Why the Order of Residualization Matters
For a gene, let \(y_i\) denote its measured expression in cell \(i\), and let \(c_i\) indicate the biological condition of that cell. A common workflow first removes fitted batch means from expression, then compares the resulting residuals across conditions. The motivation is understandable: if batch shifts expression, subtracting batch means seems to leave a cleaner biological comparison.
The difficulty is that batch and condition can be associated. When they are, part of the condition pattern is also part of the batch pattern. Regressing expression on batch indicators alone cannot distinguish those shared components; it removes them together. A subsequent residual comparison is not generally the same as estimating condition while adjusting for batch.
This tutorial focuses on a single gene and a linear contrast. In a single-cell study, cells from the same biological sample are not automatically independent biological replicates; the appropriate unit for inference must reflect the study design. The calculations below explain the geometry of residualization, not a license to treat every cell as an independent replicate.
The Model and the Projection
Suppose an additive linear model is appropriate:
Here \(\gamma\) represents batch effects, \(\beta\) is the condition effect of interest, and \(\varepsilon\) is the unexplained variation. For a common condition effect to have this interpretation, the model assumes that the condition effect is the same across batches, that the relevant condition contrast is estimable from the design, and that the mean of \(\varepsilon\) is zero conditional on the design. Standard model-based tests additionally require suitable variance assumptions; for example, the classical Gaussian linear-model test uses independent errors with constant variance and a normal distribution. If effects differ by batch, a single common \(\beta\) may not describe the data.
The joint model estimates the condition effect after accounting for batch. Blind residualization instead forms \(y^*=M_By\) and compares \(y^*\) with the original condition indicator \(c\). To see the distinction, write the ordinary least-squares condition slope from that residual comparison, including an intercept, as
The corresponding condition coefficient in the joint batch-and-condition regression is
when \(c'M_Bc>0\). The numerator is shared, but the denominators generally are not. The joint model uses the condition variation left after projecting out batch; the blind comparison uses the total centered condition variation, including its between-batch component.
A Short Derivation: What Is Removed?
Because \(M_B\) is a symmetric projection matrix, \(M_B'=M_B\), and because it removes every vector in the batch-design space, \(M_BB=0\). Applying it to the model gives
Thus the part of the condition pattern that remains in the batch-residualized expression is \(M_Bc\beta\), not \(c\beta\). If condition is perfectly determined by batch, then \(c\) lies in the column space of \(B\), so \(M_Bc=0\). The condition term disappears from the residualized expression. More generally, only the within-batch condition variation remains.
For the joint estimate, multiply the residualized model by \(c'M_B\). Since \(c'M_Bc=(M_Bc)'(M_Bc)\), solving for the coefficient gives the joint formula above. In the blind comparison, the denominator instead measures all centered condition variation. By the projection decomposition,
The second term is condition variation explained by batch. It is nonnegative, so the blind denominator is at least as large as the joint denominator. Under the additive model, this makes the blind slope attenuated toward zero when the shared numerator is nonzero. If \(M_Bc=0\), the numerator is zero too: there is no residual condition contrast to estimate.
Worked Example: Perfect Confounding Removes the Group Difference
Consider hypothetical data for one gene. All three control cells are in batch A, with expression values 8, 10, and 12. All three treated cells are in batch B, with values 18, 20, and 22. The batch means are 10 and 20.
After subtracting the batch mean, the control residuals are \(8-10=-2\), \(10-10=0\), and \(12-10=2\). The treated residuals are \(18-20=-2\), \(20-20=0\), and \(22-20=2\). Each condition has residual mean 0, so the residualized treated-minus-control contrast is \(0-0=0\).
This zero does not show that the biological condition has no effect. Condition and batch coincide exactly here: the observed 10-unit difference between group means is also the difference between batch means. The data cannot separate those explanations. Residualization has removed the shared difference, and no calculation on these observations can recover a distinct condition effect without additional assumptions or data containing condition overlap across batches.
Partial Confounding: A Contrast Can Be Shrunk
When both conditions occur in at least some of the same batches, within-batch condition comparisons can remain estimable. But blind residualization still uses a different denominator from joint adjustment. Consequently, a residualized difference may understate the condition coefficient from the joint model.
Worked Example: Partial Overlap and Attenuation
Suppose batch A has four control cells with expression 10 and two treated cells with expression 13. Batch B has two control cells with expression 20 and four treated cells with expression 23. These hypothetical values follow an additive pattern with a 10-unit batch difference and a 3-unit condition difference. Both conditions occur in both batches, so the condition coefficient is estimable under the common-effect model.
The batch A mean is \((4(10)+2(13))/6=11\), and the batch B mean is \((2(20)+4(23))/6=22\). Subtracting these means gives residuals of \(-1\) for batch A controls, \(2\) for batch A treated cells, \(-2\) for batch B controls, and \(1\) for batch B treated cells.
The pooled treated residual mean is \((2(2)+4(1))/6=1.333\). The pooled control residual mean is \((4(-1)+2(-2))/6=-1.333\). Therefore, the blind residualized contrast is \(1.333-(-1.333)=2.667\).
The within-batch contrasts are both 3: \(2-(-1)=3\) in batch A and \(1-(-2)=3\) in batch B. The joint additive-model coefficient is therefore 3. The blind residual contrast is smaller because it divides the shared within-batch signal by the total condition variation, including the condition imbalance between batches. These numbers illustrate coefficient differences only; they do not provide a valid uncertainty estimate or p-value.
Residual Contrasts Can Answer a Different Question
Even when the residual comparison is calculated correctly, its meaning depends on the design and on how effects vary across batches. Batch residuals express deviations from each batch’s mean. A pooled residual contrast therefore summarizes a particular mixture of within-batch deviations; it is not automatically the raw condition difference, an equal-weight average of batch-specific effects, or a common condition coefficient.
Worked Example: Imbalance and Batch-Specific Effects
Consider two hypothetical batches. Batch A has nine controls at expression 100 and one treated cell at 104, giving a within-batch condition contrast of \(+4\). Batch B has one control at 0 and nine treated cells at \(-1\), giving a within-batch contrast of \(-1\). There is strong condition imbalance, and the condition effect differs by batch.
The raw control mean is \((9(100)+0)/10=90\); the raw treated mean is \((104+9(-1))/10=9.5\). The raw treated-minus-control contrast is \(9.5-90=-80.5\), largely reflecting the different batch compositions.
Now subtract each batch mean. Batch A’s mean is \((9(100)+104)/10=100.4\), so its residuals are \(-0.4\) for controls and \(3.6\) for the treated cell. Batch B’s mean is \((0+9(-1))/10=-0.9\), so its residuals are \(0.9\) for the control and \(-0.1\) for treated cells. The pooled residual control mean is \((9(-0.4)+0.9)/10=-0.27\); the pooled treated residual mean is \((3.6+9(-0.1))/10=0.27\). The residual contrast is therefore \(0.27-(-0.27)=0.54\).
The positive residual contrast is not a contradiction of the negative raw contrast: the two summarize different mixtures. Nor is \(0.54\) the equal-batch average of the within-batch effects, which is \((4+(-1))/2=1.5\). With batch-specific effects, there is no single common effect unless a target and weighting scheme are specified. Calling the pooled residual contrast “the condition effect” without stating that target is misleading.
Choosing an Analysis and Checking Its Assumptions
A direct joint model makes the intended adjustment explicit: fit expression as a function of batch and condition, and estimate the condition coefficient only if the design supports it. An equivalent residualization calculation must transform both expression and condition and retain the correct uncertainty calculation. If the scientific target is a within-batch comparison, estimate those comparisons from batches containing both conditions and combine them using a stated rule. If condition effects may differ across batches, report or model those differences rather than silently imposing a common effect.
Before interpreting any residual contrast, check the batch-by-condition counts and whether the condition contrast is estimable. A batch containing only one condition contributes no within-batch comparison. If all batches contain only one condition, batch and condition are perfectly confounded. If overlap exists only in a few cells or batches, the estimate may be unstable and depend heavily on those observations. These are design limitations, not problems that residualization itself can repair.
Also distinguish a transformed data vector from a completed statistical analysis. Residuals from a fitted batch model are dependent because they share estimated batch means, and the projection uses degrees of freedom. Treating those residuals as if they were original independent observations, then applying a standard test without accounting for the fit and design, can give an invalid standard error or p-value. Inspect batch-specific residual distributions and contrasts, assess whether a common-effect model is plausible, and use an inferential method that matches the experimental unit and variance structure.
| Approach | What it estimates or preserves | Main caution |
|---|---|---|
| Residualize expression on batch only, then compare condition | A pooled contrast among batch-centered expression values | Not generally the joint-model condition coefficient; may erase shared signal or obscure the target |
| Fit batch and condition jointly | The condition coefficient conditional on batch, if estimable | Requires an appropriate model and sufficient condition overlap |
| Residualize both expression and condition | The same condition coefficient as the joint linear model | Inference must use the residualized design and correct degrees of freedom |
| Estimate within-batch contrasts | Condition comparisons supported by within-batch overlap | Requires a stated method for combining contrasts; batches without overlap add no direct comparison |
Common Mistakes and Practical Consequences
- Assuming batch removal preserves condition. If condition is associated with batch, the projection removes the component of condition aligned with batch. Check overlap before interpreting residuals.
- Residualizing only the response and calling the slope adjusted. The joint-model coefficient uses \(M_Bc\) as well as \(M_By\). Leaving \(c\) unchanged generally changes the denominator and the coefficient.
- Interpreting a zero residual contrast as evidence of no biology. Under perfect confounding, zero can result because the shared batch-condition pattern was removed. The distinct condition effect is not identified by that design.
- Assuming every residual contrast has one universal biological meaning. With imbalance or batch-specific effects, the pooled value depends on the design and weighting. State the contrast and target explicitly.
- Testing residuals as though they were untouched observations. Batch fitting affects both residual dependence and uncertainty. Use an inferential procedure appropriate to the fitted model and the biological replicate structure.
Check Your Understanding
Answer each question, then compare with the brief solutions.
- If every control is in batch A and every treated cell is in batch B, what can batch-only residualization establish about the distinct condition effect?
- In a partial-overlap design, why can the blind residual slope be smaller in magnitude than the joint-model coefficient?
- What must be residualized, in addition to expression, to obtain the joint linear-model coefficient by the Frisch–Waugh–Lovell method?
- If batch-specific condition effects differ, why might a pooled residual contrast be misleading when described simply as “the condition effect”?
- Does subtracting fitted batch means recover a comparison for a batch that contains only one condition? Explain.
Solutions. (1) It cannot distinguish condition from batch: their effects are not separately identified, and the shared group difference is removed. (2) The blind slope divides the within-batch signal by total centered condition variation, including between-batch variation; the joint coefficient divides by the remaining within-batch condition variation. (3) Residualize the condition vector as well as expression, then regress the two residualized vectors. (4) The pooled contrast uses a design-dependent mixture of within-batch deviations; it need not equal a common effect or an equal-batch average. (5) No. Without both conditions in that batch there is no within-batch condition comparison, so residualization cannot create the missing comparison.