DNA methylation in cocaine use disorder–An epigenome-wide approach in the human prefrontal cortex

Background: Cocaine use disorder (CUD) is characterized by a loss of control over cocaine intake and is associated with structural, functional, and molecular alterations in the human brain. At the molecular level, epigenetic alterations are hypothesized to contribute to the higher-level functional and structural brain changes observed in CUD. Most evidence of cocaine-associated epigenetic changes comes from animal studies while only a few studies have been performed using human tissue.

Methods: We investigated epigenome-wide DNA methylation (DNAm) signatures of CUD in human post-mortem brain tissue of Brodmann area 9 (BA9). A total of N = 42 BA9 brain samples were obtained from N = 21 individuals with CUD and N = 21 individuals without a CUD diagnosis. We performed an epigenome-wide association study (EWAS) and analyzed CUD-associated differentially methylated regions (DMRs). To assess the functional role of CUD-associated differential methylation, we performed Gene Ontology (GO) enrichment analyses and characterized co-methylation networks using a weighted correlation network analysis. We further investigated epigenetic age in CUD using epigenetic clocks for the assessment of biological age.

Results: While no cytosine-phosphate-guanine (CpG) site was associated with CUD at epigenome-wide significance in BA9, we detected a total of 20 CUD-associated DMRs. After annotation of DMRs to genes, we identified Neuropeptide FF Receptor 2 (NPFFR2) and Kalirin RhoGEF Kinase (KALRN) for which a previous role in the behavioral response to cocaine in rodents is known. Three of the four identified CUD-associated co-methylation modules were functionally related to neurotransmission and neuroplasticity. Protein-protein interaction (PPI) networks derived from module hub genes revealed several addiction-related genes as highly connected nodes such as Calcium Voltage-Gated Channel Subunit Alpha1 C (CACNA1C), Nuclear Receptor Subfamily 3 Group C Member 1 (NR3C1), and Jun Proto-Oncogene, AP-1 Transcription Factor Subunit (JUN). In BA9, we observed a trend toward epigenetic age acceleration (EAA) in individuals with CUD remaining stable even after adjustment for covariates.

Conclusion: Results from our study highlight that CUD is associated with epigenome-wide differences in DNAm levels in BA9 particularly related to synaptic signaling and neuroplasticity. This supports findings from previous studies that report on the strong impact of cocaine on neurocircuits in the human prefrontal cortex (PFC). Further studies are needed to follow up on the role of epigenetic alterations in CUD focusing on the integration of epigenetic signatures with transcriptomic and proteomic data.

1. Introduction

Cocaine is one of the most frequently consumed psychostimulants in the world with around 21 million individuals having used cocaine in 2020 (1). During the last decade, a stable increase in cocaine consumption and cocaine-associated hospitalization rates has been observed (1, 2). While a substantial fraction of individuals retain control over cocaine intake, up to 21% of regular cocaine users develop cocaine dependence characterized by repeated cycles of binge intoxication, withdrawal, and compulsive cocaine-seeking (3). Chronic use of cocaine, especially crack cocaine, increases the risk for the development of cardiovascular, cerebrovascular, and respiratory diseases (4). Also, it is associated with substance abuse of other drugs such as alcohol, cannabis, and opioids (5). Thus, chronic cocaine consumption as observed in cocaine use disorder (CUD) enhances the overall disease burden, at the individual level significantly reduces life quality, and thereby contributes to an increased risk for premature death. Investigating the pathophysiology of CUD enables the characterization of the underlying disease processes and at the same time provides a basis for the development of novel therapeutic approaches.

Substance use disorders (SUDs) including CUD are understood as brain disorders with profound alterations at the structural, functional, and molecular levels in the human brain (6–8). Drug-induced effects on epigenetics are hypothesized to contribute to structural and functional changes in the brain via the establishment of altered transcriptional programs [for review see (9–11)]. Thus, profiling of epigenetic alterations contributes to a better understanding of molecular mechanisms underlying the neuroplastic changes associated with CUD. DNA methylation (DNAm) represents the most prominent epigenetic mechanism involved in gene expression regulation and DNAm changes in the brain following acute and chronic cocaine exposure are well described (12). So far, the effects of cocaine on DNAm were most extensively studied in the rodent brain, which has led to the identification of cocaine-associated differentially methylated genes. Interestingly, several of these genes were previously related to the neurobiology of addiction, such as Bdnf (13) and Pp1c (14) as well as the immediate early genes c-Fos (15) and fosB (14). Cocaine was also shown to affect methylation and expression levels of DNA methyltransferases (Dnmt1, Dnmt3a) (15–18) and Tet demethylases (19) thereby directly altering DNAm homeostasis.

In contrast, much less is known about cocaine-associated DNAm changes in humans. Due to an easy and minimally invasive sampling procedure, epigenome-wide association studies (EWASs) of different cocaine use phenotypes have been performed in human blood samples (20, 21). These studies revealed an enrichment of cocaine-associated differential methylation within biological pathways involved in synaptic transmission, neurogenesis, cell division, and chromatin organization. However, as DNAm levels are known to be tissue-specific, the predictive value of differential methylation in blood for methylation signatures of CUD in the brain might be limited. Thus, profiling DNAm levels in brain tissue is of critical importance. So far, two epigenome-wide characterizations in CUD have been performed in human post-mortem brain tissue from the caudate nucleus (CN) (22) and the nucleus accumbens (NAc) (23). By analyzing post-mortem tissue of the CN from N = 25 CUD cases and N = 25 controls using a reduced representation bisulfite sequencing approach, a total of 173 differentially methylated regions (DMRs) were found that were significantly associated with CUD status. The IRX2 gene was among the top findings of CUD-associated differential methylation in the CN. Using a similar approach in the NAc, the key finding was a hypermethylated region within the tyrosine hydroxylase (TH) gene encoding an enzyme involved in dopamine metabolism. Interestingly, the identified DMR in TH was shared between the CN and the NAc. Both studies provide some first mechanistic insights into the biological relevance of CUD-associated DNAm changes in the human brain.

While epigenome-wide alterations associated with CUD have been studied in the human striatum, they have not been reported for cortical brain regions so far. The prefrontal cortex (PFC) has a key role in the addiction cycle due to its relation to dysregulated cognitive processes in SUDs such as inhibitory control and response to drug-related cues (24). In CUD, resting-state connectivity alterations have been found in frontostriatal circuits representing functional connections between the PFC and the striatum (25). As DNAm alterations in the PFC have been described for the Brodmann area 9 (BA9) subregion in other SUDs such as alcohol use disorder (AUD) (26, 27) and opioid use disorder (OUD) (28, 29), BA9 depicts a promising brain region for the epigenome-wide profiling of DNAm in CUD.

Based on epigenome-wide DNAm data, epigenetic clocks provide the opportunity to estimate biological age, a measure that increases naturally with chronological age, but is strongly influenced by disease processes such as chronic inflammation (30). Two commonly used methodologies are Horvath’s (31) and Levine’s et al. (30) (PhenoAge) epigenetic clocks. Epigenetic age acceleration (EAA), defined as the difference between the estimated epigenetic and chronological age, is of particular interest in the context of diseases, with its interpretation as faster or slower aging considering chronological age as the reference. Existing literature reports on a trend toward positive EAA in the PFC of individuals with SUDs such as AUD and OUD compared to a control population, especially for Levine’s clock (28, 32). This implies that the SUD brain appears to be epigenetically older compared to a brain from an individual without SUD at the same chronological age.

In the present study of CUD in human post-mortem brain tissue of BA9, we aimed to expand the epigenome-wide characterization of CUD to the PFC as another key brain region in addiction. Next to performing an EWAS of CUD in BA9, we aimed to identify functional enrichment of CUD-associated differential methylation within biological pathways and utilized epigenetic clocks to investigate the relationship between chronological and biological age in CUD.

2. Materials and methods

2.1. Post-mortem human brain samples

Post-mortem human brain tissue of BA9 from N = 42 male individuals was received from the Douglas Bell Canada Brain Bank (DBCBB, Montréal, QC, Canada). Brain samples were obtained as fresh-frozen tissue and were stored at −80°C until processing. All sample donors were of European ancestry. At the brain bank, phenotypic information was routinely obtained using a psychological autopsy including next-of-kin interviews complemented by the coroner’s notes and the medical record of the deceased donor. Psychiatric diagnoses related to substance use and comorbidities were based on the Diagnostic and Statistical Manual of Mental Disorders (DSM)-IV classification system, however, we used the more recent nomenclature of DSM-5, i.e., SUD instead of substance dependence throughout this study. We included subjects with donor age > 18 and information on CUD status (yes = CUD, no = no CUD). Exclusion criteria were severe neurodevelopmental and neuropsychiatric disorders, SUDs other than CUD or AUD, and psychiatric diagnoses except for major depressive disorder diagnoses (MDD, depressive disorder NOS). Matching was performed by age and comorbidities resulting in N = 21 CUD case/control pairs.

2.2. DNA extraction and methylation data analysis

Extraction of DNA from N = 42 BA9 post-mortem brain samples was performed using the DNeasy Blood and Tissue Kit (Qiagen, Hilden, Germany). DNA was stored at −20°C and epigenome-wide DNAm was determined using the Illumina Infinium MethylationEPIC BeadChip (Illumina, San Diego, CA, USA). Samples were randomized based on CUD case/control status and comorbidities such as depressive disorders and AUD prior to the array processing steps.

Quality control (QC) was performed on raw data using a modified version of the CPACOR-pipeline (33) as described in Zillich et al. (27). All statistical analyses were performed using the statistical computing environment R v4.2.1 (34). In brief, sample exclusion criteria were missing rate > 0.10 and discordance between phenotypic and predicted sex based on DNAm. cytosine-phosphate-guanine (CpG) sites were excluded if (1) the call rate was below 0.95, (2) they are located on sex chromosomes, and (3) their genomic position is related to a genetic variant with MAF > 0.10.

Raw intensities were quantile-normalized followed by a conversion to beta-values of methylation. Beta-values were logit-transformed to M-values as recommended by Du et al. (35). To assess DNAm signatures of CUD in post-mortem brain tissue, we applied a linear regression model using M-values as the dependent variable and CUD status as the predictor. Donor age, post-mortem interval (PMI), brain pH, diagnosis of a depressive disorder, and AUD status were included as covariates in the regression model. We performed a variance partition analysis using the R package variancePartition v1.26.0 (36) to confirm the covariates (see also Supplementary Figure 1). As variance explanation by cause of death was only minimal and would have made the statistical model significantly more complex, we did not include it as a covariate. Based on the methylation data, we estimated the fraction of neuronal and non-neuronal cells with the Houseman (37) approach using the dorsolateral PFC (dlPFC) dataset (38) as reference. We included the fraction of neuronal cells as a covariate to correct for cell-type heterogeneity in the samples. To adjust for technical effects from microarray processing, we derived principal components from the internal control probes of the EPIC array and included the first 10 principal components (cpPC1–10) into the statistical model.

This resulted in the following linear model:

$$DNAm\sim CUDstatus+age+PMI+pH+depressivedisorderstatus+AUDstatus+neuronalfraction+cpPC1+cpPC2\mathrm{\dots }+cpPC10(1)$$

We performed Benjamini–Hochberg false discovery rate (FDR) correction (39) to adjust for multiple testing. Annotation of CpG sites to genes and genomic background was based on the Illumina manifest as downloaded on 10th August 2018:

http://webdata.illumina.com.s3-website-us-east-1.amazonaws.com/downloads/productfiles/methylationEPIC/infinium-methylationepic-v-1-0-b4-manifest-file-csv.zip.

To test for enrichment of CUD-associated CpG sites within categories of gene regulatory features (TSS1500, TSS200, 5′-UTR, 1st Exon, ExonBnd, Body, 3′-UTR, and intergenic regions) and the CpG island (CGI) background (N-Shelf, N-Shore, Island, S-Shore, S-Shelf, and open sea regions), we performed chi-squared tests against the distribution of features among all probes on the EPIC array. The results of the enrichment analysis were adjusted for multiple testing using Bonferroni correction.

2.3. Differentially methylated regions

Following the assessment of CUD-associated methylation differences at single CpG sites, we investigated DMRs associated with CUD status. We used the dmrff algorithm v1.0.0 (40) for the identification of DMRs which was recently shown to outperform several other DMR identification tools when applied to methylation array data (41). We applied dmrff on the EWAS summary statistics (N = 654,448 CpG sites). The maximum distance between individual CpG sites (maxgap) was set to 500 bp and an association p-value threshold (p. cutoff) of 0.01 was used. The implemented Bonferroni correction was used to adjust for multiple testing during the identification of DMRs.

2.4. Replication of CpG-sites identified in the striatum and in whole blood

We further tested whether DNAm signatures, which had previously been associated with CUD in whole blood (20), the CN (22), and the NAc (23) were also associated with CUD in the present sample of BA9. For DMRs overlapping between the two studies in the striatum, we tested CpG sites in close proximity (500 bp up- and downstream) for their association with CUD in our sample. The comparison with CUD-associated CpG sites derived from whole blood was restricted to the subset of CpG sites that was available in our analysis in BA9.

2.5. Gene-Ontology enrichment analysis

To investigate the potential biological relevance of CUD-associated differential methylation, we analyzed the overrepresentation of differentially methylated CpG sites within biological pathways. Under the assumption of hypomethylation leading to increased gene expression and vice versa, we applied missMethyl v1.30.0 (42) separately to hypermethylated and hypomethylated CpG sites after selection for a CUD association p-value < 0.001. The set of N = 654,448 CpG sites remaining after QC was used as the reference set during the enrichment analysis. Further, we performed an enrichment analysis on DMRs resulting from dmrff using the “goregion” function in missMethyl. Correction for multiple testing in results from missMethyl was performed using FDR.

2.6. WGCNA

To further characterize the biological and functional relevance of CUD-associated methylation signatures in the brain, we performed a weighted correlation network analysis using WGCNA v1.71 (43). WGCNA identifies co-methylation modules and provides a correlation estimate between modules and traits such as CUD status. We prioritized CpG sites based on their annotation to promoter/transcription start sites (TSS 200 and TSS1500) under the assumption of a cis-regulatory role of their methylation levels on gene expression. Quantile normalized beta values of the resulting 117,876 CpG sites were used as input data for WGCNA. Soft power threshold picking was performed according to the R_signed² > 0.90 criterion resulting in a power threshold of five for the methylation data in BA9. Parameters for constructing the signed network in a block-wise approach were minModuleSize = 30, maxBlockSize = 36,000, and mergeCutHeight = 0.25. For the CUD-associated co-methylation modules, enrichment analysis was performed in missMethyl using the N = 117,876 input CpG sites as background. Module hub genes were derived by annotation of CpG sites to genes and were defined as the top 0.5% of annotated genes in a co-methylation module ranked by the product of module membership (correlation of module eigengene with methylation profile) and gene significance (correlation between CpG site and trait). Gene significance was calculated for each CpG site. Protein-protein interaction (PPI) networks of hub genes from the CUD-associated co-methylation modules were generated in Cytoscape v3.9.1 (44) using the stringApp v1.7.0 (45). Networks of the “full STRING network” type were generated by applying the STRING protein query option on the hub gene list. An interaction score threshold of 0.7 was used to include only high-confidence interactions. PPI networks were visualized using the “Prefuse Force Directed Layout” option while unconnected nodes (“Singletons”) and hidden edges were removed from the visualization. Information on stringdb identifiers corresponding to the nodes of the final PPI network for each module is provided in Supplementary Tables 7A–D.

2.7. Epigenetic age estimates

Using quantile-normalized beta-values of the BA9 methylation data, we estimated EAA using epigenetic clocks as implemented in the R package methylclock v1.2.1 (46). We applied Horvath’s (31) (Horvath) and Levines’s (30) (PhenoAge) epigenetic clocks to assess epigenetic age in the context of CUD directly in the brain. We first investigated the relationship between epigenetic age and chronological age using correlation analysis. Next, we estimated EAA (ageAccel), defined as the residuals of regressing chronological age on epigenetic age. We performed t-tests to assess the group difference between individuals with and without CUD and linear regression models to evaluate the association with CUD status while adjusting for the covariates neuronal cell fraction, pH, PMI, depressive disorder, and AUD status.

2.8. Ethics approval statement

Post-mortem brain tissue sampling at the DBCBB was performed according to their established ethical standards including written informed consent from the next-of-kin for each subject. Our study design was approved by the Ethics Committee II of the University of Heidelberg, Medical Faculty Mannheim, Germany, under the register number 2021-681.

3. Results

Demographic data of sample donors are summarized in Table 1, while a detailed overview of toxicology and cause of death is provided in Supplementary Table 1. After QC of methylation data, all N = 42 subjects, and a total of 654,448 CpG sites remained for further analysis.

TABLE 1

Table 1. Demographic data of sample donors included in the present study on cocaine use disorder (CUD).

3.1. Epigenome-wide association study of CUD

None of the CpG sites passed epigenome-wide significance (p < 7.64e–08). The strongest association with CUD in BA9 was detected for the CpG site cg04659218 (p = 7.98e–07, FDR = 0.52) on chromosome 11. A Manhattan plot summarizing the results is displayed in Figure 1A. Descriptively, several of the CpG sites showing stronger associations with CUD status were annotated to genes encoding transcription factors and DNA binding domain-containing proteins: ETS Transcription Factor ERG (ERG, p = 4.15e–06), Zinc Finger and BTB Domain Containing 4 (ZBTB4, p = 1.01e–05), RNA Polymerase II (POL-II) Subunit A (POLR2A, p = 1.01e–05), and Zinc Finger and BTB Domain Containing 22 (ZBTB22, p = 1.97e–05, see Figure 1B). EWAS summary statistics of the 1,000 CpG sites displaying the strongest association with CUD status are provided in Supplementary Table 2. A Q–Q plot (λ = 0.969) derived from the EWAS p-values is shown in Supplementary Figure 2.

FIGURE 1

Figure 1. Results of the EWAS of CUD assessing N = 654,448 cytosine-phosphate-guanine (CpG) sites in N = 42 post-mortem human brain samples from Brodmann area 9 (BA9). (A) In the Manhattan plot, the 10 CpG sites displaying the strongest association with CUD status were highlighted in blue and were annotated to gene name if available. The 20 CUD-associated DMRs are highlighted in orange. The solid blue line depicts epigenome-wide significance (p < 7.64e–08). (B) EWAS summary statistics for the 10 CpG sites showing the strongest association with CUD status.

3.2. Differentially methylated regions associated with CUD status

A total of 20 DMRs were significantly associated with CUD status after multiple testing correction. We detected more hypermethylated (N = 17, 85%) than hypomethylated (N = 3, 15%) DMRs. The DMR consisting of six hypermethylated CpG sites annotated to the genes LOC101927196 and Fibrous Sheath Interacting Protein 2 (FSIP2) displayed the strongest association with CUD status in BA9 (p = 1.74e–16). Results of the DMR analysis including summary statistics for all 20 identified regions are visualized in Table 2.

TABLE 2

Table 2. Cocaine use disorder (CUD)-associated differentially methylated regions (DMRs) (N = 20).

3.3. Replication of CpG sites identified in whole blood and striatal post-mortem brain tissue

Of the 186 CpG sites associated with cocaine and crack cocaine use in whole blood, 138 were available in our dataset after QC. Six CpG sites annotated to the genes PDSS2, EHMT2, C2, PDE12, ARF4, ARF4-AS1, RMDN2, RMDN2-AS1, and CDKN1C were nominally associated in the present analysis (see also Supplementary Table 3A) but did not remain significant after multiple testing correction. The DMRs previously identified as being conserved between the CN and NAc annotated to the genes SPEG, TH, PDXDC1, LDHD, and SYT5 were not replicated in our sample of BA9. Similarly, CpG sites located in a genomic window 500 bp up- and downstream of these overlapping striatal DMRs were not significantly associated with CUD in BA9 (Supplementary Table 3B). However, when we compared the DMRs in BA9 with DMRs detected in the CN and the NAc individually, we found one hypermethylated CUD-associated DMR annotated to the Aldehyde Dehydrogenase 3 Family Member B1 (ALDH3B1) gene on chromosome 11 that was conserved between BA9 and the CN.

3.4. Functional annotation and enrichment analysis of CUD-associated CpG sites

To analyze enrichment within genomic features (Figure 2A) and the CGI background (Figure 2B), we prioritized our findings based on CpG sites associated with CUD at nominal significance (p < 0.05). This resulted in a subset of N = 29,176 CpG sites of which N = 14,737 were hypermethylated and N = 14,439 were hypomethylated. We detected a balanced ratio between hyper- and hypomethylation in this subset (Figure 2C). For the CUD-associated CpG sites, a statistically significant enrichment within 1st Exon, N-Shore, and CGI regions was detected, whereas a significant depletion from intergenic and open sea regions was found.

FIGURE 2

Figure 2. Results from the epigenome-wide association study (EWAS) of CUD were prioritized based on association p-values resulting in a subset of N = 29,176 CpG sites with p_assoc < 0.05 that was further analyzed by annotation to panel (A) gene regulatory features and (B) the CpG island (CGI) background. The proportion of CUD-associated CpG sites (blue) within different groups was compared to the EPIC array background (gray). Differences remaining statistically significant after Bonferroni multiple testing correction are highlighted using asterisks (* = p_adj < 0.05, ** = p_adj < 0.01, *** = p_adj < 0.001). (C) Hypermethylation (light-blue) and hypermethylation (dark blue) among the N = 29,176 CUD-associated CpG sites. The top 20 CUD-associated GO terms resulting from missMethyl (42) are shown for (D) hypermethylated and (E) hypomethylated CpG sites associated with CUD (p < 0.001).

The most significantly enriched Gene Ontology (GO) term was response to exogenous dsRNA (p = 3.71e–04) for hypermethylated CpG sites and regulation of calcium ion transmembrane transporter activity (p = 2.39e–04) for the subset of hypomethylated CpG sites. No enrichment remained statistically significant after multiple testing correction (Supplementary Table 4A). Descriptively, genes harboring hypermethylated CpG sites were overrepresented within biological pathways involved in cation and especially mechanosensitive calcium ion channel activity, cyclic adenosine monophosphate (cAMP)-dependent protein kinase activity (PKA), and immunological signaling mainly related to the cyclic GMP-AMP synthase (CGAS)-mediated response to nucleic acids (Figure 2D). In contrast, the gene set derived from hypomethylated CpG sites was enriched within pathways related to calcium transmembrane transport, negative regulation of muscle relaxation, and cAMP binding (Figure 2E). Also, positive regulation of inositol-1,4,5-trisphosphate-sensitive calcium release channel activity was among the enriched GO terms detected for hypomethylated CpG sites. The strongest overrepresentation resulting from the enrichment analysis using DMRs was observed for “thiosulfate-thiol sulfurtransferase activity” (Supplementary Table 4B).

3.5. Weighted correlation network analysis

Due to the known cis-regulatory role of DNAm on gene expression, we extracted TSS1500- and TSS200-associated CpG sites from the methylation matrix to identify transcriptionally relevant co-methylation networks. After hierarchical clustering in WGCNA, the sample tree was pruned at the height of 20 resulting in the removal of one sample from the analysis that was considered as an outlier (N = 41 remaining). The sample dendrogram was related to traits resulting in the clustering of CUD cases and non-CUD individuals within subgroups (Figure 3A). Network construction led to the identification of a total of 69 modules ranging in size between 38 and 35,758 CpG sites. A module-trait correlation heatmap for all resulting modules is shown in Supplementary Figure 3.

FIGURE 3

Figure 3. Functional analysis of the CUD-associated co-methylation modules identified in WGCNA. Quantile-normalized beta values of methylation for N = 117,876 promoter-associated cytosine-phosphate-guanine (CpG) sites (TSS200 and TSS1500) were extracted and used as an input dataset for WGCNA. (A) Dendrogram of hierarchical sample clustering with trait annotation. (B) Correlation plots of module membership and CUD gene significance for the four CUD-associated co-methylation modules (brown, brown4, steelblue, blue). (C) Module-trait relationships for the CUD-associated modules. Module-trait correlation coefficients are color-coded, and correlation p-values are shown in brackets. (D) Results of the Gene Ontology (GO) enrichment analysis. The 10 GO terms with the strongest enrichment within the respective co-methylation module are displayed (* = FDR < 0.05). CUD, cocaine use disorder; PMI, post-mortem interval; neuronal fraction, estimated neuronal cell type proportion based on Houseman et al. (37); AUD, alcohol use disorder.

Of the 69 modules, four were significantly associated with CUD status: brown (N = 6,309 CpG sites, p = 9.15e–03), brown4 (N = 187 CpG sites, p = 0.04), blue (N = 14,050 CpG sites, p = 0.02), and steelblue (N = 340 CpG sites, p = 0.02). A positive, highly significant correlation between module membership and gene significance for CUD status was detected within these modules (Figure 3B). Whereas module brown was positively associated with CUD status (r = 0.40), a negative correlation was identified for modules brown4, blue, and steelblue (−0.37 ≤ r ≤ −0.32, see Figure 3C). In addition, for the brown module, significant negative correlations were detected with age (r = −0.5, p = 8.11e–04) and the fraction of neuronal cells (r = −0.91, p = 7.69e–17). The blue module also depicts a cell type-specific network displaying a highly significant positive correlation with neuronal fraction (r = 0.96, p = 1.78e–22). To further assess the biological processes related to CUD-associated modules, we performed GO enrichment analyses (Figure 3D). After multiple testing correction, GO enrichment in modules brown and blue remained statistically significant. Descriptively, module brown was enriched for synaptic signaling and neurotransmitter transport, whereas module brown4 revealed GO terms related to transcription factor binding and neuroblast as well as synapse maturation. The blue module displayed enrichment of large GO terms involved in basic cellular processes and compartments, such as extracellular region, plasma membrane, and cell adhesion. Within module steelblue, nominally significant pathway enrichment was related to ion channel activity, membrane potential homeostasis, and vesicle transport. Full enrichment statistics for the top enriched GO terms in each module are provided in Supplementary Table 5.

To further characterize the molecular processes enriched in the CUD-associated modules, we constructed PPI networks based on module hub genes (listed in Supplementary Table 6). The resulting PPI network plots are displayed in Supplementary Figures 4A–D. Among PPI network hub nodes of modules brown and blue (Supplementary Table 8), we found members of the Transforming Growth Factor Beta (TGF-β) signaling cascade such as TGF Alpha (TGFA), TGFB1, and Endoglin (ENG). Reflecting the pathway enrichment for neuron-related biological processes in modules brown, brown4, and steelblue, we detected multiple hub nodes with an important role in neuronal activity and function within these modules. Among them were the Ca_v1.2 subunit (CACNA1C, module steelblue), the glucocorticoid receptor (NR3C1, module brown), the acetylcholine transporter Solute Carrier Family 18 Member A3 (SLC18A3, module brown4), and Vesicle Associated Membrane Protein 2 (VAMP2, module steelblue). Further, several tyrosine kinases and their downstream targets were found among hub nodes of PPI networks derived from modules brown4 and blue, such as Fibroblast Growth Factor Receptor 1 (FGFR1), Erb-B2 Receptor Tyrosine Kinase 3 (ERBB3), KRAS Proto-Oncogene GTPase (KRAS), and Signal Transducer And Activator Of Transcription 5A (STAT5A). The JUN was the node displaying the strongest connectivity in module steelblue.

3.6. Epigenetic age in CUD

A moderate (R_Levine² = 0.44) to strong (R_Horvath² = 0.80) correlation was found between estimated epigenetic age and chronological age (Figures 4A, B). Levine’s epigenetic clock produced strongly negative PhenoAge estimates (Figure 4B). When assessing ageAccel we observed a trend toward a more positive ageAccel in CUD cases compared to individuals without CUD, however, not statistically significant (p_Horvath = 0.39, p_Levine = 0.15, Figures 4C, D). For both epigenetic clocks, we consistently found positive associations for CUD status while adjusting for covariates in the regression models for ageAccel, although not statistically significant (Supplementary Table 9).

FIGURE 4

Figure 4. Evaluation of epigenetic clocks in cocaine use disorder (CUD) using post-mortem human brain samples of Brodmann Area 9 (BA9). Correlation of chronological age with DNA methylation (DNAm) based (A) Horvath’s and (B) Levine’s epigenetic age estimates as implemented in the R package methylclock (46). Comparison of ageAccel between individuals with and without CUD using (C) Horvath’s and (D) Levine’s epigenetic clocks.

4. Discussion

We performed an epigenome-wide analysis of DNAm differences in BA9 post-mortem human brain tissue from N = 21 CUD cases and N = 21 individuals without CUD. In the BA9 subregion of the human PFC, none of the associations of CpG site methylation levels with CUD status passed the level of epigenome-wide significance. However, we identified 20 DMRs significantly associated with CUD.

Among the strongest single-site associations with CUD status, there were multiple CpG sites annotated to genes involved in transcriptional regulation either via intrinsic transcription factor activity or due to the interaction with RNA POL-II. The transcription factor ERG is involved in the regulation of vasculature integrity (47). For ZBTB4, binding to methylated CpG sites was reported (48) and ZBTB22 is assumed to be involved in transcription regulation via interaction with POL-II (49). In our EWAS, we also detected CUD-associated differential methylation within the largest subunit of POL-II (POL2RA), suggesting a role of CUD-associated methylation signatures in gene expression regulation. Results of the enrichment analysis for CUD-associated CpG sites within genomic features are in line with this: the significant enrichment of CUD-associated differential methylation within features containing transcription factor binding sites and the depletion from intergenic and non-CGI-associated (open sea) regions points toward a regulatory potential of CUD-associated DNAm patterns on gene expression.

For several of the genes harboring DMRs significantly associated with CUD in BA9, there is evidence for their relevance for neuronal function. The DMR displaying the strongest association with CUD status was related to the LOC101927196 gene. Overexpression of the long non-coding RNA LOC101927196 was shown to reduce cell proliferation and support apoptotic processes in the brain in a rat model of autism (50). For NPFFR2 and KALRN, we found direct links to their importance for the biological response to cocaine. NPFFR2 interferes with the opioid system in the brain and its activation contributes to anxiety-like behavior in mice (51). Intraventricular application of NPFF reduced cocaine-induced conditioned place preference in rats pointing toward an important role of NPFF signaling in cocaine reward processing (52). The KALRN gene encodes several isoforms and Kalirin-7 was shown to be involved in the development and function of dendritic spines (53). After repeated cocaine exposure in a mouse model, expression levels of Kalirin-7 were elevated in the NAc contributing to the formation of new dendritic spines. Kalirin-7 knock-out mice displayed increased cocaine-induced locomotor sensitization and a decreased conditioned place preference (54).

Looking at multiple brain regions allows a more complete understanding of the role of DNAm alterations in CUD. In line with results from previous studies on epigenome-wide alterations of CUD in the CN (22) and in the NAc (23), we observed more hypermethylation than hypomethylation at the DMR level in BA9. In our systematic comparison of the DMRs that were consistently associated with CUD in both the CN and NAc, none were detected in our analysis. We detected one DMR annotated to the ALDH3B1 that was conserved between BA9 and the CN. Hypermethylation of this DMR was also consistent between studies. It has to be noted that the EPIC array covers only a finite number of CpG sites which may limit the comparison with results from reduced representation bisulfite sequencing. In addition, we also compared our results with DNAm signatures of CUD in whole blood. While some nominally significant overlap was observed, there was little evidence for a systematic relationship between the CUD-associated CpG sites in BA9 and whole blood, as the overlap was not significant after multiple testing correction. DNAm is highly tissue-specific, which could partly explain the missing overlap. At the same time, the samples were phenotypically different in age and substance use.

In our analysis in BA9, no enrichment of hypo- and hypermethylated CUD-associated CpG sites among GO terms remained statistically significant after multiple testing correction. Further, pathway enrichment was driven by relatively few genes leading to the emergence of multiple biologically related pathways among the results. The functional relevance of CUD-associated differential methylation should thus be evaluated together with the enriched biological processes and hub genes derived from WGCNA co-methylation modules. Three of the four CUD-associated modules were enriched for neuron-specific processes such as synaptic signaling, ion channel activity, and neurotransmitter transport. Next to the association with CUD, modules brown and blue displayed a strong correlation with neuronal fraction, which might also contribute to the enrichment of neuronal processes within these co-methylation modules. PPI networks based on module hub genes revealed CACNA1C and NR3C1 as highly connected genes for which a previous relation to cocaine use is known. CACNA1C is a well-described risk gene for psychiatric disorders such as bipolar disorder, MDD, and schizophrenia (55). By encoding a subunit of the voltage-gated Ca_v1.2 calcium channel, it depicts a key regulator of neuronal excitability. Further, the expression of cocaine-induced locomotor sensitization was shown to be dependent on Ca_v1.2 signaling in mice (56). The NR3C1 gene encodes a glucocorticoid receptor for which downregulation in blood was observed in chronic cocaine users (57) and its DNAm levels were shown to predict substance use phenotypes in adolescents (58). Further, brain-specific conditional knock-out of NR3C1 flattened the dose-response curve for cocaine self-administration and reduced behavioral sensitization in an intravenous cocaine self-administration model in mice (59). Thus, CACNA1C and NR3C1 depict important regulators of the reinforcing effects of cocaine.

We further identified JUN as the hub node with the strongest connectivity in module steelblue. By forming heterodimers with FOS family proteins, JUN is involved in the formation of the AP-1 transcription factor complex, a key regulator of neuroplasticity in addiction (60). In a mouse model, decreased cocaine-induced conditioned place preference was found during expression of a dominant-negative mutant of c-Jun in the brain (61). Epigenetic dysregulation of the JUN interaction network might therefore contribute to neurocircuit alterations in CUD.

In our analysis of epigenetic age in the brain, we observed a trend toward a more positive ageAccel in CUD cases compared to individuals without CUD that was more pronounced with Levine’s clock compared to Horvath’s clock. This is in line with results from previous studies that report on positive ageAccel, particularly for Levine’s epigenetic clock, in the PFC of individuals with opioid intoxication, AUD, and OUD (28, 32). As cocaine intake is associated with neuroinflammatory processes in the brain (62, 63), ageAccel in CUD could be a consequence of an inflammatory response that is particularly reflected by Levine’s (30) epigenetic clock. However, as comorbid diseases are frequently observed in CUD, an increase in ageAccel could at least in part be driven by systemic diseases that are themselves associated with an accelerated epigenetic age. Due to the lack of information on somatic comorbidities, we were not able to address this point in the present study. Further, for the interpretation of epigenetic age in the brain, it must be noted that epigenetic clocks are at least to some extent tissue-specific and might not be applicable to every tissue. We only assessed epigenetic clocks that were either trained on non-blood tissues (Horvath) or at least showed a high correlation of findings across tissues (Levine). The negative values of Levine’s PhenoAge we observed in some of the brain samples might represent a post-mortem brain tissue artifact, as Levine’s clock was initially trained on blood methylation data. Similarly, in a recent study in post-mortem human brain tissue, negative PhenoAge estimates were observed for Levine’s clock (32). Within-subject comparisons using blood and brain tissue from the same person should be performed to investigate such potential tissue artifacts. Further, analyses with larger sample sizes, careful adjustment for comorbidities, and appropriate estimation methods are required to evaluate if accelerated epigenetic aging in the brain is a true phenomenon in CUD.

Several limitations apply to our study on DNAm signatures of CUD in post-mortem human brain tissue. First, our sample size of N = 42 is comparatively low for an epigenome-wide analysis and is associated with limitations in statistical power. Another limitation is the frequent comorbidity of CUD with other SUDs and mood disorders, which is well reflected in our sample. Drug intoxication at death and medication prior to death could have affected methylation levels and therefore might confound the results for CUD. Next to the CUD phenotype, assessment of used substances and comorbidities including their treatment is crucial to adjust for confounding effects on DNAm signatures. Due to the small sample size, we were not able to adjust for drug intoxication at death and medication prior to death. Thus, we cannot rule out a potential confounding of the identified CUD-associated methylation signatures by epigenetic effects of other substances. Third, the case-control design of the EWAS does not allow the assessment of temporal changes in DNAm levels. Thus, within the present study, we are not able to separate differential methylation acquired during the disease course of CUD from an epigenetic predisposition to CUD.

In the present study, we identified CUD-associated DNAm signatures in the BA9 subregion of the human PFC. We detected CUD-associated differential methylation within several transcription factor genes, an enrichment within genomic features of transcription regulation, and an association with gene regulatory networks involved in neurotransmission and neuroplasticity. Specifically, epigenetic alterations related to calcium and cAMP signaling as well as the CACNA1C and JUN PPI networks might contribute to the development of the altered neurocircuits observed in the CUD brain.

To follow up on these hypotheses, further studies are required. These should be based on larger sample sizes with detailed phenotyping to minimize the confounding of results in association studies. Further, integration of epigenome-wide methylation data with other omics data such as transcriptomics or proteomics is necessary to prioritize and substantiate findings from the EWAS. Such multi-omics studies depict promising approaches for deeper profiling of altered biological processes in SUDs as recently shown for AUD (64) and OUD (65). To address the temporal dynamics of molecular changes in the disease course of CUD, a longitudinal study assessing different stages of CUD within the same individuals should be performed. While a longitudinal design is nearly impossible using post-mortem brain tissue, it is a feasible approach for blood samples. The application of such innovative methods could pave the way for the development of precision medicine approaches thereby addressing the urgent need for novel therapeutic options in CUD.

Data availability statement

Methylation data presented in this study have been deposited at the European Genome-phenome Archive (EGA), which is hosted by the EBI and the CRG, under accession number: EGAS00001006826 (https://ega-archive.org/studies/EGAS00001006826).

Ethics statement

The studies involving human participants were reviewed and approved by the Ethics Committee II of the Medical Faculty Mannheim, Heidelberg University, Mannheim, Germany under the register number: 2021-681. Written informed consent was not provided because post-mortem brain tissue was obtained from individuals that died by suicide or sudden death. Following the established ethical standards of the DBCBB, written informed consent was obtained from the next-of-kin for each subject.

Author contributions

SW, RS, AH, and MR: conceptualization and funding acquisition. EP, FS, JCF, and LZ: data curation. EP and LZ: formal analysis and writing—original draft. EP, LZ, JCF, MF, AH, RS, and SW: investigation. LZ, EP, FS, and JF: methodology. SW, MN, AH, and RS: project administration. GT, NM, AH, and MN: resources. SW, JF, MR, and RS: supervision. All authors contributed to the article writing—review and editing and approved the submitted version.

Funding

Funding was provided by the German Federal Ministry of Education and Research (BMBF) through the e:Med research program “A systems-medicine approach toward distinct and shared resilience and pathological mechanisms of substance use disorders” (01ZX01909 to RS, MR, SW, and AH). Further, by “Toward Targeted Oxytocin Treatment in Alcohol Addiction (Target-OXY)” (031L0190A to JCF), the ERA-NET program: Psi-Alc (FKZ: 01EW1908), the Hetzler Foundation for Addiction Research (to AH), and the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) Project-ID 402170461–TRR 265 (66) to RS, AH, and MR.

Acknowledgments

We thank Claudia Schäfer-Arnold and Elisabeth Röbel for technical assistance.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpsyt.2023.1075250/full#supplementary-material

References

1. UNODC. World Drug Report 2022. New York, NY: United Nations publication (2022).

Google Scholar

2. European Monitoring Centre for Drugs and Drug Addiction. European Drug Report 2022: Trends and Developments. Luxembourg: Publications Office of the European Union (2022).

Google Scholar

3. Lopez-Quintero C, Pérez de los Cobos J, Hasin D, Okuda M, Wang S, Grant B, et al. Probability and predictors of transition from first use to dependence on nicotine, alcohol, cannabis, and cocaine: results of the National Epidemiologic Survey on Alcohol and Related Conditions (NESARC). Drug Alcohol Depend. (2011) 115:120–30. doi: 10.1016/j.drugalcdep.2010.11.004

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Butler AJ, Rehm J, Fischer B. Health outcomes associated with crack-cocaine use: systematic review and meta-analyses. Drug Alcohol Depend. (2017) 180:401–16. doi: 10.1016/j.drugalcdep.2017.08.036

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Bierut LJ, Strickland J, Thompson J, Afful S, Cottler L. Drug use and dependence in cocaine dependent subjects, community-based individuals, and their siblings. Drug Alcohol Depend. (2008) 95:14–22. doi: 10.1016/j.drugalcdep.2007.11.023

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Volkow ND, Michaelides M, Baler R. The neuroscience of drug reward and addiction. Physiol Rev. (2019) 99:2115–40. doi: 10.1152/physrev.00014.2018

PubMed Abstract | CrossRef Full Text | Google Scholar

7. Volkow ND, Fowler JS, Wang GJ. The addicted human brain: insights from imaging studies. J Clin Invest. (2003) 111:1444–51. doi: 10.1172/JCI18533

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Koob GF, Volkow ND. Neurocircuitry of addiction. Neuropsychopharmacology. (2010) 35:217–38. doi: 10.1038/npp.2009.110

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Nielsen DA, Utrankar A, Reyes J, Simons D, Kosten T. Epigenetics of drug abuse: predisposition or response. Pharmacogenomics. (2012) 13:1149–60. doi: 10.2217/pgs.12.94

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Nestler EJ, Lüscher C. The molecular basis of drug addiction: linking epigenetic to synaptic and circuit mechanisms. Neuron. (2019) 102:48–59. doi: 10.1016/j.neuron.2019.01.016

PubMed Abstract | CrossRef Full Text | Google Scholar

11. Stewart AF, Fulton SL, Maze I. Epigenetics of drug addiction. Cold Spring Harb Perspect Med. (2021) 11:a040253. doi: 10.1101/cshperspect.a040253

PubMed Abstract | CrossRef Full Text | Google Scholar

12. Vaillancourt K, Ernst C, Mash D, Turecki G. DNA methylation dynamics and cocaine in the brain: progress and prospects. Genes. (2017) 8:138. doi: 10.3390/genes8050138

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Tian W, Wang J, Zhang K, Teng H, Li C, Szyf M, et al. Demethylation of c-MYB binding site mediates upregulation of Bdnf IV in cocaine-conditioned place preference. Sci Rep. (2016) 6:22087. doi: 10.1038/srep22087

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Anier K, Malinovskaja K, Aonurm-Helm A, Zharkovsky A, Kalda A. DNA methylation regulates cocaine-induced behavioral sensitization in mice. Neuropsychopharmacology. (2010) 35:2450–61. doi: 10.1038/npp.2010.128

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Wright KN, Hollis F, Duclot F, Dossat AM, Strong CE, Francis TC, et al. Methyl supplementation attenuates cocaine-seeking behaviors and cocaine-induced c-Fos activation in a DNA methylation-dependent manner. J Neurosci. (2015) 35:8948–58. doi: 10.1523/JNEUROSCI.5227-14.2015

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Anier K, Zharkovsky A, Kalda A. S-adenosylmethionine modifies cocaine-induced DNA methylation and increases locomotor sensitization in mice. Int J Neuropsychopharmacol. (2013) 16:2053–66. doi: 10.1017/S1461145713000394

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Tian W, Zhao M, Li M, Song T, Zhang M, Quan L, et al. Reversal of cocaine-conditioned place preference through methyl supplementation in mice: altering global DNA methylation in the prefrontal cortex. PLoS One. (2012) 7:e33435. doi: 10.1371/journal.pone.0033435

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Cannella N, Oliveira AM, Hemstedt T, Lissek T, Buechler E, Bading H, et al. Dnmt3a2 in the nucleus accumbens shell is required for reinstatement of cocaine seeking. J Neurosci. (2018) 38:7516–28. doi: 10.1523/JNEUROSCI.0600-18.2018

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Feng J, Shao N, Szulwach KE, Vialou V, Huynh J, Zhong C, et al. Role of Tet1 and 5-hydroxymethylcytosine in cocaine action. Nat Neurosci. (2015) 18:536–44. doi: 10.1038/nn.3976

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Camilo C, Maschietto M, Vieira HC, Tahira AC, Gouveia GR, Santos AC, et al. Genome-wide DNA methylation profile in the peripheral blood of cocaine and crack dependents. Braz J Psychiatry. (2019) 41:485–93. doi: 10.1590/1516-4446-2018-0092

PubMed Abstract | CrossRef Full Text | Google Scholar

21. Shu C, Jaffe AE, Sabunciyan S, Ji H, Astemborski J, Sun J, et al. Epigenome-wide association analyses of active injection drug use. Drug Alcohol Depend. (2022) 235:109431. doi: 10.1016/j.drugalcdep.2022.109431

PubMed Abstract | CrossRef Full Text | Google Scholar

22. Vaillancourt K, Yang J, Chen GG, Yerko V, Théroux J, Aouabed Z, et al. Cocaine-related DNA methylation in caudate neurons alters 3D chromatin structure of the IRXA gene cluster. Mol Psychiatry. (2021) 26:3134–51. doi: 10.1038/s41380-020-00909-x

PubMed Abstract | CrossRef Full Text | Google Scholar

23. Vaillancourt K, Chen GG, Fiori L, Maussion G, Yerko V, Théroux JF, et al. Methylation of the tyrosine hydroxylase gene is dysregulated by cocaine dependence in the human striatum. iScience. (2021) 24:103169. doi: 10.1016/j.isci.2021.103169

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Goldstein RZ, Volkow ND. Dysfunction of the prefrontal cortex in addiction: neuroimaging findings and clinical implications. Nat Rev Neurosci. (2011) 12:652–69. doi: 10.1038/nrn3119

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Hu Y, Salmeron BJ, Gu H, Stein EA, Yang Y. Impaired functional connectivity within and between frontostriatal circuits and its association with compulsive drug use and trait impulsivity in cocaine addiction. JAMA Psychiatry. (2015) 72:584–92. doi: 10.1001/jamapsychiatry.2015.1

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Lohoff FW, Roy A, Jung J, Longley M, Rosoff DB, Luo A, et al. Epigenome-wide association study and multi-tissue replication of individuals with alcohol use disorder: evidence for abnormal glucocorticoid signaling pathway gene regulation. Mol Psychiatry. (2021) 26:2224–37. doi: 10.1038/s41380-020-0734-4

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Zillich L, Frank J, Streit F, Friske M, Foo J, Sirignano L, et al. Epigenome-wide association study of alcohol use disorder in five brain regions. Neuropsychopharmacology. (2021) 47:832–9. doi: 10.1038/s41386-021-01228-7

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Shu C, Sosnowski DW, Tao R, Deep-Soboslay A, Kleinman JE, Hyde TM, et al. Epigenome-wide study of brain DNA methylation following acute opioid intoxication. Drug Alcohol Depend. (2021) 221:108658. doi: 10.1016/j.drugalcdep.2021.108658

PubMed Abstract | CrossRef Full Text | Google Scholar

29. Liu A, Dai Y, Mendez EF, Hu R, Fries GR, Najera KE, et al. Genome-wide correlation of DNA methylation and gene expression in postmortem brain tissues of opioid use disorder patients. Int J Neuropsychopharmacol. (2021) 24:879–91. doi: 10.1093/ijnp/pyab043

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Levine ME, Lu AT, Quach A, Chen BH, Assimes TL, Bandinelli S, et al. An epigenetic biomarker of aging for lifespan and healthspan. Aging. (2018) 10:573–91. doi: 10.18632/aging.101414

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Horvath S. DNA methylation age of human tissues and cell types. Genome Biol. (2013) 14:R115. doi: 10.1186/gb-2013-14-10-r115

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Cabrera-Mendoza B, Stertz L, Najera K, Selvaraj S, Teixeira A, Meyer T, et al. Within subject cross-tissue analyzes of epigenetic clocks in substance use disorder postmortem brain and blood. Am J Med Genet B Neuropsychiatr Genet. (2022) 192:13–27. doi: 10.1002/ajmg.b.32920

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Lehne B, Drong AW, Loh M, Zhang W, Scott WR, Tan S, et al. A coherent approach for analysis of the Illumina HumanMethylation450 BeadChip improves data quality and performance in epigenome-wide association studies. Genome Biol. (2015) 16:37.

Google Scholar

34. R Core Team. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing. Vienna: R Foundation (2022).

Google Scholar

35. Du P, Zhang X, Huang C, Jafari N, Kibbe WA, Hou L, et al. Comparison of Beta-value and M-value methods for quantifying methylation levels by microarray analysis. BMC Bioinformatics. (2010) 11:587. doi: 10.1186/1471-2105-11-587

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Hoffman GE, Schadt EE. VariancePartition: interpreting drivers of variation in complex gene expression studies. BMC Bioinformatics. (2016) 17:483. doi: 10.1186/s12859-016-1323-z

PubMed Abstract | CrossRef Full Text | Google Scholar

37. Houseman EA, Accomando WP, Koestler DC, Christensen BC, Marsit CJ, Nelson HH, et al. DNA methylation arrays as surrogate measures of cell mixture distribution. BMC Bioinformatics. (2012) 13:86. doi: 10.1186/1471-2105-13-86

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Jaffe A, Kaminsky Z. FlowSorted.DLPFC.450k: Illumina HumanMethylation Data on Sorted Frontal Cortex Cell Populations. R Package Version 1.24.0. Vienna: R Foundation (2020).

Google Scholar

39. Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J R Stat Soc Ser B. (1995) 57:289–300. doi: 10.1111/j.2517-6161.1995.tb02031.x

CrossRef Full Text | Google Scholar

40. Suderman M, Staley J, French R, Arathimos R, Simpkin A, Tilling K, et al. DMRFF: identifying differentially methylated regions efficiently with power and control. bioRxiv [Preprint]. (2018) 508556. doi: 10.1101/508556

CrossRef Full Text | Google Scholar

41. Lent S, Cardenas A, Rifas-Shiman SL, Perron P, Bouchard L, Liu C, et al. Detecting differentially methylated regions with multiple distinct associations. Epigenomics. (2021) 13:451–64. doi: 10.2217/epi-2020-0344

PubMed Abstract | CrossRef Full Text | Google Scholar

42. Phipson B, Maksimovic J, Oshlack A. missMethyl: an R package for analyzing data from Illumina’s HumanMethylation450 platform. Bioinformatics. (2016) 32:286–8. doi: 10.1093/bioinformatics/btv560

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. (2008) 9:559. doi: 10.1186/1471-2105-9-559

PubMed Abstract | CrossRef Full Text | Google Scholar

44. Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. (2003) 13:2498–504. doi: 10.1101/gr.1239303

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Doncheva NT, Morris JH, Gorodkin J, Jensen LJ. Cytoscape StringApp: network analysis and visualization of proteomics data. J Proteome Res. (2019) 18:623–32. doi: 10.1021/acs.jproteome.8b00702

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Pelegí-Sis DÓ, Prado PD, Ronkainen J, Bustamante M, González JR. Methylclock: a bioconductor package to estimate DNA methylation age. Bioinformatics. (2021) 37:1759–60. doi: 10.1093/bioinformatics/btaa825

PubMed Abstract | CrossRef Full Text | Google Scholar

47. Birdsey G, Shah A, Dufton N, Reynolds L, Osuna Almagro L, Yang Y, et al. The endothelial transcription factor ERG promotes vascular stability and growth through Wnt/β-Catenin signaling. Dev Cell. (2015) 32:82–96. doi: 10.1016/j.devcel.2014.11.016

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Roussel-Gervais A, Naciri I, Kirsh O, Kasprzyk L, Velasco G, Grillo G, et al. Loss of the Methyl-CpG–binding protein ZBTB4 alters mitotic checkpoint, increases aneuploidy, and promotes tumorigenesis. Cancer Res. (2017) 77:62–73. doi: 10.1158/0008-5472.CAN-16-1181

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Safran M, Rosen N, Twik M, BarShir R, Stein T, Dahary D, et al. The GeneCards suite. In: Abugessaisa I, Kasukawa T editors. Practical Guide to Life Science Databases. Singapore: Springer Nature Singapore (2021). p. 27–56. doi: 10.1007/978-981-16-5812-9_2

CrossRef Full Text | Google Scholar

50. Yao W, Huang J, He H. Over-expressed LOC101927196 suppressed oxidative stress levels and neuron cell proliferation in a rat model of autism through disrupting the Wnt signaling pathway by targeting FZD3. Cell Signal. (2019) 62:109328. doi: 10.1016/j.cellsig.2019.05.013

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Lin Y-T, Yu Y, Hong W, Yeh T, Chen T, Chen J, et al. NPFFR2 activates the HPA axis and induces anxiogenic effects in rodents. Int J Mol Sci. (2017) 18:1810. doi: 10.3390/ijms18081810

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Kotlinska J, Pachuta A, Silberring J. Neuropeptide FF (NPFF) reduces the expression of cocaine-induced conditioned place preference and cocaine-induced sensitization in animals. Peptides. (2008) 29:933–9. doi: 10.1016/j.peptides.2008.01.008

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Ma X, Kiraly DD, Gaier ED, Wang Y, Kim E, Levine ES, et al. Kalirin-7 is required for synaptic structure and function. J Neurosci. (2008) 28:12368–82. doi: 10.1523/JNEUROSCI.4269-08.2008

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Kiraly DD, Nemirovsky NE, LaRese TP, Tomek SE, Yahn SL, Olive MF, et al. Constitutive knockout of kalirin-7 leads to increased rates of cocaine self-administration. Mol Pharmacol. (2013) 84:582–90. doi: 10.1124/mol.113.087106

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Dedic N, Pöhlmann ML, Richter JS, Mehta D, Czamara D, Metzger MW, et al. Cross-disorder risk gene CACNA1C differentially modulates susceptibility to psychiatric disorders during development and adulthood. Mol Psychiatry. (2018) 23:533–43. doi: 10.1038/mp.2017.133

PubMed Abstract | CrossRef Full Text | Google Scholar

56. Giordano TP, Tropea TF, Satpute SS, Sinnegger-Brauns MJ, Striessnig J, Kosofsky BE, et al. Molecular switch from L-type Ca v 1.3 to Ca v 1.2 Ca2+ channel signaling underlies long-term psychostimulant-induced behavioral and molecular plasticity. J Neurosci. (2010) 30:17051–62. doi: 10.1523/JNEUROSCI.2255-10.2010

PubMed Abstract | CrossRef Full Text | Google Scholar

57. Schote AB, Jäger K, Kroll SL, Vonmoos M, Hulka LM, Preller KH, et al. Glucocorticoid receptor gene variants and lower expression of NR3C1 are associated with cocaine use. Addict Biol. (2019) 24:730–42. doi: 10.1111/adb.12632

PubMed Abstract | CrossRef Full Text | Google Scholar

58. Raffetti E, et al. DNA methylation of the glucocorticoid receptor gene predicts substance use in adolescence: longitudinal data from over 1000 young individuals. Transl Psychiatry. (2021) 11:477. doi: 10.1038/s41398-021-01601-6

PubMed Abstract | CrossRef Full Text | Google Scholar

59. Deroche-Gamonet V, Sillaber I, Aouizerate B, Izawa R, Jaber M, Ghozland S, et al. The glucocorticoid receptor as a potential target to reduce cocaine abuse. J Neurosci. (2003) 23:4785–90. doi: 10.1523/JNEUROSCI.23-11-04785.2003

PubMed Abstract | CrossRef Full Text | Google Scholar

60. McClung CA, Nestler EJ. Neuroplasticity mediated by altered gene expression. Neuropsychopharmacology. (2008) 33:3–17. doi: 10.1038/sj.npp.1301544

PubMed Abstract | CrossRef Full Text | Google Scholar

61. Peakman M, Colby C, Perrotti LI, Tekumalla P, Carle T, Ulery P, et al. Inducible, brain region-specific expression of a dominant negative mutant of c-Jun in transgenic mice decreases sensitivity to cocaine. Brain Res. (2003) 970:73–86. doi: 10.1016/s0006-8993(03)02230-3

PubMed Abstract | CrossRef Full Text | Google Scholar

62. Sil S, Niu F, Tom E, Liao K, Periyasamy P, Buch S. Cocaine mediated neuroinflammation: role of dysregulated autophagy in pericytes. Mol Neurobiol. (2019) 56:3576–90. doi: 10.1007/s12035-018-1325-0

PubMed Abstract | CrossRef Full Text | Google Scholar

63. López-Pedrajas R, Ramírez-Lamelas D, Muriach B, Sánchez-Villarejo M, Almansa I, Vidal-Gil L, et al. Cocaine promotes oxidative stress and microglial-macrophage activation in rat cerebellum. Front Cell Neurosci. (2015) 9:279. doi: 10.3389/fncel.2015.00279

PubMed Abstract | CrossRef Full Text | Google Scholar

64. Zillich L, Poisel E, Frank J, Foo JC, Friske MM, Streit F, et al. Multi-omics signatures of alcohol use disorder in the dorsal and ventral striatum. Transl Psychiatry. (2022) 12:190. doi: 10.1038/s41398-022-01959-1

PubMed Abstract | CrossRef Full Text | Google Scholar

65. Mendez EF, Wei H, Hu R, Stertz L, Fries GR, Wu X, et al. Angiogenic gene networks are dysregulated in opioid use disorder: evidence from multi-omics and imaging of postmortem human brain. Mol Psychiatry. (2021) 26:7803–12. doi: 10.1038/s41380-021-01259-y

PubMed Abstract | CrossRef Full Text | Google Scholar

66. Heinz A, Kiefer F, Smolka MN, Endrass T, Beste C, Beck A, et al. Addiction research consortium: losing and regaining control over drug intake (ReCoDe)-from trajectories to mechanisms and interventions. Addict Biol. (2020) 25:e12866. doi: 10.1111/adb.12866

PubMed Abstract | CrossRef Full Text | Google Scholar