Sampling Uncertainty In Genome-Wide Correlations

Why Correlation Uncertainty Matters

A co-expression analysis may estimate a correlation for every pair of genes, but each estimate is based on a finite set of biological samples. Even when the underlying association is weak or absent, sampling variation produces nonzero correlations. With thousands of genes, the number of pairs can be enormous, so ordinary-looking estimation noise can create many large-looking entries in a correlation matrix. The issue is not only whether one pair's correlation is uncertain; it is whether the matrix as a whole can be trusted as a stable description of relationships.

This tutorial focuses on uncertainty for the Pearson correlation. Let \(X_{ig}\) and \(X_{ih}\) denote the expression measurements for genes \(g\) and \(h\) in sample \(i\), for \(i=1,\ldots,n\). Write \(\rho_{gh}\) for their population correlation in the population of biological units of interest, and \(r_{gh}\) for the sample Pearson correlation. These are measures of marginal association, not evidence that one gene causes changes in the other.

Assumptions and the Confidence-Interval Procedure

For a standard Fisher-transformation interval, the \(n\) pairs of measurements for the two genes should come from independent biological sampling units. Within a sampled unit, the two measurements form a pair. The usual approximation is justified most directly when the pairs follow a bivariate normal distribution; in practice it is an approximation whose accuracy improves with sample size. The measurements should be sufficiently well behaved that a few extreme observations do not dominate the Pearson correlation. The usual formula also assumes the same \(n\) complete pairs are used for both genes.

Under these conditions, transform the observed correlation, calculate an approximately normal interval on the transformed scale, then apply the hyperbolic tangent to its endpoints. For a two-sided confidence level of \(100(1-\alpha)\%\), the interval is:

Here \(z_{1-\alpha/2}\) is the corresponding standard-normal quantile. For a 95% interval it is approximately 1.96. The reported interval is for the population correlation \(\rho\), not for a transformed quantity. Under repeated samples from the stated model, the procedure produces intervals that cover the same \(\rho\) about 95% of the time. It does not mean that a particular interval has a 95% probability of containing a fixed \(\rho\).

How Sample Size Controls Uncertainty

On the Fisher scale, the standard error depends on \(n\), not on the observed \(r\): it is \(1/\sqrt{n-3}\). Increasing the number of independent samples narrows the transformed interval, and therefore generally narrows the interval after back-transformation as well. Conversely, a correlation estimated from a small study may be quite uncertain even when its point estimate looks large. The transformation is especially useful because the raw correlation scale is bounded between -1 and 1 and has a sampling distribution that can be skewed, particularly when the true correlation is far from zero.

Why a Genome-Wide Matrix Can Be Unstable

For \(G\) genes, the number of distinct off-diagonal gene pairs is \(G(G-1)/2\). Each pair has its own sampling uncertainty, and the uncertainty does not disappear merely because the matrix contains many entries. In fact, a large number of pairs makes it likely that some estimated correlations will be far from their population values. The entries are also dependent: correlations that share a gene are computed from overlapping measurements, so one cannot generally treat all pairwise estimates as independent.

A second warning comes from matrix rank. If all genes are measured on the same \(n\) samples and the data are centered by gene, the sample expression matrix has rank at most \(n-1\). Its sample correlation matrix therefore also has rank at most \(n-1\). When \(G>n-1\), the matrix is necessarily singular: some directions in the \(G\)-dimensional space cannot be estimated from the available samples. This is a structural limitation of the sample information, not a numerical software defect.

Diagnostics, Alternatives, and Limits

Fisher intervals address uncertainty for an individual Pearson correlation under a stated sampling model. Before interpreting them, check that observations correspond to independent biological units, that both genes use the same complete samples, and that the scatterplot is not dominated by outliers or strong nonlinear structure. A batch effect, unmodeled group structure, or repeated measurements can make the ordinary independent-pair model inappropriate. In those cases, an interval calculated from the nominal \(n\) may be too narrow or may describe a different estimand than intended.

A bootstrap can provide an alternative uncertainty estimate: resample biological units, recompute the correlation for each resample, and use the resulting distribution to form an interval. This can be useful when the Fisher approximation is questionable, but it still requires a defensible sampling unit and enough independent units to represent the population. Resampling genes instead of biological samples does not estimate sampling uncertainty in the observed subjects. For rank-based correlations, the Pearson Fisher formula is not automatically an exact interval procedure; use an interval method matched to the chosen correlation measure.

Fisher intervals are typically pointwise: their nominal coverage applies to one specified gene pair. They do not imply that all intervals across a genome-wide analysis jointly cover their correlations at the same confidence level. Nor do they make a noisy matrix stable. A small sample size, many genes, influential samples, or changes in preprocessing can cause estimated entries and downstream network summaries to vary substantially. Stability should be examined by checking sensitivity to reasonable analysis choices and, where possible, by comparing estimates across independent samples. A pointwise interval is a measure of uncertainty, not a guarantee of reproducibility.

A common mistake is to interpret an estimated correlation of 0.30 as compelling evidence of a gene link simply because it exceeds a chosen numerical threshold. That reasoning ignores both the sample size and the number of pairs examined. The correction is to quantify uncertainty for the relevant pair, recognize that thousands of noisy estimates are being considered, and avoid treating a sample correlation as a stable biological relationship without supporting evidence. A large point estimate can have a wide interval, and a matrix with many entries can be unstable even when each individual calculation is correct.