The epigenome is thought to mediate between genes and the environment, particularly in response to adverse life experiences. Similar to other psychiatric diseases, the suicide liability of an individual appears to be influenced by many genetic factors of small effect size as well as by environmental stressors. To identify epigenetic marks associated with suicide, which is considered the endpoint of complex gene-environment interactions, we compared the cortex DNA methylation patterns of 6 suicide completers versus 6 non-psychiatric sudden-death controls, using Illumina 450K methylation arrays. Consistent with a multifactorial disease model, we found DNA methylation changes in a large number of genes, but no changes with large effects reaching genome-wide significance. Global methylation of all analyzed CpG sites was significantly (0.25 percentage point) lower in suicide than in control brains, whereas the vast majority (97%) of the top 1,000 differentially methylated regions (DMRs) were higher methylated (0.6 percentage point) in suicide brains. Annotation analysis of the top 1,000 DMRs revealed an enrichment of differentially methylated promoters in functional categories associated with transcription and expression in the brain. In addition, we performed a comprehensive literature research to identify suicide genes that have been replicated in independent genetic association, brain methylation and/or expression studies. Although, in general, there was no significant overlap between different published data sets or between our top 1,000 DMRs and published data sets, our methylation screen strengthens a number of candidate genes (APLP2, BDNF, HTR1A, NUAK1, PHACTR3, MSMP, SLC6A4, SYN2, and SYNE2) and supports a role for epigenetics in the pathophysiology of suicide.
Suicide is a major preventable health problem accounting for ∼0.1% of deaths in developed countries, with men's suicide rates being 3-4 times that of women [OECD, 2011]. In the vast majority of cases, suicide is a symptom (or sequelae) of psychiatric disease, such as major depression or bipolar disorder. In addition, its etiology includes complex interactions between genetic dispositions, personality styles and psychosocial factors [Hawton and van Heeringen, 2009]. For example, the MARCKS gene shows an association with suicide and violent behavior [Punzi et al., 2014] and FKBP5 with suicide and emotional responses to stress [Pérez-Ortiz et al., 2013]. Apart from mental disease-specific mechanisms, there are further independent neurobiological factors predisposing to suicide [Mann et al., 2001; Yen et al., 2009]. Although the heritability of suicidal behavior is estimated to be around 50% [Mann et al., 2001], so far genome-wide association studies did not find genetic variants of large effect size [Perlis et al., 2010; Schild et al., 2013]. In light of this ‘missing heritability', it is plausible to assume that the predisposition to suicide is mainly due to changes in gene regulation, which can be influenced by environmental factors, rather than changes in the DNA sequence itself, such as single nucleotide polymorphisms (SNPs) and copy number variations.
Epigenetic mechanisms controlling the temporally and spatially highly coordinated gene expression patterns may underlie brain pathophysiology associated with suicide. Epigenetic information is not encoded by the DNA sequence itself but by reversible modifications of DNA and/or histones, that can be transmitted from cells to daughter cells. Promoter DNA methylation during development, differentiation or disease processes is associated with posttranslational histone modifications that lead to a locally condensed inactive chromatin structure and gene silencing [Jaenisch and Bird, 2003; Hatchwell and Greally, 2007; Feil and Fraga, 2012]. It is generally assumed that epigenetics provides the molecular mechanisms for gene-environment interactions. Accumulating experimental evidence in both animal models and humans suggests that adverse, in particular stress-related experiences, can lead to persistent epigenetic changes in the brain that can influence metabolism, behavior, cognition, mood, and neuroendocrine functions [Murgatroyd et al., 2009; Szyf, 2011; Hoffmann and Spengler, 2012; Lutz and Turecki, 2014].
Candidate gene studies [Poulter et al., 2008; Ernst et al., 2009; Guipponi et al., 2009; Klempan et al., 2009; Keller et al., 2010; Fiori and Turecki, 2011; Gross et al., 2013; Guintivano et al., 2014] have associated cortex DNA methylation changes in several genes, i.e. AMD1, ARG2, BDNF, GABRA1, NTRK2, QKI, SKA2, and SSAT, with suicidal behavior. Genome-wide analyses of hippocampal DNA methylation patterns in 2 groups of suicide completers showed dysregulation of several hundred loci each, but without much overlap between groups [Labonté et al., 2012, 2013]. This may be explained by different early life experiences between groups (with and without early childhood trauma). On the other hand, such types of studies always rely on a limited number of brain samples, which makes it difficult to reach genome-wide significance for individual genes with small effect sizes. Genome-wide expression analyses also revealed several hundred suicide-specific changes in the prefrontal cortex [Sequeira et al., 2012; Zhurov et al., 2012]. Gene ontology analysis can serve for the interpretation and organization of dysregulated candidate genes into known pathways/processes.
The epigenome is the sum of the epigenetic modifications which bring a cellular phenotype into being, and each of the many different cell types and tissues of the body is endowed with a specific combination of silenced and expressed genes. Functional impairment of the prefrontal cortex, in particular dysfunction of serotonergic neurotransmission, is thought to play an important role in the development of suicidal behavior [van Heeringen, 2001; Raust et al., 2007; Sequeira et al., 2012]. Since genome-wide methylation patterns of suicide completers so far were only reported for hippocampal DNA [Labonté et al., 2012, 2013], we here compared prefrontal cortical DNA methylation in male suicide completers and controls.
Materials and Methods
This study was approved by the ethics committee of the Landesärztekammer Rheinland-Pfalz. Brain samples (online suppl. table 1; for all online suppl. material, see www. karger.com/doi/10.1159/000435778) were obtained within 1 day postmortem from 6 suicide completers and 6 non-psychiatric sudden-death comparison subjects. Only 1 individual in the suicide group (S2 had suffered from depression) and none of the controls had a medical history of psychiatric disease. All samples were from Caucasian males. There was no between-group difference in age or post-mortem time (Mann-Whitney U test). Because previous studies [Labonte et al., 2012, 2013] suggested that suicide-associated DNA methylation patterns are independent of the underlying psychiatric disorders and because of the small sample size, we did no further subgrouping of the suicide completers, as to the presence of comorbidities, medication, early-life adversity, or other possible confounding factors.
Cortex tissue was dissected from the frontal pole (Brodmann area BA10) and immediately frozen at -80°C until further use. Genomic DNA was extracted using the Precellys Tissue DNA Kit (PEQLAB, Erlangen, Germany). Amount and quality of cortex DNA were determined with a Nanodrop spectrophotometer (NanoDrop, Wilmington, Del., USA). The ratio of absorbance at 260 versus 280 nm was around 1.8 for all samples, indicative of pure DNA. Bisulfite conversion was performed with 1 µg DNA each using the EZ-96 DNA methylation kit (Zymo Research, Irvine, Calif., USA).
After bisulfite conversion, the 12 samples (of 6 suicide completers and 6 controls) were whole-genome amplified, enzymatically fragmented, and hybridized to 1 Illumina HumanMethylation450 BeadChip according to the manufacturer's protocol (Illumina, San Diego, Calif., USA). The arrays were scanned with an Illumina iScan. This assay allows quantitative measurements at >485,000 CpG sites per sample at single-nucleotide resolution, covering 99% of RefSeq genes (promoter, first exon, gene body, 5′ and 3′ UTRs) and 96% of CpG islands.
Microarray data were exported as idat files and phenotype data as tab delimited csv files into the R-based RnBeads package (http://rnbeads.mpi-inf.mpg.de). Prefiltering was performed to remove probes containing more than 2 SNPs (n = 18,998). The Greedycut algorithm then iteratively removed additional probes (n = 61) showing the highest impurity (online suppl. fig. 1). Background noise was subtracted using ‘Noob' [Triche et al., 2013], and signal intensity was normalized, as recommended by the Illumina protocol. Tools for background subtraction and data normalization are implemented in the Lumi package [Du et al., 2008, 2010]. Data quality was inferred from positive and negative control probes on the array. The generated control plots (online suppl. fig. 2) report information on (1) the signal intensity distribution for both channels, (2) the intensity for each sample measured with the different probe types, and (3) the intensity distribution for ~600 negative control probes. All samples passed the quality thresholds for staining efficiency (online suppl. fig. 2A), hybridization efficiency (online suppl. fig. 2B), bisulfite conversion rate for Infinium I and II probes (online suppl. fig. 2C, D), allele-specific extension of Infinium I and II probes (online suppl. fig. 2E, F), extension efficiency of A, T, G, and C nucleotides (online suppl. fig. 2G), and overall efficiency and specificity of all steps (online suppl. fig. 2H) by querying non-polymorphic bases in the genome, signal monitoring of probes without CpG dinucleotides, and calculating the scaling factor for normalization. In a post-filtering step, 3,144 probes in a non-CpG context (CC, CAG, CAH, CTG, and CTH with H = A, C, or T) were also removed.
The β values were controlled for batch effects and analyzed for CpG site-specific as well as region-specific (i.e. promoter, gene, CpG island, and tiling) methylation. A promoter is defined as the region spanning 1,500 bp upstream and 500 bp downstream of the transcription start site. In total, 30,145 Ensembl (version 69) gene promoters could be analyzed. Genes (including promoter, gene body, 5′ and 3′ UTRs) represent 30,039 Ensembl genes (version 69), downloaded with BioMart. CpG islands (n = 26,576) were downloaded from the UCSC Genome Browser (http://genome.ucsc.edu). Tiling region pairs (n = 131,992) represent contiguous DNA segments of 5,000 bp length.
Between-group methylation differences were identified using p values adjusted for multiple testing by the false discovery rate and/or combined ranks. For single CpG sites, the 3 combined ranks consisted of (1) differences in the methylation means, (2) log2 of the methylation quotient, and (3) p values calculated by a 2-sided Welch t test. For region-specific methylation differences, the combined ranks were (1) the mean of mean methylation levels across all sites in a region for the suicide group and controls, (2) the mean difference in means across all sites in a region, and (3) the log2 of the mean quotient of means across all sites in a region.
Functional relevance of our results was analyzed using the Database for Annotation, Visualization, and Integrated Discovery (DAVID), version 6.7 (http://david.abcc.ncifcrf.gov). Fisher's exact test was used to calculate the p values for enrichment in the DAVID annotation system. All p values were corrected for multiple testing with the Benjamini-Hochberg procedure. Network analysis was performed with STRING version 9.1 [Franceschini et al., 2013], a database of known and predicted protein interactions, including direct (physical) and indirect (functional) associations derived from high-throughput experiments, genomic context, conserved coexpression and previous knowledge.
To define a list of suicide candidate genes, we performed a PubMed (March 2015) search using the terms ‘gene, expression, suicide, variant, methylation'. The 205 studies, which met the search criteria, were classified into 3 categories: genetic association, gene expression and gene methylation studies.
Genome-Wide Methylation Analysis
Cortex DNA samples of 6 male suicide completers and 6 sudden-death comparison subjects (online suppl. table 1) were analyzed with Infinium Human Methylation 450 arrays which quantify methylation of >485,000 CpG sites covering essentially all genes. To avoid batch effects, all samples were hybridized to 1 BeadChip and processed together. After prefiltering and data preprocessing, 463,374 probes were evaluated for each of the 12 samples, yielding a data set with 5,560,488 measurements. Firstly, we performed an unsupervised cluster analysis of site-specific and regional methylation β values to detect major effects in the entire data set (online suppl. fig. 3). The 2 identified clusters did not separate suicide completers and controls (Rand index of 0.04). One cluster consisted of 5 suicide and 4 control samples and the second smaller cluster of 1 suicide and 2 control samples. The between-group variability of site-specific and regional β values ranged from 0.02 to 0.04, which is low. Similarly, in a principal component analysis, we found 2 clusters on the first axis, which were not related to the 2 studied groups (online suppl. fig. 4).
Then, we compared methylation of single CpG sites between groups. The mean methylation of all analyzed sites was 50.69% in suicide completers and 50.94% in sudden-death comparison subjects, indicating a minor (0.25 percentage point) but statistically significant (Mann-Whitney U test; p = 0.001) between-group difference in global DNA methylation. Using p values and combined ranks, none of the differentially methylated CpG sites reached genome-wide significance. In addition, we looked for region-specific methylation differences in promoter regions, genes, CpG islands, and tiling regions (5,000-bp windows). Similar to the single CpG site analysis, none of the identified differentially methylated regions (DMRs) reached genome-wide significance. Online supplementary table 2 shows the top 25 promoter DMRs. Considering the small group size with 6 samples each and the generally small effect sizes of epigenetic variants [El Hajj et al., 2013], this result is not unexpected. Indeed, post-hoc power analysis (effect size 20%; α = 0.05; 2-sided) revealed a power of only ∼10%.
The identified DMRs provide a list of potential candidates for epigenetic programming of suicidal behavior. Functional annotation enrichment analysis (with the DAVID database) was used to test whether the identified DMRs belong to predefined categories more than by chance alone. Because conceptually related genome-wide expression studies [Sequeira et al., 2012; Zhurov et al., 2012] found several percent of the analyzed genes to be differentially regulated in prefrontal cortex (BA10) of suicide completers, we estimated that 1,000 of the ∼21,000 genes assessed by the HumanMethylation450 array may be epigenetically regulated. Therefore, we generated a list with the top 1,000 genes and the top 1,000 promoter regions showing (non-significant) between-group methylation differences in our study. There was only a small overlap (n = 92) between these 2 lists. We compared our gene and promoter lists with reference lists of genes for (1) disease associations, (2) functional categories, (3) biological processes (Gene Ontology terms, GO), (4) molecular functions (GO), (5) pathways, and (6) enrichment in tissues. The top 1,000 gene DMRs were enriched in some categories; however, after multiple-testing correction none of these enrichments was significant. The top 1,000 promoter DMRs were significantly (∼1.3-fold) enriched in 7 functional categories (3 of them related to transcription) and 2 tissue expressions (brain and epithelium; table 1). Overall, we assessed 25,033 promoters which were assigned to an official gene symbol and represented by at least 2 CpG sites on the array. Of these, 14,716 showed higher methylation in suicide cortices and 10,317 in controls. However, when only considering the top 1,000 promoter DMRs, 967 (97%) were higher methylated (on average by 0.6 percentage point) in suicide completers and only 33 (3%) in controls. Considering that mean methylation was even slightly (0.25 percentage point) lower in suicide brains, this generally increased methylation in potentially dysregulated genes (top 1,000 DMRs) is all the more suspicious.
Because CpG island methylation in promoter regions plays a most prominent role in the regulation of transcription [Deaton and Bird, 2011], we explored the functional category ‘transcription regulation' in more detail. A network analysis [Franceschini et al., 2013] was performed with the 127 promoter DMRs enriched in this category. The main network was formed by 28 promoters, including 7 genes that have been associated with depression (fig. 1). One of these 7 genes, NOTCH3, and the NOTCH signaling pathway, are associated with suicide behavior [Monsalve et al., 2014] and 3 genes, CRY1, PER3 and RORB, with circadian rhythm, which is frequently dysregulated in mood and anxiety disorders [Partonen, 2012; Robillard et al., 2013]. With 1 notable exception, all promoters in this network showed higher methylation in suicide brains.
Based on 112 genetic association studies in PubMed (March 2015), we identified 145 candidate genes (online suppl. table 3) which have been associated with the development of suicide [for reviews, see Turecki et al., 2012; Coon et al., 2013; Mandelli and Seretti, 2013; Sadkowski et al., 2013]. Fifteen (10%) of these 145 genes have been replicated in independent studies (table 2); BDNF, HTR1A, HTR1B, SLC6A4 (5-HTTLPR), SSAT (SAT1), TPH1, and TPH2 even showed association in multiple studies. Eleven (8%) of the 145 suicide-associated genes, BDNF, CDH13, DISC1, DRD2, EVC, GDNF, HTR1A, NXPH1, PAH, SKA2, and SLC6A4, belong to the group of 1,000 promoter DMRs, which generally showed a somewhat higher methylation in suicide brains. However, this overlap is not greater than expected by chance.
Similarly, we identified 18 PubMed studies on brain DNA methylation with 594 candidate genes (online suppl. table 3) showing methylation differences between suicide and control brains. Comparison of these 594 candidate genes with our top 1,000 promoter DMRs yielded 31 matches, which is not a significant overlap. Eighty-eight (15%) of the 594 genes have been replicated in independent methylation assays (table 2). Five of the 88 genes, BDNF, MSMP, NUAK1, PHACTR3, and SYNE2, belong to the top 1,000 promoter DMRs.
Based on 75 expression studies in PubMed, we defined a list of 906 suicide candidate genes (online suppl. table 3) exhibiting expression differences between suicide completers and controls, predominantly in prefrontal cortex and/or other brain regions. When we compared our top 1,000 promoter DMRs with the list of 906 differentially expressed genes, we found 46 matches, which again is no significant overlap. Sixty-five (7%) of the 906 genes have been replicated in independent expression assays (table 2). BDNF, CTNNB1, GLUL, HTR2A, LEPR, NR3C1 (GR), NTRK2, SLC6A4, SSAT, TPH2, and UBC were differentially expressed in 3 or more studies. Four of the 65 genes, APLP2, BDNF, SLC6A4, and SYN2, belong to the top 1,000 promoter DMRs.
Since most expression studies did not report the direction of change (up- or downregulation) in suicide brains, for only 10 DMRs the methylation and expression changes can be correlated. For 6 genes, ATXN2L, BDNF, LMNB2, NOTCH3, SYMPK, and TRRAP, increased promoter methylation was associated with reduced expression. Four genes, CDH22, PTPN1, SYN2, and SLC25A23, showed increased methylation and expression. Although in general there is an inverse correlation between methylation and expression, it is well known that some genes can be activated by methylation, i.e. of a repressor element in the promoter region [Farcas et al., 2009].
Table 2 presents the 153 suicide candidate genes that have been replicated in independent genetic association, methylation and/or expression studies. It is noteworthy that only 2 genes, BDNF and SSAT, were previously linked with suicide in all 3 categories. Nine of the 153 top candidates, APLP2, BDNF, HTR1A, MSMP, NUAK1, PHACTR3, SLC6A4, SYN2, and SYNE2, belong to the group of top 1,000 DMRs.
Human behavior, including suicide, results from complex interactions of genes and environment. Although a large number of genetic variants, most importantly in BDNF [Dwivedi, 2010; Pregelj et al., 2011] and serotonergic genes [Antypa et al., 2013; Sadkowski et al., 2013], have been associated with suicidal behavior, major contributions of single genes have not been found [Perlis et al., 2010; Schild et al., 2013]. This is consistent with polygenic/multifactorial models with small effect sizes of individual genes and gene-environment interactions influencing the vulnerability to suicide. Environmentally induced epigenetic changes are the most likely mechanism mediating such complex interactions. Accumulating adverse environmental factors such as childhood maltreatment, acute or chronic stress, social stress, loss/separation and inappropriate family support are thought to increase the risk for suicidal behavior [Paykel, 1976]. It has been shown in both animal models and humans that such social stressors can lead to persistent epigenetic changes in the brain [Weaver et al., 2004; McGowan et al., 2009; Lutz and Turecki, 2014; Zannas and West, 2014].
Previous genetic association studies as well as expression and methylation studies in different brain regions (online suppl. table 3) suggested that regulation of multiple genes and pathways may be altered in the brain of suicide completers; however, the overlap between different genome-wide data sets was generally small. This may be partially due to differences in the brain regions studied, the patient and/or control cohorts analyzed, and the technical platforms used. On the other hand, existing suicide candidate gene lists are likely polluted with false positives and false negatives. Therefore, it is not unexpected that the top 1,000 DMRs in our study did not significantly overlap with candidate genes from previous studies. Only genes that have been replicated in independent genetic association, expression and/or methylation studies should be considered further as primary candidate genes (table 2).
Here, we compared for the first time genome-wide cortex (BA10) methylation patterns of suicide completers versus sudden-death comparison subjects. The small size of the analyzed groups (n = 6) explains why we did not find individual genes showing significant methylation differences after genome-wide multiple testing correction. However, using functional annotation (DAVID) analysis of the top 1,000 promoter DMRs (of >25,000 analyzed promoters), we observed several significantly enriched associations. It is interesting to note that the vast majority (97%) of the top 1,000 promoter DMRs showed somewhat (on average 0.6 percentage point) higher methylation in suicide brains than in controls. This is conspicuous, compared to only 46% of all analyzed CpG sites (n = 465,374) and 59% of all analyzed promoters (n = 25,033) showing higher methylation in suicide completers. Three of 7 significant terms (with false discovery rate-adjusted p values <0.05) in the functional category (alternative splicing, transcription and transcription regulation) were directly related to transcription. Of the identified network of 28 transcriptional regulators, 7 genes have been associated with depression, including EZH2, a transcriptional repressor of the polycomb group family, NOTCH3, a key regulator of neural development, and the circadian rhythm genes CRY1, PER3 and RORB. Significant enrichment of the promoter DMRs was also observed in the category tissue expression (brain and epithelium). Collectively, these results suggest that transcription may be dysregulated by increased promoter methylation in frontal cortex of suicide completers, maybe particularly in cells derived from the neuroepithelium.
It is reassuring that the identified promoter DMRs, most prominently the neurotrophic factor BDNF, the serotonin transporter SLC6A4, the serotonin receptor HTR1A, and the synapsin SYN2 have been implicated with suicidal behavior and neuropsychiatric disorders in the literature (online suppl. table 3). Consistent with an earlier study demonstrating increased BDNF methylation in the Wernicke area [Keller et al., 2010], we found a somewhat higher methylation of the BDNF promoter in frontal cortex of suicide completers (6.9%), compared to controls (6.3%). BDNF mRNA and protein expression appears to be reduced in cortex, hippocampus and blood of suicide subjects regardless of psychiatric co-morbidities including depression, post-traumatic stress disorders and schizophrenia [Dwivedi, 2010; Sher, 2011]. Interestingly, BDNF methylation was also increased in blood of patients with suicidal behavior [Kang et al., 2013], suggesting that it may be exploited as an epigenetic biomarker for this trait. Gathering evidence suggests that despite enormous between-tissue differences of DNA methylation patterns, the blood epigenome can reflect at least to some extent the methylation variation in the brain [Gladkevich et al., 2004; Davies et al., 2012].
Epigenetic variation is regulated in a tissue-specific manner by stochastic, genetic and environmental factors [Kerkel et al., 2008; Grundberg et al., 2013]. Because epigenomes largely differ between tissues, it is crucial to study the appropriate target tissue in the first place. In addition, it is important to consider that cortex tissue is composed of different cell types. It seems plausible to assume that suicidal behavior is associated with changes in the neuronal epigenome(s). Future methylation analyses should be performed on isolated neuronal versus glial cell populations. Although it is desirable to increase sample size in the next few years, it will not be possible to reach genetic study standards using human brain material. Bisulfite conversion as used in this study does not discriminate between DNA methylation and hydroxymethylation. Distinct regulatory epigenetic processes, including histone modifications, can be expected to contribute to cortex gene regulation in complex behavioral phenotypes.
Despite the small scale of this study, our results shed new light on the etiology and pathophysiology of suicidal behavior, which has been associated with small but numerous changes in cortex DNA methylation. At the individual gene level, it is hardly possible to appraise the phenotypic consequences of methylation changes in the order of 0.5-1 percentage point. However, the accumulation of small changes acting in the same direction in multiple genes interacting with each other in a pathway(s)/network(s) may well play a role in the etiopathogenesis of complex disease. Similar to genome-wide association studies with genetic markers, even minor methylation differences between suicide completers and controls may uncover genes and pathways underlying suicidal behavior. We propose that an accumulation of unfavorable epigenetic events, which are programmed by adverse life experiences, works in concordance with other factors to increase the risk of suicidal behavior.
The authors thank Prof. Andreas Reif for his helpful comments on the paper.