DNA methylation in adolescents with anxiety disorder: a longitudinal study

Anxiety disorders (AD) typically manifest in children and adolescents and might persist into adulthood. However, there are still few data concerning epigenetic mechanisms associated with onset, persistence or remission of AD over time. We investigated a cohort of adolescents and young adults at baseline (age; 13.19 ± 2.38) and after 5 years and classified them according to the AD diagnosis and their longitudinal trajectories into 4 groups: (1) Typically Developing Comparisons (TDC; control group, n = 14); (2) Incident (AD in the second evaluation only, n = 11); (3) Persistent (AD in both evaluations, n = 14) and (4) Remittent (AD in the first evaluation only, n = 8). DNA methylation was evaluated with the Infinium HumanMethylation450 BeadChip from saliva samples collected at both evaluations. Gene set enrichment analysis was applied to consider biological pathways. We found decreased DNA methylation in TDC group while the chronic cases of AD presented hypermethylation in central nervous system development pathways. Moreover, we showed that this persistent group also presented hypermethylation while the other three groups were associated with hypomethylation in nervous system development pathway. Incidence and remission groups were associated with increased and decreased methylation in neuron development pathways, respectively. Larger studies are likely to detect specific genes relevant to AD.

Anxiety disorders are a group of syndromes characterized by maladaptive responses to threats, that develop as a result of the interplay between biological factors, psychological mechanisms and environmental influences [1][2][3][4] . These disorders usually have their onset during childhood and adolescence and frequently persist into adulthood [5][6][7] . Understanding the genetically mediated factors that lead to persistence (or remission) of anxiety over time might reveal new ways of preventing and treating those disorders.
Epigenetic mechanisms, mainly DNA methylation, are the result of the interface between genes and environmental factors 4,8,9 . Some studies have suggested that these mechanisms might be able to assure enduring environmental influences on the body, thus mediating the increased susceptibility for psychiatric disorders 10,11 . Studies in rodents described that the negative impact of early stress on behavioral responses as well as on brain changes may influence the regulation of DNA methylation across generations 12,13 . Furthermore, a study in humans found that changes in the offspring's perceptions of maternal bonding were also related to DNA hypermethylation in cell signaling process over time 14 . Therefore, these findings are consistent with the idea that epigenetic mechanisms may be associated to the behavioral and emotional response to environmental factors in long term.
Several candidate genes (e.g: BDNF; FKBP5; SLC6A4; AVP; NR3C1; CRH; COMT; MAOA; OXTR and APOE) to psychiatric disorders have been quoted in epigenetic studies using locus-specific assays 11,[15][16][17][18] . In general, these studies reported that the effects of adverse environment might generate epigenetic modifications that, in turn, can alter physiological processes important to the development of mental disorders. However, no hub genes have consistently been associated to diagnosis or the course of anxiety disorders. Given the current state of knowledge on the genetics of anxiety disorders, unbiased hypothesis-free methodology might be more fruitful for investigating the DNA methylation fingerprints related to the course of these disorders [19][20][21] . To our knowledge there is a paucity of studies involving genome-wide DNA methylation and anxiety in humans [22][23][24][25] and in animals 26 studies.
This is an unbiased hypothesis-free study with no predictions with respect to specific genes aiming to investigate which biological pathways are associated with the course of anxiety disorders based on gene ontology. This analysis allows for the investigation of sets of genes instead of looking for isolated genes 27 . Our aim is to explore biological pathways that may confer an epigenetic signature according to favorable and unfavorable trajectories of anxiety disorders.

Methods
Sample selection and Participants. This is a longitudinal 5-year follow-up study that involves a subsample of adolescents and young adults from a community cohort enriched for participants with anxiety disorders originated in 2008 that were re-evaluated 5 years later, in 2013.
In 2008, 240 non-medicated adolescents from a total of 2.457 students that answered a screening scale for anxiety disorders (Screen for Child and Anxiety Related Emotional Disorders -SCARED) 28 in their schools underwent an extensive psychiatric diagnostic assessment (Schedule for Affective Disorder and Schizophrenia for School-Age Children-Present and Lifetime Version -K-SADS-PL) 29,30 and evaluation in the anxiety outpatient research program at Hospital de Clínicas de Porto Alegre. Exclusion criteria included: (1) a significant organic illness; (2) a history of bipolar disorder, pervasive developmental disorder, or any psychotic disorder; (3) a history of alcohol or drug dependence or abuse; or (4) a clinical suspicion of intellectual disability. Of these selected individuals, six were excluded due to intellectual disability remaining a total of 234 adolescents: 134 with anxiety disorders and 100 not anxious individuals classified as controls to anxiety disorders. Details about initial selection procedures are described elsewhere 31 .
The initial 2008 sample (n = 234) was contacted again in 2013 and of these, 76 adolescents agreed to participate in the 5-year follow-up survey. Adolescents were re-evaluated throughout the K-SADS-PL or the Mini International Neuropsychiatric Interview (MINI) 32 depending on their ages. After the second evaluation, 47 subjects (Table 1) were selected to undergo genome-wide DNA methylation analysis due to logistical strategies and financial limitations. These anxious and non-anxious individuals were categorized into four groups of subjects according to their anxiety disorder trajectories from 2008 to 2013, carefully paired by age and gender: (1) Typically Developing Comparisons (TDC; n = 14; mean age = 17.96; SD = 2.38; 57.14% females). This group is considered control to anxiety disorders because they were not diagnosed with anxiety disorder in both psychiatric evaluations (2008 and 2013), according to the instruments used (described in section 2.2). (2) Incident (n = 11; mean age = 17.27; SD = 1.95; 45.45% females). This group is composed by those who were controls for anxiety disorders in 2008 but were considered cases for any anxiety disorder in 2013, which means that these participants were not previously anxious (2008), but have the onset of anxiety disorder in the second evaluation (2013). (4) Remittent (n = 8; mean age = 19.73; SD = 2.71; 50% females). This group is composed by youths that had anxiety disorder in 2008, but not in 2013. They remitted from the disorder.
This study was approved by the Research Ethics Committee of Hospital de Clínicas de Porto Alegre (HCPA; protocol number 12-0254) and all participants and the caregivers who were participant's legal guardians provided written informed consent in order to participate in the study. All methods were performed in accordance with the relevant guidelines and regulations.
Anxiety Disorder Diagnosis. The psychiatric diagnosis was assessed in 2008 using the Brazilian version of K-SADS-PL 29 and in 2013 with the same instrument or with the MINI 32 depending on the age of the individual (K-SADS for subjects with age lower than 18 years old and MINI for those with 19 years old or higher).
The K-SADS-PL is a semi-structured diagnostic interview that ascertains both lifetime and current psychiatric diagnostic status in children and adolescents based on Diagnostic and Statistical Manual of Mental Disorders, 4 th (DSM-IV) criteria 30 . The Brazilian version of this instrument presents adequate psychometric properties (reliability and validity) and can be used in both clinical practice and research in order to evaluate children and adolescents mental health 29 . Our Inter-rater reliability in this sample was 0.932 (kappa value) for anxiety diagnoses 33 .
The MINI is a structured clinical diagnostic psychiatric interview that evaluates the main diagnoses according DSM-IV criteria in individuals above 18 years old 32 . This instrument also has good psychometric properties with kappa coefficients ranging between 0.65 and 0.85 34 .
Array-Based Genome-Wide DNA Methylation Assays. DNA was extracted from saliva sample, collected in both evaluations, 2008 and 2013, using Oragene™ kits (DNA Genotek, Ottawa, Ontario, Canada). We treated the extracted genomic DNA (500 ng) with sodium bisulphite using the EZ-96 DNA Methylation-Gold™ Kit (Zymo Research, Orange, CA, USA) according to the manufacturer's protocol. Bisulfite-converted DNA was subsequently assessed for DNA methylation status at 485.577 CpG loci with the Infinium HumanMethylation450 (IHM450) BeadChip 35,36 . The IHM450 BeadChip covers 99% of Ref Seq genes regions and involves targeted gene regions with sites in the promoter region, 3′and 5′untranslated regions (3′UTR and 5′UTR), first exon, and gene body in order to explore the genome-wide DNA methylation. This multiple-site approach is extended to CpG islands/CpG island regions for which 96% of islands were covered overall, with multiple sites within islands and island shores, as well as those regions flanking island shores (island shelves) 35,37,38 . Pre-Processing of Raw Data of IHM450 BeadChip. Raw data was evaluated according to the quality control of samples and probes; background correction; normalization; type 1 and 2 probe scaling and batch/plate/ chip/confounder adjustment 36 following the processing workflows described in the R packages lumi 39 methyAnalysis 40 and sva 41 . Areas with single nucleotide polymorphisms and sexual chromosomes were removed. After this cleaning process, we had 382.264 probes of IHM450 BeadChip. All data including DNA methylation values/ subject, and methylated vs. unmethylated probes are deposited in GEO and are accessible via the GEO identifier GSE78975.

Differential Methylation in Signatures.
We used limma R package to evaluate differential methylation signatures comparing the longitudinal and the case-control contrasts 42 . In the longitudinal study, we compared times t1 and t2 to (A) controls, (B) persistent (cases), (C) incident and (D) remittent groups resulting in four contrasts: A t2 -A t1 , B t2 -B t1 , C t2 -C t1 and D t2 -D t1 . In the case-control study, we compared two contrasts, (A) controls and (B) persistent (cases) groups at each time: B t1 -A t1 and B t2 -A t2 . Therefore, each probe has been evaluated by 6 different contrasts. Results from each contrast, for each probe, are summarized in log2-fold change values. The probe with the highest variance across samples was selected to represent a gene when we had multiples probes mapped to the promoter region. A "differential methylation signature" represents the differential methylation of a given contrast mapped to 19.556 unique genes. The differential methylation signatures have been ranked by log2-fold change values for the enrichment analysis described bellow.

Gene Set Enrichment and statistical analysis. Package HTSanalyzeR was used to perform Gene Set
Collection Analysis (GSCA) on the differential methylation signatures 43 . We used the Sub catalog C5 of Molecular Signature Database (MSigBD) composed by Gene Ontology (GO; http://geneontology.org/) that describes gene products in terms of their associated biological pathways in a species-independent manner 27,44 .
Gene set enrichment analysis (GSEA) was applied to GO. GSEA is a method that groups genes that share common biological functions or regulation. It has three keys elements: calculation of an Enrichment Score (ES); estimation of significance level of ES and adjustment for multiple hypothesis testing 45 .
We used the high scoring gene sets that were grouped according to the basis of leading-edge subsets of shared genes. Leading-edge subset is used to extract the core members of high scoring gene sets that contribute to the enrichment signal and thus can reveal a biologically important subset within a gene set 45 . The enrichment analysis was performed with differential methylation signatures mapped to 19.556 genes. The GSEA scores reflected the difference in DNA methylation pattern between two groups (contrasts) and did not aim to identify hub genes. Moreover, it avoids any conclusion at a single gene level. One advantage of this approach is that it does not rely on arbitrary statistical thresholds to assign significance for individual genes, but it uses all genes as a differential methylation signature.
We performed two sets of analysis. First, we compared differences in DNA methylation patterns from GSEA analysis considering time (baseline versus 5-year follow-up) in each of the 4 groups (TDC, incident, persistent and remittent). Second, we performed case-control analysis comparing anxious vs. non-anxious cases across both time points. We analyzed 825 biological pathways. We considered 30 as the minimum gene set size in the gene set enrichment analysis 43 . The significant gene set cutoff p-value (adjusted) was set to <0.001. To estimate the p value we considered 1000 permutations in GSEA and P-values were corrected by BH method 46 . can be considered representative of the whole sample from 2008. There were also no differences regarding age, sex and ethnicity in the 4 groups with the different anxiety trajectories (TDC, incident, remittent and persistent). Furthermore, all analyses were controlled for age and sex. These analyses are available in the supplemental material.

Gene Signatures in longitudinal approach
From the 825 biological pathways analyzed, we observed 50 differentially methylated pathways in TDC group, 13 in the incident group, 27 in the persistent group and 11 in the remittent group (Table 2). Overall we can see a pattern of DNA hypomethylation from several cellular processes (e.g., intracellular signaling cascade, regulation of signal transduction and regulation of metabolic process) in the TDC group that were not seen in the other groups. On the other hand, we found patterns of DNA hypomethylation in other biological pathways in the other three groups. These results are depicted in Fig. 1.
Of note, several biological pathways related to nervous systems were differently methylated according to the different trajectory of anxiety disorders. There are two main findings that should be considered worth noting. First, we could see a pattern of DNA hypomethylation in biological pathways associated with neurogenesis, neuron differentiation and generation of neurons for TDC, incident and remittent groups that were not significant in the persistent group. In addition, the same pattern of DNA hypomethylation could be seen in nervous system development process in TDC, incident and remittent groups, whereas a pattern of DNA hypermethylation was observed in the persistent group. Furthermore, DNA hypermethylation was also observed in biological pathways related to the development of the central nervous system in individuals with chronic anxiety whereas in the TDC group we found DNA hypomethylation in this pathway (Table 2). Second, neuron development pathway was hypermethylated in the incident group and hypomethylated in the remittent group, but no significant results were found for TDC and persistent groups ( Fig. 1 and Table 2). The genes that comprised each biological pathway in all four groups can be observed in supplementary information (Supplementary Appendix 1). Gene Signatures in cross-sectional approach: "Case -Control contrasts". Considering the cross-sectional approach in which the contrast is between cases as compared to control groups, we found 54 biological pathways significant at false discovery rate (FDR) correction in 2008 and 21 biological pathways in 2013 (Table 3). We observed DNA methylation in opposite directions in several biological pathways comparing both evaluations (2008 and 2013), as can be seen in Fig. 2.
When we considered the time of the evaluation, we found that many biological pathways involved with nervous system (e.g, nervous system development, central nervous system development and generation of neurons) were DNA hypomethylated in 2008 and hypermethylated in 2013. It is noteworthy that other biological pathways (e.g, regulation of development process, regulation of gene expression and regulation of metabolic process) showed a pattern of DNA hypomethylation in 2008, which did not reach significance in 2013 ( Fig. 2 and Table 3). The genes that compose each biological pathway in control and case groups can be observed in supplementary information (Supplementary Appendix 1).

Discussion
We were able to show different DNA methylation patterns in many biological pathways depending on the longitudinal trajectory of anxiety disorders. Furthermore, we explore cross-sectional data in two important periods of life. Our results demonstrated that anxiety disorder trajectories in adolescence might be reflected on different DNA methylation patterns evaluated from saliva samples.
Even though some biological pathways showed DNA hypomethylation patterns either in the TDC (controls) or in the persistent (cases) groups (e.g, signal transduction, regulation of cellular metabolic process and regulation of metabolic process), when we considered only nervous system development and central nervous system development processes, we found opposite directions of DNA methylation (hypomethylation in controls as opposed to hypermethylation in cases) within the longitudinal contrast. Moreover, we were able to observe that the neuron development process was associated to DNA hypermethylation in the incident group while this same biological pathway was hypomethylated in the remittent group, but without significance in TDC or persistent groups. Up to now, we have no way to infer the consequence of DNA methylation in transcriptional levels but maybe the severity of anxiety symptoms over time could act silencing genes involved within nervous system. Surprisingly, when we investigated cases and controls to anxiety disorders in cross-sectional approach, we found DNA methylation patterns with opposite directions with the "cases -controls" contrast. Moreover, the beginning of adolescence was characterized by DNA hypomethylation patterns whereas young adults presented DNA hypermethylation patterns. Although the majority of biological pathways were observed at both evaluations, some of them could be observed only in one or in another period. These results suggest that the time factor could be more important than anxiety diagnosis per se in changing epigenetic development biological patterns.
There is no consensus in the literature regarding the increase or decrease of DNA methylation in a longitudinal lifetime course. There are studies in healthy populations that suggested DNA hypermethylation patterns as people get older 47,48 , as well as to DNA hypomethylation over the time 49 . Talens et al. (2012) suggested small or even no differences in mean DNA methylation between the young versus the old age healthy twins, while increases in variation were generally more substantial in older individuals. It is suggested that epigenetic changes might either accumulate with age being generally nondirectional or they are the outcome of many smaller changes with opposite directions 50 Table 2. GSEA of Biological pathways in groups defined by trajectories of anxiety disorders. Note: GSEA, Gene set enrichment analysis; TDC, Typically Developing Controls; NS, No significant (p > 0.001). All biological pathways significant at false discovery rate (FDR) correction. Signal negative represents DNA hypomethylated biological pathways.

Regulation processes
Negative  According to our knowledge this is the first study that addresses the direct association between anxiety disorders and genome-wide DNA methylation, with a longitudinal approach in adolescents. There are very few data evaluating anxiety (trait, state or disorder) and genome-wide DNA methylation in humans [21][22][23][24]51 . Shimada et al. 23 studied genome-wide DNA methylation in patients with Panic Disorder (PD) and healthy control subjects. They depicted forty CpG sites significant associated with PD at 5% FDR correction, but with small different rates of the DNA methylation. A study that compared medication-free PD patients with healthy controls 25 suggested a possible sex-specific methylation process (hypermethylation in female patients) in the HECA gene, but no global changes in DNA methylation. Schartner et al. 52 suggested that patients with PD presented significantly DNA hypomethylation in promotor region of CRHR1 gene when compared with a sample of healthy controls. Unfortunately, we cannot add information on these specific genes or any other gene because our data is limited to genome-wide DNA methylation.
Recent studies that explored the genome-wide DNA methylation in order to understand the etiology of psychiatric disorders did not consider anxiety disorders [53][54][55][56][57][58] . Most studies evaluated the contribution of environmental adversity in DNA methylation differences associated to depression 53 ; post-traumatic stress disorder 16,59-62 , schizophrenia 63 ; borderline personality disorder, attention deficit hyperactivity disorder and bipolar disorders 64,65 . Also, some cross-sectional studies have investigated DNA methylation in specific genes (e,g; BDNF and MB-COMT) in individuals with different psychiatric disorders 17,18,66 and their data suggested higher DNA methylation in individuals with internalizing disorders and decreased DNA methylation in patients diagnosed with schizophrenia and bipolar in comparison to subjects without psychiatric disorders.
The few studies that explored a longitudinal approach involving the genome-wide DNA methylation in early development suggested that adversities in childhood are associated to differences in DNA methylation and risk to mental health 67,68 . In parallel to the influence of early adversities, DNA methylation changes may be observed as a result of therapeutic treatment using the same concept of environmental interference in DNA regulation 69 . An interesting study compared changes in SLC6A4 methylation after Cognitive Behavior Treatment (CBT) in anxious children, considering those that remitted as compared to those that did not remit from their primary anxiety diagnoses. While responders showed an increase percentage of DNA methylation in SLC6A4, non-responders had decreased SLC6A4 DNA methylation percentage 70 . These results demonstrated the complexity of the epigenetic mechanisms underlying anxiety as well as psychiatric disorders in general.
Our study has some limitations, such as the small sample size, which meant that we did not have sufficient power to explore individual genes associated with anxiety. In addition, we did not have information regarding the anti-anxiety medication use over time, or data on childhood trauma. We exclude participants with drug abuse or dependence, but we did not evaluate eventual smoking. The use of drugs or smoking can influence genome-wide DNA methylation [71][72][73] and the presence of trauma or early life stressors is known to increase the risk to psychopathology during lifetime and are able to reprogram the epigenetic mechanisms 74 . Moreover, we did not evaluate the effects of genotype on the DNA methylation differences because we did not consider specific genes. In addition, even though we know that saliva samples contain a heterogeneous mixture of different cell-types: epithelial cells and leukocytes 75 , we were not able to estimate the proportion of epithelial cells in comparison to leukocytes in our sample. Although this issue could be a limitation of our study, the study of Smith et al. 75 verified that comparisons of DNA methylation between saliva and four brain regions (cerebellum, frontal cortex, entorhinal cortex, and superior temporal gyrus) seems to be more similar than comparisons between blood and these same brain regions validating our measure. There are also some limitations in relation to gene-set analyses. Gene-gene correlations and the contributions of some genes to multiple Gene Ontology terms mean that some biological pathways are spuriously represented, when the association with the phenotype of interest could be explained by the overlap between genes within the biological pathways.
However, our study has important strengths that should be acknowledged. We used DNA from buccal cells present in saliva sample in our genome-wide DNA methylation study that were collected in two different evaluations in long term. This approach is considered a gold pattern design to evaluate epigenetics mechanisms. The DNA extracted from saliva sample is more representative than blood DNA due its consistence with typical  Table 3. GSEA of Biological pathways in groups defined by "case-control" contrasts. Note: GSEA, Gene set enrichment analysis; NS, No significant (p > 0.001). All biological pathways significant at false discovery rate (FDR) correction. Signal negative represents DNA hypomethylated biological pathways.
SCIeNtIFIC RepoRtS | (2018) 8:13800 | DOI:10.1038/s41598-018-32090-1 methylation patterns in the brain regions once it is originated from the same ectodermal layer during development of brain tissue 57,75,76 . Although there are tissue-specific epigenetic variation across brain and blood 57,77 , the use of peripheral tissue is the only feasible method concerning biological investigation of the central nervous system 76 . We also used the GSEA method, which is the most appropriate way to generate hypotheses in a large-scale experiment in order to identifying biological pathways, avoiding the challenge of indicating single genes. In this way, we used our rich dataset to study the genome-wide DNA methylation in a more exploratory approach considering the different longitudinal trajectories of anxiety disorders. Our data may help to further focus on specific sets of genes in order to delineate more a priori hypothesis studies 45 .

Conclusions
We used a rich dataset to study the genome-wide DNA methylation in a more exploratory approach considering the presence of either hypo or hyper DNA methylation in the distinct biological pathways as well as the participants'age along the different longitudinal trajectories of anxiety disorders. In the longitudinal analysis we observed that the incident and remittent youth cases of anxiety presented biological pathways associated to DNA hypomethylation patterns. A pattern of DNA hypermethylation within the nervous system development was observed in the persistent group over the five years. On the other hand, individuals presented mainly DNA hypomethylation when they were younger, and DNA hypermethylation patterns with increased age. Further studies are needed to deepen the understanding of the biological pathways involved with anxiety disorders, their longitudinal course and epigenetic changes. Our data may further help looking for more specific genes associated to this complex disorder in the future.