Two Outputs, Two Questions
A negative-binomial (NB) depth model can produce more than one useful normalized matrix, but those matrices do not represent the same quantity. Corrected counts rescale observed counts to a chosen reference depth. Pearson residuals measure how far an observed count is from its fitted mean, in units of its model-based standard deviation. The first is intended to retain a count-like expression signal while adjusting for depth; the second expresses departure from the fitted model.
This distinction matters when summarizing data. A mean of corrected counts can be described as average expression on a reference-depth scale. A mean of residuals instead describes average signed deviation from the model. Neither output is universally preferable: the appropriate choice follows from the scientific question and the downstream calculation.
Model and Definitions
Let \(Y_{ij}\) be the observed count for gene \(j\) in cell \(i\), and let \(s_i>0\) be the cell’s depth factor. In the depth-only NB model considered here, the count has conditional mean \(\mu_{ij}\) and dispersion \(\theta_j>0\), with \(Y_{ij}\mid s_i \sim \operatorname{NB}(\mu_{ij},\theta_j)\). We use the parameterization \(\operatorname{Var}(Y_{ij}\mid s_i)=\mu_{ij}+\mu_{ij}^2/\theta_j\). The mean is modeled as \(\mu_{ij}=\exp(\beta_{0j})s_i^{\beta_{1j}}\). Here \(\beta_{0j}\) is a gene-specific intercept and \(\beta_{1j}\) is its fitted depth coefficient. In practice these parameters are estimated; formulas below use their fitted values, with hats omitted for readability.
This model requires counts to be nonnegative integers, positive depth factors, and a mean and dispersion specification suitable for the data. It also assumes that the depth relationship is adequately represented by the stated mean model for the purpose at hand. In particular, a depth adjustment alone does not account for every biological or technical difference among cells.
The corrected count is an observed count multiplied by a model-derived adjustment factor. It is not the fitted mean itself. The residual, in contrast, subtracts the fitted mean from the observation and divides by the model-based standard deviation. Thus the corrected count is nonnegative, while the residual can be negative, zero, or positive.
Why the Formulas Differ
The mean model gives the expected count at the observed depth and the expected count at the chosen reference depth. Their ratio isolates the depth adjustment:
Multiplying the observed count by that ratio expresses it on the reference-depth scale. The intercept cancels; the correction depends on the depth factors and the gene’s fitted depth coefficient. This operation retains the observed count’s sampling variation, rather than replacing the count with its fitted expectation.
For the residual, the NB variance assumption supplies the denominator. Under the stated parameterization, \(\operatorname{Var}(Y_{ij}\mid s_i)=\mu_{ij}+\mu_{ij}^2/\theta_j\); therefore the square root of this variance is the model-based standard deviation. Subtracting the mean and dividing by that standard deviation puts departures on a scale that accounts for the mean-dependent NB variability. This standardization does not guarantee that fitted residuals are exactly independent or have exactly unit variance: parameters are estimated, counts are discrete, and the model may be imperfect.
Worked Example: Correcting Counts Can Change Their Order
Suppose a gene has fitted depth coefficient \(\beta_{1j}=1\) and reference depth \(s_{\mathrm{ref}}=1\). Cell A has observed count \(Y_{Aj}=1\) and depth \(s_A=1\); cell B has \(Y_{Bj}=2\) and depth \(s_B=4\). The counts are hypothetical.
For cell A, \(C_{Aj}=1(1/1)^1=1\). For cell B, \(C_{Bj}=2(1/4)^1=0.5\). Although the observed count in B is larger, its corrected count is smaller. Each adjustment factor is positive, so the correction preserves nonnegativity, but it does not generally preserve the ordering of counts across cells. Ordering is preserved when the compared cells have the same adjustment factor, or when their particular counts and factors happen to retain that ordering.
This is why a shared fitted depth coefficient does not establish shared correction factors: the factors also depend on each cell’s depth. Corrected counts should be compared on their reference-depth scale, not assumed to reproduce the raw ranking of cells.
What Happens to Zeros?
Since the corrected-count factor \((s_{\mathrm{ref}}/s_i)^{\beta_{1j}}\) is positive, a zero observed count remains zero: \(Y_{ij}=0\) implies \(C_{ij}=0\). A positive observed count remains positive, although its corrected value may be fractional. Consequently, corrected counts retain the zero pattern of the observed matrix, but they are not generally integer-valued counts and should not be treated as if they followed the original NB sampling model.
For a Pearson residual, a zero observed count has a different representation. If \(\mu_{ij}>0\), then \(R_{ij}=-\mu_{ij}/\sqrt{\mu_{ij}+\mu_{ij}^2/\theta_j}<0\). So an observed zero becomes a negative residual, not a zero in the residual matrix. Its magnitude depends on the fitted mean and dispersion: a zero where the model predicts little has a smaller-magnitude residual than a zero where the model predicts more.
This change in representation is useful when the goal is to compare observations with expectations, but it means a residual matrix is not a sparse count matrix in the same sense as the original data. The presence of a nonzero residual does not mean that a transcript was observed.
Worked Example: A Zero Count and a Positive Residual
For a hypothetical gene and cell, suppose the fitted mean is \(\mu_{ij}=2\) and the dispersion is \(\theta_j=4\). The model variance is \(2+2^2/4=3\), so the model-based standard deviation is \(\sqrt{3}\approx1.732\).
If the observed count is zero, then \(R_{ij}=(0-2)/\sqrt{3}\approx-1.155\). The zero is represented by a negative residual. If instead the observed count is five, then \(R_{ij}=(5-2)/\sqrt{3}\approx1.732\), a positive residual indicating an observation above the fitted mean.
These residuals describe model-relative departures, not transcript counts. In particular, \(R_{ij}=1.732\) does not mean that the cell contains 1.732 transcripts; the observed count in this example is five.
Choosing Summaries for Each Output
Corrected counts and residuals support different kinds of summaries. The table gives a practical distinction; the model and the intended interpretation should still be checked for the analysis being performed.
| Output | What a value represents | Suitable summaries and uses | Interpretive caution |
|---|---|---|---|
| Corrected counts \(C_{ij}\) | Observed count rescaled to the chosen reference-depth scale | Gene-wise means or sums across cells; relative expression comparisons on that scale; summaries that require nonnegative, count-like values | Values can be fractional and are not raw counts from the NB sampling model; sums do not restore the original library totals |
| Pearson residuals \(R_{ij}\) | Signed difference from the fitted mean, divided by the fitted standard deviation | Summaries of model departures; feature matrices for methods that use centered, signed inputs, such as many dimension-reduction procedures | Negative values are meaningful; sums or means describe deviations, not abundance or transcript totals |
The earlier material on log-normalized values introduces another expression-scale output. A log transformation is useful when a downstream method expects transformed expression, but it is not the same operation as NB residual standardization: it does not directly express each observation’s discrepancy from a fitted NB mean in model-based standard-deviation units. Conversely, residuals are not a replacement for corrected counts when an abundance-like, nonnegative summary is required.
Worked Example: Matching a Summary to Its Interpretation
Suppose a gene has fitted coefficient \(\beta_{1j}=1\), reference depth \(s_{\mathrm{ref}}=2\), and two hypothetical observations. Cell A has \(s_A=1\), \(Y_{Aj}=2\); cell B has \(s_B=2\), \(Y_{Bj}=3\). For the residual calculations, suppose the fitted intercept gives means \(\mu_{Aj}=2\) and \(\mu_{Bj}=4\), and \(\theta_j=4\).
The corrected counts are \(C_{Aj}=2(2/1)=4\) and \(C_{Bj}=3(2/2)=3\). Their mean is \((4+3)/2=3.5\) corrected counts per cell on the reference-depth scale.
For cell A, the residual variance is \(2+2^2/4=3\), and \(R_{Aj}=(2-2)/\sqrt{3}=0\). For cell B, the variance is \(4+4^2/4=8\), and \(R_{Bj}=(3-4)/\sqrt{8}\approx-0.354\). The mean residual is \((0-0.354)/2\approx-0.177\). That mean says the two observations, on average, fall slightly below their fitted means in standardized units; it is not an average expression level of \(-0.177\).
If the question is average expression after depth correction, the corrected-count mean of 3.5 is the relevant summary. If the question is whether these observations tend to depart from the fitted model, the residuals are relevant. Reporting the residual mean as an expression quantity would give it the wrong interpretation.
Assumptions, Diagnostics, and Limitations
Both outputs inherit the limitations of the fitted mean model. If the depth factors are unsuitable, the depth relationship is misspecified, or important structure is omitted, the adjustment can leave systematic patterns or create misleading comparisons. A correction factor is model-based; it cannot make an inappropriate model appropriate. The chosen reference depth also sets the scale of corrected counts. Changing that reference multiplies the corrected values by a common positive factor for a given gene under this model, so absolute corrected values depend on the reference convention.
Useful basic checks include verifying that all \(s_i\) are positive, that corrected values are nonnegative, and that cells with zero observed counts still have zero corrected counts. For residuals, verify that the sign matches \(Y_{ij}-\mu_{ij}\) and that the denominator is the square root of the NB variance, not the variance itself. Patterns of residuals across depth or unusually large residuals can indicate that the fitted model has not captured the data adequately; residuals are not a guarantee that depth effects have been removed.
An additional practical limitation is the treatment of extreme residuals. Some workflows cap or clip residuals before using them as features. Clipping changes their values and therefore changes their influence on summaries and downstream analyses; clipped values are no longer the unmodified Pearson residuals defined above. Such a choice must be stated rather than silently treating the clipped matrix as the original residual matrix.
Common Mistake
A frequent mistake is to call corrected counts and Pearson residuals interchangeable “normalized expression.” That reasoning overlooks what each number measures. Corrected counts are scaled observations and remain nonnegative; residuals are signed deviations from model expectations. The correction does not necessarily preserve rankings across cells, and a residual’s magnitude is not a transcript count. Select the output based on whether the task needs reference-scale expression or standardized model departure, then describe summaries in those terms.
Check Your Understanding
Answer each question, then compare with the solutions below.
- If \(Y_{ij}=0\) and \(s_i>0\), what are the signs or values of the corrected count and Pearson residual, assuming \(\mu_{ij}>0\)?
- For \(Y=3\), \(\beta_1=1\), \(s=2\), and \(s_{\mathrm{ref}}=1\), calculate the corrected count.
- With \(\mu=4\) and \(\theta=4\), calculate the Pearson residual when \(Y=2\).
- Which output should be used to report a mean on the reference-depth expression scale, and why is the other output unsuitable for that interpretation?
Solutions. (1) The corrected count is zero because multiplying zero by a positive factor leaves zero. The residual is negative: \((0-\mu)/\sqrt{\mu+\mu^2/\theta}<0\). (2) \(C=3(1/2)^1=1.5\) corrected counts on the chosen reference scale. (3) The variance is \(4+4^2/4=8\), so \(R=(2-4)/\sqrt{8}\approx-0.7071\). (4) Use corrected counts: their mean is an average of observed counts rescaled to the reference depth. A residual mean is an average signed, standardized departure from fitted means, not an expression abundance.