A methylome-wide association study of major depression with out-of-sample case-control classification and trans-ancestry comparison

Abstract

Major depression (MD) is a leading cause of global disease burden, and both experimental and population-based studies suggest that differences in DNA methylation may be associated with the condition. However, previous DNA methylation studies have, so far, not been widely replicated, suggesting a need for larger meta-analysis studies. Here we conducted a meta-analysis of methylome-wide association analysis for lifetime MD across 18 studies of 24,754 European-ancestry participants (5,443 MD cases) and an East Asian sample (243 cases, 1,846 controls). We identified 15 CpG sites associated with lifetime MD with methylome-wide significance. The methylation score created using the methylome-wide association analysis summary statistics was significantly associated with MD status in an out-of-sample classification analysis (area under the curve 0.53). Methylation score was also associated with five inflammatory markers, with the strongest association found with tumor necrosis factor beta. Mendelian randomization analysis revealed 23 CpG sites potentially causally linked to MD, with 7 replicated in an independent dataset. Our study provides evidence that variations in DNA methylation are associated with MD, and further evidence supporting involvement of the immune system.

Keywords: DNA methylation; Depression; Epigenomics; Predictive markers.

Fig. 1. Meta-analysis of MWAS for MD.

a, Workflow of meta-analysis. b, Manhattan plot for the meta-analysis of MWAS for MD. Each dot represents a CpG probe. The x axis represents the relative positions of the probes in the genome. The y axis represents −log₁₀-transformed, two-sided P values. The red and yellow dashed lines represent the significance threshold for Bonferroni and FDR correction, respectively.

Fig. 2. Comparison between basic and complex models.

a,b, In both comparison of effect sizes indicated by standardized coefficients of linear regression (a) and comparison of double-sided P values (b), each dot represents a CpG probe. In a, the x axis represents effect size, and in b, the x axis represents −log₁₀-transformed, double-sided P value. Error bars represent ± one standard error of the mean. Differences in P values reflect both the model used and the sample sizes. The y axis represents individual CpG sites. The yellow dashed line in b represents the significance threshold for nominal significance (double-sided P < 0.05).

Fig. 3. Out-of-sample classification of MD using MS.

a, Box plot for case–control classification using MS created using P-value threshold at 5 × 10⁻⁸. MS was standardized and residualized against age, sex and values for aryl hydrocarbon receptor repressor (AHRR) probes. Residuals of MSs were used for the plot. The x axis represents MS, and the y axis represents MD case and control groups. Box boundaries are determined by the top and bottom quartiles. Box whiskers are defined by 1.5 interquartile range. Any values outside of the box whiskers are presented as individual data points. Median values per group are indicated as vertical lines within the boxes. b, Statistics for MS association tests. Cohen’s d was calculated by comparing MD case–control differences in MS. Double-sided P values are reported. P_T, P threshold for creating the MS.

Fig. 4. MR analysis of the causal effect of DNAm on MD.

a, P plot for discovery MR analysis in GS. Each dot represents a CpG site. The x axis represents chromosomes and base pair position. The y axis represents −log₁₀-tranformed P value of Wald ratio MR analysis. The gray dashed line shows the FDR-corrected significance threshold. Discovery MR analysis was performed on CpG sites available in the EPIC array. b, Replication analysis in GoDMC for the significant CpGs found in the discovery analysis. Replication analysis was performed on CpG sites available on the 450K array. Of the 23 significant CpG sites in the discovery analysis, 19 were available on both EPIC and 450K arrays and therefore were included in the replication analysis. The x axes represent effect size and −log₁₀-transformed P value on the left and right panels. The y axis represents individual CpG sites. The yellow dashed line in the right panel represents FDR-significance threshold.