Schizophrenia (SZ), which affects 20 million people worldwide , is one of the most common types of serious mental illness. To date, attempts to discover genes that are directly associated with SZ have been largely unsuccessful [2,3]. A recent genome-wide association study (GWAS) of SZ revealed 108 loci associated with the disease . However, the proportion of the variance in liability explained by single nucleotide polymorphisms (SNPs) is small; significant genome-wide loci explain only 3.4% of the variance in liability, and the cumulative effect of common loci expressed as a polygenic risk score (PRS) was estimated to explain only 7%.
Therefore, both genetic and environmental factors must be considered in research on the etiology of SZ. Epidemiological studies have suggested that multiple environmental factors, including prenatal infection/immune activation, paternal age, malnutrition, hypoxia-related obstetric complications, childhood/adolescent trauma and cannabis abuse are associated with an increased risk of SZ . The interaction of risk genes with the environment can result in aberrant epigenetic mechanisms in SZ. After first being proposed by Petronis (2004) , the epigenetic component of the pathophysiology of SZ has been widely recognized in recent years. DNA methylation, one of the main epigenetic mechanisms, can be examined by utilizing two complementary approaches. The first involves hypothesis-driven analysis of DNA methylation in candidate genes. Using this approach, abnormal methylation of some genes, such as DRD2 , DLGAP2 , and COMT , has been implicated in the onset and progression of SZ. The second approach involves microarray-based epigenomic profiling. With this microarray-based technique, assessment of DNA methylation in about 450,000 or 853,307 CpG sites can be performed. However, this represents only a small fraction of all CpG sites (∼28 million) . Several methylome-wide association studies (MWAS) have been performed using postmortem brain tissues [11-14] and the peripheral blood cells of patients with SZ [15-20]. These MWAS have revealed differences genome-wide methylation patterns between SZ patients and controls, but the number of CpGs analyzed, and the sample sizes, have been limited. Two large-scale methylome studies linked SZ to hypoxia and infection , and some differentially methylated positions were located in a top-ranked SZ region based on GWAS analyses . Few MWAS have been performed in first-episode SZ [19,21-23]. Moreover, with a few exceptions [18,24,25], most previous studies did not control for potential confounders of methylation, such as lifestyle  and dietary [26-28] factors.
The present study aimed to identify methylation changes in CpG sites associated with recent-onset (RO) psychosis using the Infinium Human Methylation BeadChip (Illumina, San Diego, CA, USA), which covers more than 850 K CpG sites. In the methylation analyses, lifestyle factors such as smoking, alcohol, exercise, and diet were controlled for. Furthermore, we evaluated the association of altered methylation patterns with symptoms and early trauma in patients and healthy controls.
As inclusion criteria, it was required that subjects be aged between 19 and 58 years and meet the Diagnostic and Statistical Manual of Mental Disorders, Fifth Edition  criteria for SZ spectrum disorders (SZ, schizophreniform disorder, psychotic disorder not otherwise specified [NOS]), brief psychotic disorder, or delusional disorder. RO was defined as duration of illness ≤ 2 years. The exclusion criteria were an intelligence quotient < 70, history of head trauma, serious neurological disorder (epilepsy, stroke, Parkinson’s disease, and/or dementia), and significant medical illness. Participants were recruited from the Korean Early Psychosis Cohort Study (KEPS), a nationwide, multicenter, prospective and naturalistic observational study. Age- and sex-matched healthy individuals were recruited for the control group via advertisements. All participants provided written informed consent in accordance with a protocol approved by the Ethics Committee of Jeonbuk National University Hospital (approval no. CUH 2014-11-002).
The severity of symptoms was evaluated within a week of blood sampling using the Positive and Negative Syndrome Scale (PANSS) [30,31]. Childhood trauma was evaluated with the Early Trauma Inventory Self Report-Short Form (ETISR-SF) . Data on factors related to lifestyle were obtained using the Fagerstrom Test for Nicotine Dependence (FTND) , Alcohol Use Disorders Identification Test (AUDIT) , Dietary Habits Questionnaire (DHQ) , and Physical Activity Rating (PA-R) . For simplicity, we used only the fourth item of the FTND (0, non-smoker; 1, ≤ 10 cigarettes/day; 2, 11−20 cigarettes/day; 3, 21−30 cigarettes/day; and 4, ≥ 31 cigarettes/day) and the mean of items 1 and 2 of the AUDIT ([sum of the scores on items 1 and 2] / 2). The DHQ is a 20-item self-administered questionnaire consisting of three subcategories: there are 5 items on diet regularity, 6 on balanced diet, and 9 on unhealthy diet and eating habits. This scale was developed based on dietary guidance published by the Korean Ministry for Health, Welfare and Family Affairs (2010) . Each item is scored on a three-point scale (1, 3, and 5 points) according to the frequency of the dietary habit. Higher scores indicate better dietary habits in the respective categories. The PA-R is a questionnaire for rating physical activity, with scores ranging from 0 (avoids walking or exercise) to 7 (runs more than 10 miles per week or spends more than 3 hours per week in comparable physical activity).
DNA methylation refers to a biochemical process in which a methyl group is bound to the cytosine of a CpG dinucleotide. In turn, this alters the levels of gene transcription, without altering the underlying DNA sequence. Genomic DNA was extracted from peripheral blood mononuclear cell samples using a QIAamp DNA Mini Kit (QIAGEN, Hilden, Germany) according to the manufacturer’s instructions. All prepared samples were bisulfite converted according to the EZ DNA Methylation Kit protocols (Zymo Research, Orange, CA, USA). Bisulfite-converted DNA was assessed using the Infinium MethylationEPIC BeadChip (Illumina), to determine the methylation levels of 853,307 CpG sites. After staining, the BeadChips were scanned, and images were extracted using GenomeStudio software (Illumina). The raw intensity data files (idat files) were preprocessed. First, we excluded 486 probes with detection levels of a p value ≥ 0.01 in 25% of the sample, because they were considered unreliable. Then, we removed 19,377 CpG probes located on the sex chromosome. During this process, 59 probes identified as SNPs were removed. The data were background-corrected and adjusted for dye bias using the methylumi and lumi R packages. To reduce bias in the Infinium I and II assays, the corrected signal value was normalized. The beta-mixture quantile method was used to decompose the density profiles of the Infinium I and II assays into three methylation states: unmethylated (close to 0), partially methylated (close to 0.5), and fully methylated (close to 1). Quantile normalization was then performed to fit each β-distribution of the Infinium II profile to the corresponding β-distribution of the Infinium I profile, leaving 846,055 CpGs for subsequent analysis (Supplementary Fig. 1; available online). Raw data were extracted as beta values ranging between 0 and 1 (β value = methylated intensity / [methylated intensity + unmethylated intensity]). All data analyses, and the visualization of differentially methylated CpGs, were conducted using R software (version 3.3.3; www.r-project.org).
To identify 2,912 differentially methylated CpG sites associated with RO psychosis, we performed linear regression adjusted for confounding factors: β = α0 + α1 age + α2 sex + α3 education + α4 alcohol + α5 smoking + α6 exercise + α7 diet + ε. To adjust for multiple testing, a false discovery rate (FDR) cutoff p value of < 0.001 was used. Correlations between the methylation values of 2,912 selected CpG sites and total and subscale scores on the PANSS and ETISR-SF were calculated. p values < 0.05 were considered significant.
To identify functions associated with the methylation-altered genes (2,535 genes corresponding to 2,912 CpG sites), we performed Gene Ontology (GO) analysis and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis using the Database for Annotation, Visualization, and Integrated Discovery (DAVID; http://david.abcc.ncifcrf.gov/). We performed KEGG enrichment analysis to identify the biological processes, pathways, and networks shared by these genes. GO terms and KEGG pathways that contained 10 genes and had an Expression Analysis Systematic Explorer score < 0.05 were considered significant. The top 15 GO terms and KEGG pathways were selected.
The study was powered to detect mean differences in DNA methylation between cases and controls at a p value < 3.44E-06, following a previously described method . The power thresholds (ranging between 0 and 1) of the probes, based on the numbers of cases and controls, are shown in Supplementary Fig. 2 (available online). Our study was adequately powered to detect a 0.05% mean difference at ∼40.6% of all sites.
Demographic and clinical characteristics of the participants (61 RO psychosis patients and 47 controls) are described in Table 1. Diagnoses included SZ (29.51%), schizophreniform disorder (31.15%), other specified SZ spectrum and psychotic disorders (22.95%), brief psychotic disorder (9.84%), and delusional disorder (6.56%). We identified 2,912 differentially methylated CpG sites associated with RO psychosis (FDR < 0.001) after controlling for age, sex, education level, alcohol use, smoking status, exercise, and diet. Of the 2,912 CpGs, 1,509 sites (51.82%) were hypomethylated, and 1,403 (48.18%) were hypermethylated. Of the 1,509 hypomethylated CpGs, 951 (63.02%), 308 (20.41%), and 11 (0.73%) sites were located in the promoter regions, gene bodies, and 3'-UTRs, respectively (Fig. 1A). When classified according to the CpG contents in the genes, 970 (64.28%), 193 (12.79%), and 49 sites (3.25%) were located in the CpG islands (CGIs), CGI shores, and CGI shelves, respectively (Fig. 1B). Of the 1,403 hypermethylated CpGs, 505 (35.99%), 644 (45.90%), and 51 (3.64%) sites were located in the promoter regions, gene bodies, and 3'-UTRs, respectively (Fig. 1C). In addition, 526 (37.49%), 375 (26.73%), and 148 sites (10.55%) were located in the CGIs, CGI shores, and CGI shelves, respectively (Fig. 1D). The top 20 differentially methylated sites with the lowest FDR values are presented in Table 2. Corresponding gene information of top 20 sites are provided in the online Supplementary Table 1. Ten sites from the ACIN1, PDCD10, SERPINI1, GNAS, HIST1H4J, ZNF668, SSRP1, ANKRD13B, HIST1H4K, CDKN2A, CDKN2B-AS1, and FKBP2 genes, and ten from the ETFA, TYRO3P, RNF31, PSME2, NUMBL, MARCH1, ANP32C, Autophagy-related 4B cysteine peptidase (ATG4B), IGF2BP2, and METTL16 genes, showed that patients were more likely to have a hypo- or hypermethylated status than controls. The full names of these genes, and the top 100 differentially methylated sites, are provided in the online Supplementary Table 2.
The top 15 significantly enriched KEGG pathways included the oxytocin signaling pathway, long-term depression pathway, axon guidance, endometrial cancer, long-term potentiation, MAPK signaling pathway, and glutamatergic pathway, among others (Table 3). In total, 341 genes were annotated with 15 KEGG pathways. The top 15 significantly enriched GO terms are provided in the online Supplementary Table 3.
Significant correlations are shown in Figure 2. In the patient group, positive correlations were found between the beta value of cg13562874 and the ETISR-SF total score, while negative correlations were identified between the beta values of cg13810931 and cg18128437 and the PANSS total score. In the control group, the ETISR-SF total score exhibited negative correlations with the beta values of cg23207361 and cg23408615, and a positive correlation with the beta value of cg25548986.
Epigenetics provides insight into the molecular link between the genetic and environmental factors that contribute to the complex etiology of diseases. The goal of the present MWAS was to identify methylation biomarkers in the peripheral blood of patients with RO psychosis. We found associations of several novel methylated genes with psychosis and early trauma. The results provide new insight into the epigenetic causes and treatment of psychotic disorders.
Analysis of the distribution of differentially methylated CpG sites within genes demonstrated that hypomethylation was more common in the promoter region, and hypermethylation in the gene body. Although the role of DNA methylation in gene bodies remains unclear, gene body DNA methylation in highly expressed genes is often seen in human cells . CGIs in gene bodies reportedly act as alternative promoters [40,41], and tissue- or cell type-specific CGI methylation is prevalent in gene bodies [41,42]. Similarly, in this study, hypomethylation was more common in CGIs, while hypermethylation was more abundant outside CGIs. This indicates that psychosis-associated CpG sites often do not reside in CGIs. Accumulating evidence suggests that phenotypically relevant variation in DNA methylation often occurs outside CGIs .
Several of the 22 genes associated with the top 20 differentially methylated CpG sites have been implicated in the pathogenesis of SZ. Two genes, apoptotic chromatin condensation inducer 1 (ACIN1) and programmed cell death 10 (PDCD10), are involved in apoptosis. Our finding that hypomethylation of cg23712970 in ACIN1 and PDCD10 may result in increased expression of these genes is partially consistent with the recent finding that expression levels of the BCL2 and CASP3 apoptosis genes were significantly higher in the peripheral blood lymphocytes of schizophrenic patients . Analysis of mRNA expression in the post-mortem brain tissue of patients with SZ also revealed altered mRNA expression in the IFITM2, TIAL1, PDCD6, and MCL genes, which are involved in programmed cell death . ATG4B is another gene involved in autophagy. Furthermore, the hypomethylation of cg15028047 seen in the present study is in line with the decreased expression of autophagy-related genes in Brodmann area 22 in previously reported in SZ patients [46,47]. Importantly, genome-wide profiling of DNA methylation in induced pluripotent stem cells demonstrated that 23 differentially methylated genes between SZ patients and controls were functionally annotated as ‘Autophagy’ . It should be noted that, as antipsychotic drugs are known to induce autophagy , our finding could be related either to disease processes or uncontrolled environmental factors. Hypomethylation of the serpin family I member 1 (SERPINI1), ring finger protein 31 (RNF31), and FKBP prolyl isomerase 2 (FKBP2) genes, which appear to enhance axonal growth, prevent neuroinflammation, and accelerate protein folding and trafficking, respectively, seems counterintuitive given that increased expression of these genes could be beneficial. A possible explanation could be that these methylation changes may reflect either compensatory mechanisms in response to disease processes or the effects of antipsychotic drugs. Hypermethylation of methyltransferase-like 16 (METTL16) is notable, in that a methionine-induced animal model of SZ , and the methylation hypothesis of SZ, provide solid support for a role of epigenetics in understanding the pathophysiology of SZ . Histones play a central role in transcription regulation, DNA repair, DNA replication, and chromosomal stability. Our findings regarding the hypomethylation of H4 clustered histone 11 (HIST1H4J) and H4 clustered histone 12 (HIST1H4K) are novel; most previous studies of SZ have considered histone methylation at H3K9 [52,53] and global histone H4 acetylation . Hypermethylation of the membrane-associated ring-CH-type finger 1 (MARCH1) gene, which is involved in regulating vesicular transport and endocytosis, merits further investigation given the novel hypothesis that altered clathrin-mediated endocytosis (CME) is responsible for the pathophysiology of SZ . When we compared our top 100 CpG sites (with 120 genes) to the 172 differentially methylated positions reported by a large-scale methylome study of SZ , only four genes, caspase recruitment domain family member 14 (CARD14), nuclear factor of activated T cell 1 (NFATC1), phosphofurin acidic cluster sorting protein 2 (PACS2), and SKI proto-oncogene (SKI), overlapped. Moreover, comparison with the 100 top-ranked SZ-associated genes in the study by Dempster et al.  revealed no overlap. This lack of overlap may be attributable to differences in disease stage and confounding factors, such as the study subjects, total number of CpG sites examined, and sample size. Overall, it should be noted that fold- change and △β values for the top 20 CpG sites were low. These low values were also reported in previous studies [18,24]. This suggests that marked differences in methylation status between patients and controls may be difficult to detect in peripheral blood samples, or by microarrays.
Some of the top 15 KEGG pathways have important implications for SZ. Oxytocin plays a central role in birthing, mother-infant bonding, and a wide range of related social behaviors, and is considered a candidate susceptibility gene for SZ. Altered oxytocin receptor (OXTR) methylation and OXTR polymorphisms have been reported in SZ, and in individuals at ultra-high risk for psychosis [22,56]. As for the axon guidance pathway, a recent analysis of European and American GWAS of SZ revealed that 15 axon guidance pathway related genes were associated with the disease . The mitogen-activated protein kinase (MAPK)-associated pathway activates transcription factors related to learning, memory, cell proliferation, and apoptosis. Hypermethylated genes related to MAPK signaling pathways were reported in SZ-affected monozygotic twins . In the endocytosis pathway, 38 genes were enriched, including MARCH1. Interestingly, DNA methylation-mediated modulation of endocytosis has been proposed as a mechanism underlying synaptic function .
Higher trauma scores in our patient group were associated with hypermethylation of the zinc finger SWIM- type containing 8 (ZSWIM8) gene. The C2H2-type zinc finger proteins (C2H2-ZNFs) are known to play significant roles in the regulation of neural stem cell activity and subsequent brain development [60,61]. Accumulating evidence suggests a role for C2H2-ZNFs in the development of SZ [62-64]. Therefore, this finding suggests that earlier trauma may decrease the expression of ZNFs, in turn leading to aberrant brain development. However, in the control group, higher trauma scores were associated with hypomethylation of the branched-chain amino acid transaminase (BAT2) gene and the intergenic region, as well as hypermethylation of the one cut homeobox 1 (ONECUT1) gene. The BAT2 and ONECUT1 genes are involved in glucose metabolism and diabetes mellitus. It is difficult to determine the implications of these findings given that the controls were all healthy. Another interesting finding was that, in the patient group, higher PANSS scores were correlated with hypomethylation of the phosphatidylinositol transfer protein membrane associated 1 (PITPNM1), PDCD10, and SERPINI1 genes. PITPNM1 plays a role in maintaining normal diacylglycerol levels in the Golgi apparatus, while the PDCD10 and SERPINI1 genes are involved in apoptosis and axonal growth. Hence, these results suggest that more severe psychopathology may be associated with more active changes in gene expression at both the cellular and axonal levels.
Some limitations should be considered when interpreting the findings of this study. First, the small sample size and diagnostic heterogeneity among the participants may limit the generalizability of the findings. However, it should be noted that the present study was adequately powered to detect a 0.05% mean difference at ∼40.6% of all sites. Second, as some of the patients were on antipsychotics (19 were antipsychotic-naïve, 19 were antipsychotic-free, and 23 were medicated), potential effects of medication on the epigenome cannot be ruled out. Considering several reports on the effects of antipsychotics on methylation profiling in animals  and patients with SZ , this factor should be controlled for in future studies. Third, significant CpG sites were not validated with independent samples. Fourth, as we did not have RNA samples, the effects of methylation changes on gene expression were not explored. Fifth, cellular heterogeneity of blood samples was not controlled for. Given that cellular composition explains much of the observed variability in DNA methylation , the cellular heterogeneity of our samples may have affected the results. Sixth, the SNPs were not all removed. As sequence polymorphisms may affect the readout of methylation signals in Infinium arrays , this should be considered in future studies.
As strengths of the present study, the patients were all in a relatively early disease stage, and various lifestyle factors were controlled for. Most importantly, this was the first study to investigate the association between early trauma and DNA methylation changes in RO psychosis. In summary, using a methylome-wide array, we identified differentially methylated CpG sites in patients with RO psychosis compared to controls. The genes associated with the top 20 differentially methylated sites were largely novel, had not been reported in previous methylation studies, and were involved in apoptosis, autophagy, axonal growth, neuroinflammation, protein folding, etc. In the patient group, significant associations of these novel methylated genes with early trauma and psychopathology were observed. These findings support the notion that altered DNA methylation could be involved in early trauma and the pathophysiology of psychosis.
The data set analyzed in this study is available in the Gene Expression Omnibus (GEO) repository (GSE 157252, https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE157252).
The biospecimens and data used in this study were provided by the Biobank of Jeonbuk National University Hospital, a member of the Korea Biobank Network supported by the Ministry of Health, Welfare and Family Affairs. All samples derived from the Korea Biobank Network were obtained with informed consent under institutional review board-approved protocols.
Korean Mental Health Technology R&D Project, Ministry of Health and Welfare, Republic of Korea (HL19C0015) Korea Health Technology R&D Project through the Korea Health Industry Development Institute (KHIDI), funded by the Ministry of Health and Welfare, Republic of Korea (HI18C2383).
No potential conflict of interest relevant to this article was reported.
Conceptualization: Young-Chul Chung. Data acquisition: Yan-Hong Piao, Yin Cui. Formal analysis: Yan-Hong Piao, Yin Cui, Fatima Zahra Rami, Ling Li, Maryam Karamikheirabad, Shi Hyun Kang, Sung-Wan Kim, Jung Jin Kim, Bong Ju Lee. Funding: Young-Chul Chung. Supervision: Young-Chul Chung. Writing—original draft: Young-Chul Chung. Writing—review & editing: Young-Chul Chung.