ComBat Empirical Bayes Batch Adjustment

Why ComBat Uses Empirical Bayes Shrinkage

A batch can shift expression values for many genes, and it can also change their spread. Removing a separate location and scale estimate for every gene and batch is possible, but estimates from small batches can be unstable. ComBat addresses this by estimating gene-specific batch effects and shrinking them toward batch-wide distributions learned across genes. The aim is to reduce technical differences while retaining modeled biological variation.

ComBat is intended for expression measurements on a suitable transformed, approximately continuous scale, such as log-scale expression, rather than raw sequencing counts. The earlier material on library-size normalization and variance-stabilizing transformations establishes why counts commonly need appropriate normalization and transformation before this kind of location-scale adjustment. ComBat does not make confounded biology and batch separately identifiable: if batch and condition are perfectly confounded, the data cannot determine which explains their shared differences.

The Location-Scale Model and Assumptions

A common ComBat model is

Here, \(\alpha_i+X_j\beta_i\) represents the gene-specific baseline and modeled biological effects, while \(\gamma_{i,b(j)}\) and \(\delta_{i,b(j)}\) represent batch location and scale effects. The model assumes \(\delta_{i,b}>0\). For the usual exact normal likelihood and conjugate posterior calculations, assume the errors are independent and normally distributed, \(\epsilon_{ij}\sim N(0,1)\), conditional on the modeled effects. Mean zero and variance one alone do not imply normality, so they are not sufficient for those exact calculations.

ComBat first estimates and removes the modeled baseline and biological covariates, then standardizes each gene using a pooled residual scale \(s_i\). Write the standardized residual as

In this standardized scale, a batch’s gene-specific observations are modeled as \(Z_{ij}=\gamma_{i,b}+\delta_{i,b}\epsilon_{ij}\). Standardization gives genes comparable scales for estimating the across-gene priors. In practice, the fitted baseline and standardization are estimated from the data; their adequacy and the resulting residual batch structure require careful assessment.

For a batch \(b\), ComBat uses distributions across genes such as \(\gamma_{i,b}\sim N(\bar{\gamma}_b,\tau_b^2)\) and \(\delta_{i,b}^2\sim \operatorname{InvGamma}(\lambda_b,\theta_b)\). The inverse-gamma parameterization here has density proportional to \(v^{-\lambda_b-1}\exp(-\theta_b/v)\) for \(v>0\). The prior hyperparameters are estimated across genes within a batch. The assumption is that genes in that batch provide a useful shared distribution for stabilizing one another’s estimates; it does not require all genes to have the same batch effect.

How Shrinkage Is Calculated

For one gene and batch, let \(n_b\) be the number of samples in the batch and \(\bar{Z}_{i,b}\) the mean standardized residual. If the batch scale variance \(\delta_{i,b}^2\) is treated as known and the standardized errors are independent normal, the batch mean has the exact likelihood

Combining this likelihood with the normal prior on \(\gamma_{i,b}\) gives a normal posterior. Its mean is a weighted average of the gene’s observed batch mean and the batch-wide prior mean:

The observed mean receives more weight when the batch has more samples or the prior variance \(\tau_b^2\) is large. The prior mean receives more weight when the batch is small or the gene-level residual variance is large. This is partial pooling: estimates are drawn toward a shared center, not forced to equal it.

For scale, conditional on a fixed location \(\gamma_{i,b}\), the standardized residual sum of squares is \(R_{i,b}=\sum_{j:b(j)=b}(Z_{ij}-\gamma_{i,b})^2\). With the inverse-gamma prior above and normal errors, the conditional posterior is inverse-gamma with shape \(\lambda_b+n_b/2\) and scale \(\theta_b+R_{i,b}/2\). When the shape exceeds one, its mean is

In the empirical Bayes procedure, location and scale estimates are updated using one another, with batch-wide hyperparameters estimated across genes. A calculation that plugs a shrunk location estimate into the residual sum of squares is a conditional empirical Bayes update, not a claim that the fully joint posterior has this simple independent form. The normal likelihood is essential for the stated exact conjugate updates.

Adjusting Expression Values

After estimating batch location and scale effects, ComBat removes them from the standardized residuals and restores the gene’s original scale and modeled biological component. In the notation above, the adjusted value is

The location estimate is subtracted and the residual is divided by the estimated batch scale. The biological component \(\hat{\alpha}_i+X_j\hat{\beta}_i\) is retained. The adjustment therefore depends on the fitted covariate model and on the assumption that batch effects are separable from the biology represented by \(X\).

What ComBat Adds, and Where It Can Fail

A standard linear adjustment with batch indicators estimates batch coefficients without borrowing strength across genes. That can be appropriate when its model is adequate and the estimates are sufficiently precise. ComBat adds empirical Bayes shrinkage, which can stabilize gene-specific estimates when many genes share information and batch sample sizes are limited. It is not a substitute for a paired or blocked design: those designs control batch through the study structure, whereas ComBat estimates and adjusts batch effects after data collection.

ComBat’s usefulness depends on design and scale. If every sample in one condition was processed in one batch and every sample in another condition in a different batch, condition and batch are confounded. No adjustment can reliably preserve the condition effect while removing the batch effect because the observed data do not distinguish them. The method also does not justify applying a normal location-scale model directly to raw counts, nor does it guarantee that all technical differences have been removed.

A key diagnostic is whether the model’s assumptions and estimated effects are plausible for the data: check the expression scale, sample sizes and condition-by-batch overlap, and whether the fitted covariates reflect the biology to retain. A later examination of ComBat fit and residual batch structure addresses diagnostics in more detail. For this tutorial, the essential caution is that shrinkage controls estimation variability under its model; it cannot repair an unidentifiable design or make a misspecified model correct.