sva R Package: Install, Worked Example and Output

The sva package (Surrogate Variable Analysis) is a Bioconductor package for R that removes batch effects and other unwanted variation from high-throughput data such as gene expression arrays and RNA-seq. This page is for analysts who need to install it, choose the right function for their data, run a worked example, and read the output. As of the Bioconductor 3.23 release (checked 2026-11-08), the package is version 3.60.0, licensed Artistic-2.0, and maintained by Jeffrey T. Leek, John D. Storey and W. Evan Johnson.

What sva is used for in research

The package vignette describes two uses: (1) estimating surrogate variables that capture unknown sources of variation, such as unrecorded processing differences, and (2) removing a known batch with ComBat. The Bioconductor description adds that it can also correct using control probes. Surrogate variables are then included as covariates in a downstream model, so that differential-expression results are not driven by the hidden structure.

Which function for which situation

Situation	Function	What it does (per the vignette)
Unknown sources of variation, continuous data (e.g. microarray)	num.sv(), sva()	Estimates the number of latent factors, then the surrogate variables themselves
Unknown sources of variation, count data (RNA-seq)	svaseq()	Surrogate variable estimation for counts; a supervised form accepts a controls argument
Known batch, continuous data	ComBat()	Adjusts for a batch variable you supply; returns a matrix of the same dimensions
Known batch, count data	ComBat_seq()	Negative binomial adjustment that keeps integer counts
Prediction on new samples	fsva()	Adjusts the training data with surrogate variables and the test data with a "frozen" algorithm

For a deeper treatment of the known-batch route, see the CASRAI guide to ComBat batch correction.

Install or access (official source only)

sva is distributed through Bioconductor, not CRAN. The package page at bioconductor.org/packages/release/bioc/html/sva.html gives the installation commands:

if (!require("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install("sva")

The package declares R (>= 3.2) as its minimum, while the Bioconductor installation instructions for release 3.23 call for R 4.6. Match your R version to the Bioconductor release you install from. The vignette is available inside R with browseVignettes("sva").

Step-by-step task: remove hidden variation from an expression matrix

The steps below follow the bladderbatch example in the official vignette. The input is a matrix with features (genes) in rows and samples in columns.

1. Load the packages and data.

   library(sva)
   library(bladderbatch)
   data(bladderdata)
   pheno = pData(bladderEset)
   edata = exprs(bladderEset)

2. Build the full and null models. The full model contains the adjustment variables plus the variable of interest; the null model contains only the adjustment variables. Both are made with model.matrix.

   mod  = model.matrix(~as.factor(cancer), data=pheno)
   mod0 = model.matrix(~1, data=pheno)

3. Estimate the number of surrogate variables.

   n.sv = num.sv(edata, mod, method="leek")

   The vignette notes that sva() will also estimate the number itself if you omit n.sv. The estimate is based on a permutation procedure (Buja and Eyuboglu, 1992) as described in Leek and Storey (2007).

4. Estimate the surrogate variables.

   svobj = sva(edata, mod, mod0, n.sv=n.sv)

5. Test with and without adjustment. f.pvalue() runs a parametric F-test per row and requires nested models. Convert to Benjamini-Hochberg adjusted values with p.adjust.

   pValues = f.pvalue(edata, mod, mod0)
   qValues = p.adjust(pValues, method="BH")
   
   modSv  = cbind(mod, svobj$sv)
   mod0Sv = cbind(mod0, svobj$sv)
   pValuesSv = f.pvalue(edata, modSv, mod0Sv)
   qValuesSv = p.adjust(pValuesSv, method="BH")

6. Or pass the surrogate variables to limma. Include them in the design matrix; they adjust the fit but are not part of the contrasts.

   fit = lmFit(edata, modSv)

Known batch with ComBat

If the batch is recorded, the vignette uses ComBat. The batch variable is passed separately and is not placed in the model matrix; the model matrix holds only the other covariates (here, an intercept-only model).

batch = pheno$batch
modcombat = model.matrix(~1, data=pheno)
combat_edata = ComBat(dat=edata, batch=batch, mod=modcombat,
                      par.prior=TRUE, prior.plots=FALSE)

A categorical covariate with p levels needs p-1 indicator columns. par.prior=TRUE uses the parametric empirical Bayes adjustment; par.prior=FALSE uses the nonparametric version, which is slower; prior.plots=TRUE draws the prior plots.

RNA-seq counts

For counts, the vignette uses svaseq() (its zebrafish example has six samples, so it sets n.sv=1), and ComBat_seq() when the batch is known:

svseq = svaseq(dat0, mod1, mod0, n.sv=1)$sv

adjusted = ComBat_seq(count_matrix, batch=batch, group=NULL)

The group argument of ComBat_seq() preserves a biological variable of interest. For supervised svaseq(), controls holds values between 0 and 1 giving the probability that a gene is affected by batch but not by the group; negative-control genes are marked TRUE. For differential expression on counts, see the CASRAI guide to differential gene expression analysis.

Output and how to read it

sva() returns a list of four elements:

• sv: a matrix whose columns are the surrogate variables, one row per sample.
• pprob.gam: the posterior probability that each gene is associated with one or more latent variables.
• pprob.b: the posterior probability that each gene is associated with the variables of interest.
• n.sv: the number of surrogate variables used.

The practical check is the before-and-after comparison in step 5. In the vignette’s unadjusted analysis, nearly 70% of genes appear strongly differentially expressed at FDR below 5%, which the authors call artificially high. Refitting with the surrogate variables added to both models is the adjustment. ComBat() instead returns a matrix with the same dimensions as the input, with the known batch effect removed.

Troubleshooting

• Very large feature sets. With more than 100,000 features, num.sv() and sva() are slow. The vfilter argument (an integer between 100 and the number of features) restricts the computation to the most variable features, for example num.sv(edata, mod, vfilter=2000, method="leek") and sva(edata, mod, mod0, n.sv=n.sv, vfilter=2000). The vignette does not recommend fewer than 1,000 features.
• Surrogate variables that track real biology. If the latent variables reflect true unknown subgroups, sva may not be appropriate, because surrogate variables can correlate with the subgroup. Direct adjustment removes only the batch variables you name.
• ComBat model matrix errors. Keep the batch variable out of mod, and use p-1 indicators for a covariate with p levels.
• Nested models. f.pvalue() needs the null model to be nested in the full model; build both from the same adjustment variables.
• Few samples. The vignette’s six-sample RNA-seq example uses a single surrogate variable (n.sv=1).

Free alternatives and companions

sva is itself open source (Artistic-2.0). Within the same ecosystem, the vignette documents these related options:

Approach	When the vignette uses it
Direct adjustment with limma (lmFit) using a known batch as a covariate	Batch is recorded and you model it directly; removes only the named variables
ComBat() / ComBat_seq() (in sva)	Known batch; data cleaned before downstream analysis
sva() / svaseq() with surrogate variables as covariates	Unknown or incompletely recorded sources of variation

Python users with a general statistical-modelling need can look at the CASRAI guide to the statsmodels API; sva itself is an R/Bioconductor package. Whatever route you choose, record the package version, function arguments and random seed so the analysis can be repeated, in line with reproducibility practice.

FAQ

What does sva stand for?

Surrogate Variable Analysis. The Bioconductor package of that name estimates surrogate variables for unwanted variation in high-throughput data.

Is sva on CRAN?

No. The documented installation route is Bioconductor, via BiocManager::install("sva").

Should I use sva or ComBat?

Use ComBat when the batch variable is known and recorded; use sva or svaseq when the sources of variation are unknown. Both are in the same package.

How many surrogate variables should I use?

Estimate with num.sv(), or let sva() estimate it when n.sv is omitted. With very small sample sizes the vignette example uses one.

How do I cite sva?

The Bioconductor page gives: Leek JT, Johnson WE, Parker HS, Fertig EJ, Jaffe AE, Zhang Y, Storey JD, Collado Torres L. sva: Surrogate Variable Analysis. doi:10.18129/B9.bioc.sva. The vignette also points to Leek and Storey (2007, 2008), Johnson et al. (2007) for ComBat, and Leek (2014) for svaseq.

Sources (official documentation)

• Bioconductor, sva package page, release 3.23, version 3.60.0, read 2026-11-08.
• Bioconductor, sva vignette (browseVignettes("sva")), read 2026-11-08.