Why an Extreme Count Needs Context
A very large RNA-seq count can have several explanations. The sample may have a larger library, the gene may be genuinely more highly expressed in one condition, or a single sample may be unusually high relative to the fitted model. These possibilities are not interchangeable: a condition effect is supported by a pattern across samples, whereas a single-sample outlier is a departure concentrated in one observation. A large number by itself does not identify which explanation is appropriate.
The first comparison should account for library size. The next should be against the gene’s fitted count distribution, which allows for both the mean and overdispersion. Finally, ask whether the estimated gene-level condition effect changes substantially when the unusual sample is omitted. These checks answer related but distinct questions: normalization addresses scale, residuals address model fit, and influence addresses sensitivity of the gene’s fitted effect.
Model and Diagnostic Quantities
For gene \(i\) in sample \(j\), let \(y_{ij}\) be the observed count and \(s_j\) the sample’s size factor. A standard negative-binomial generalized linear model (GLM) assumes, conditional on the design and fitted parameters, that \(y_{ij}\sim NB(\mu_{ij},\alpha_i)\), with mean \(\mu_{ij}\) and variance \(\mu_{ij}+\alpha_i\mu_{ij}^2\). The dispersion \(\alpha_i\) is nonnegative and describes extra variation beyond a Poisson model. The log-link model is \(\log(\mu_{ij})=\log(s_j)+\mathbf{x}_j^\mathsf{T}\boldsymbol{\beta}_i\), where \(\mathbf{x}_j\) is the design row and \(\boldsymbol{\beta}_i\) contains gene \(i\)’s coefficients. This account assumes an appropriate design, positive size factors, adequately estimated dispersion, and no important unmodelled structure such as a batch effect confounded with condition.
A normalized count \(y_{ij}/s_j\) is useful for comparing relative expression across samples with different library sizes. It is not itself a replacement for the count-model fit: in particular, it does not account for the gene’s mean-dependent variability. For a fitted mean \(\widehat{\mu}_{ij}\) and dispersion estimate \(\widehat{\alpha}_i\), a Pearson residual is the observed-minus-fitted count divided by its fitted standard deviation:
A large positive residual flags a count that is high relative to the fitted mean and variance; a large negative residual flags a low one. Residual size is a diagnostic clue, not a definitive outlier test. Counts are discrete, parameters have been estimated, and residuals need not follow a standard normal distribution exactly. Interpret them alongside the other samples, the design, and the fit.
Worked Example: Normalizing Before Judging a Count
Suppose one gene has hypothetical counts 80, 160, and 240 in three samples, with size factors 0.5, 1.0, and 1.5. Calculate the normalized counts and decide whether the largest raw count is evidence of an outlier.
Normalize each count. The normalized counts are \(80/0.5=160\), \(160/1.0=160\), and \(240/1.5=160\). Although 240 is the largest raw count, all three normalized counts are equal.
Interpretation. The raw-count difference is explained by the sample size factors in this hypothetical example. It does not, by itself, indicate unusual expression in the third sample. This calculation does not establish that the gene follows the fitted negative-binomial model; residuals and the sample pattern must still be checked.
Use the Sample Pattern, Not Just the Maximum
A genuine condition effect should generally be visible as a group-level pattern after accounting for sample size and the model. A count that is extreme in one sample but not supported by other samples in that condition raises a different concern. Replicates matter: with few samples, one observation can strongly affect a group mean and its estimated condition coefficient. The count should be examined in the context of the design, not compared with an arbitrary universal threshold.
Worked Example: A Consistent Condition Difference
Consider hypothetical normalized counts for one gene: 20, 22, and 19 in controls, and 80, 75, and 85 in treated samples. Assume equal size factors for this illustration and a two-condition model with separate fitted means. For a simple numerical residual check, use dispersion \(\alpha_i=0.1\).
Compare the group pattern. The control mean is \((20+22+19)/3=20.33\). The treated mean is \((80+75+85)/3=80. The ratio is \(80/20.33=3.934\), so the estimated treated-to-control log base 2 ratio is \(\log_2(3.934)=1.976\). The three treated counts are all much higher than the three controls; the difference is not confined to one treated sample.
Check the treated residuals. For a fitted treated mean of 80, the negative-binomial variance is \(80+0.1(80^2)=720\), and the standard deviation is \(\sqrt{720}=26.83\). The Pearson residuals for counts 80, 75, and 85 are \(0/26.83=0\), \(-5/26.83=-0.186\), and \(5/26.83=0.186\). None is a conspicuous departure from the fitted treated mean.
Conclusion. These hypothetical data show a consistent condition-associated difference, not a single unusually high treated observation. This pattern is compatible with differential expression under the stated model assumptions; the calculation is a diagnostic illustration, not a complete inferential test.
Influence: Does One Sample Drive the Gene Effect?
A residual measures the discrepancy between one observation and its fitted mean relative to the model’s variance. Influence asks a different question: how much would the fitted coefficients change if that observation were left out? A point can have a large residual but little influence if the design gives it little leverage. Conversely, a point with substantial leverage may materially affect the fitted coefficients even if its residual is not the largest.
For a GLM, a common gene-level influence summary is Cook’s distance, an approximate measure of how much the coefficient fit changes when one sample is omitted. One standard GLM approximation uses the Pearson residual \(r_{P,ij}\), leverage \(h_{ij}\), and \(p\), the number of fitted coefficients for the gene:
Leverage describes how unusual a sample’s design position is for estimating the coefficients. In a weighted GLM, it is obtained from the diagonal of the weighted hat matrix. Cook’s distance is a diagnostic scale for comparing samples within a gene and identifying fits worth investigating; its approximate value depends on the model and implementation. There is no universal cutoff that proves an observation erroneous or mandates deletion.
Worked Example: Reading an Approximate Cook’s Distance
Suppose a hypothetical gene has \(p=2\) fitted coefficients, and one sample has Pearson residual \(r_P=4.0\) and leverage \(h=0.60\). Calculate the approximate Cook’s distance.
Substitute into the GLM approximation: \(D\approx 4.0^2(0.60)/[2(1-0.60)^2]\). The numerator is \(16(0.60)=9.6\), and the denominator is \(2(0.40)^2=0.32\). Therefore \(D\approx9.6/0.32=30\).
This value flags the observation for investigation because its residual and leverage combine to suggest potential influence. It does not establish that the count is a technical error. The next step is to inspect the data and compare the full fit with a fit omitting that sample, while retaining the original result unless there is a justified reason to handle the observation differently.
Why Leave-One-Sample Refitting Helps
A direct way to assess influence is to refit the gene model after removing one sample and compare the relevant coefficient or test statistic with the full-data result. The comparison should be made for the same design and gene. A large shift warns that the reported gene effect is sensitive to that observation. It does not tell you whether the sample is faulty: a genuine biological extreme can also be influential.
The connection between omission and coefficient change can be justified by a first-order score expansion. Let \(U(\boldsymbol{\beta})\) be the full-data score, \(u_j\) the score contribution from sample \(j\), and \(\mathbf{I}\) the full-data observed information at the full-data estimate \(\widehat{\boldsymbol{\beta}}\). The full score is approximately zero, so the reduced-data score evaluated at that same estimate is \(U(\widehat{\boldsymbol{\beta}})-u_j\approx-u_j\). Write \(\widehat{\boldsymbol{\beta}}_{(-j)}\) for the estimate after omitting sample \(j\). Using the full information as a local approximation, expand the reduced-data score around the full estimate. Since the score’s derivative is approximately \(-\mathbf{I}\), the reduced-data score at the reduced estimate is approximately \(-u_j-\mathbf{I}(\widehat{\boldsymbol{\beta}}_{(-j)}-\widehat{\boldsymbol{\beta}})\). Setting it to zero gives the approximate coefficient change below.
This is a local approximation, not a guarantee that a one-step calculation equals a full refit. It explains why an observation’s score contribution, uncertainty, and design position can all matter to influence. When the result is consequential, the actual leave-one-sample refit is more informative than relying on the approximation alone.
Worked Example: One Count Changes the Apparent Group Effect
Suppose hypothetical normalized counts are 10, 12, and 11 in controls, and 10, 11, and 100 in treated samples. With equal size factors and a model using separate means for the two groups, compare the estimated treated-to-control ratio before and after omitting the count of 100. This calculation illustrates sensitivity; it does not establish why that count is high.
Full data. The control mean is \((10+12+11)/3=11\). The treated mean is \((10+11+100)/3=121/3=40.33\). The ratio is \(40.33/11=3.667\), corresponding to \(\log_2(3.667)=1.874\).
Omitting the high count. The remaining treated mean is \((10+11)/2=10.5\). The treated-to-control ratio is \(10.5/11=0.9545\), corresponding to \(\log_2(0.9545)=-0.067\). The estimated direction and size of the apparent effect change markedly.
Conclusion. The full-data group difference is highly sensitive to the count of 100, so that sample deserves careful review. The leave-one-out calculation is a diagnostic comparison, not justification for removing the sample or reporting the reduced-data estimate as the primary result. Check sample quality, design, and biological context, and document any decision about how to handle the observation.
A Practical Diagnostic Sequence
Compare counts after applying the sample size factors. A large raw count that is ordinary after normalization is not, on that basis alone, a count outlier.
For the gene and design, inspect fitted means and residuals using the negative-binomial variance. Ask whether the observation is unusual relative to its expected mean and dispersion.
Compare the sample with other samples in its condition and inspect leverage or Cook’s distance. If needed, refit without the sample to quantify sensitivity.
Review sample identity, sequencing or processing records, and plausible biological explanations. Do not delete a count solely because it is extreme or influential.
These steps complement, rather than replace, residual diagnostics for negative-binomial GLMs. A count can look extreme on a normalized-count plot yet be plausible under a large dispersion. Alternatively, a count can appear modest in raw units but be surprising relative to a small fitted mean and its model variance. Always interpret diagnostics relative to the fitted model and experimental design.
Failure Modes and Limits
- Using raw counts alone. A larger library can produce larger counts broadly. Normalize for sample size before comparing relative expression.
- Calling every large residual an error. A large residual identifies a model discrepancy, not its cause. Check the sample pattern and data quality; a valid biological observation can be extreme.
- Treating Cook’s distance as a deletion rule. Influence means the fitted result is sensitive to a sample. It does not show that the sample is invalid, nor does a universal cutoff resolve the question.
- Ignoring the design. A missing batch term, confounding, or incorrect sample labels can make residual and influence diagnostics misleading. Revisit the design and data before interpreting a flagged point.
- Confusing a group effect with one unusual replicate. A condition difference supported across replicates is different from a fitted difference that disappears when one sample is omitted. Report that sensitivity accurately.
Small sample sizes limit the reliability of any outlier diagnosis: fitted means, dispersions, and influence measures are themselves uncertain. In particular, a leave-one-out fit has fewer observations and can be unstable. The purpose of these diagnostics is to identify and communicate sensitivity, then support a defensible investigation—not to make biological or technical judgments from a statistic alone.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions below.
- A count of 300 comes from a sample with size factor 3.0. What is its normalized count, and why is the raw count alone insufficient to label it an outlier?
- For a fitted mean of 25, dispersion 0.2, and observed count 55, calculate the Pearson residual.
- In one gene, all four treated samples are higher than all four controls. What pattern would distinguish this from a single-sample outlier?
- Why can a high Pearson residual and a high Cook’s distance answer different diagnostic questions?
- A gene’s estimated condition coefficient changes greatly after one sample is omitted. Does that prove the sample should be removed? Explain.
Solutions. (1) The normalized count is \(300/3.0=100\). Without the other normalized counts and the fitted mean and dispersion, the raw count cannot show whether the observation is unusual. (2) The fitted variance is \(25+0.2(25^2)=150\), so the residual is \((55-25)/\sqrt{150}=30/12.25=2.449\). It flags a positive departure for investigation, not a confirmed error. (3) A condition effect is supported when the separation is shared across treated samples rather than being driven by one unusually high treated count; model fit and design still matter. (4) A Pearson residual measures discrepancy relative to fitted variability, while Cook’s distance also accounts for leverage and summarizes potential change in the fitted coefficients. (5) No. It establishes sensitivity to that sample, not invalidity. Inspect data quality and biological context, and document any justified handling decision.