Skip to main content

Effect of short-term prescription opioids on DNA methylation of the OPRM1 promoter

Abstract

Background

A long-term opioid use has been associated with hypermethylation of the opioid receptor mu 1 (OPRM1) promoter. Very little is currently known about the early epigenetic response to therapeutic opioids. Here, we examine whether we can detect DNA methylation changes associated with a few days’ use of prescribed opioids. Genome-wide DNA methylation was assayed in a cohort of 33 opioid-naïve participants who underwent standard dental surgery followed by opioid self-administration. Saliva samples were collected before surgery (visit 1), and at two postsurgery visits at 2.7 ± 1.5 days (visit 2), and 39 ± 10 days (visit 3) after the discontinuation of opioid analgesics.

Results

The perioperative methylome underwent significant changes over the three visits that were primarily due to postoperative inflammatory response and cell heterogeneity. To specifically examine the effect of opioids, we started with a candidate gene approach and evaluated 10 CpGs located in the OPRM1 promoter. There was a significant cross-sectional variability in opioid use, and for participants who self-administered the prescribed drugs, the total dosage ranged from 5–210 morphine milligram equivalent (MME). Participants were categorized by cumulative dosage into three groups: < 25 MME, 25–90 MME, and ≥ 90 MME. Using mixed-effects modeling, 4 CpGs had significant positive associations with opioid dose at two-tailed p value < 0.05, and overall, 9 of the 10 OPRM1 promoter CpGs showed the predicted higher methylation in the higher dose groups relative to the lowest dose group. After adjustment for age, cellular heterogeneity, and past tobacco use, the promoter mean methylation also had positive associations with cumulative MME (regression coefficient = 0.0002, one-tailed p value = 0.02) and duration of opioid use (regression coefficient = 0.003, one-tailed p value = 0.001), but this effect was significant only for visit 3. A preliminary epigenome-wide association study identified a significant CpG in the promoter of the RAS-related signaling gene, RASL10A, that may be predictive of opioid dosage.

Conclusion

The present study provides evidence that the hypermethylation of the OPRM1 promoter is in response to opioid use and that epigenetic differences in OPRM1 and other sites are associated with a short-term use of therapeutic opioids.

Background

Prescription opioids were once considered as a relatively benign treatment for pain management [1, 2]. However, over the past decade, prescribed analgesics have emerged as a major socio-environmental factor that has contributed to the opioid epidemic [3, 4]. For many individuals who develop opioid use disorder (OUD), the initiation phase may begin with treatment for acute pain or minor surgery, with primary care physicians and dentists accounting for a large fraction of prescribed opioids [5,6,7,8,9,10,11]. Even short-term use (e.g., up to 3 days) is a risk factor for some individuals, and the risk for addiction increases proportionally with dosage and duration of use [8, 9, 12,13,14,15].

Drug addiction is a chronic disease that is triggered by an exposure to an environmental agent. Following the initial exposure, the addictive substance continues to have a persistent effect, and this suggests a form of cellular memory. There is strong evidence that epigenetic processes, including DNA methylation, play a key role in maintaining the long-term effects of the additive substance [16, 17]. Studies particularly in model organisms have shown that drugs of abuse trigger intracellular signaling cascades that alter gene transcription, repeated exposure to the drug then results in remodeling of the epigenome that persists over time, and these epigenetic processes maintain the long-term changes in steady-state gene expression that underlie addiction [16, 18,19,20]. Work in humans generally relies on postmortem tissue from long-term drug users, and studies have found significant epigenetic differences in the brains of former addicts compared to non-addicts [21, 22]. While the brain is the most relevant tissue in terms of neuroadaptation and drug-seeking behavior, epigenetic markers of addiction have also been detected in peripheral tissues such as the blood and sperm [23,24,25,26,27,28]. Easily accessible peripheral tissues are clearly the practical choice when it comes to defining biomarkers of drug use and/or predictors of individual risk for addiction.

The μ-opioid receptor gene (OPRM1) encodes the primary target for both endogenous and exogenous opioids and plays a central role in mediating the rewarding and therapeutic effects. The CpG island located in the promoter of this gene is a potential sensor for drug use, and multiple studies in leukocytes and sperm have found higher DNA methylation among long-term opioid users compared to control samples [23, 29,30,31,32,33]. Hypermethylation of the promoter region has also been found among people with alcohol dependence [34]. However, as all these studies are cross-sectional comparisons between opioid-exposed individuals and controls, there is no definite way to discern whether the epigenetic differences are the cause, or effect, of drug use. Since genetic variants both within and near the OPRM1 gene have also been associated with susceptibility to addiction and drug sensitivity [35,36,37], it is plausible that such epigenetic markers represent genetic effects that preceded drug use. Another lingering question is, if the epigenetic changes are induced by drug use, does the hypermethylation of the promoter CpGs occur only after repeated and sustained exposure, or are these indicators of the early epigenomic, and potentially transcriptomic, responses to drugs? In the case of potent drugs such as opioids, the initial exposure is a crucial phase in the pathway to drug dependence and addiction, and it is reasonable to expect that some of the modification to the epigenome occurs within the first few exposures.

To address these questions, we applied a longitudinal design and collected saliva samples and self-reports of the opioid use from a group of opioid naïve dental patients before oral surgery and at two follow-up visits after surgery. We assayed genome-wide DNA methylation and explored (1) the methylome during the perioperative period, (2) how demographic variables such as age and race/ethnicity relate to methylome changes and immune response, and (3) whether we can discern opioid-associated CpGs from the highly heterogeneous methylome data. As the site of surgery and postsurgery inflammatory response, the saliva presents particular challenges due to immune-related cellular heterogeneity. To overcome this, we applied in silico approaches to deconvolute the underlying cellular heterogeneity and demonstrate the utility of the methylome-based cell estimates as proxies for the immune changes induced by surgery. For the effect of opioids, we specifically focused on the OPRM1 promoter CpGs and evaluated whether the data replicates the CpG hypermethylation. Overall, our results show a dose-dependent increase in methylation at the OPRM1 promoter that can be discerned despite extensive heterogeneity in the methylome data, and this indicates that the epigenetic response to opioids occurs within the first few days to weeks following drug exposure. Additionally, we also performed an epigenome-wide association study (EWAS), and this identified a few CpGs that may be predictive of opioid dosing.

Results

The number of enrolled participants (N = 41) and the timeline of sample collection are shown in Fig. 1. Only 33 patients (19 females) received prescription opioids after an oral procedure. The baseline characteristics, other diagnosed diseases, casual use of other drugs (specifically tobacco and marijuana; no participant reported the use of cocaine, psychedelics, and other hard drugs), and prescription opioid self-administration are reported only for these 33 participants (Table 1). Following the pre-surgery visit (visit 1 or v1), the second visit (visit 2 or v2) occurred after surgery and within a week of the last opioid dose (average number of days between last opioid dose and visit 2 was 2.7 ± 1.5 days). The last sample collection (visit 3 or v3) occurred between 32 and 88 days from surgery, and the number of days between the last opioid dose and visit 3 was 39 ± 10 days. In total, 26 participants provided saliva samples at all three visits, 6 participants provided saliva at two visits, and one provided saliva only at v1 (Table 1). The mean age was 33.61 ± 13.84 years and ranged from 19 to 61 years (Table 1). Based on self-reported race/ethnicity, there were 13 Caucasians (mean age = 31.69 ± 14.11 years), 13 African-Americans (mean age = 39.92 ± 13.91), and the remaining 7 were of “other” racial/ethnic group (mostly Hispanic/Latino; mean age = 25.43 ± 8.02). The African-American group was slightly older, but there was no statistically significant difference in age between the groups (p value = 0.06). Sex distribution was not significantly different between the race/ethnic groups. Individual-level information, including comorbidities, is provided as Additional file 1: Table S1.

Fig. 1
figure 1

Timeline of sample collection. Saliva samples were collected before surgery and the two follow-up visits after surgery and end of opioid self-administration. The notations above the arrows show the range of days between events

Table 1 Participant characteristics and postoperative opioid use

Postoperative opioid dosing data was based on self-reported pill counts converted to morphine milligram equivalent (MME). With the exception of one individual who used no opioids (and we considered this individual to represent a dose of 0 MME with 0 days of use), all patients started opioid treatment generally within 24 h of surgery and continued use for an average of 6 ± 4 days for up to 17 days (Additional file 1: Table S1). As expected, cumulative dosage correlated with length of use (r = 0.67, p value < 0.0001). For the 32 participants that self-administered opioids, the total cumulative dosage over the course of treatment ranged from 5–210 MME. Based on the quantile distribution of the cumulative MME, participants were classified into three groups: < 25 MME (those below the 25th percentile or quartile 1 for opioid dosage), 25–90 MME (those within the interquartile range), and ≥ 90 MME (those above the 75th percentile or quartile 3) (Table 1). Opioid dosage showed no significant association with age, sex, and self-reported race/ethnicity. There was no significant association between opioid dosage and the presence or absence of other comorbidities. Dosage was also not associated with past marijuana use. However, the group that reported using tobacco within the past 12 months had significantly higher self-administered opioid dosage (mean of 95.63 ± 57.78 MME among tobacco users, and 49.55 ± 36.19 MME among non-users; p value = 0.008).

Global shift in postoperative methylome

For an overview of the methylome and the variance structure, we started with a principal component analysis (PCA) using the full set of high-quality probes (736,432 probes passed QC criteria). The top PC (PC1) captured a vast portion of the variance at 63.5%, and following that, PC2 and PC3 captured only 2.5% and 1.6% of variance, respectively (PCs for each methylome data in Additional file 1: Table S2). PC1 was not significantly associated with the demographic variables (sex, age, self-reported race/ethnicity) or with comorbidities and past use of marijuana or tobacco. Instead, the visit was the most significant explanatory variable for PC1 (F2,86 = 5.94, p value = 0.004), and the pattern indicated a significant change in the methylome with the strongest contrast between v3 and v2 (Tukey-Kramer post hoc p value = 0.003) (Fig. 2a). To deduce whether the longitudinal variance capture by PC1 could be explained by the length of time from surgery or opioid self-administration, we performed bivariate analyses between PC1 and the following variables: opioid dose, days from surgery to sample collection, and days from last opioid self-administration to sample collection. This analysis was done for the three visits separately, and at v2, PC1 had a modest but significant correlation with days from surgery to v2 (r = 0.40, p value = 0.03, n = 31 participants with methylome data at v2; Fig. 2b). Similarly, at v3, PC1 was correlated with days from surgery to v3 (r = 0.41, p value = 0.04, n = 27 participants with methylome data at v3). PC1 was not correlated with opioid dose or the number of days from the last opioid use. From this, we can infer that the longitudinal shift in the methylome is primarily due to surgery.

Fig. 2
figure 2

Global patterns in DNA methylation across visits. a The top principal component (PC1) extracted from the methylome-wide data explained 63.5% of the variance, and the ANOVA plot shows significant differences between the three visits (F2,86 = 5.94, p value = 0.004). b At visit 2, PC1 is correlated with the number of days from surgery to the second visit (r = 0.40, p value = 0.03, n = 31). c The epigenome-wide association plot depicts the location of each CpG (autosomal chromosomes 1 to 22 on the x-axis) and the –log10(p value) for the effect of visit (y-axis). The genome-wide significant threshold was set at p value = 5 × 10-8 (upper red horizontal line); the suggestive threshold was set at p value = 10-5 (lower blue horizontal line). d Distribution of p values for the effect of visit shows a significant deviation from the null hypothesis. e For the CpGs above the suggestive threshold, comparison of mean differences between visit 1 and visit 2 (x-axis), and visit 3 and visit 2 (y-axis) indicates a reversal in methylation patterns from visit 2 to visit 3, with the majority of sites showing lower methylation at visit 2, and then increasing in methylation by visit 3

To profile the CpGs that changed longitudinally over the three visits, we performed a mixed-effects ANOVA with a visit as a fixed variable and the person ID as a random effect (Fig. 2c). The p values for visit showed a significant deviation from the null hypothesis (Fig. 2d histogram). However, only 2 intergenic CpGs (cg05639411 and cg24904009) were above the genome-wide significant threshold of 5.0e−8 (Fig. 2c), and overall, the pattern indicated the modest shift in the methylome across several CpGs. At a genome-wide suggestive threshold of p value = 1.5e−5, there were 1701 CpGs that underwent change over the visits (Additional file 1: Table S3). The majority of these CpGs (> 65%) decreased in methylation between v1 and v2 and regained methylation by v3 such that these sites showed significantly higher levels of methylation at v3 compared to both v1 and v2 (Fig. 2e). Similarly, for the ~ 35% of CpGs that gained methylation between v1 and v2, these sites generally declined in methylation by v3 resulting in significantly lower methylation compared to both v1 and v2 (Fig. 2e). Gene set enrichment analysis (GSEA) of the 1133 annotated genes represented by the CpGs conveyed mostly an innate immune inflammatory response (Additional file 1: Table S4). The most overrepresented pathway was natural killer cell-mediated cytotoxicity (KEGG ID hsa04650; normalized enrichment score = − 1.93, FDR = 0.03), and the most overrepresented function was for genes involved in cellular defense response (GO ID 0006968; normalized enrichment score = − 1.83, p = 0.001, FDR = 0.3), and these immunity-related categories were enriched among the CpGs that decreased in methylation at v2. The opioid receptors were not represented in the list of visit-associated CpGs. Based on these observations, a possible explanation for the shift in the methylome is that it is the result of surgery-induced immune response and changes in the oral cell composition. Opioid use, if it had an impact, is likely to exert a weaker signal and, given the limited sample size, more suitable for a focused candidate gene study.

Deconvolution of cellular heterogeneity

To decompose cell types from the composite DNA methylation signal, we applied a reference-free approach [38]. The bootstrapping method described in Houseman et al. [38] determined K = 4 cell types (Additional file 1: Table S2). Cell 1, which represented the most abundant cell type, showed an increase at v2 right after surgery followed by a decline by v3 (Table 2). Aside from cell 1, no other cell showed significant change over the visits (Table 2). To deduce what cell types are represented by the 4 groups, we also estimated blood leukocyte proportions (mainly lymphocytes and granulocytes/neutrophils) using a reference-based approach [39] and compared correlations between the 4 cell types to the reference-based cell estimates (Table 2; Additional file 1: Table S2). Cell 1 had a strong positive correlation with granulocytes, and cell 4 had a strong positive correlation with lymphocytes indicating that cells 1 and 4 are chiefly representative of the leukocyte population in saliva and serve as a proxy for the increase in granulocyte proportions after surgery. Cells 2 and 3 had only modest correlations with leukocyte estimates (|r| of 0.4–0.5) and may be more representative of the epithelial cells.

Table 2 Reference-free and reference-based estimates of cellular proportions

The cell estimates were not associated with opioid dose. To evaluate if the baseline characteristics were related to cellular composition, we tested associations with age, sex, and race/ethnicity. Cell 1 had a significant negative correlation with age (r = − 0.48, p value = 0.006) only at v2 that suggests an age-dependent immune response in the days immediately after surgery (Fig. 3a). Cell 3 had the strongest association with age at all three visits (Fig. 3b). Cells 2 and 3 showed extensive cross-sectional variability without longitudinal change, and both were significantly associated with race/ethnicity at all three visits, indicating that these could serve as proxies for the cellular composition differences between populations (Fig. 3c, d). Cell 4 was not associated with any of the baseline variables, and sex was not a factor for any of the cell types.

Fig. 3
figure 3

Estimated cell-type proportions and associated variables. a Cell 1 shows both longitudinal and cross-sectional variability. The proportion of cell 1 is negatively correlated with age at visit 2 (r = − 0.48, p value = 0.006, n = 31; black squares and dashed line), but not at visit 1 (r = − 0.13, p value = 0.49, n = 31; red × markers and dotted line), and only slightly at visit 3 (r = − 0.32, p value = 0.10, n = 27; gray circles and solid line). b Cell 3 is associated with cross-sectional variability but no significant longitudinal change. The estimated proportion has a strong positive correlation with age at all three visits. At visit 1, r = 0.36 (p value = 0.05); visit 2, r = 0.57 (p value = 0.0009); visit 3, r = 0.48 (p value = 0.01). c Cell 3 also shows a significantly higher proportion in African-Americans at all three visits (F2,28 = 15.66, p value < 0.0001 at visit 1). d Cell 2 is also ethnicity-specific and associated with a lower proportion in African-Americans at all three visits (F2,28 = 4.77, p value < 0.02 at visit 1)

Effect of opioid dose on OPRM1 promoter methylation

To examine if the higher opioid dose is related to higher promoter methylation, we started with a candidate gene approach and focused on the CpGs located in the OPRM1 promoter. In total, 10 promoter CpGs were targeted by the Illumina probes and these encompassed the CpG island described by Nielsen et al. and replicated by Chorbov et al. [29, 30] (Fig. 4a, Table 3; individual level β values in Additional file 1: Table S2). We first applied a mixed regression model with an opioid dosage group and visit as fixed categorical variables, and each participant ID as a random intercept. With the exception of the last CpG, the regression estimates for all the OPRM1 promoter CpGs were positive, with higher methylation levels for the two higher dosage groups (i.e., 25–90 MME and ≥ 90 MME) relative to the lowest dosage group (< 25 MME) (Table 3). At a nominal p value of 0.05, 4 CpGs were significantly associated with opioid dosage groups. The ANOVA plots for these CpGs showed that the difference between dosage groups was pronounced at v3 (for CpG1, CpG2, CpG6) and v2 (for CpG7) but not at v1 (Fig. 4a). As tobacco use was associated with a higher self-administered dosage of opioids, we considered it as a potential contributing factor. However, including the past tobacco use in the regression model did not alter the results, and this indicated that the higher methylation at the OPRM1 CpGs is a specific effect of opioids.

Fig. 4
figure 4

OPRM1 promoter CpG methylation. a The OPRM1 promoter and the CpG island (green block) are depicted along with the base pair coordinates (black line, GRCh37/hg19) and location of the 10 CpGs (filled circles). Residual β values were extracted after fitting participant ID as a random intercept, and the plots show the methylation patterns across the three visits for CpG1, CpG2, CpG6, and CpG7 (panels with ANOVA line plots; error bars are standard error). The difference between the dosage groups appears at visit 3 (for CpG1, CpG2, CpG6) and visit 2 (for CpG7). The lowest cumulative dose group (< 25 MME: blue dotted line) has lower average methylation compared to the two higher cumulative dose groups (25–90 MME, yellow dashed line; ≥ 90 MME, red solid line). b The promoter mean methylation was taken as the average β values for CpG1 to CpG9. After fitting a regression model with adjustment for age and cell proportions, the leverage plots show a significantly higher average promoter methylation (y-axes) associated with higher MME (x-axis, left panel) and longer duration of use (x-axis, right panel). MME is morphine milligram equivalent

Table 3 Dose-dependent methylation of individual OPRM1 promoter CpGs

To check whether the association with opioid dosage can be robustly detected, we summarized the overall methylation pattern in the promoter by taking the mean DNA methylation β values for the nine CpGs that were positively associated with opioid dosage (CpG1 to CpG9). We applied a linear regression model and tested whether higher mean methylation was associated with either higher MME or longer length of opioid use. This analysis was done for the three visits separately and adjusted for age, tobacco use, and cellular heterogeneity. Both MME and length of opioid use were associated with higher mean methylation, but this effect was significant only at v3, further indicating that the hypermethylation of the OPRM1 receptor is more likely a response rather than a predisposing factor (Fig. 4b, Table 4). Our results are consistent with the opioid-associated hypermethylation and indicate that even a relatively short-term opioid use may induce an increase in methylation that is proportional to dosage at the OPRM1 promoter.

Table 4 Mean methylation in the OPRM1 promoter and association with opioid dose and days of use

Preliminary epigenome-wide association study for the opioid dose

Since the candidate gene approach indicated that the short-term use of prescribed opioids can have an impact on CpG methylation, we expanded the analysis to an EWAS using the same mixed model to test association with opioid dosage. A CpG (cg08105965) located in the promoter CpG island of the GTPase signaling gene, RASL10A (RAS like family 10 member A), was genome-wide significant (p value of 5.0e−8, Fig. 5a). Unlike the pattern for the OPRM1 promoter, the lowest dose group had a significantly lower methylation level at both visits 1 and 3, indicating that the difference preceded opioid use (Fig. 5b). In addition to the RASL10A promoter CpG, 5 other CpG sites were associated with opioid dosage at the suggestive threshold (p value of 1.0e−5, Fig. 5c–g). For most of these, the methylation differences were apparent at v1 and preceded opioid use. For these top CpGs, adjusting for past tobacco use did not alter the results, and none of these sites were significantly associated with tobacco use. To explore if any of the CpGs that were above the suggestive threshold have been previously implicated in opioid use or dependence, we referred to recent human EWAS for opioid dependence [24] and methadone treatment dosage [28]. Based on the comparison of probe IDs and genes, none of the CpGs we report here have been previously linked to opioid-related traits.

Fig. 5
figure 5

Epigenome-wide test for prescription opioid dosage. a The epigenome-wide Manhattan plot depicts the location of each CpG (autosomal chromosomes 1 to 22 on the x-axis) and the –log10(p value) for the effect of doage (y-axis). The genome-wide significant threshold was set at p value = 5 × 10-8 (upper red horizontal line); the suggestive threshold was set at p value = 10-5 (lower blue horizontal line). For the six CpG sites that were above the genome-wide suggestive threshold, residual β values were extracted after fitting participant ID as random intercept, and the plots show the methylation patterns across the three visits (error bars are standard error) for b cg08105965 (RASL10A), c cg18500286 (AFF1), d cg04718304 (intergenic), e cg16719801 (VSNL1), f cg08163714 (ANXA2), and g cg25124965 (PAIP2). For most of these, the difference between the dosage groups appears at visits 1 and 3

While this is preliminary results from a small study cohort and is yet to be replicated, we provide the list of 64 CpGs that were associated with opioid dosage at a nominal uncorrected p value of 1.0e−4, along with the gene ontology IDs and KEGG pathways for the corresponding genes in Additional file 1: Table S5.

Discussion

Here, we report results from a longitudinal study of DNA methylation in a cohort of opioid naïve dental patients who received prescription opioids following oral surgery. To summarize the main result, we found increased methylation at the OPRM1 promoter associated with a higher cumulative opioid dose. This replicates the hypermethylated profile among long-term opioid users and alcohol-dependent individuals [23, 29,30,31,32,33,34]. The pattern of methylation we observed indicates that the increase in methylation is more likely the response to, rather than the cause of, opioid use [20]. The present study provides evidence that such epigenetic modifications are induced within the early days of drug use and may represent early epigenomic responses to an addictive substance.

A peculiar challenge we faced was that the site of sample collection was also the site of surgery. Saliva has a highly heterogeneous cellular makeup and is estimated to constitute about ~ 45% epithelial cells and about ~ 55% leukocytes from the circulating blood [40]. The main goal of the study was to detect the effect of short-term and comparatively low-dose opioids, while accounting for the larger perturbation caused by surgery. Although we do not have details on the severity of the oral surgery, most were third molar extractions and were relatively minor and non-invasive. Nonetheless, the patients would have experienced an injury-induced inflammatory response that can result in changes in numbers of circulating immune cells [41], and consequently, changes in oral cell composition. As DNA methylation is highly cell-type specific, the heterogeneity in cells will be a major source of “noise” in the methylome data [39, 42,43,44]. The longitudinal variability in DNA methylation that was captured by the top PC can therefore be attributed to cell composition rather than opioid use. We could deduce this by the significant correlation between PC1, and the number of days from surgery to the follow-up visits. We were able to partly resolve the cell heterogeneity by applying reference-free and reference-based estimates of cell proportions. The reference-free method estimated four major cell types. Although saliva is highly heterogeneous, and certainly has more than just 4 types of cells [40], the classification into 4 broad groups likely reflects the limitation in the in silico approach to resolve finer differences between cellular subtypes. Cell 1 most likely represented the granulocyte population (chiefly neutrophils), which constitutes the most abundant leukocyte subtype in circulating blood, and is responsible for innate immunity and acute inflammatory response. Consistent with the known increase in granulocyte-to-lymphocyte ratio in the few days following surgery [45], we also found an increase in cell 1 and in the relative abundance of granulocytes compared to lymphocytes at visit 2. This was followed by a compensatory decrease in granulocyte proportions by visit 3. Cell 2 and cell 3 are presumed to represent a portion of the epithelial cell population, and these showed no significant within-individual changes over the visits. However, these cells exhibited a significant association with age and self-reported race/ethnicity. Although cell type decomposition was not the primary objective of the study, our analyses demonstrated that the saliva methylome can be highly informative of individual differences in perioperative immune profiles.

For the effect of postsurgical opioid use, we first focused on the OPRM1 promoter region as an epigenetic sensor of opioid dose. The CpG-rich promoter harbors a CpG island, and several studies in different populations have demonstrated higher DNA methylation at this site among opioid users and methadone-maintained heroin addicts [23, 29,30,31,32,33]. The increased methylation of the OPRM1 promoter is not only limited to OUD but has also been detected among individuals with alcohol dependence, suggesting that the hypermethylation is generally associated with substance use disorder and addiction [34]. A question has been whether such epigenetic differences are the result of drug use or the cause of increased vulnerability to addiction [20]. To address this, we interrogated 10 CpGs in a 550-bp region that encompassed the promoter CpG island investigated by Nielson et al. and Chorbov et al. (the CpG island is depicted in Fig. 4) [29, 30]. With the exception of the last CpG, the remaining 9 CpGs showed higher methylation in the two higher-dose groups relative to the low-dose group, and four of these CpGs were significantly associated with dosage at a nominal alpha of 0.05. Comparison of mean methylation differences between the dosage groups across the three visits indicated that higher methylation in the higher dose groups is more apparent at the postsurgery visits, particularly visit 3. The positive association between the mean promoter methylation and cumulative MME, and mean promoter methylation and days of opioid use, were also significant only at visit 3. The heightened inflammatory state at visit 2, which occurred within a few days of surgery, may have been the reason why more subtle effect of opioids was not significant at visit 2, and the positive association emerged only at visit 3.

The OPRM1 locus presents a prime site for gene x environment interaction, a critical aspect of addiction since the addictive substance is an environmental agent that has a long-lasting biological effect. The OPRM1 gene has been the subject of several candidate gene studies for addiction. Much attention has been paid to the missense SNP that alters the OPRM1 protein function, although its impact on addiction traits and OUD is somewhat ambiguous [46, 47]. Several studies have also identified non-coding variants in the OPRM1 locus that alters DNA methylation and gene expression [33, 35, 48]. At least one genome-wide association study has also identified a genome-wide significant association between a SNP upstream of OPRM1 and methadone-maintenance dosing [37]. These studies collectively provide evidence that common genetic variants in the proximal region of OPRM1 affect DNA methylation and gene expression and could have a downstream impact on an opioid response that could potentially influence vulnerability to addiction. Our present work was carried out in a small sample size, and our primary goal was to track the within-individual trajectory across visits. If there were genetically modulated small cross-sectional differences at baseline, this sample size would be underpowered to detect the differences, and the significant association with the opioid dose that we found may have been the result of opioid-induced augmentation of differential methylation at the postoperative visits.

The hypermethylation of the OPRM1 promoter is likely only a small part of a larger network of genes involved in the cellular response to drug exposure. We therefore followed-up with a preliminary EWAS exploration to identify other CpGs that may be associated with opioid self-administration. To our surprise, despite the small sample size, one CpG, located in the promoter region of RASL10A, a RAS-related GTPase signaling gene, was genome-wide significant. Perhaps this is due to the power of the longitudinal design in capturing differentially methylated sites that are significantly different between dosage groups at more than one visit. For instance, the differential methylation of cg08105965 at RASL10A is apparent at both visits 1 and 3. Although RASL10A has not been previously implicated in opioid response or addiction, it is notable that the OPRM1 protein is a G-protein-coupled receptor, and its activation results in a cellular signaling cascade that also involve RAS GTPase activity [49,50,51]. In addition to the RASL10A CpG, five other CpGs were at or above the suggestive threshold, including sites located in AFF1, VSNL1, ANXA2, and PAIP2. To our knowledge, DNA methylation at these genes has not been previously linked to opioid use or dependence. However, one recent study of gene expression in the rat model has shown an upregulation of RASL10A and VSNL1 in the brain following acute morphine treatment [52]. Similar to RASL10A, VSNL1 (visinin-like 1) also codes for an intracellular signaling molecule with high expression in the brain [53]. The list of CpGs that were associated with opioid dosage at a nominal uncorrected p value < 0.0001 included few other cellular signaling genes (e.g., ANXA2, RET, ADRB1 ). Taken together, the EWAS results hint that epigenetic modulation of genes involved in intracellular signal transduction may play a role during the early phase of opioid use.

We must emphasize that the small sample size and the heterogeneity in methylome signal, partly due to cell composition and partly due to the heterogeneous population group, are major limitations, and the EWAS results await replication in an independent cohort. The CpGs identified by the present EWAS were differentially methylated even at v1, prior to opioid use, and this suggests that there may be genetic variants underlying these epigenetic differences. However, such potential effects of genetic variation are not addressed in the present study due to the lack of genotype data. Another weakness that we should note is that the main variable of interest, therapeutic opioid dosing, was based on patient self-reports rather than objective measures of drug use [54]. A future strategy would be to use existing technologies such as wearable devices that can provide additional means of tracking the physiological responses to opioids [55]. The present study also does not address whether these epigenetic changes linger or diminish over time in the absence of continued drug use. A more comprehensive longitudinal epigenomic study of the early effects of prescription opioids that also integrates genetic effects and with a longer follow-up period would be the next phase of the study.

Regarding the potential for epigenetic persistence, we must point out that any peripheral tissue serves only as a proxy for the possible epigenetic changes in the brain. A distinction is that the blood and epithelium are mitotically active tissues, and cells are renewed within a few days to a few weeks, with the exception of long-lived memory T cells. For methylation signals to persist, it will require either continued presence of the perturbation (i.e., continued exposure to opioids) or methylation changes in mitotically active stem cells that can be faithfully transmitted to the daughter cells. The brain, on the other hand, is mitotically inactive and consists of mostly terminally differentiated cells that last a lifetime. If the relatively modest dose and short-term use of prescription opioids have a similar impact on the brain cells, the effects may not readily decay and may be long-lasting in the central nervous system.

Conclusion

In conclusion, our study replicates the hypermethylation of the OPRM1 promoter with the opioid use. Previous studies reported on the effects of chronic opioid use; here, we provide evidence that the epigenetic restructuring begins within the initial stage of opioid exposure. The present findings on the acute effects of prescription opioids, as well as the CpGs that may be predictive of opioid dosing, require further replication with a well-powered and more comprehensive study in a larger cohort.

Methods

Participants

Eligible participants were scheduled for tooth extractions, mostly third molar extractions, at an oral and maxillofacial surgery clinic that were typically followed by postoperative prescriptions of hydrocodone/acetaminophen (7.5/325 mg q4–6h prn pain) or oxycodone/acetaminophen (5 mg/325 mg q6h prn pain). For inclusion in the study, individuals were required to be 18 years of age or older, opioid naïve, able to consent, able to understand and speak English, and willing to provide saliva samples. Individuals were excluded if they reported previous use of opioids, had current substance use dependence, were pregnant, were incarcerated, had other causes of pain, were unable to consent, or had a developmental disability that prevented participation. The study received approval by the university Institutional Review Board. Eligible participants were provided a summary of the consent form by the study coordinator and allowed to read and ask questions before enrollment. All participants provided written informed consent.

Forty-one individuals consented to the study and provided contact information and responded to a demographic questionnaire. The enrolled participants also provided information on existing diagnosed diseases and were assessed for casual substance use within the past 12 months (tobacco, marijuana, cocaine, psychedelics, other hard drugs; details in Additional file 1: Table S1). The clinical staff provided routine opioid medication and recovery instructions for all participants right after surgery. The opioids prescribed to participants were hydrocodone, oxycodone, and codeine in doses that varied between 5 mg and 30 mg (Table 1). The study coordinator also provided opioid medication logs to participants to record self-administration including the date, time, individual dose per opioid pill, and number of opioid pills taken. For the 33 participants who received opioid medication, only one participant (person ID 142) reported no opioid usage.

Sample processing and DNA methylation assay

The saliva was collected using the Oragene DNA sample collection kit (OGR 500) by DNA Genotek (http://www.dnagenotek.com). The first set of samples was collected before surgery. The second saliva sample was collected a few days after opioid discontinuation, and the third was collected on a follow-up visit (Fig. 1, Table 1). DNA was purified using the DNA Genotek PrepIT L2P kit according to the manufacturer’s instruction. Genome-wide DNA methylation was assayed on the Illumina Infinium Human MethylationEPIC BeadChips following the manufacturer’s standard protocol at the HudsonAlpha Genomic Services Lab (https://gsl.hudsonalpha.org).

Data processing

Raw intensity IDAT files were loaded to R, and all quality checks, data preprocessing, and normalization were carried out using the R package, minfi (v.1.31) [56]. Methylation levels were estimated as β values (ratio of methylated by unmethylated probes) and quantile normalized. The initial QC involved a comparison between the log median intensities of methylated and unmethylated channels, and the density plots for β values (Additional file 2: Fig. S1a). All samples passed these checks. Sex estimated from the DNA methylation data also matched the self-reported sex. To retain only high-quality data, probes with detection p value > 0.01 (14,676 probes) were excluded. Probes that target CpGs on the sex chromosomes were also removed (18,605 probes). Finally, a total of 96,146 probes that overlapped annotated SNPs and/or were flagged for poor mapping quality (MASK.general list from [57]) were also filtered out. A total of 736,432 high-quality probes were retained and used for downstream analysis.

As further QC, we performed unsupervised hierarchical clustering using the full set of high-quality probes (Additional file 2: Fig. S1b). While samples longitudinally collected from the same individual tended to cluster together, there were also several samples that did not cluster with self. To check for possible errors in sample labeling, we repeated the cluster analysis using a subset of 30,435 probes that had been filtered out due to overlap with common SNPs. While these were deemed poor-quality probes and unfit for differential methylation analysis, in terms of sample identity check, these probes can serve as proxy genotype markers that can help verify if the samples came from the same person. Using this set, almost all samples collected from the same participant clustered with self, and for the most part, the clusters also aligned with self-reported race/ethnicity groups (Additional file 2: Fig. S1c). Only one of the 33 participants who received prescription opioids (person ID 108) did not pair with self, and data from this person were excluded from all downstream analyses.

Estimation of cellular proportions

To infer the relative proportions of the major cell types, we first implemented a reference-free deconvolution of the methylome data using the R package RefFreeEWAS (v2.2) [58]. The RefFreeEWAS algorithm applies a non-negative matrix factorization to decompose a matrix Y = , where M represents an m × K matrix with m as CpG specific methylation for an unknown number of K cell types, and Ω as the cell-type proportion constrained to sum to a value ≤ 1. For computational efficiency, the K cell types have to be first specified, and as described in Houseman et al. [58], we set the K to vary from 2 to 10 cell types and decomposed the Y = MΩ. Following this, we applied bootstrapping to estimate the optimal K value. For this estimation, we applied 10 iterations with replacement every 1000 times. The optimal K = 4 was selected based on the minimum value of the average of bootstrapped deviances for each putative cell type. While this method provides the relative proportions of cell types, the identity of the four cells is unknown. Since a significant proportion of saliva consists of leukocytes, we also applied a reference-based approach to estimate the relative proportions of lymphocytes and neutrophils [39, 42, 43]. To infer the putative identities of cells, we performed Pearson correlations between the K cells and the proportions of leukocyte types.

Statistical analyses

For the global analysis, PCA was done on the full set of 736,432 probes using the prcomp function in R. In order to evaluate which variable had the most significant association with PC1, we examined the association between PC1 and the following variables: sex, age, self-reported race/ethnicity, presence or absence of other diseases, tobacco or marijuana usage in the past 12 months, and opioid dose, days from surgery to sample collection, and days from last opioid dose to sample collection. We used ANOVA for categorical variables and Pearson correlation for continuous variables, and these tests were conducted separately for the three visits. We also performed similar analyses for estimated cell proportions to examine whether the variables were significantly related to the cell proportions. PC1 and the cell proportions were also related to visit using ANOVA. Since the visit was the most significant explanatory variable for PC1, we identified the CpGs that showed longitudinal change over the three visits by applying a mixed-effects ANOVA: aov(β value ~ visit + Error(ID/visit)). This epigenome-wide analysis was done for the 26 participants with data from all 3 visits. For the set of genome-wide suggestive CpGs that changed over the visits (uncorrected p value ≤ 10-5), GSEA was implemented on the WebGestalt platform (http://www.webgestalt.org) with each CpG ranked by the mean β value difference between v2 and v1.

For candidate gene analysis, we surveyed the promoter region of OPRM1. The CpG island that was interrogated by Nielson et al. and Chorbov et al. is located at 154360587–154360922 bp of chromosome 6 (GRCh37/hg19) [29, 30]. Within that exact coordinate, our data only had 4 CpG probes. We therefore considered a slightly wider region (550 bp), and in total, the array data contained 10 probes that targeted promoter CpGs at chr6:154360344-154360894 bp. To evaluate methylation at individual CpGs, we applied a linear mixed-effects model with dosage group and visit as fixed categorical variables, and person ID as random intercept: lmer(β value ~ dosage + visit + (1|ID)). To test if the history of tobacco use could account for some of the effects, we repeated the test with the model: lmer(β value ~ MME + tobacco-use + visit + (1|ID)). This was done using the “lmertest” R package, and to get the p values for the main effect of dosage groups, the degrees of freedom were computed by the Satterthwaite’s method [59, 60]. Following the CpG level analysis, we estimated the general methylation trend for the promoter by averaging the β values for the 9 CpGs that had a positive regression coefficient with the dosage groups. We then tested the association between the promoter mean methylation score, and two opioid-related continuous variables: length of opioid use in days and cumulative MME. This analysis was done for the three visits separately and adjusted for age and cellular heterogeneity using the equations lm(mean-β ~ MME + tobacco-use + age + cell1 + cell2 + cell3), and lm(mean-β ~ days-of-use + tobacco-use + age + cell1 + cell2 + cell3).

Following the candidate gene study, we then performed an EWAS for opioid dosage using the same mixed-effect model: lmer(β value ~ MME + visit + (1|ID)), and for the CpGs identified by the EWAS at above the genome-wide suggestive threshold, we checked for the effect of past tobacco use using the model: lmer(β value ~ MME + tobacco-use + visit + (1|ID)).

Availability of data and materials

The full de-identified raw DNA methylation data are available from the NCBI NIH Gene Expression Omnibus repository (GEO accession ID: GSE151485).

Abbreviations

EWAS:

Epigenome-wide association study

FDR:

False discovery rate

GO:

Gene ontology

GSEA:

Gene set enrichment analysis

KEGG:

Kyoto encyclopedia of genes and genomes

MME:

Morphine milligram equivalent

OPRM1:

Opioid receptor mu 1

OUD:

Opioid use disorder

PC:

Principal component

PCA:

Principal component analysis

QC:

Quality control

SNP:

Single-nucleotide polymorphism

References

  1. Porter J, Jick H. Addiction rare in patients treated with narcotics. N Engl J Med. 1980;302(2):123.

    CAS  PubMed  Google Scholar 

  2. Rummans TA, Burton MC, Dawson NL. How good intentions contributed to bad outcomes: the opioid crisis. Mayo Clin Proc. 2018;93(3):344–50.

    Article  PubMed  Google Scholar 

  3. Administration SAaMHS: Results from the 2017 National Survey on Drug Use and Health: Detailed Tables. 2017.

  4. Understanding the epidemic [https://www.cdc.gov/drugoverdose/epidemic/index.html].

  5. Edlund MJ, Martin BC, Russo JE, DeVries A, Braden JB, Sullivan MD. The role of opioid prescription in incident opioid abuse and dependence among individuals with chronic noncancer pain: the role of opioid prescription. Clin J Pain. 2014;30(7):557–64.

    PubMed  PubMed Central  Google Scholar 

  6. National Academies of Sciences E, and Medicine: Pain management and the opioid epidemic: balancing societal and individual benefits and risks of prescription opioid use. Washington, DC: The National Academies Press; 2017.

  7. Schroeder AR, Dehghan M, Newman TB, Bentley JP, Park KT. Association of opioid prescriptions from dental clinicians for US adolescents and young adults with subsequent opioid use and abuse. JAMA Intern Med. 2019;179(2):145–52.

    Article  PubMed  Google Scholar 

  8. Shah A, Hayes CJ, Martin BC. Characteristics of initial prescription episodes and likelihood of long-term opioid use — United States, 2006–2015. Mmwr-Morbid Mortal W. 2017;66(10):265–9.

    Article  Google Scholar 

  9. Denisco RC, Kenna GA, O'Neil MG, Kulich RJ, Moore PA, Kane WT, Mehta NR, Hersh EV, Katz NP. Prevention of prescription opioid abuse. The role of the dentist. J Am Dent Assoc. 2011;142(7):800–10.

    Article  PubMed  Google Scholar 

  10. Rasubala L, Pernapati L, Velasquez X, Burk J, Ren YF. Impact of a mandatory prescription drug monitoring program on prescription of opioid analgesics by dentists. PLoS One. 2015;10(8):e0135957.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  11. Tong ST, Hochheimer CJ, Brooks EM, Sabo RT, Jiang V, Day T, Rozman JS, Kashiri PL, Krist AH. Chronic opioid prescribing in primary care: factors and perspectives. Ann Fam Med. 2019;17(3):200–6.

    Article  PubMed  PubMed Central  Google Scholar 

  12. Alam A, Gomes T, Zheng H, Mamdani MM, Juurlink DN, Bell CM. Long-term analgesic use after low-risk surgery: a retrospective cohort study. Arch Intern Med. 2012;172(5):425–30.

    Article  PubMed  Google Scholar 

  13. Brummett CM, Waljee JF, Goesling J, Moser S, Lin P, Englesbe MJ, Bohnert ASB, Kheterpal S, Nallamothu BK. New persistent opioid use after minor and major surgical procedures in US adults. JAMA Surg. 2017;152(6):e170504.

    Article  PubMed  PubMed Central  Google Scholar 

  14. Cicero TJ, Ellis MS, Surratt HL, Kurtz SP. The changing face of heroin use in the United States: a retrospective analysis of the past 50 years. JAMA Psychiatry. 2014;71(7):821–6.

    Article  PubMed  Google Scholar 

  15. Sun EC, Darnall BD, Baker LC, Mackey S. Incidence of and risk factors for chronic opioid use among opioid-naive patients in the postoperative period. JAMA Intern Med. 2016;176(9):1286–93.

    Article  PubMed  PubMed Central  Google Scholar 

  16. Nestler EJ: Epigenetic mechanisms of drug addiction. Neuropharmacology 2014, 76 Pt B:259-268.

  17. Browne CJ, Godino A, Salery M, Nestler EJ. Epigenetic mechanisms of opioid addiction. Biological psychiatry. 2020;87(1):22–33.

    Article  CAS  PubMed  Google Scholar 

  18. Flagel SB, Chaudhury S, Waselus M, Kelly R, Sewani S, Clinton SM, Thompson RC, Watson SJ Jr, Akil H. Genetic background and epigenetic modifications in the core of the nucleus accumbens predict addiction-like behavior in a rat model. Proceedings of the National Academy of Sciences of the United States of America. 2016;113(20):E2861–70.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Koo JW, Mazei-Robison MS, LaPlant Q, Egervari G, Braunscheidel KM, Adank DN, Ferguson D, Feng J, Sun H, Scobie KN, et al. Epigenetic basis of opiate suppression of Bdnf gene expression in the ventral tegmental area. Nature neuroscience. 2015;18(3):415–22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Nielsen DA, Utrankar A, Reyes JA, Simons DD, Kosten TR. Epigenetics of drug abuse: predisposition or response. Pharmacogenomics. 2012;13(10):1149–60.

    Article  CAS  PubMed  Google Scholar 

  21. Kozlenkov A, Jaffe AE, Timashpolsky A, Apontes P, Rudchenko S, Barbu M, Byne W, Hurd YL, Horvath S, Dracheva S: DNA methylation profiling of human prefrontal cortex neurons in heroin users shows significant difference between genomic contexts of hyper- and hypomethylation and a younger epigenetic age. Genes 2017, 8(6).

  22. Egervari G, Landry J, Callens J, Fullard JF, Roussos P, Keller E, Hurd YL. Striatal H3K27 acetylation linked to glutamatergic gene dysregulation in human heroin abusers holds promise as therapeutic target. Biological psychiatry. 2017;81(7):585–94.

    Article  CAS  PubMed  Google Scholar 

  23. Doehring A, Oertel BG, Sittl R, Lotsch J. Chronic opioid use is associated with increased DNA methylation correlating with increased clinical pain. Pain. 2013;154(1):15–23.

    Article  CAS  PubMed  Google Scholar 

  24. Montalvo-Ortiz JL, Cheng Z, Kranzler HR, Zhang H, Gelernter J. Genomewide study of epigenetic biomarkers of opioid dependence in European-American women. Sci Rep. 2019;9(1):4660.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  25. Bleich S, Lenz B, Ziegenbein M, Beutler S, Frieling H, Kornhuber J, Bonsch D. Epigenetic DNA hypermethylation of the HERP gene promoter induces down-regulation of its mRNA expression in patients with alcohol dependence. Alcohol Clin Exp Res. 2006;30(4):587–91.

    Article  CAS  PubMed  Google Scholar 

  26. Bonsch D, Lenz B, Reulbach U, Kornhuber J, Bleich S. Homocysteine associated genomic DNA hypermethylation in patients with chronic alcoholism. J Neural Transm (Vienna). 2004;111(12):1611–6.

    Article  CAS  Google Scholar 

  27. Knothe C, Doehring A, Ultsch A, Lotsch J. Methadone induces hypermethylation of human DNA. Epigenomics. 2016;8(2):167–79.

    Article  CAS  PubMed  Google Scholar 

  28. Marie-Claire C, Crettol S, Cagnard N, Bloch V, Mouly S, Laplanche JL, Bellivier F, Lepine JP, Eap C, Vorspan F. Variability of response to methadone: genome-wide DNA methylation analysis in two independent cohorts. Epigenomics. 2016;8(2):181–95.

    Article  CAS  PubMed  Google Scholar 

  29. Chorbov VM, Todorov AA, Lynskey MT, Cicero TJ. Elevated levels of DNA methylation at the OPRM1 promoter in blood and sperm from male opioid addicts. Journal of opioid management. 2011;7(4):258–64.

    Article  PubMed  PubMed Central  Google Scholar 

  30. Nielsen DA, Yuferov V, Hamon S, Jackson C, Ho A, Ott J, Kreek MJ. Increased OPRM1 DNA methylation in lymphocytes of methadone-maintained former heroin addicts. Neuropsychopharmacology: official publication of the American College of Neuropsychopharmacology. 2009;34(4):867–73.

    Article  CAS  Google Scholar 

  31. Ebrahimi G, Asadikaram G, Akbari H, Nematollahi MH, Abolhassani M, Shahabinejad G, Khodadadnejad L, Hashemi M. Elevated levels of DNA methylation at the OPRM1 promoter region in men with opioid use disorder. The American journal of drug and alcohol abuse. 2018;44(2):193–9.

    Article  PubMed  Google Scholar 

  32. Nielsen DA, Hamon S, Yuferov V, Jackson C, Ho A, Ott J, Kreek MJ. Ethnic diversity of DNA methylation in the OPRM1 promoter region in lymphocytes of heroin addicts. Hum Genet. 2010;127(6):639–49.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Yuferov V, Levran O, Proudnikov D, Nielsen DA, Kreek MJ. Search for genetic markers and functional variants involved in the development of opiate and cocaine addiction and treatment. Ann N Y Acad Sci. 2010;1187:184–207.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Zhang H, Herman AI, Kranzler HR, Anton RF, Simen AA, Gelernter J. Hypermethylation of OPRM1 promoter region in European Americans with alcohol dependence. J Hum Genet. 2012;57(10):670–5.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Hancock DB, Levy JL, Gaddis NC, Glasheen C, Saccone NL, Page GP, Hulse GK, Wildenauer D, Kelty EA, Schwab SG, et al. Cis-expression quantitative trait loci mapping reveals replicable associations with heroin addiction in OPRM1. Biological psychiatry. 2015;78(7):474–84.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. LaForge KS, Yuferov V, Kreek MJ. Opioid receptor and peptide gene polymorphisms: potential implications for addictions. Eur J Pharmacol. 2000;410(2-3):249–68.

    Article  CAS  PubMed  Google Scholar 

  37. Smith AH, Jensen KP, Li J, Nunez Y, Farrer LA, Hakonarson H, Cook-Sather SD, Kranzler HR, Gelernter J. Genome-wide association study of therapeutic opioid dosing identifies a novel locus upstream of OPRM1. Mol Psychiatry. 2017;22(3):346–52.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Houseman EA, Kile ML, Christiani DC, Ince TA, Kelsey KT, Marsit CJ. Reference-free deconvolution of DNA methylation data and mediation by cell composition effects. BMC Bioinformatics. 2016;17:259.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  39. Houseman EA, Accomando WP, Koestler DC, Christensen BC, Marsit CJ, Nelson HH, Wiencke JK, Kelsey KT. DNA methylation arrays as surrogate measures of cell mixture distribution. BMC Bioinformatics. 2012;13:86.

    Article  PubMed  PubMed Central  Google Scholar 

  40. Theda C, Hwang SH, Czajko A, Loke YJ, Leong P, Craig JM. Quantitation of the cellular content of saliva and buccal swab samples. Sci Rep. 2018;8(1):6944.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  41. Ni Choileain N, Redmond HP. Cell response to surgery. Arch Surg. 2006;141(11):1132–40.

    Article  PubMed  Google Scholar 

  42. Koestler DC, Christensen B, Karagas MR, Marsit CJ, Langevin SM, Kelsey KT, Wiencke JK, Houseman EA. Blood-based profiles of DNA methylation predict the underlying distribution of cell types: a validation analysis. Epigenetics. 2013;8(8):816–26.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Reinius LE, Acevedo N, Joerink M, Pershagen G, Dahlen SE, Greco D, Soderhall C, Scheynius A, Kere J. Differential DNA methylation in purified human blood cells: implications for cell lineage and studies on disease susceptibility. PLoS One. 2012;7(7):e41361.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  44. Jaffe AE, Irizarry RA. Accounting for cellular heterogeneity is critical in epigenome-wide association studies. Genome Biol. 2014;15(2):R31.

    Article  PubMed  PubMed Central  Google Scholar 

  45. Tabuchi T, Shimazaki J, Satani T, Nakachi T, Watanabe Y, Tabuchi T. The perioperative granulocyte/lymphocyte ratio is a clinically relevant marker of surgical stress in patients with colorectal cancer. Cytokine. 2011;53(2):243–8.

    Article  CAS  PubMed  Google Scholar 

  46. Hancock DB, Markunas CA, Bierut LJ, Johnson EO. Human genetics of addiction: new insights and future directions. Curr Psychiatry Rep. 2018;20(2):8.

    Article  PubMed  PubMed Central  Google Scholar 

  47. Kalsi G, Euesden J, Coleman JR, Ducci F, Aliev F, Newhouse SJ, Liu X, Ma X, Wang Y, Collier DA, et al. Genome-wide association of heroin dependence in Han Chinese. PLoS One. 2016;11(12):e0167388.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  48. Oertel BG, Doehring A, Roskam B, Kettner M, Hackmann N, Ferreiros N, Schmidt PH, Lotsch J. Genetic-epigenetic interaction modulates mu-opioid receptor regulation. Hum Mol Genet. 2012;21(21):4751–60.

    Article  CAS  PubMed  Google Scholar 

  49. Elam C, Hesson L, Vos MD, Eckfeld K, Ellis CA, Bell A, Krex D, Birrer MJ, Latif F, Clark GJ. RRP22 is a farnesylated, nucleolar, RAS-related protein with tumor suppressor potential. Cancer Res. 2005;65(8):3117–25.

    Article  CAS  PubMed  Google Scholar 

  50. Belcheva MM, Vogel Z, Ignatova E, Avidor-Reiss T, Zippel R, Levy R, Young EC, Barg J, Coscia CJ. Opioid modulation of extracellular signal-regulated protein kinase activity is RAS-dependent and involves Gbetagamma subunits. J Neurochem. 1998;70(2):635–45.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  51. Bian JM, Wu N, Su RB, Li J. Opioid receptor trafficking and signaling: what happens after opioid receptor activation? Cell Mol Neurobiol. 2012;32(2):167–84.

    Article  CAS  PubMed  Google Scholar 

  52. Valderrama-Carvajal A, Irizar H, Gago B, Jimenez-Urbieta H, Fuxe K, Rodriguez-Oroz MC, Otaegui D, Rivera A. Transcriptomic integration of D4R and MOR signaling in the rat caudate putamen. Sci Rep. 2018;8(1):7337.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  53. Braunewell KH, Dwary AD, Richter F, Trappe K, Zhao C, Giegling I, Schonrath K, Rujescu D. Association of VSNL1 with schizophrenia, frontal cortical function, and biological significance for its gene product as a modulator of cAMP levels and neuronal morphology. Transl Psychiatry. 2011;1:e22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Ready LB, Sarkis E, Turner JA. Self-reported vs. actual use of medications in chronic pain patients. Pain. 1982;12(3):285–94.

    Article  CAS  PubMed  Google Scholar 

  55. Ferreri F, Bourla A, Mouchabac S. Karila L: e-Addictology: an overview of new technologies for assessing and intervening in addictive behaviors. Front Psychiatry. 2018;9:51.

    Article  PubMed  PubMed Central  Google Scholar 

  56. Aryee MJ, Jaffe AE, Corrada-Bravo H, Ladd-Acosta C, Feinberg AP, Hansen KD, Irizarry RA. Minfi: a flexible and comprehensive bioconductor package for the analysis of infinium DNA methylation microarrays. Bioinformatics. 2014;30(10):1363–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  57. Zhou W, Laird PW, Shen H. Comprehensive characterization, annotation and innovative use of Infinium DNA methylation BeadChip probes. Nucleic Acids Res. 2017;45(4):e22.

    PubMed  Google Scholar 

  58. Houseman EA, Molitor J, Marsit CJ. Reference-free cell mixture adjustments in analysis of DNA methylation data. Bioinformatics. 2014;30(10):1431–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  59. Bates D, Machler M, Bolker BM, Walker SC. Fitting linear mixed-effects models using lme4. J Stat Softw. 2015;67(1):1–48.

    Article  Google Scholar 

  60. Kuznetsova A, Brockhoff PB. Christensen RHB: lmerTest Package: tests in linear mixed effects models. J Stat Softw. 2017;82(13):1–26.

    Article  Google Scholar 

Download references

Acknowledgements

We thank the UTHSC Office of Research and the CORNET Awards for supporting this project.

Funding

Funded by the University of Tennessee Health Science Center CORNET Clinical Awards

Author information

Authors and Affiliations

Authors

Contributions

JVSS performed the lab work and data analysis and contributed to the manuscript. FISG contributed to the data analysis. JHB identified suitable patients and facilitated participant recruitment at the dental clinic. KJD contributed to the study conception and design. KM contributed to the study conception, design and data analysis, and wrote the manuscript. All authors contributed to and approved the final version of the manuscript.

Corresponding author

Correspondence to Khyobeni Mozhui.

Ethics declarations

Ethics approval and consent to participate

All participants provided written informed consent and the study received IRB approval.

Consent for publication

Not applicable

Competing interests

We have no financial or non-financial conflicts of interest.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary information

Additional file1: Table S1.

Participant-level characteristics and opioid use data. Table S2. Longitudinally collected sample-level data on DNA methylation chracteristics and cell estimates. Table S3. Top 1701 CpGs considering a genome-wide suggestive threshold of p-value= 1.5e-5. Table S4. Gene set enrichment analysis (GSEA) result using 1133 annotated genes obtained from 1701 CpGs. Table S5. CpG sites associated with opioid dosage at p ≤ 0.0001.

Additional file 2: Fig. S1.

Microarray data quality checks. a The density plots for the β-values (left: raw values; right: quantile normalized) using the full set of >850K probes show the expected bimodal distribution. b Unsupervised hierarchical clustering using 736,432 high quality probes shows that while many of the longitudinally collected data cluster appropriately with self, several samples do not cluster with self. This may be due to within-individual heterogeneity over time; however, a concerning possibility is errors in sample identity. c Repeating the unsupervised hierarchical clustering using probes that were flagged due to overlap with SNPs shows that all samples collected longitudinally from the same participant cluster with self, with the exception of only one participant who received prescription opioids (greyed: person ID 108). Person 108 paired with a participant who received no prescription opioids (person ID 122) and had no follow-up data (both 108 and 122 were excluded from the main statistical analyses). Note that this QC checks were done using data from the full set of recruited participants; however, participants who received no prescription opioids (and therefore no follow-up visits) were excluded from the main analyses.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Sandoval-Sierra, J.V., Salgado García, F.I., Brooks, J.H. et al. Effect of short-term prescription opioids on DNA methylation of the OPRM1 promoter. Clin Epigenet 12, 76 (2020). https://doi.org/10.1186/s13148-020-00868-8

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s13148-020-00868-8

Keywords