Mean-Variance Trends in Microarray Models

Why Let the Variance Prior Depend on Expression?

In expression microarrays, genes can differ not only in their average measured intensity but also in how variable their measurements are. In some data sets, the residual variance tends to change systematically across the range of average expression. A single variance prior shared by every gene can then be a poor reference for genes at one end of that range. A mean-variance trend addresses this by allowing the prior variance to depend on average expression.

The trend does not say that a gene’s average expression determines its biological variability, nor does it replace the gene-wise linear model. It describes a pattern across genes that informs the prior used to stabilize noisy residual-variance estimates. The earlier material on estimating the variance prior across genes supplies the empirical Bayes framework; here the specific change is that the prior scale varies with average expression.

A typical gene-wise model is \(y_g=X\beta_g+\varepsilon_g\), where \(y_g\) is the vector of observed expression values for gene \(g\), \(X\) is the design matrix, and \(\beta_g\) is the vector of gene-specific coefficients. Under the standard linear-model assumptions, the errors have mean zero, variance \(\sigma_g^2\), and are independent with a normal distribution conditional on the model. The design must be of adequate rank, and the residual variance estimate must have positive residual degrees of freedom. The average expression \(A_g\) is calculated from the normalized expression measurements; it is a predictor of the prior scale, not the residual variance estimate itself.

The trend is useful when the expected residual variance changes with average expression. For example, measurement noise may be larger among dim probes, or variability may follow another smooth pattern. The relationship is data-dependent: neither its direction nor its shape should be assumed in advance.

Prior and Moderated Variance

Conditional on the average expression, a common formulation gives each gene a scaled inverse-chi-squared prior for its residual variance. The prior scale is the trend value at that gene’s average expression. The gene-wise residual variance estimate supplies the likelihood information, and empirical Bayes combines it with the prior.

The formula shows exactly how the trend changes moderation. A gene with a small residual degrees of freedom receives more weight from its prior scale. A gene with a large residual degrees of freedom receives relatively more weight from its own residual variance estimate. Since the prior scale is local to average expression, two genes with equal \(s_g^2\) can receive different moderated variances.

For a contrast \(c^\mathsf{T}\beta_g\), let \(\hat\theta_g=c^\mathsf{T}\hat\beta_g\), and let \(v_g\) be the design-based variance factor for that contrast. If the coefficient estimate has variance \(\sigma_g^2v_g\), its moderated standard error and t-statistic are

Under the model and prior assumptions, the moderated statistic is compared with a t distribution having \(d_g+d_0\) degrees of freedom. The variance trend changes the estimated scale in the denominator; it does not change the observed coefficient estimate or the contrast’s design-based variance factor.

Why the Weighted Formula Holds

The gene-wise sampling model contributes a likelihood for \(\sigma_g^2\) proportional to \((\sigma_g^2)^{-d_g/2}\exp[-d_gs_g^2/(2\sigma_g^2)]\). The scaled inverse-chi-squared prior contributes a factor proportional to \((\sigma_g^2)^{-(d_0/2+1)}\exp[-d_0s_{0g}^2/(2\sigma_g^2)]\). Multiplication adds the powers and the exponential terms:

This has the scaled inverse-chi-squared form, with degrees of freedom \(d_0+d_g\). Matching its scale term to the usual form \((d_0+d_g)s_{g,\mathrm{mod}}^2\) gives \(s_{g,\mathrm{mod}}^2=(d_0s_{0g}^2+d_gs_g^2)/(d_0+d_g)\). Thus the weighted estimate follows directly from combining the prior and sampling information. It is the posterior scale used for moderated inference, not in general the posterior mean of \(\sigma_g^2\).

Fitting, Diagnostics, and Alternatives

In practice, average expression and gene-wise residual variance estimates are calculated from the normalized data and fitted gene-wise models. A smooth relationship between average expression and variance information is then estimated across genes. Because individual residual variances are noisy, the trend should represent the broad relationship rather than follow each gene’s estimate. Estimation methods may account for the sampling variability of log residual variances; simply drawing a curve through raw variances can give excessive influence to noisy genes.

The resulting fitted scale \(s_0^2(A_g)\) must be positive and defined across the expression range where it is used. The model also relies on the underlying gene-wise linear models being appropriate, the residual variance estimates having their stated sampling behavior, and the trend being a reasonable summary of how prior variance changes with average expression. The empirical Bayes borrowing further assumes that genes provide useful information about the prior after conditioning on the modeled trend. Genes need not have identical variances, but a smooth common relationship must be informative.

Useful diagnostics include a plot of residual variance information against average expression, with the fitted trend shown; the distribution of departures from that trend; and checks for systematic residual patterns that the smooth curve misses. Inspect whether the fitted curve is being driven by a small number of genes or behaves erratically near the ends of the expression range. A flat or unstable fit, a sharp bend supported by few genes, or a broad pattern of residuals above the curve in one region and below it in another suggests the trend is not capturing the structure adequately.

The main alternatives answer related but distinct needs. A common variance prior is simpler and can be appropriate when there is little evidence of a mean-variance relationship. No moderation uses each gene’s residual variance directly, avoiding prior misspecification but retaining its small-sample instability. A trend prior can improve the local reference for variance moderation when the pattern is real, but it adds a modeling choice and can mislead if the smooth fit is poor. Modeling a mean-dependent variance through observation-level weights is another possible strategy; it changes the variance representation in the linear model rather than merely choosing a gene-specific prior scale. These approaches are not interchangeable without checking their assumptions and the data structure.

Common Pitfalls and Limitations

• Confusing mean expression with residual variance. Average expression predicts the prior scale across genes; it is not a substitute for the gene-wise residual variance.
• Assuming the trend must slope downward. The direction and shape must be supported by the data. A trend that is forced to decrease can impose a false prior pattern.
• Assuming a high or low prior scale determines the gene’s result. The moderated estimate also depends on \(s_g^2\), \(d_g\), and \(d_0\). Genes with more residual degrees of freedom are less influenced by the prior.
• Treating the moderated variance as the posterior mean. The weighted formula is the posterior scale used for moderated inference. The posterior mean, when it exists, has a different denominator.
• Ignoring model or fit failures. A smooth trend cannot repair an incorrect design, unmodeled dependence, poor normalization, or a relationship that changes in ways the fitted curve cannot represent.

Trend-based moderation is valuable when there is a reproducible, sufficiently smooth relationship between average expression and residual variance. It is limited by the quality of the gene-wise models and the trend estimate. A fitted pattern is an empirical prior assumption, not a causal explanation for why expression measurements vary.