Empirical-Bayes Batch Adjustment with ComBat

Why ComBat Moderates Batch Estimates

A batch can shift many measured features in a similar way, while the number of samples in any one batch may be too small to estimate each feature’s batch effect precisely. Estimating every shift and scale separately can therefore produce noisy corrections. ComBat addresses this by combining feature-specific estimates with information pooled across features: estimates that are imprecise are pulled toward a batch-level distribution, while estimates supported by more data are changed less.

ComBat is a location-and-scale adjustment for continuous measurements. It is useful when samples were processed in identifiable batches and a biological covariate, such as treatment group, should remain represented after adjustment. It does not make a confounded study identifiable: if batch and biology cannot be separated in the observed design, shrinkage cannot recover the missing information.

Model and Assumptions

Definition: Let $Y_{ij}$ be the measurement for feature $i$ in sample $j$, and let $b(j)$ denote that sample’s batch. A standard ComBat model is $$ Y_{ij}=\alpha_i+X_j\beta_i+\gamma_{i,b(j)}+\delta_{i,b(j)}\epsilon_{ij}, \qquad \epsilon_{ij}\sim N(0,\sigma_i^2). $$ Here, $\alpha_i$ is a feature-specific intercept, $X_j\beta_i$ represents specified biological or other protected covariates, $\gamma_{i,b}$ is the batch-specific location effect, and $\delta_{i,b}$ is a positive batch-specific scale multiplier. The errors are assumed independent across samples within a feature, conditional on the model, and approximately normal after any appropriate transformation.

ComBat first estimates the covariate contribution and feature-specific scale, then standardizes the measurements. Write the standardized residual as $Z_{ij}$. In the idealized standardized model, a sample from batch $b$ has $Z_{ij}$ centered at $\gamma_{ib}$ with variance $\delta_{ib}^{2}$. The raw batch-specific location estimate is $\hat{\gamma}_{ib}$, and the raw scale estimate is $\hat{\delta}_{ib}^{2}$. The notation $\delta^2$ denotes a variance multiplier; $\delta$ itself is a standard-deviation multiplier.

The empirical-Bayes step treats the feature-specific effects within a batch as draws from shared distributions. In the parametric version, a common model is

Model: For features $i=1,\ldots,m$ in batch $b$, assume $\gamma_{ib}\sim N(\bar{\gamma}_b,\tau_b^2)$ and $\delta_{ib}^{2}\sim \operatorname{InverseGamma}(a_b,b_b)$. The hyperparameters are estimated from the collection of features in that batch. This requires features pooled together to have sufficiently comparable standardized batch effects for a shared distribution to be meaningful.

These are modeling assumptions, not consequences of having many features. ComBat also requires enough samples in each batch to estimate the within-batch variation. The biological covariates to preserve must be included in the design, and the design must contain information that distinguishes those covariates from batch. If, for example, every treated sample is in one batch and every control sample is in another, treatment and batch are perfectly confounded; their effects cannot be separated by this model.

How Location Shrinkage Works

Under the standardized model, the sampling variance of the raw location estimate is approximately $\delta_{ib}^{2}/n_b$, where $n_b$ is the number of samples in batch $b$. Combining this sampling distribution with the normal prior gives the empirical-Bayes posterior mean

\tilde{\gamma}_{ib} = \frac{n_b\tau_b^2\hat{\gamma}_{ib}+\hat{\delta}_{ib}^{2}\bar{\gamma}_b} {n_b\tau_b^2+\hat{\delta}_{ib}^{2}}.

The raw estimate receives more weight when it is based on more samples or when the prior distribution across features is broad. The batch mean $\bar{\gamma}_b$ receives more weight when the raw estimate is noisy or the across-feature prior is concentrated. The expression is a weighted average: its weights are proportional to inverse variances.

Estimate feature-specific effects.
Fit the covariate model and calculate standardized batch-specific location and scale estimates.

Estimate batch-level distributions.
Use the collection of feature estimates in each batch to estimate the location and scale prior hyperparameters.

Moderate the estimates.
Combine each raw feature estimate with its batch-level prior, producing moderated location and scale estimates.

Adjust and restore the modeled signal.
Remove the moderated batch effects from standardized measurements, rescale, and add back the fitted covariate contribution.

Worked Example: Moderating a Batch Location Estimate

Suppose a hypothetical batch contains $n_b=4$ samples. For one standardized feature, the raw location estimate is $\hat{\gamma}_{ib}=0.80$, and its raw scale estimate is $\hat{\delta}_{ib}^{2}=0.36$. Across features in this batch, suppose the estimated location prior has mean $\bar{\gamma}_b=0.20$ and variance $\tau_b^2=0.25$. The location estimate is moderated as follows:

\tilde{\gamma}_{ib} = \frac{(4)(0.25)(0.80)+(0.36)(0.20)} {(4)(0.25)+0.36} = \frac{0.800+0.072}{1.000+0.360} = \frac{0.872}{1.360} \approx 0.641.

The moderated estimate is about $0.641$ standardized units, between the raw estimate $0.80$ and the prior mean $0.20$. Because this batch has only four samples and the feature’s estimated variance is not negligible relative to the prior variance, the estimate is pulled toward the batch-level mean. It is not replaced by that mean.

Scale Moderation and the Adjustment

ComBat also moderates batch-specific variances. A common parametric empirical-Bayes update uses an inverse-gamma prior and, conditional on the moderated location estimate, updates the scale variance by

\tilde{\delta}_{ib}^{2} = \frac{b_b+\frac{1}{2}\sum_{j:b(j)=b}(Z_{ij}-\tilde{\gamma}_{ib})^2} {a_b+\frac{n_b}{2}-1}.

This is a plug-in empirical-Bayes update: the location estimate is treated as fixed in the variance update, and the prior hyperparameters $a_b$ and $b_b$ are estimated across features. Thus, the formula describes the standard moderated estimation procedure, not a claim that all uncertainty in the location estimate has been integrated out. Positive variance estimates and a suitable prior are required.

Worked Example: Moderating a Batch Scale Estimate

Consider a hypothetical feature with four standardized residuals in one batch: $-1.2,-0.4,0.4,1.2$. Their mean is $0$, and their sum of squared deviations from zero is $1.44+0.16+0.16+1.44=3.20$. Suppose the moderated location is $\tilde{\gamma}_{ib}=0$, and the batch’s estimated inverse-gamma hyperparameters are $a_b=4$ and $b_b=2$. Then

\tilde{\delta}_{ib}^{2} = \frac{2+\frac{1}{2}(3.20)} {4+\frac{4}{2}-1} = \frac{2+1.60}{5} = 0.720.

The moderated scale variance is $0.720$, in squared standardized units. The raw sample variance of these four values, using denominator $n_b-1$, is $3.20/3\approx1.067$. The prior moderates the estimate toward its implied mean variance $b_b/(a_b-1)=2/3\approx0.667$. The adjusted standard-deviation multiplier is $\sqrt{0.720}\approx0.849$, not $0.720$.

After estimating the moderated effects, ComBat adjusts the standardized measurement and restores the fitted non-batch component. In the notation above, the adjustment is

Y_{ij}^{\mathrm{adj}} = \left(\frac{Z_{ij}-\tilde{\gamma}_{i,b(j)}}{\tilde{\delta}_{i,b(j)}}\right)\hat{\sigma}_i +\hat{\alpha}_i+X_j\hat{\beta}_i.

Subtracting the moderated location removes the estimated batch shift; dividing by the moderated standard-deviation multiplier removes the estimated batch-specific scale difference. Restoring $\hat{\alpha}_i+X_j\hat{\beta}_i$ retains the fitted covariate signal included in the model. ComBat therefore adjusts the specified batch effects, not every difference between samples.

Worked Example: Applying the Adjustment

Suppose a hypothetical sample has standardized value $Z_{ij}=1.50$. Its batch’s moderated location is $0.40$, its moderated variance multiplier is $0.64$, and the feature’s estimated standard deviation is $\hat{\sigma}_i=2.0$ measurement units. The fitted covariate contribution, including the intercept, is $3.0$ measurement units. The adjusted value is

Y_{ij}^{\mathrm{adj}} = \left(\frac{1.50-0.40}{\sqrt{0.64}}\right)(2.0)+3.0 = \left(\frac{1.10}{0.80}\right)(2.0)+3.0 = 2.75+3.0 = 5.75.

The adjusted measurement is $5.75$ units. The square root is necessary because $0.64$ is the moderated variance multiplier, while the adjustment divides by its standard-deviation multiplier, $0.80$. This calculation removes the estimated batch contribution while retaining the stated covariate contribution.

Checks, Alternatives, and Limitations

Before using the adjusted matrix, check whether its assumptions and design are plausible. These checks do not prove the model is correct, but they can reveal serious mismatches:

Inspect the design. Tabulate batch against the covariates to be preserved. Empty or nearly empty batch-by-covariate combinations indicate weak separation; perfect confounding makes the desired decomposition unidentifiable.
Check measurement suitability. The model is intended for continuous measurements whose residual distributions are reasonably represented by the normal model after suitable preprocessing. It is not automatically appropriate for raw counts, binary outcomes, or highly censored values.
Check pooled-feature comparability. The shared prior is estimated across features. A mixture of very different feature classes, scales, or batch-effect distributions can make one prior a poor summary. Consider analyzing defensible feature groups separately rather than pooling indiscriminately.
Check sample sizes and estimates. Very small batches yield uncertain feature-specific variances, making results sensitive to prior estimation. Zero or near-zero within-batch variation can also make scale estimation unstable.
Compare raw and adjusted summaries in context. Large changes concentrated in one batch, unexpected compression of within-batch variation, or altered covariate contrasts warrant investigation. An adjusted matrix is not evidence that every remaining difference is biological.

ComBat differs from the covariate-preserving linear batch adjustment in that it explicitly models and moderates both location and scale effects across features. An unmoderated feature-by-feature adjustment estimates batch effects separately, which may be reasonable with ample data but can be noisy in small batches. ComBat’s gain is stabilization through pooling; its cost is reliance on a shared across-feature distribution. Parametric ComBat uses specified normal and inverse-gamma priors, whereas nonparametric empirical-Bayes variants estimate the across-feature distributions more flexibly. Greater flexibility does not solve confounding, poor measurement distributions, or an inappropriate choice of features to pool.

ComBat is not a substitute for study design, and it does not guarantee that technical effects have been removed without affecting biology. Its correction is conditional on the covariates and batches supplied to the model. If an important biological variable is omitted, its signal may be treated as unexplained variation; if it is inseparable from batch, no adjustment can reliably preserve it while removing batch.

Common Mistake

Pitfall: “Empirical Bayes borrows information across features, so ComBat can separate a biological effect from batch even when all cases are in one batch and all controls are in another.” This is incorrect. Shrinkage stabilizes estimates under the model; it does not create overlap or identify effects that the design confounds. Include the biological covariate and inspect its overlap with batch before adjustment. If the effects are perfectly confounded, their separate contributions are not estimable from these data.

Key takeaway: ComBat combines feature-specific batch estimates with batch-level empirical-Bayes distributions, moderating both location and scale. Its usefulness depends on a suitable measurement model, comparable features for pooling, adequate batch information, and a design that distinguishes batch from protected biology.

Check Your Understanding

Answer each question, then compare with the solutions below.

In the location update, what happens to shrinkage toward $\bar{\gamma}_b$ when $n_b$ increases, with the other quantities held fixed?
A moderated variance multiplier is $0.49$. What standard-deviation multiplier is used in the adjustment?
Why does including a treatment covariate not solve perfect confounding between treatment and batch?
For $n_b=2$, $\hat{\gamma}_{ib}=1.0$, $\hat{\delta}_{ib}^{2}=0.25$, $\bar{\gamma}_b=0$, and $\tau_b^2=0.25$, calculate $\tilde{\gamma}_{ib}$.

Solutions. (1) The raw estimate receives more weight because $n_b\tau_b^2$ increases relative to $\hat{\delta}_{ib}^{2}$, so shrinkage toward the prior mean decreases. (2) The standard-deviation multiplier is $\sqrt{0.49}=0.70$. (3) If treatment and batch have exactly the same sample assignments, their effects cannot be distinguished from the observed data; adding a covariate cannot supply missing overlap. (4) $\tilde{\gamma}_{ib}=[(2)(0.25)(1.0)+(0.25)(0)]/[(2)(0.25)+0.25]=0.50/0.75\approx0.667$ standardized units.