Why link methylation and expression?
Promoter methylation and gene expression are measurements of different molecular features. A disease cohort may show higher methylation near a gene's promoter and lower expression of that gene, but the two findings alone do not show that the samples with higher methylation are the samples with lower expression. Nor does an inverse association establish that methylation caused the expression change. The analysis must preserve which measurements came from the same subject and make clear whether the question concerns a disease-related difference, a within-disease association, or both.
For a focused analysis, define the promoter region and its probes before examining the results. Summarize methylation for that promoter in each sample using a prespecified method, such as the mean of quality-controlled probe beta values. Beta values are proportions on a scale from 0 to 1. Pair each subject's promoter summary with expression for the corresponding gene from the same biological sample. Expression values are commonly analyzed on a normalized log scale; the scale and units must be stated.
Separate paired disease differences from within-disease association
If each disease sample has a matched reference sample from the same subject, calculate the within-subject difference for each measurement. For methylation, \(D_{Mi}=M_{i,\mathrm{disease}}-M_{i,\mathrm{reference}}\); define \(D_{Ei}\) analogously. A paired t procedure estimates the mean of these differences. Its assumptions are that pairs are independent across subjects and that the differences are approximately normally distributed, especially when the number of pairs is small. It does not require disease and reference values within a pair to be independent.
The reason pairing matters is visible from the variance of a difference:
When measurements from the same subject are positively correlated, the covariance term reduces the variability of the paired differences. Treating the samples as unrelated discards that information and answers a different question. Conversely, paired disease differences do not establish that methylation and expression covary across diseased subjects. That requires an analysis of the matched measurements within the disease samples.
Worked Example: Paired Disease and Reference Differences
Suppose eight hypothetical subjects each contribute a disease sample and a matched reference sample. The promoter beta-value differences (disease minus reference) are 0.08, 0.10, 0.06, 0.12, 0.09, 0.05, 0.11, and 0.07. The corresponding expression differences, in log2 expression units, are -0.6, -0.8, -0.4, -1.0, -0.7, -0.3, -0.9, and -0.5. Assume the pairs are independent across subjects and the difference distributions are reasonably compatible with a t analysis.
For methylation, the mean difference is \(\bar{D}_M=0.68/8=0.0850\) beta-value units. The sum of squared deviations from this mean is 0.00420, so \(s_{D_M}=\sqrt{0.00420/7}=0.02449\), and \(SE=s_{D_M}/\sqrt{8}=0.008660\). Thus \(t=0.0850/0.008660=9.815\), with \(df=7\). A 95% confidence interval, using \(t_{0.975,7}=2.365\), is \(0.0850\pm2.365(0.008660)=(0.0645,0.1055)\). In this hypothetical cohort, the mean promoter methylation difference is estimated as an increase of 0.0850 beta-value units; the interval describes uncertainty in the population mean paired difference under the model.
For expression, \(\bar{D}_E=-5.2/8=-0.6500\) log2 units. The squared deviations sum to 0.420, giving \(s_{D_E}=\sqrt{0.420/7}=0.2449\) and \(SE=0.2449/\sqrt{8}=0.08660\). Therefore \(t=-0.6500/0.08660=-7.506\), with \(df=7\). The 95% interval is \(-0.6500\pm2.365(0.08660)=(-0.8548,-0.4452)\) log2 units. Both paired tests give evidence of a mean disease-versus-reference difference, but these two results alone do not quantify the methylation–expression relationship among disease samples.
Estimate the association in matched disease samples
For the within-disease question, a simple linear regression uses one matched pair of values per subject:
Here \(\beta_1\) is the expected change in expression, in its stated units, for a one-unit increase in promoter beta value. The model assumes a linear conditional mean, independent errors across subjects, constant error variance, and normally distributed errors for exact small-sample t-based inference. The methylation measurement must be aligned to the expression measurement from the same subject. A slope is an association conditional on this model, not a causal effect.
The ordinary least-squares slope is \(\hat{\beta}_1=S_{ME}/S_{MM}\), where \(S_{ME}=\sum_i(M_i-\bar{M})(E_i-\bar{E})\) and \(S_{MM}=\sum_i(M_i-\bar{M})^2\). A useful derivation of its uncertainty follows from the model: the fitted residuals estimate the error variance by \(s^2=\sum_i(E_i-\hat{E}_i)^2/(n-2)\). Conditional on the observed methylation values, the slope has variance \(\sigma^2/S_{MM}\), because the slope is a weighted sum of independent errors with weights \((M_i-\bar{M})/S_{MM}\). Replacing \(\sigma^2\) by \(s^2\) gives \(SE(\hat{\beta}_1)=\sqrt{s^2/S_{MM}}\). The test of \(H_0:\beta_1=0\) uses \(t=\hat{\beta}_1/SE(\hat{\beta}_1)\) with \(df=n-2\), under the stated assumptions.
Worked Example: Inverse Association Within the Disease Group
Consider six hypothetical disease samples with promoter beta values \(M=(0.20,0.30,0.40,0.50,0.60,0.70)\) and log2 expression values \(E=(7.9,7.7,7.2,6.8,6.3,6.1)\). The measurements are matched by subject. Assume independent subjects, a linear relationship, and no severe outlier or variance pattern; with only six subjects, those assumptions require cautious assessment.
For convenient arithmetic, set \(z=(M-0.45)/0.10\), giving \(z=(-2.5,-1.5,-0.5,0.5,1.5,2.5)\). The expression mean is 7.0, and direct calculation gives \(S_{zz}=17.5\), \(S_{zE}=-6.8\), and \(S_{EE}=2.68\). The fitted slope per one-unit increase in \(z\) is \(-6.8/17.5=-0.3886\) log2 units. Since one unit of \(z\) is 0.10 beta-value units, the slope per one beta-value unit is \(-0.3886/0.10=-3.886\) log2 units. The fitted line is approximately \(\hat{E}=8.749-3.886M\).
The residual sum of squares is \(SSE=S_{EE}-S_{zE}^2/S_{zz}=2.68-46.24/17.5=0.03771\). Thus \(s^2=0.03771/(6-2)=0.009429\). The slope standard error on the \(z\) scale is \(\sqrt{0.009429/17.5}=0.02321\), and the t statistic is \(-0.3886/0.02321=-16.74\), with \(df=4\). The two-sided p-value is less than 0.001 under the regression assumptions. A 95% interval uses \(t_{0.975,4}=2.776\): on the beta-value scale it is \(-3.886\pm2.776(10\times0.02321)=(-4.530,-3.242)\) log2 units per one-unit increase in beta value. This small hypothetical data set illustrates a strong inverse association, but the estimate is sensitive to individual samples and should not be interpreted as causal.
Multiple genes, diagnostics, and alternatives
An analysis may test methylation–expression association for many genes. Testing each at level \(\alpha\) without adjustment raises the chance of false discoveries across the collection. The Benjamini–Hochberg procedure controls the false discovery rate under independence and also under certain positive dependence conditions. Sort the \(m\) p-values as \(p_{(1)}\le\cdots\le p_{(m)}\); find the largest \(k\) for which \(p_{(k)}\le kq/m\), where \(q\) is the target false discovery rate, and reject the hypotheses corresponding to ranks 1 through \(k\). Dependence among gene-level tests should be considered; for arbitrary dependence, a more conservative procedure such as Benjamini–Yekutieli may be appropriate.
Worked Example: Applying Benjamini–Hochberg to Gene-Level Tests
Suppose five hypothetical genes have valid, unadjusted association-test p-values 0.001, 0.008, 0.021, 0.040, and 0.200. With target \(q=0.05\), the ordered Benjamini–Hochberg cutoffs are \(0.05(1/5)=0.010\), 0.020, 0.030, 0.040, and 0.050. The first four p-values meet their corresponding cutoffs, including \(0.040\le0.040\); the fifth does not meet its cutoff. The largest qualifying rank is \(k=4\), so the first four hypotheses are rejected by this procedure. This controls the expected false discovery proportion under its conditions; it does not mean that each selected gene has a 5% probability of being a false finding.
Diagnostics should address both the data construction and the fitted relationship. Confirm that promoter probes map to the intended promoter and that each expression value belongs to the same subject and biological sample as its methylation summary. Inspect scatterplots with fitted lines, residuals versus fitted values, and residual distributions. Curvature suggests that a straight-line slope may be inadequate; a single high-leverage sample can dominate the slope; and increasing residual spread across fitted values challenges constant variance. With small cohorts, plots can reveal concerns but cannot reliably verify normality or rule out influential observations.
Also consider batch, disease severity, age, treatment, and tissue composition when relevant to the design. For example, differing proportions of tumor and non-tumor cells can shift both methylation and expression and create an association that does not represent a within-cell relationship. A regression that includes prespecified covariates can estimate an adjusted association, but requires adequate sample size, appropriate model form, and no severe collinearity. Report the covariates and the interpretation as conditional on them. Avoid choosing probes, genes, covariates, or transformations because they produce the most favorable result without accounting for that selection.
| Approach | Question answered | Main consideration |
|---|---|---|
| Paired disease-reference tests | Does the mean methylation or expression differ between matched conditions? | Analyze within-subject differences; does not itself test their association. |
| Regression in disease samples | How does expression vary with methylation across matched disease samples? | Provides a slope and model-based inference; check linearity, residuals, and confounding. |
| Spearman rank association | Do higher methylation ranks tend to accompany higher or lower expression ranks? | Can capture monotone patterns without a linear-scale assumption, but does not estimate a change in expression units per beta-value unit. |
| Unpaired group comparison | Do two unrelated groups differ? | Not a substitute for a paired analysis when the design contains matched subjects. |
A common mistake is to observe a significant increase in disease-sample promoter methylation and a significant decrease in disease-sample expression, then conclude that methylation predicts expression. Those are marginal condition differences, not a test of subject-level association. The correction is to analyze the matched \(M_i,E_i\) pairs within the relevant group, report the estimated association and uncertainty, and distinguish it from paired disease effects. A further limitation remains: observational matched measurements cannot, by themselves, establish whether methylation changes expression, expression changes methylation, or other factors influence both.
Check Your Understanding
Answer each question, then compare with the solutions below.
- A study has disease and reference samples from the same eight subjects. What values should be analyzed for a paired test of mean methylation change?
- In the regression example, what does the slope of \(-3.886\) mean, and what does it not establish?
- Why can two significant disease-versus-reference differences fail to demonstrate a methylation–expression association?
- When might Spearman association be preferred to simple linear regression, and what information is then not directly provided?
- For five tests at \(q=0.05\), the fourth ordered p-value is 0.040 and the fifth is 0.200. What does the Benjamini–Hochberg rule conclude if the first three also meet their rank cutoffs?
Solutions. (1) Calculate one disease-minus-reference methylation difference per subject and analyze those eight independent pair differences. (2) Under the fitted linear model, a one-unit increase in beta value is associated with an estimated 3.886-log2-unit decrease in expression; the slope is not a causal effect. (3) The two tests address marginal mean changes, not whether subjects with higher methylation tend to have lower expression. (4) Spearman can suit a monotone but nonlinear relationship or data where ranks are more defensible; it does not directly give a slope in expression units per beta-value unit. (5) Rank four is the largest qualifying rank because \(0.040\le0.05(4/5)=0.040\), while rank five fails because \(0.200>0.050\); reject the first four hypotheses under the procedure's dependence conditions.