Tutorials › High-Dimensional Genomic Data Analysis › Shrinkage Of Differential-Expression Effects

Count-Based Differential Expression for RNA-seq · Tutorial 449 of 795

Shrinkage Of Differential-Expression Effects

You will derive and apply shrinkage of RNA-seq log fold changes, assess its assumptions and limitations, and distinguish stabilized effect estimates from tests of differential expression.

Intermediate 10 min read

What You'll Learn

  • Define log fold-change shrinkage and its target.
  • Derive the normal-prior posterior estimate from a normal likelihood approximation.
  • Calculate shrunk estimates for genes with different standard errors.
  • Explain why shrinkage changes effect estimates but does not automatically change a Wald test or its p-value.
  • Recognize when a prior or likelihood approximation may produce misleading shrinkage.
  • Compare shrinkage with unshrunk estimation and heavy-tailed prior approaches.

Why Shrink Differential-Expression Effects?

In RNA-seq analysis, a fitted count model can estimate how expression changes between conditions. For a gene with little information, that estimate may be highly variable: a small change in observed counts can produce a large change in the estimated log fold change. Shrinkage stabilizes such estimates by pulling uncertain effects toward zero, while allowing well-estimated effects to move less.

The target is important. A log fold change of zero means no difference for the specified contrast. Shrinking toward zero is not a claim that every gene has no effect; it is a way to temper noisy estimates using an assumption about the distribution of effects across genes. The result is an effect estimate, not by itself a test of whether a gene is differentially expressed.

Model and Definition

For gene \(g\) in sample \(i\), a common starting point is a negative-binomial model \(Y_{gi}\sim NB(\mu_{gi},\alpha_g)\), with variance \(\mu_{gi}+\alpha_g\mu_{gi}^2\). A log-link model writes \(\log(\mu_{gi})=o_{gi}+x_i^\mathsf{T}\gamma_g\), where \(o_{gi}\) is an offset accounting for library size or normalization, \(x_i\) contains design variables, and \(\gamma_g\) contains gene-specific coefficients. The model must represent the design and contrast of interest; the count observations are assumed independent between samples conditional on the model, with dispersion handled appropriately.

For a particular contrast, let \(\beta_g\) denote the true log fold change in base 2, and let \(\hat{\beta}_g\) be its unshrunk estimate with standard error \(s_g\). If a model reports a coefficient on the natural-log scale, convert it to base 2 by dividing by \(\log(2)\), and convert its standard error by the same factor. On the log2 scale, an effect of \(1\) corresponds to a twofold increase; an effect of \(-1\) corresponds to half the expression, for that contrast.

Definition: Effect shrinkage replaces an unshrunk estimate \(\hat{\beta}_g\) with an estimate pulled toward a target, commonly zero. The amount of shrinkage depends on the estimate's uncertainty and on an assumed or estimated distribution of effects across genes. Shrinkage stabilizes effect estimates; it does not, by definition, perform a significance test.

A simple derivation uses a normal approximation to the likelihood and a normal prior:

$$ \hat{\beta}_g\mid\beta_g \ \dot{\sim}\ N(\beta_g,s_g^2), \qquad \beta_g\sim N(0,\tau^2). $$

Here \(s_g\) is treated as the standard error from the fitted count model, and \(\tau\) describes the spread of effects across genes. The prior is centered at zero and assigns less prior weight to very large effects than to moderate ones. Combining the likelihood and prior gives a posterior distribution proportional to their product. Its exponent, apart from a constant, is

$$ -\frac{(\hat{\beta}_g-\beta_g)^2}{2s_g^2} -\frac{\beta_g^2}{2\tau^2}. $$

Collecting the terms involving \(\beta_g\) gives a quadratic with precision \(1/s_g^2+1/\tau^2\). Completing the square therefore yields posterior variance \(v_g\) and posterior mean \(\tilde{\beta}_g\):

$$ v_g=\left(\frac{1}{s_g^2}+\frac{1}{\tau^2}\right)^{-1} =\frac{s_g^2\tau^2}{s_g^2+\tau^2}, \qquad \tilde{\beta}_g =v_g\frac{\hat{\beta}_g}{s_g^2} =\frac{\tau^2}{\tau^2+s_g^2}\hat{\beta}_g. $$

Thus the estimate is a weighted version of the unshrunk estimate. The weight \(w_g=\tau^2/(\tau^2+s_g^2)\) lies between zero and one: when \(s_g\) is large relative to \(\tau\), the weight is small; when \(s_g\) is small, the weight approaches one. This calculation is exact for the stated normal likelihood and prior. In practice, the normal approximation is applied to a fitted count-model coefficient, and \(\tau\) is generally estimated from the collection of gene-wise effects. That makes the procedure empirical Bayes: information across genes helps set the prior spread, while each gene's own uncertainty determines how strongly its estimate is pulled.

Worked Estimates

Worked Example: A Noisy Estimate Is Pulled Substantially Toward Zero

Suppose a hypothetical negative-binomial model for a specified condition contrast gives a gene an unshrunk log2 fold-change estimate of \(1.20\), with standard error \(0.80\). Assume the fitted design and dispersion adequately describe the samples, and use a normal prior with \(\tau=0.50\) log2 units. The data are not precise relative to the prior spread.

The shrinkage weight is

$$ w=\frac{0.50^2}{0.50^2+0.80^2} =\frac{0.25}{0.89} =0.2809. $$

The shrunk estimate is

$$ \tilde{\beta}=0.2809(1.20)=0.3371\text{ log2 units}. $$

The unshrunk estimate corresponds to a fold change of \(2^{1.20}\approx2.30\), while the shrunk estimate corresponds to \(2^{0.3371}\approx1.263\). Under this model and prior, the gene's estimated expression is about 1.26 times as high in the numerator condition as in the denominator condition. The estimate is pulled closer to no change because its standard error is large; shrinkage does not prove that the true effect is small.

Worked Example: A Precise Estimate Is Changed Less

For a second hypothetical gene, suppose the same contrast gives \(\hat{\beta}=1.20\) log2 units and \(s=0.20\) log2 units. Use the same prior spread, \(\tau=0.50\). Assume the count-model fit is appropriate and the relatively small standard error reflects substantial information for this gene.

The weight and shrunk estimate are

$$ w=\frac{0.25}{0.25+0.20^2} =\frac{0.25}{0.29} =0.8621, \qquad \tilde{\beta}=0.8621(1.20)=1.0345\text{ log2 units}. $$

The corresponding fold change is \(2^{1.0345}\approx2.05\), compared with \(2^{1.20}\approx2.30\) before shrinkage. Both estimates indicate higher expression in the numerator condition, and the precise estimate is changed less than the noisy estimate in the preceding example. The amount of shrinkage is tied to uncertainty, not simply to whether the estimated effect is positive or negative.

Worked Example: A Large Effect and the Choice of Prior

Consider a third hypothetical gene with \(\hat{\beta}=3.00\) log2 units and \(s=0.40\) log2 units. These values come from a hypothetical count-model fit that is assumed to be valid, including an adequate treatment of dispersion and sample design. Under the same normal prior with \(\tau=0.50\), the estimate becomes

$$ w=\frac{0.25}{0.25+0.16}=0.6098, \qquad \tilde{\beta}=0.6098(3.00)=1.829\text{ log2 units}. $$

The unshrunk estimate corresponds to an eightfold change, \(2^3=8\), while the shrunk estimate corresponds to about \(2^{1.829}\approx3.55\)-fold. Although the estimate is precise relative to the first example, this zero-centered normal prior still pulls a very large effect appreciably inward. That is a consequence of the prior's relatively thin tails, not evidence that the gene's true effect must be 1.829. In applications where very large effects are plausible, a heavy-tailed prior can preserve them more readily while still shrinking weakly informed estimates. Prior choice and its fit therefore matter, especially for extreme estimated effects.

Effect Estimation Is Not Significance Testing

A Wald test for the null hypothesis \(H_0:\beta_g=0\) uses the unshrunk estimate and its model-based standard error. The test statistic is \(z=\hat{\beta}_g/s_g\), with a two-sided p-value calculated from the standard normal distribution when the Wald approximation is appropriate. Shrinking \(\hat{\beta}_g\) does not automatically replace this statistic, change the p-value, or establish that \(H_0\) is true.

For the first worked estimate, the unshrunk Wald statistic is \(1.20/0.80=1.50\), giving a two-sided normal-approximation p-value of about \(0.134\). The shrunk estimate, \(0.3371\), is useful as a stabilized description of effect size; it is not a new Wald statistic to divide by the original standard error and present as the same test. Some methods define separate posterior or threshold-based tests, but those require their own stated model and inferential procedure. They should not be silently substituted for the original test.

This distinction is useful in reporting. A significance test addresses compatibility of the data with a specified null under a testing model. A shrunken log fold change describes an effect after regularization. A gene may have a modest p-value but a substantially shrunken estimate, or an estimated effect may be large while its uncertainty prevents strong evidence against the null. Neither quantity should be used as a synonym for the other.

Assumptions, Diagnostics, and Alternatives

The normal-prior calculation relies on a reasonably calibrated standard error and a likelihood that is adequately approximated as normal on the effect scale. The negative-binomial count model must also be appropriate for the data: the design and offsets should represent the comparison, and dispersion should be handled adequately. Shrinkage cannot repair a misspecified design, an unmodeled batch effect, an unsuitable dispersion estimate, or a problematic influential observation.

The prior adds another assumption: gene-wise effects for the coefficient being shrunk are reasonably described by a distribution centered at zero with a common spread, or by an appropriate extension of that idea. This is a population-level regularization assumption, not a statement that all genes have the same effect. A fitted prior that is too narrow can over-shrink genuine large effects; one that is too wide may do little to stabilize noisy estimates.

Useful checks include comparing unshrunk and shrunken estimates across genes, examining how the change relates to standard error, checking the fitted distribution of effects and the prior scale, and investigating striking outliers in the count-model diagnostics. Under the simple normal-prior formula, estimates with larger \(s_g\) should generally receive smaller weights. If large effects are all pulled implausibly close to zero, examine whether a thin-tailed prior is a poor match. If the standard errors are systematically unreliable, investigate the count model before interpreting the shrunken results.

ApproachWhat it providesMain caution
Unshrunk model estimateThe fitted coefficient without a cross-gene prior adjustmentCan be highly variable for genes with little information
Normal-prior shrinkageA simple, transparent estimate pulled toward zero according to uncertaintyMay pull very large effects strongly toward zero
Heavy-tailed-prior shrinkageCan stabilize noisy estimates while preserving large effects more readilyRequires an appropriate prior and does not remove model-fitting concerns

No shrinkage is a useful reference, but it can leave unstable low-information estimates that are difficult to compare across genes. A heavy-tailed prior is an alternative when a normal prior's treatment of large effects is unsuitable. The choice is about effect estimation and its assumptions; it does not replace the need to select and interpret a valid test for differential expression.

Common Pitfall

Common mistake: “The shrunken estimate is closer to zero, so the gene is not significant,” or “I can divide the shrunken estimate by the original standard error to get the original Wald test.” Correction: Shrinkage modifies the effect estimate under a prior model. The original Wald test is based on the unshrunk estimate and its standard error. Use a separately defined inferential procedure if the goal is to test a posterior or threshold-based claim.

A second pitfall is interpreting a shrunk log2 fold change as a change in raw counts. It describes a relative expression contrast on the log2 scale; converting it to a fold change requires exponentiation with base 2. Neither the log2 estimate nor its fold-change version establishes causation.

Key takeaway: Shrinkage stabilizes uncertain gene-wise log fold changes by combining each gene's estimate and uncertainty with information about the distribution of effects across genes. Interpret the resulting effect estimate separately from the p-value and check both the count model and the prior assumptions.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions below.

  1. With \(\tau=0.60\), which estimate is shrunk more: a gene with \(s=0.15\) or one with \(s=0.90\)? Explain using the shrinkage weight.
  2. For \(\hat{\beta}=0.80\), \(s=0.60\), and \(\tau=0.40\), calculate the normal-prior shrunk estimate.
  3. Why does a shrunken estimate not by itself provide the original Wald-test p-value?
  4. What might it suggest if a normal prior pulls nearly all very large effects far toward zero?

Solutions. 1. The gene with \(s=0.90\) is shrunk more. Its weight is \(0.60^2/(0.60^2+0.90^2)=0.36/1.17\approx0.3077\); for \(s=0.15\), it is \(0.36/(0.36+0.0225)\approx0.9412\). 2. The weight is \(0.40^2/(0.40^2+0.60^2)=0.16/0.52\approx0.3077\), so \(\tilde{\beta}=0.3077(0.80)\approx0.2462\) log2 units. 3. The original Wald statistic is defined using the unshrunk estimate divided by its standard error; shrinkage changes the estimate and does not automatically define a new test. 4. The prior may be too thin-tailed or too narrow for the effect distribution, so a heavy-tailed prior or a reassessment of prior fitting may be warranted.

Clinical Biostats