Abstract
The increased prevalence of opioid use disorder (OUD) makes it imperative to disentangle the biological mechanisms contributing to individual differences in OUD vulnerability. OUD shows strong heritability, however genetic variants contributing to vulnerability remain poorly defined. We performed a genome-wide association study using over 850 male and female heterogeneous stock (HS) rats to identify genes underlying behaviors associated with OUD such as nociception, as well as heroin-taking, extinction and seeking behaviors. By using an animal model of OUD, we were able to identify genetic variants associated with distinct OUD behaviors while maintaining a uniform environment, an experimental design not easily achieved in humans. Furthermore, we used a novel non-linear network-based clustering approach to characterize rats based on OUD vulnerability to assess genetic variants associated with OUD susceptibility. Our findings confirm the heritability of several OUD-like behaviors, including OUD susceptibility. Additionally, several genetic variants associated with nociceptive threshold prior to heroin experience, heroin consumption, escalation of intake, and motivation to obtain heroin were identified. Tom1, a microglial component, was implicated for nociception. Several genes involved in dopaminergic signaling, neuroplasticity and substance use disorders, including Brwd1, Pcp4, Phb1l2 and Mmp15 were implicated for the heroin traits. Additionally, an OUD vulnerable phenotype was associated with genetic variants for consumption and break point, suggesting a specific genetic contribution for OUD-like traits contributing to vulnerability. Together, these findings identify novel genetic markers related to the susceptibility to OUD-relevant behaviors in HS rats.
Similar content being viewed by others
Introduction
Opioid use disorder (OUD) is characterized by compulsive opioid use and relapse propensity, despite efforts to maintain abstinence. OUD diagnosis has increased substantially in the past decade [1], affecting approximately 16 million people per year worldwide [2]. The rise in diagnosis is not demographically constrained [1], but is prevalent in individuals across socio-economic class, age and sex [3], increasing the urgency to better understand the biological mechanisms contributing to OUD vulnerability. In twin studies, OUD shows high heritability (~50%) [4, 5], indicating that specific genetic variants contribute to biological basis of OUD susceptibility.
Genome-wide association studies (GWAS) are a useful tool for identifying the genetic variants that influence heritable traits contributing to OUD. GWAS of OUD, though limited, have already yielded valuable insights into genetic contributions to diagnostic vulnerability [6,7,8,9,10,11,12,13]. However, human GWAS of OUD are limited to observational data and are thus confounded by numerous poorly defined environmental factors. Understandably, due to a lack of empirical data from study participants, human GWAS have largely focused on OUD diagnosis, which is comprised of poorly quantified amounts of drug-taking, extinction and seeking behaviors. Using rodent models of OUD, we can quantify specific traits (e.g., consumption, escalation, reinstatement) associated with addiction-related behaviors to directly interrogate the relationship between these metrics and genetic variants. By doing so, we can discover novel genetic variants that may contribute to the different phases of OUD, which can then be compared to data from human studies. Crucially, such genes may also provide targets for the development of novel OUD treatments.
In the current study, we assessed several measures of heroin-taking, extinction and seeking behaviors in outbred male and female heterogeneous stock (HS) rats at two different testing locations. Additionally, anxiety- and stress-like behaviors, as well as nociceptive threshold, were assessed prior to and following heroin experience. HS rats are a diverse outbred population that was created by interbreeding 8 inbred strains at NIH in 1984 to capture genetic and behavioral heterogeneity reflective of the human population [14,15,16]. Compared to inbred and commercially available outbred rat lines, HS rats exhibit comparatively greater diversity [16] and have been used in several studies assessing addiction-related behaviors [17,18,19,20,21,22,23], including those focusing on genetic analyses [24,25,26]. In addition, HS rats are supported by complementary expression quantitative trait loci datasets [27], furthering the translational utility of our findings. In this study, we found genetic variants specific to different OUD behavioral traits, many of which had known involvement in neurobiological processes associated with substance use disorder (SUD). Additionally, we employed a novel non-linear network-based clustering approach to identify subpopulations of OUD vulnerable and resilient rats [23] and thereby identify genetic variants associated with OUD susceptibility.
Methods
Behavioral testing occurred at two geographically distinct locations: Medical University of South Carolina (MUSC, USA) and the University of Camerino (UCAM, Italy). Two testing locations rendered the advantage of accounting for environmental differences on behavior, akin to the human experience, while maintaining standardized experimental procedures and equipment used between the sites. All experimental procedures were approved by the MUSC Institutional Animal Care and Use Committee and the Italian Ministry of Health (UCAM). All testing procedures were compliant with the National Institute of Health Guide for the Care and Use of Laboratory Animals and the Assessment and Accreditation of Laboratory Animals Care, as well as the European Community Council Directive for Care and Use of Laboratory Animals.
Subjects
Male and female HS rats (NMcwiWFsm #13673907, RRID:RGD13673907) bred at Wake Forest University were shipped to MUSC and UCAM in batches of 40 (20 females and 20 males/site) at 5 weeks of age. Each site had an equal representation of littermates thereby creating genetically comparable testing cohorts. Whenever possible, each shipment contained only one male and one female per litter, to reduce genetic relationships, which can reduce power in GWAS. Upon arrival, rats were left undisturbed for 3 weeks in a climate-controlled vivarium room with a standard 12-hr light/dark cycle prior to testing. Rats were pair-housed and received ad libitum access to food and water throughout the duration of testing. Behavioral tests for stress, anxiety and nociception occurred during the light cycle, whereas all heroin-associated testing occurred during the dark phase. A total of 1160 rats underwent testing from February of 2019 through December of 2022, with the following animals excluded from all analyses: 147 for behavioral control studies (e.g., did not undergo heroin self-administration procedures), 107 due to death (post-surgical illness, n = 84; surgery, n = 23); 29 for behavioral data collection technical issues, and 3 for spleen identification issues. We sought to collect data from approximately 1000 rats as this size is sufficient to identify alleles accounting for ~2.5–3% of the variance, which are commonly observed in HS rats. Ongoing analysis is working to increase this sample size for increased genetic variant identification. Final behavioral and genetic analyses were conducted on 874 rats (MUSC, n = 479; UCAM, n = 395).
Drugs
Heroin hydrochloride was supplied to both sites by the National Institute on Drug Abuse (Bethesda, MD) and dissolved in 0.9% sterile saline.
Behavioral testing
The experimental timeline is shown in Fig. 1a. Experimental testing occurred as previously described [17, 23]. Extensive behavioral analyses for these animals, including sex and site differences, have been conducted and reported in prior publications [17, 23], and we refer the reader to these manuscripts for further edification. The focus of this manuscript is on identifying the genetic variants that contribute to individual differences in these behaviors. A brief description of behavioral testing procedures is provided below with greater detail in the Supplemental Information.
a Rats underwent testing for stress- and anxiety-related behaviors (OFT: open field test; EPM: elevated-plus maze) and nociceptive threshold (TF: tail flick) prior to jugular catheterization surgery, with tests repeated following heroin experience. Heroin self-administration occurred for 3 weeks (12-h/session, 4 days a week) followed by a progressive ratio test and an additional 3 days of self-administration training. Next, rats underwent a within-session extinction-heroin prime reinstatement test (6-h). Extinction training sessions (2-h/session) then occurred followed by a cue-induced reinstatement test (2-h). Spleen sample were collected for genetic analyses at the conclusion of training. b Correlation matrix for behavioral trait covariance. Behavioral testing that occurred prior to heroin experience is designated with a 1, and testing following heroin training as 2. Traits in bold are those with at least one significant QTL identified in GWAS analysis (n = 874).
Stress, anxiety and nociceptive-like behavioral measures
Stress and anxiety-related behaviors were assessed in the elevated-plus maze (EPM) and the open field test (OFT). Nociception under baseline (administration of 1 mg/kg saline, sc) and test (0.75 mg/kg heroin, sc) conditions were evaluated using the tail flick (TF) test. All tests were performed prior to and following heroin self-administration (see Supplemental Information and Kuhn et al. [17, 23]).
Heroin taking, extinction and seeking measures
After EPM, OFT and TF tests, rats were outfitted with an indwelling jugular catheter prior to undergoing heroin self-administration training procedures. During training, presses on the active lever resulted in presentation of a tone and cue light stimulus, as well an infusion of heroin (20 µg/kg/100 µl infusion over 3-s). Training occurred for 3 weeks with four 12-h sessions/week. Next, motivation to take heroin was assessed using a progressive ratio test whereupon the number of lever presses needed to earn an infusion exponentially increased according to the following formula: (5*e0.2n)-5 [28]. The maximum number of lever presses exerted to receive an infusion was considered the break point. Heroin self-administration training was then re-established for 3 days, followed by a within-session extinction-prime test lasting 6 h. Rats were in extinction training conditions (i.e., no cue presentation or heroin delivery) for the entire 6-h session. However, with two hours left in the session, animals received a heroin prime injection (0.25 mg/kg, sc; heroin primed reinstatement test). Animals next underwent 6 daily 2-h extinction training sessions prior to a 2-h test for cue-induced reinstatement (i.e., active lever presses resulted in cue presentation but no heroin infusion). The EPM, OFT and TF tests were then repeated. See Supplemental Information for further details.
Tissue collection and genotyping
Following all behavioral testing, animals were sacrificed and spleen samples were collected and shipped to University of California San Diego for genotyping. Genotyping-by-sequencing was employed as previously described to assess genotypes [29, 30]. A total of 7,955,513 single nucleotide polymorphisms (SNPs) were identified with an error rate <1%. However, prior to GWAS analysis, the SNPs were filtered for missingness not above 10%, minor allele frequency (MAF) of no less than 0.5% and the Hardy-Weinberg equilibrium (HWE) deviation p-value no less than 1e−10, resulting in a final set of 5,308,726 SNPs for analysis.
Statistical analyses
Behavioral traits assessed in GWAS
Behaviors included in GWAS analysis included behavioral testing traits (i.e., EPM, OFT and TF) both prior to (time point 1) and following (time point 2) heroin experience, as well as the difference in individual testing measures between the time points (Supplemental Table 1). Traits representative of heroin taking, extinction and seeking measures were also analyzed (Supplemental Table 1). Heroin-taking behaviors included: escalation of intake in the first hour of the session across training (ug/kg; average consumption hour 1 days 1–3 subtracted from hour 1 days 10–12; Escalation of intake: Hour 1), escalation of intake across the entire session (ug/kg; average consumption days 1–3 subtracted from days 10–12; Escalation of intake: Hour 1–12), total consumption (ug/kg) across training sessions 1–12, and break point achieved during the progressive ratio test. Extinction traits include active lever presses made during the following: hours 1–2 (Extinction: Burst) and 3–4 (Extinction: Prime) of the 6-h within-session extinction-heroin prime reinstatement test, the last day of extinction training prior to the test for cued reinstatement (Extinction: Day 6) and the difference in lever presses between the first and last day of extinction training (Extinction de-escalation). Active lever presses made during the tests for heroin prime and cued reinstatement were assessed for seeking measures.
Behavioral clustering
In addition to specific traits, we wanted to determine if any genetic variants were associated with OUD susceptibility. To do so, we employed a non-linear network-based clustering approach [31] using an R package ‘mlsbm’ [32] and characterized rats into OUD vulnerable, intermediate or resilient clusters. Behaviors included for analysis are indicated in Supplemental Table 1. Briefly, behavioral data was standardized within sex and site and combined to create a rat-rat similarity network, whereupon a Bayesian stochastic block model was employed to define clusters assignments [31]. Behavioral and neurobiological validation of this model are extensively detailed in prior publications [23]. Cluster assignment was used as phenotypes (Vulnerable vs intermediate/resilient; Intermediate vs vulnerable/resilient; Resilient vs vulnerable/intermediate) for GWAS. Differences in cluster composition between the top SNP allelic variants were assessed using a Chi-square test.
Phenotypic correlations
To address difference in data variance between sites (MUSC vs UCAM), covariates (including sex) that explained more than 2% of trait variance (Supplemental Tables 2 and 3) were regressed out and data was quantile normalized within site. This resulted in trait variances being equal between sites, thus weighing data from MUSC and UCAM equally. Next, data from both MUSC and UCAM were merged and quantile normalized again to maintain normality. Behavioral co-variance between traits was assessed using Spearman correlations with Bonferroni adjustment to correct for multiple comparisons using GraphPad Prism version 10. SNP heritability for behavioral traits was determined using GCTA-GREML analysis with significance determined by the p-value of rejecting the null hypothesis that the heritability is zero (p < 0.05).
Quantitative trait loci (QTL) and SNP identification
To perform the GWAS we used GCTA software employing a genetic relatedness matrix and the Leave One Chromosome Out method to account for proximal contamination [33, 34]. Because all data was quantile normalized, only one significance threshold, which was determined by permutation, was necessary [35]. We report QTLs that surpassed the permutation-derived threshold of –log10 p > 5.58 (p < 0.05). Significant QTLs are only reported when there was at least one additional SNP located within 0.5 Mb with a significance within 2 –log10(p) units of the first SNP [36]. Independence of QTLs on the same chromosome was established by making the top SNP in the most significant QTL a covariate and then repeating the GWAS analysis for the chromosome in question. This process was repeated until no more significant QTLs were identified [36]. Linkage disequilibrium (LD; r2 ≥ 0.6) was used to define the boundaries of each QTL. SNPEff was used to identify coding SNPs within each QTL interval [37]. To identify associations from other GWAS that have been performed in HS rats that map to the same loci, a phenome wide association study (PheWAS) was conducted. To explore functional mechanisms whereby SNPs may affect OUD-like traits, we assessed the co-___location of the peak SNP in each QTL with known SNPs that can alter gene expression (eQTL), slicing patterns (sQTL) or directly alter the amino acid sequence of genes within a 3 Mb window using RatGTEx.org [27]. The reference database for eQTLs was comprised of animals not in this study, therefore, simple correlations between expression and behavior were not possible. Importantly, the absence of coding variants, eQTLs or sQTLs could reflect the fact that RNA sequencing for the eQTL/sQTL analysis used a polyA pulldown protocol, limiting our findings to RNAs that are polyadenylated.
Principal component analysis
To assess whether variance in behavioral traits within significant QTLs could be explained by common factors, a principal component analysis was applied to the data. Number of components was determined using a minimum threshold of Eigenvalue of 1 and Varimax method with a minimum factor loading of 0.3 was used to rotate identified components. Analyses used SPSS version 27.
Results
Behavioral correlations
Behavioral traits showed weak to moderate correlations overall (Fig. 1b, Supplemental Table 4). However, behavior during the progressive ratio test, cued reinstatement and total heroin consumption were positively correlated (r2 ≥ 0.48), suggesting overlapping behavioral mechanisms guiding these behaviors.
Heritability estimates
SNP heritability (h2) estimates the proportion of the phenotypic variance explained by the genetic variance and is dependent on the genetic diversity of the population and environmental sources of variability. Aspects of heroin taking (Consumption, h2 = 0.17) and the inability to initially refrain from seeking (Extinction: Burst, h2 = 0.10) were heritable (Fig. 2a, Supplemental Table 5). While some measures of heroin taking, extinction and seeking were not heritable, OUD susceptibility cluster assignment [23] was heritable (Vulnerable vs resilient/intermediate, h2 = 0.17; Intermediate vs vulnerable/resilient, h2 = 0.09; Resilient vs vulnerable/intermediate, h2 = 0.15), suggesting a heritable genetic contribution to OUD susceptibility (Fig. 2a, Supplemental Table 5).
a SNP heritability (h2) measure ± SEM for all behavioral traits assessed. Significant h2 measures (p < 0.05) are indicated in blue, and non-significant traits are in orange. b Porcupine plot for selected behavioral traits. Top SNPs exceeding the genome-wide significance threshold are highlighted and identified (red: -log(p) = 5.36 and blue: -log(p) = 5.58, which correspond to permutation derived p values of 0.10 and 0.05, respectively (n = 874).
Several behavioral measures prior to and following heroin administration experience showed significant heritability (Fig. 2a, Supplemental Table 5). Nociceptive threshold under basal conditions both prior to (TF baseline 1, h2 = 0.11) and following (TF baseline 2, h2 = 0.15) heroin experience were both heritable. Heroin-induced nociception after prolonged heroin exposure was also heritable (TF test 2, h2 = 0.12), reflecting the genetic contribution to the long-term effects of heroin exposure on nociceptive threshold. Aspects of anxiety- and stress-like behaviors were also heritable, including high anxiety prior to heroin experience (EPM closed 1, h2 = 0.08), but lower levels following heroin self-administration (EPM open 2, h2 = 0.07). Our results support previous findings that locomotor behavior is heritable [24, 25] (OFT distance 1, h2 = 0.14; OFT time 1, h2 = 0.10) and show that this heritability persists even after exposure to heroin administration (OFT distance 2, h2 = 0.23; OFT time 2, h2 = 0.16).
GWAS
We identified genome-wide significant associations for several traits: escalation of intake across the entire session on Chromosome 10, total heroin consumption on Chromosome 11, break point achieved during the progressive ratio test on Chromosome 19, and baseline nociception prior to heroin experience (TF baseline 1) on Chromosomes 2 and 19 (Fig. 2b). Shared behavioral variance between these traits was assessed using PCA, with heroin taking traits loading onto one factor and baseline nociception onto its own factor, reflecting the independence of these measures in our data (Supplemental Fig. 1). While both nociception and break point were associated with loci on Chromosome 19, the associated loci are separated by approximately 4 Mb and are not in strong LD with one another, indicating that they are due to separate causal loci. See Supplemental Table 6 for MAF and HWE values for each significant loci and the Supplementary Information for a detailed report of all GWAS findings. Both sites contributed comparably to significant loci (Supplemental Table 7).
Escalation of intake
We identified a significant locus on Chromosome 10 for escalation of heroin intake (Fig. 3). The peak SNP was located at 10:35,034,818 and spans from positions 34,961,844 to 35,103,415. This locus did not contain any coding variants that are in strong LD with the top SNP. Though eQTLs and sQTLs were identified, none are in strong LD (i.e., r2) with the top SNP, suggesting they are unlikely to be causally related to the escalation of heroin intake behavior. Identified eQTLs include metabotropic glutamate receptor 6 (Grm6) in whole brain (r2 = 0.63), collagen type XXIII alpha 1 chain (Col23a1) in the orbitofrontal cortex (OFC; r2 = 0.63) and zinc finger protein (Zfp2) in the prelimbic cortex (PrL; r2 = 0.64). Two sQTLs were also identified, jade family PHD finger 2 (Jade2) in the nucleus accumbens core (NAcc) and basolateral amygdala (BLA; both r2 = 0.63), and Col23a1 in the PrL (r2 = 0.63). The lack of any putatively causal variants in this interval suggests that the causal variant remains unknown, possibly because there is an eQTL or sQTL in a tissue that is not represented in the RatGTEx database, or an eQTL that was not detected, for example, because it is significant at a different developmental time point or under different environmental conditions.
Regional association plot for chromosome 10, with the top SNP in purple and chromosome ___location indicated on the plot. Other SNPs are colored based on their linkage disequilibrium with the top SNP. Known genes of SNPs are indicated below the plots. Significance threshold: red line, p < 0.10; blue line, p < 0.05 (n = 874).
Total consumption
Total heroin consumed across training (Supplemental Fig. 2) was associated with a locus on Chromosome 11 (Fig. 4a) that spanned from positions 35,033,985 to 36,454,941 with the peak SNP at 35,034,818. Several significant causal coding variants were identified that are in strong LD with this locus, including bromodomain and WD repeat ___domain containing 1 (Brwd1: Leu2187Pro, r2 = 0.99), Purkinje cell protein 4 (Pcp4: Lys12Arg, r2 = 0.99) and two different amino acid coding variants for SH3 ___domain binding glutamate-rich protein (Sh3bgr: Glu169Lys, r2 = 0.99; Glu177Lys, r2 = 0.99). There were also eQTLs for Sh3bgr that was in high LD with this locus (r2 = 0.99 for all) in whole brain, BLA, NAcc, PrL, lateral habenula (LHb) and OFC. Two additional genes with eQTLs in high LD with the locus were identified: Pcp4 (PrL, r2 = 0.88) and a long non-coding RNA (ENSRNOG00000068065; infralimbic cortex (IL), r2 = 0.99). We also identified an sQTL for ENSRNOG00000068065 in several brain regions (whole brain, BLA, LHb, OFC, NAcc, PrL, and PrL2; r2 = 0.99–1.00). Lastly, there was an sQTL for the gene guided entry of tail-anchored proteins factor 1 (Get1; r2 = 0.99) in the IL. Heroin consumption is likely mediated, in part, by these coding differences and or by heritable differences in expression of one or more of these genes. Because several plausible variants are in LD with one another, mechanistic or other complementary approaches will be needed to identify the causal gene. Additionally, PheWAS analysis showed that OUD vulnerability (Vulnerable vs Intermediate/resilient cluster) was associated with the top SNP at a suggestive level (−log(p) = 4.5), with cluster composition differing between the allelic variants (x2 = 18.12, p < 0.0001; Fig. 4b, c).
a Regional association plot for chromosome 11, with the purple dot indicating the top SNP and ___location on chromosome indicated on the plot. The color of additional SNPS refers to extent of linkage disequilibrium with the top SNP. Known genes of SNPs are indicated below the plots. b Effect plots for total heroin consumption during training (left: normalized data; right: raw data) for OUD phenotypes (Vulnerable: pink; Intermediate/resilient: yellow) according to peak SNP allelic variant. Half of a violin plot for each phenotype is shown, with the solid line indicating the median value of the data and the dash lines representing the upper and lower interquartile range. c Top SNP allelic variant OUD phenotype cluster composition differ (x2(1, 853) = 18.12, p < 0.0001), with the minor allele associated with greater OUD resiliency. Significance threshold: red line, p < 0.10; blue line, p < 0.05. (n = 853; Vulnerable (pink), n = 361; Intermediate/resilient (yellow), n = 492).
Break point
The break point achieved during the progressive ratio test was associated with a locus on Chromosome 19 (Fig. 5a) spanning from positions 8,647,094 to 8,801,084 with the peak SNP at 8,648,697. Two potentially causal coding variants in prohibitin 1 like 2 (Phb1l2; Lys70Asn and Ile231Asn), which were in perfect LD with the top SNP (r2 = 1.0), were identified. An eQTL with perfect LD with the locus for the gene matrix metallopeptidase 15 (Mmp15; r2 = 1.0) was also identified in the whole brain. An additional eQTL (solute carrier family 38 member 7, Slc38a7) and sQTL (casein kinase 1 alpha 2, Csnk2a2) were identified in the BLA, however, the relationship between these genes and the locus are modest (r2 = 0.65 for both) rendering them much less likely to be the causes of the association with break point behavior. Similar to the Chromosome 11 association for total heroin consumption, PheWAS analysis showed a suggestive association between the top SNP and the OUD vulnerable cluster (Vulnerable vs intermediate/resilient clusters; −log(p) = 5.16; Fig. 5b). OUD cluster composition differed between the alleles at the peak SNP (x2 = 18.91, p < 0.0001; Fig. 5c). While there is a small bias toward the AA variant conferring OUD resiliency, the minor variant, AC, strongly coincides with OUD vulnerability.
a Regional association plot for chromosome 19, with the top SNP in purple and the color of other significant SNPs referring to level of linkage disequilibrium with the top SNPs. The ___location of the top SNP is indicated on the plot. Known genes of SNPs are indicated below the plots. b Effect plots showing break point achieved during the progressive ratio test (left: normalized data, right: raw data) for OUD vulnerable (pink) and intermediate/resilient (yellow) rats in relation to the top SNP allelic variant. Half a violin plot is shown for each OUD phenotype, with the solid line indicating the median value of the data set, whereas the dashed lines represent the upper and lower interquartile range. c OUD phenotypic cluster composition for both peak SNP allelic variants. Cluster composition between variants differed (x2(1, 864) = 18.91, p < 0.0001), with the minor allele heavily biased toward an OUD vulnerable phenotype. Significance threshold: red line, p < 0.10; blue line, p < 0.05. (n = 864; Vulnerable (pink), n = 345; Intermediate/resilient (yellow), n = 519).
Baseline nociception prior to heroin experience
A locus on Chromosome 2 spanning positions 51,450,191 to 52,982,078 with the peak SNP at 51,588,876 was identified for basal nociception prior to heroin taking procedures (Fig. 6a). Several coding variants in modest LD with the peak SNP were identified in genes within this locus, including coiled-coil ___domain containing 152 (Ccdc152: Ala242Ser and Thr14Asn, both r2 = 0.77), growth hormone receptor (Ghr: Ala546Val and Lys541Asn, both r2 = 0.77), zinc finger protein 131 (Zfp131: Gly506Ser, r2 = 0.82) and serine/threonine-protein kinase NIM1 (Nim1k; r2 = 0.64). The gene selenoprotein P (Selenop) had an sQTL in both BLA and whole brain (r2 = 0.78–0.79), and an eQTL for both whole brain and the NAcc (both r2 = 0.66). These data provide support of the involvement of these genes in the Chromosome 2 association with baseline nociception.
Regional association plots for (a) chromosome 2 and (b) chromosome 19. The top SNP identified for each chromosome is in purple with chromosome ___location indicated on the plot. The color of other SNPs corresponds to extent of linkage disequilibrium with the top SNP. Known genes of SNPs are indicated below the plots. Significance threshold: red line, p < 0.10; blue line, p < 0.05 (n = 874).
An additional significant association for this behavior was identified on Chromosome 19 (Fig. 6b). The peak SNP was located at 19:12,470,318 with the QTL spanning from positions 11,508,674 to 12,474,455. A potentially causal coding variant in modest LD with the top SNP in the gene LARGE xylosyl- and glucuronyltransferase 1 (Large1: Ala92Val, r2 = 0.75) was present. There were also several eQTLs for Large1 (whole brain, r2 = 0.69; BLA, r2 = 0.78; IL, r2 = 0.68; NAcc, r2 = 0.70; OFC, r2 = 0.69; PrL, r2 = 0.70; and PrL2, r2 = 0.75). In addition, there was an eQTL for ribosomal protein S18 (Rps18I2) within the LHb (r2 = 0.98). An sQTL for the gene myb1 trafficking protein (Tom1) was in strong LD with this locus (r2 = 1.0) within the PrL. Finally, there was a significant sQTL for metallothionein 3 (Mt3) in the BLA (r2 = 0.75). These coding variants, eQTLs and sQTLs could point to heritable differences and underlying biological mechanisms for basal nociception prior to drug experience at the Chromosome 19 locus.
Discussion
The economic toll on society for OUD totals over one trillion dollars annually [38]. Furthermore, over one million individuals seek treatment for OUD annually [39] and more than 80,000 people die from opioid overdoses every year in the US alone [40]. Thus, there is an urgent need to understand the factors that influence vulnerability to OUD. Using a rat line emulating within species genetic and phenotypic variability [14] akin to humans, we assessed the genetic basis of several OUD-like behaviors and further extended analysis to account for OUD susceptibility. While heroin was used, we expect these findings will be relevant for other opioids of abuse whose primary site of action is mu opioid receptors. We found several traits related to OUD are heritable, including overall OUD susceptibility. Genetic variants associated with nociception prior to heroin experience as well as heroin consummatory behaviors were identified, several of which are either biomarkers for other SUDs or are involved in neuroplasticity. Using this model, we are able to identify genetic variants possibly associated with an OUD diagnosis, but also specific behavioral traits that comprise vulnerability or resiliency thereby further contributing to our knowledge of the genetic basis of OUD.
Heritability
OUD is heritable in humans and we found several OUD-like behaviors were also heritable in HS rats. Human twin studies have yielded heritability estimates between 23–54% for OUD [4, 5] and human GWAS of OUD have produced SNP heritability estimates, which are expected to be lower, ranging from 6–12.7% [9,10,11,12,13]. However, given the complex and multi-symptomatic nature of OUD [41], assessing heritability of distinct behavioral traits associated with human OUD vulnerability is difficult. In a large sample of HS rats, we identified several heritable heroin addiction-like traits that model different aspects of human OUD. Both heroin consumption and initial extinction behavior were significantly heritable, with SNP heritability estimates for these traits similar to SNP heritability obtained from human studies of OUD. However, behaviors following prolonged heroin use and forced abstinence, such as cued reinstatement, did not exhibit significant heritability, suggesting weaker genetics influences and/or limitations in our procedure’s ability to detect the cascade of neuroplastic adaptations that occur following prolonged heroin self-administration [42, 43].
Novelty-induced locomotion, a heritable behavior [24, 25], has been associated with vulnerability in opioid addiction-related behaviors in both rodent models [17, 44,45,46,47,48] and humans [49]. We confirmed the heritability of locomotor behavior both prior to and following heroin experience in HS rats, suggesting an enduring genetic contribution to this behavior that remains relatively unaffected by factors (e.g., epigenetic changes) associated with prolonged heroin experience. Lastly, in agreement with prior work in humans, we showed that nociceptive threshold to heat stimuli under baseline conditions is heritable [50]. Following opioid use, both humans [51, 52] and rats [53] demonstrate opioid-induced hyperalgesia in response to a painful stimuli [17], a finding replicated in our study. However, the heritability of nociception in humans is unknown. Given the inherent behavioral and genetic heterogeneity present in HS rats, and the translational relevance of the findings reported here, HS rats offer an opportunity to disentangle the genetic factors contributing to opioid-induced nociception.
GWAS
Genetic variants associated with OUD-like behaviors in rats can contribute novel insights that improve our understanding of the neurobiological mechanisms underlying distinct behavioral traits associated with OUD in humans. This study implicated several loci and genes within those loci that had coding variants, eQTLs and/or sQTLs associated with affective or behavioral traits associated with SUD (Supplemental Table 8).
Escalation of heroin intake across the entire session
While no causal coding variants, eQTLs or sQTLs were in strong LD with the QTL on Chromosome 10 for escalation of heroin-taking behavior across the entire session, several genes in this interval are of potential interest. For example, Adamts2, an extracellular matrix proteinase [54] localized within the mesolimbic and mesocortical regions within the brain [55], is a biomarker for heavy tobacco smoking and is downregulated following long-term morphine exposure and withdrawal [56]. Similarly, the gene heterogeneous nuclear ribonucleoprotein H1 (Hnrnph1) is located in this interval. Hnrnph1 was identified as the causal gene for a mouse QTL for methamphetamine response [57,58,59,60] and subsequent work has shown it also modulates the response to fentanyl [61] and alcohol [62]. Further investigation is required to determine if there are overlooked coding variants or eQTLs for Adamts2 or Hnrnph1.
Heroin consumption
The region on Chromosome 11 that was associated with total consumption of heroin across all training sessions contained a number of genes with coding and expression variants. Three genes harbored coding variants that were in strong LD with the peak SNP: Brwd1, Pcp4 and Sh3bgr. In humans, BRWD1 has been associated with potentially relevant traits including cigarette smoking [63,64,65,66,67,68], externalizing behaviors [69], ADHD [66, 70] and hyperthymia in bipolar disorder [71]. PCP4 mediates calcium signaling and through downstream mechanisms also affects dopamine and acetylcholine signaling [72]. PCP4 levels are attenuated in post-mortem tissue of individuals diagnosed with alcohol use disorder [73] as well as in the NAc of rats following 2-weeks of nicotine withdrawal [74]. This gene has also been linked to cigarette smoking [63, 75] and peripheral effects of alcohol drinking [76] in humans. We identified an eQTL for Pcp4 in the PrL, a critical component of drug-reward circuitry and plasticity [72, 77]. The relationship between PCP4 and OUD is unknown, but our results suggest it may act to alter cortical to NAc transmission as a consequence of sustained heroin use. Finally, Sh3bgr had two coding variants and numerous eQTLs in whole brain, BLA, NAcc, OFC, LHb, and PrL, however the homologous human gene has not been associated with any psychiatric traits in prior human GWAS. In addition to these genes, a long non-coding RNA (ENSRNOG00000068065) and the gene Get1 had sQTLs, but there is no supporting literature for their role in OUD-relevant behaviors.
Ets2, which contained the peak SNP for this locus, did not identify any coding variants or an e/sQTL. Ets2 encodes a ubiquitous transcription factor that is found in many cells and has be linked to molecular processes [73] and cell types [74] that modulate synaptic remodeling, including those associated with OUD [78]. Our results suggest Ets2 and its allelic variants may impact OUD susceptibility, though further investigation is necessary.
While the functional association between heroin consumption and the other identified e/sQTLs are not known, DSCAM is correlated with the top SNP and encodes a neuronal cell adhesion molecule that regulates cellular connectivity [79] and glutamatergic plasticity associated with learning [80, 81]. Glutamate-mediated neuroplasticity occurs following heroin experience and is considered a driving force behind relapse vulnerability [43]. DSCAM is necessary for tasks requiring memory formation [82, 83], including those associated with rewarded behaviors [83], and was identified in human GWAS as a biomarker for cigarette cessation success [84]. Based on these data, DSCAM warrants further investigation for OUD.
Break point
The QTL on Chromosome 19 for break point achieved during the progressive ratio test contained two coding variants in strong LD with the top variant for the poorly characterized gene Phb1l2. Phb1I2 is localized to both neuronal and non-neural cell types in the central nervous system and helps to maintain mitochondria respiration by attenuating reactive oxygen species which contribute to cellular dysfunction [85]. Chronic morphine administration disrupts mitochondrial respiration and restoring proper functioning attenuates opioid withdrawal [86]. Additionally, Phb1, an analogue of Phb1l2, has been implicated in mediating glutamatergic and dopaminergic neural transmission and becomes dysregulated in several disorders, including alcohol use disorder [85]. It is thus possible that differences in the coding sequence of Phb1l2 influence break point through glutamate and dopamine dependent processes. Additionally, an eQTL for Mmp15 in whole brain that was in strong LD with the top SNP was identified, suggesting Mmp15 expression could modulate break point. MMP15 activates MMP2 [87] which contributes to NAcc cell-specific regulation of extinction behavior following heroin self-administration [88]. It is possible that attenuated Mmp15-induced Mmp2 expression contributes to higher break point levels. Furthermore, PheWAS analysis identified this locus as being associated with OUD vulnerability with the coding variant likely impacting OUD susceptibility.
Several of the heroin taking traits are both phenotypically and genetically correlated. Despite this, the significant loci identified for these traits are not the same (Supplemental Table 9). We attribute this observation to the fact that these traits are highly polygenic, meaning that many loci make small contributions to the observed phenotype. Because our power is limited, we only expect to identify a small subset of the true positive loci, and those subsets are expected to differ from one trait to another for purely stochastic reasons.
Baseline nociceptive threshold prior to heroin experience
Genetic variants on two separate chromosomes were found for nociceptive threshold under baseline conditions, prior to heroin experience. On Chromosome 19, Large1 was notable because it contained both putatively causal coding variants and an eQTL. Large1 has been associated with various behavioral traits in humans, including educational attainment [89,90,91], smoking [63], alcohol dependence [92], neuroticism [93] and depression [93, 94]; however, no prior evidence that we are aware of link it to nociception. Heritable differences in splicing patterns (sQTLs) for both Tom1 in the PrL and Mt3 in the BLA are also candidates for nociception, though neither has been previously directly linked to pain processing. Mt3 is a brain zinc binding protein involved in several metabolic processes and has been implicated in cellular dysfunction of diseases including Alzheimer and Huntington disease [95]. Tom1 is an integral component of microglial function [96], which is involved in neuropathic pain [97], but has yet to be implicated in acute painful stimuli. This gene has also been associated with cigarette smoking in humans [75].
On chromosome 2, two genes that operate in tandem with one another were also associated with nociception. There was modest LD between the top SNP for this QTL and splicing patterns in Selenop which encodes for a selenoprotein. This protein is secreted by the liver and acts to deliver selenium, a trace element necessary for proper brain function [98], via the bloodstream to regions including the brain. Selenium has been implicated in mediating dopamine release and transmission within the substantia nigra [99, 100] and the NAc [101]. The presence of Selenop on NAc dopamine terminals is proposed to attenuate dopamine release into the synaptic cleft via direct mediation of dopamine 2 autoreceptors, protecting against neurotoxic effects of drugs of abuse such as methamphetamine [101]. Ccdc152, which also had coding variants that were in modest LD with this QTL, is encoded within the Selenop gene, and acts to attenuate Selenop levels via post-translational modification [102]. Though the periaqueductal gray is the predominant brain nuclei processing pain signaling [103], striatal dopamine signaling has recently been implicated in affecting pain processing and perception [104, 105] and warrants further investigation.
Conclusion and future directions
We used a rodent model to identify several heritable traits and genetic variants associated with OUD. Unlike GWAS performed in humans, we were able to identify genetic variants associated with distinct OUD behaviors that contribute toward diagnosis. Several of the identified genes are supported by prior relationships to SUDs and other relevant behavioral and psychiatric traits, suggesting that they warrant further investigation for OUD. Furthermore, using an established rodent model that characterizes rats using a non-linear compilation of behaviors across the SUD phases to confer OUD susceptibility, we identified genetic variants for heroin consumption and motivation that were associated specifically with OUD vulnerability. These findings strengthen the translational utility of the behavioral model and can help to assess genetic variants and biomarkers associated with OUD vulnerability versus resiliency. However, complementing human OUD [5], phenotypic variance can only be partially explained by the genetic variants discovered here. Other genetic variants as well as epigenetic [106] and uncontrolled environmental factors also contribute to OUD behaviors in HS rats. Furthermore, both preclinical and clinical GWAS analysis focusing on genetic variants associated with OUD resiliency, a neurobiological active process [107, 108], should be considered in future analyses.
Data availability
All genetic data sets generated and analyzed during this study are included in the Supplemental Information of this manuscript. The behavioral data set used for these analyses are available from the corresponding author on reasonable request.
References
Wide-ranging online data for epidemiologic research (WONDER) [Internet]. 2020. Available from: http://wonder.cdc.gov/.
Chang HY, Kharrazi H, Bodycombe D, Weiner JP, Alexander GC. Healthcare costs and utilization associated with high-risk prescription opioid use: a retrospective cohort study. BMC Med. 2018;16:69.
Jones CM, Logan J, Gladden RM, Bohm MK. Vital signs: demographic and substance use trends among heroin users - United States, 2002–2013. MMWR Morb Mortal Wkly Rep. 2015;64:719–25.
Berrettini W. A brief review of the genetics and pharmacogenetics of opioid use disorders. Dialogues Clin Neurosci. 2017;19:229–36.
Crist RC, Reiner BC, Berrettini WH. A review of opioid addiction genetics. Curr Opin Psychol. 2019;27:31–5.
Cheng Z, Zhou H, Sherva R, Farrer LA, Kranzler HR, Gelernter J. Genome-wide association study identifies a regulatory variant of RGMA associated with opioid dependence in European Americans. Biol Psychiatry. 2018;84:762–70.
Gelernter J, Kranzler HR, Sherva R, Koesterer R, Almasy L, Zhao H, et al. Genome-wide association study of opioid dependence: multiple associations mapped to calcium and potassium pathways. Biol Psychiatry. 2014;76:66–74.
Nelson EC, Agrawal A, Heath AC, Bogdan R, Sherva R, Zhang B, et al. Evidence of CNIH3 involvement in opioid dependence. Mol Psychiatry. 2016;21:608–14.
Sanchez-Roige S, Fontanillas P, Jennings MV, Bianchi SB, Huang Y, Hatoum AS, et al. Genome-wide association study of problematic opioid prescription use in 132,113 23andMe research participants of European ancestry. Mol Psychiatry. 2021;26:6209–17.
Zhou H, Rentsch CT, Cheng Z, Kember RL, Nunez YZ, Sherva RM, et al. Association of OPRM1 functional coding variant with opioid use disorder: a Genome-Wide Association study. JAMA Psychiatry. 2020;77:1072–80.
Deak JD, Zhou H, Galimberti M, Levey DF, Wendt FR, Sanchez-Roige S, et al. Genome-wide association study in individuals of European and African ancestry and multi-trait analysis of opioid use disorder identifies 19 independent genome-wide significant risk loci. Mol Psychiatry. 2022;27:3970–9.
Polimanti R, Walters RK, Johnson EC, McClintick JN, Adkins AE, Adkins DE, et al. Leveraging genome-wide data to investigate differences between opioid use vs. opioid dependence in 41,176 individuals from the Psychiatric Genomics Consortium. Mol Psychiatry. 2020;25:1673–87.
Song W, Kossowsky J, Torous J, Chen CY, Huang H, Mukamal KJ, et al. Genome-wide association analysis of opioid use disorder: a novel approach using clinical data. Drug Alcohol Depend. 2020;217:108276.
Consortium S, Saar K, Beck A, Bihoreau MT, Birney E, Brocklebank D, et al. SNP and haplotype mapping for genetic analysis in the rat. Nat Genet. 2008;40:560–6.
Johannesson M, Lopez-Aumatell R, Stridh P, Diez M, Tuncel J, Blazquez G, et al. A resource for the simultaneous high-resolution mapping of multiple quantitative trait loci in rats: the NIH heterogeneous stock. Genome Res. 2009;19:150–8.
Hansen C, Spuhler K. Development of the National Institutes of Health genetically heterogeneous rat stock. Alcohol Clin Exp Res. 1984;8:477–9.
Kuhn BN, Cannella N, Crow AD, Roberts AT, Lunerti V, Allen C, et al. Novelty-induced locomotor behavior predicts heroin addiction vulnerability in male, but not female, rats. Psychopharmacology. 2022;239:3605–20.
Kallupi M, Carrette LLG, Kononoff J, Solberg Woods LC, Palmer AA, Schweitzer P, et al. Nociceptin attenuates the escalation of oxycodone self-administration by normalizing CeA-GABA transmission in highly addicted rats. Proc Natl Acad Sci. 2020;117:2140–8.
Deal A, Cooper N, Kirse HA, Uneri A, Raab-Graham K, Weiner JL, et al. Early life stress induces hyperactivity but not increased anxiety-like behavior or ethanol drinking in outbred heterogeneous stock rats. Alcohol. 2021;91:41–51.
Wang T, Han W, Chitre AS, Polesskaya O, Solberg Woods LC, Palmer AA, et al. Social and anxiety-like behaviors contribute to nicotine self-administration in adolescent outbred rats. Sci Rep. 2018;8:18069.
King CP, Tripi JA, Hughson AR, Horvath AP, Lamparelli AC, Holl KL, et al. Sensitivity to food and cocaine cues are independent traits in a large sample of heterogeneous stock rats. Sci Rep. 2021;11:2223.
Cannella N, Tambalo S, Lunerti V, Scuppa G, de Vivo L, Abdulmalek S, et al. Long-access heroin self-administration induces region specific reduction of grey matter volume and microglia reactivity in the rat. Brain Behav Immun. 2024;118:210–20.
Kuhn BN, Cannella N, Crow AD, Lunerti V, Gupta A, Walterhouse SJ, et al. A multi-symptomatic model of heroin use disorder in rats reveals distinct behavioral profiles and neuronal correlates of heroin vulnerability versus resiliency. Am J Psychiatry. 2025;182:198–208.
Chitre AS, Hebda-Bauer EK, Blandino P, Bimschleger H, Nguyen KM, Maras P, et al. Genome-wide association study in a rat model of temperament identifies multiple loci for exploratory locomotion and anxiety-like traits. Front Genet. 2022;13:1003074.
Gunturkun MH, Wang T, Chitre AS, Garcia Martinez A, Holl K, St Pierre C, et al. Genome-wide association study on three behaviors tested in an open field in heterogeneous stock rats identifies multiple loci implicated in psychiatric disorders. Front Psychiatry. 2022;13:790566.
Zhou JL, de Guglielmo G, Ho AJ, Kallupi M, Pokhrel N, Li HR, et al. Single-nucleus genomics in outbred rats with divergent cocaine addiction-like behaviors reveals changes in gene amygdala GABAergic inhibition. Nat Neurosci. 2023.
Munro D, Wang T, Chitre AS, Polesskaya O, Ehsan N, Gao J, et al. The regulatory landscape of multiple brain regions in outbred heterogeneous stock rats. Nucleic Acids Res. 2022;50:10882–95.
Richardson NR, Roberts DC. Progressive ratio schedules in drug self-administration studies in rats: a method to evaluate reinforcing efficacy. J Neurosci Methods. 1996;66:1–11.
Parker CC, Gopalakrishnan S, Carbonetto P, Gonzales NM, Leung E, Park YJ, et al. Genome-wide association study of behavioral, physiological and gene expression traits in outbred CFW mice. Nat Genet. 2016;48:919–26.
Gileta AF, Gao J, Chitre AS, Bimschleger HV, St Pierre CL, Gopalakrishnan S, et al. Adapting genotyping-by-sequencing and variant calling for heterogeneous stock rats. G3 (Bethesda). 2020;10:2195–205.
Allen C, Kuhn BN, Cannella N, Crow AD, Roberts AT, Lunerti V, et al. Network-based discovery of opioid use vulnerability in rats using the bayesian stochastic block model. Front Psychiatry. 2021;12:745468.
Allen C, Chung D mlsbm: efficient estimation of bayesian SBMs & MLSBMs. R package version 0.99.2 ed2021
Yang J, Lee SH, Goddard ME, Visscher PM. GCTA: a tool for genome-wide complex trait analysis. Am J Hum Genet. 2011;88:76–82.
Cheng R, Parker CC, Abney M, Palmer AA. Practical considerations regarding the use of genotype and pedigree data to model relatedness in the context of Genome-Wide Association studies. G3 (Bethesda). 2013;3:1861–7.
Cheng R, Palmer AA. A simulation study of permutation, bootstrap, and gene dropping for assessing statistical significance in the case of unequal relatedness. Genetics. 2013;193:1015–8.
Chitre AS, Polesskaya O, Holl K, Gao J, Cheng R, Bimschleger H, et al. Genome-Wide Association study in 3173 outbred rats identifies multiple loci for body weight, adiposity, and fasting glucose. Obesity (Silver Spring). 2020;28:1964–73.
Cingolani P, Platts A, Wang le L, Coon M, Nguyen T, Wang L, et al. A program for annotating and predicting the effects of single nucleotide polymorphisms, SnpEff: SNPs in the genome of Drosophila melanogaster strain w1118; iso-2; iso-3. Fly (Austin). 2012;6:80–92.
Florence C, Luo F, Rice K. The economic burden of opioid use disorder and fatal opioid overdose in the United States, 2017. Drug Alcohol Depend. 2021;218:108350.
National Survey of Substance Abuse Treatment Services (N-SSATS): 2020 [Internet]. Substabce Abuse and Mental Health Services Administration. 2021. Available from: https://www.samhsa.gov/data/sites/default/files/reports/rpt35313/2020_NSSATS_FINAL.pdf.
Drug overdose death rates [Internet]. National Institute on Drug Abuse. 2023. Available from: https://nida.nih.gov/research-topics/trends-statistics/overdose-death-rates.
Association AP. Diagnostic and statistical manual of mental disorders. 5 ed2013
Mazei-Robison MS, Nestler EJ. Opiate-induced molecular and cellular plasticity of ventral tegmental area and locus coeruleus catecholamine neurons. Cold Spring Harb Perspect Med. 2012;2:a012070.
Kruyer A, Chioma VC, Kalivas PW. The opioid-addicted tetrapartite synapse. Biol Psychiatry. 2020;87:34–43.
Ambrosio E, Goldberg SR, Elmer GI. Behavior genetic investigation of the relationship between spontaneous locomotor activity and the acquisition of morphine self-administration behavior. Behav Pharmacol. 1995;6:229–37.
Lamarque S, Taghzouti K, Simon H. Chronic treatment with Delta(9)-tetrahydrocannabinol enhances the locomotor response to amphetamine and heroin. Implications for vulnerability to drug addiction. Neuropharmacology. 2001;41:118–29.
Xigeng Z, Yonghui L, Xiaojing L, Lin X, Dongmei W, Jie L, et al. Social crowding sensitizes high-responding rats to psychomotor-stimulant effects of morphine. Pharmacol Biochem Behav. 2004;79:213–8.
Swain Y, Muelken P, LeSage MG, Gewirtz JC, Harris AC. Locomotor activity does not predict individual differences in morphine self-administration in rats. Pharmacol Biochem Behav. 2018;166:48–56.
Chang SE, Krueger LD, Flagel SB. Investigating individual differences in opioid-taking and opioid-seeking behavior in male rats. Psychopharmacology. 2022.
Wingo T, Nesil T, Choi JS, Li MD. Novelty seeking and drug addiction in humans and animals: from behavior to molecules. J Neuroimmune Pharmacol. 2016;11:456–70.
Angst MS, Phillips NG, Drover DR, Tingle M, Ray A, Swan GE, et al. Pain sensitivity and opioid analgesia: a Pharmacogenomic Twin study. Pain. 2012;153:1397–409.
Ren ZY, Shi J, Epstein DH, Wang J, Lu L. Abnormal pain response in pain-sensitive opiate addicts after prolonged abstinence predicts increased drug craving. Psychopharmacology. 2009;204:423–9.
Carcoba LM, Contreras AE, Cepeda-Benito A, Meagher MW. Negative affect heightens opiate withdrawal-induced hyperalgesia in heroin dependent individuals. J Addict Dis. 2011;30:258–70.
Marchette RCN, Gregory-Flores A, Tunstall BJ, Carlson ER, Jackson SN, Sulima A, et al. kappa-opioid receptor antagonism reverses heroin withdrawal-induced hyperalgesia in male and female rats. Neurobiol Stress. 2021;14:100325.
Shiomi T, Lemaitre V, D’Armiento J, Okada Y. Matrix metalloproteinases, a disintegrin and metalloproteinases, and a disintegrin and metalloproteinases with thrombospondin motifs in non-neoplastic diseases. Pathol Int. 2010;60:477–96.
Ruso-Julve F, Pombero A, Pilar-Cuellar F, Garcia-Diaz N, Garcia-Lopez R, Juncal-Ruiz M, et al. Dopaminergic control of ADAMTS2 expression through cAMP/CREB and ERK: molecular effects of antipsychotics. Transl Psychiatry. 2019;9:306.
Thiagarajan SK, Mok SY, Ogawa S, Parhar IS, Tang PY. Receptor-mediated AKT/PI3K signalling and behavioural alterations in zebrafish larvae reveal association between schizophrenia and opioid use disorder. Int J Mol Sci. 2022;23:4715.
Yazdani N, Parker CC, Shen Y, Reed ER, Guido MA, Kole LA, et al. Hnrnph1 is a quantitative trait gene for methamphetamine sensitivity. PLoS Genet. 2015;11:e1005713.
Ruan QT, Yazdani N, Blum BC, Beierle JA, Lin W, Coelho MA, et al. A mutation in Hnrnph1 that decreases methamphetamine-induced reinforcement, reward, and dopamine release and increases synaptosomal hnRNP H and mitochondrial proteins. J Neurosci. 2020;40:107–30.
Ruan QT, Yazdani N, Reed ER, Beierle JA, Peterson LP, Luttik KP, et al. 5′ UTR variants in the quantitative trait gene Hnrnph1 support reduced 5′ UTR usage and hnRNP H protein as a molecular mechanism underlying reduced methamphetamine sensitivity. FASEB J. 2020;34:9223–44.
Borrelli KN, Langan CR, Dubinsky KR, Szumlinski KK, Carlezon WA Jr, Chartoff EH, et al. Intracranial self-stimulation and concomitant behaviors following systemic methamphetamine administration in Hnrnph1 mutant mice. Psychopharmacology. 2021;238:2031–41.
Bryant CD, Healy AF, Ruan QT, Coehlo MA, Lustig E, Yazdani N, et al. Sex-dependent effects of an Hnrnph1 mutation on fentanyl addiction-relevant behaviors but not antinociception in mice. Genes Brain Behav. 2021;20:e12711.
Fultz EK, Coelho MA, Lieberman D, Jimenez-Chavez CL, Bryant CD, Szumlinski KK. Hnrnph1 is a novel regulator of alcohol reward. Drug Alcohol Depend. 2021;220:108518.
Saunders GRB, Wang X, Chen F, Jang SK, Liu M, Wang C, et al. Genetic diversity fuels gene discovery for tobacco and alcohol use. Nature. 2022;612:720–4.
Liu M, Jiang Y, Wedow R, Li Y, Brazel DM, Chen F, et al. Association studies of up to 1.2 million individuals yield new insights into the genetic etiology of tobacco and alcohol use. Nat Genet. 2019;51:237–44.
Xu K, Li B, McGinnis KA, Vickers-Smith R, Dao C, Sun N, et al. Genome-wide association study of smoking trajectory and meta-analysis of smoking status in 842,000 individuals. Nat Commun. 2020;11:5302.
Karlsson Linner R, Biroli P, Kong E, Meddens SFW, Wedow R, Fontana MA, et al. Genome-wide association analyses of risk tolerance and risky behaviors in over 1 million individuals identify hundreds of loci and shared genetic influences. Nat Genet. 2019;51:245–57.
Kichaev G, Bhatia G, Loh PR, Gazal S, Burch K, Freund MK, et al. Leveraging polygenic functional enrichment to improve GWAS power. Am J Hum Genet. 2019;104:65–75.
Wootton RE, Richmond RC, Stuijfzand BG, Lawn RB, Sallis HM, Taylor GMJ, et al. Evidence for causal effects of lifetime smoking on risk for depression and schizophrenia: a Mendelian Randomisation study. Psychol Med. 2020;50:2435–43.
Karlsson Linner R, Mallard TT, Barr PB, Sanchez-Roige S, Madole JW, Driver MN, et al. Multivariate analysis of 1.5 million people identifies genetic associations with traits related to self-regulation and addiction. Nat Neurosci. 2021;24:1367–76.
Rao S, Baranova A, Yao Y, Wang J, Zhang F. Genetic relationships between attention-deficit/hyperactivity disorder, autism spectrum disorder, and intelligence. Neuropsychobiology. 2022;81:484–96.
Greenwood TA, Akiskal HS, Akiskal KK, Bipolar Genome S, Kelsoe JR. Genome-wide association study of temperament in bipolar disorder reveals significant associations with three novel loci. Biol Psychiatry. 2012;72:303–10.
Kalivas PW, Volkow ND. The neural basis of addiction: a pathology of motivation and choice. Am J Psychiatry. 2005;162:1403–13.
Trojanowska M. Ets factors and regulation of the extracellular matrix. Oncogene. 2000;19:6464–71.
Cornell J, Salinas S, Huang HY, Zhou M. Microglia regulation of synaptic plasticity and learning and memory. Neural Regen Res. 2022;17:705–16.
Argos M, Tong L, Pierce BL, Rakibuz-Zaman M, Ahmed A, Islam T, et al. Genome-wide association study of smoking behaviours among Bangladeshi adults. J Med Genet. 2014;51:327–33.
Simino J, Sung YJ, Kume R, Schwander K, Rao DC. Gene-alcohol interactions identify several novel blood pressure loci including a promising locus near SLC16A9. Front Genet. 2013;4:277.
Hearing M. Prefrontal-accumbens opioid plasticity: implications for relapse and dependence. Pharmacol Res. 2019;139:158–65.
Green JM, Sundman MH, Chou YH. Opioid-induced microglia reactivity modulates opioid reward, analgesia, and behavior. Neurosci Biobehav Rev. 2022;135:104544.
Ronn LC, Berezin V, Bock E. The neural cell adhesion molecule in synaptic plasticity and ageing. Int J Dev Neurosci. 2000;18:193–9.
Li HL, Huang BS, Vishwasrao H, Sutedja N, Chen W, Jin I, et al. Dscam mediates remodeling of glutamate receptors in Aplysia during de novo and learning-related synapse formation. Neuron. 2009;61:527–40.
Stachowicz K, Panczyszyn-Trzewik P, Sowa-Kucma M, Misztak P. Changes in working memory induced by lipopolysaccharide administration in mice are associated with metabotropic glutamate receptors 5 and contrast with changes induced by cyclooxygenase-2: involvement of postsynaptic density protein 95 and down syndrome cell adhesion molecule. Neuropeptides. 2023;100:102347.
Chen BE, Kondo M, Garnier A, Watson FL, Puettmann-Holgado R, Lamar DR, et al. The molecular diversity of Dscam is functionally required for neuronal wiring specificity in Drosophila. Cell. 2006;125:607–20.
Keene AC, Krashes MJ, Leung B, Bernard JA, Waddell S. Drosophila dorsal paired medial neurons provide a general mechanism for memory consolidation. Curr Biol. 2006;16:1524–30.
Uhl GR, Liu QR, Drgon T, Johnson C, Walther D, Rose JE, et al. Molecular genetics of successful smoking cessation: convergent Genome-Wide Association study results. Arch Gen Psychiatry. 2008;65:683–93.
Bernstein HG, Smalla KH, Keilhoff G, Dobrowolny H, Kreutz MR, Steiner J. The many “Neurofaces” of prohibitins 1 and 2: crucial for the healthy brain, dysregulated in numerous brain disorders. J Chem Neuroanat. 2023;132:102321.
Kowalczyk P, Sulejczak D, Kleczkowska P, Bukowska-Osko I, Kucia M, Popiel M, et al. Mitochondrial oxidative stress-a causative factor and therapeutic target in many diseases. Int J Mol Sci. 2021;22:13384.
Wright JW, Harding JW. Contributions of matrix metalloproteinases to neural plasticity, habituation, associative learning and drug addiction. Neural Plast. 2009;2009:579382.
Chioma VC, Kruyer A, Bobadilla AC, Angelis A, Ellison Z, Hodebourg R, et al. Heroin seeking and extinction from seeking activate matrix metalloproteinases at synapses on distinct subpopulations of accumbens cells. Biol Psychiatry. 2021;89:947–58.
Okbay A, Wu Y, Wang N, Jayashankar H, Bennett M, Nehzati SM, et al. Polygenic prediction of educational attainment within and between families from genome-wide association analyses in 3 million individuals. Nat Genet. 2022;54:437–49.
Schoeler T, Speed D, Porcu E, Pirastu N, Pingault JB, Kutalik Z. Participation bias in the UK Biobank distorts genetic associations and downstream analyses. Nat Hum Behav. 2023;7:1216–27.
Lee JJ, Wedow R, Okbay A, Kong E, Maghzian O, Zacher M, et al. Gene discovery and polygenic prediction from a genome-wide association study of educational attainment in 1.1 million individuals. Nat Genet. 2018;50:1112–21.
Kapoor M, Wang JC, Wetherill L, Le N, Bertelsen S, Hinrichs AL, et al. Genome-wide survival analysis of age at onset of alcohol dependence in extended high-risk COGA families. Drug Alcohol Depend. 2014;142:56–62.
Baselmans BML, Jansen R, Ip HF, van Dongen J, Abdellaoui A, van de Weijer MP, et al. Multivariate genome-wide analyses of the well-being spectrum. Nat Genet. 2019;51:445–51.
Li QS, Wajs E, Ochs-Ross R, Singh J, Drevets WC. Genome-wide association study and polygenic risk score analysis of esketamine treatment response. Sci Rep. 2020;10:12649.
Koh JY, Lee SJ. Metallothionein-3 as a multifunctional player in the control of cellular processes and diseases. Mol Brain. 2020;13:116.
Wang H, Wang X, Shen Y, Wang Y, Yang T, Sun J, et al. SENP1 modulates chronic intermittent hypoxia-induced inflammation of microglia and neuronal injury by inhibiting TOM1 pathway. Int Immunopharmacol. 2023;119:110230.
Inoue K, Tsuda M. Microglia in neuropathic pain: cellular and molecular mechanisms and therapeutic potential. Nat Rev Neurosci. 2018;19:138–52.
Pillai R, Uyehara-Lock JH, Bellinger FP. Selenium and selenoprotein function in brain disorders. IUBMB Life. 2014;66:229–39.
Bellinger FP, Raman AV, Rueli RH, Bellinger MT, Dewing AS, Seale LA, et al. Changes in selenoprotein P in substantia nigra and putamen in Parkinson’s disease. J Parkinsons Dis. 2012;2:115–26.
Romero-Ramos M, Venero JL, Cano J, Machado A. Low selenium diet induces tyrosine hydroxylase enzyme in nigrostriatal system of the rat. Brain Res Mol Brain Res. 2000;84:7–16.
Torres DJ, Yorgason JT, Mitchell CC, Hagiwara A, Andres MA, Kurokawa S, et al. Selenoprotein P modulates methamphetamine enhancement of vesicular dopamine release in mouse nucleus accumbens via dopamine D2 receptors. Front Neurosci. 2021;15:631825.
Mita Y, Uchida R, Yasuhara S, Kishi K, Hoshi T, Matsuo Y, et al. Identification of a novel endogenous long non-coding RNA that inhibits selenoprotein P translation. Nucleic Acids Res. 2021;49:6893–907.
Averitt DL, Eidson LN, Doyle HH, Murphy AZ. Neuronal and glial factors contributing to sex differences in opioid modulation of pain. Neuropsychopharmacology. 2019;44:155–65.
Taylor AMW, Becker S, Schweinhardt P, Cahill C. Mesolimbic dopamine signaling in acute and chronic pain: implications for motivation, analgesia, and addiction. Pain. 2016;157:1194–8.
Harris HN, Peng YB. Evidence and explanation for the involvement of the nucleus accumbens in pain processing. Neural Regen Res. 2020;15:597–605.
Werner CT, Altshuler RD, Shaham Y, Li X. Epigenetic mechanisms in drug relapse. Biol Psychiatry. 2021;89:331–8.
Ersche KD, Meng C, Ziauddeen H, Stochl J, Williams GB, Bullmore ET, et al. Brain networks underlying vulnerability and resilience to drug addiction. Proc Natl Acad Sci. 2020;117:15253–61.
Ersche KD, Williams GB, Robbins TW, Bullmore ET. Meta-analysis of structural brain abnormalities associated with stimulant drug dependence and neuroimaging of addiction vulnerability and resilience. Curr Opin Neurobiol. 2013;23:615–24.
Acknowledgements
The authors would like to thank the numerous technicians within the Kalivas and Ciccocioppo lab for their help performing the behavioral assays on the rats used in these analyses, as well as members of the Kalivas lab for their helpful feedback on the manuscript. This work was supported by the National Institute on Drug Abuse U01DA45300 (PWK), P50DA037844 (AAP), T32DA007288 (BNK) and K99DA057390 (BNK).
Funding
Open access funding provided by the Carolinas Consortium.
Author information
Authors and Affiliations
Contributions
Experimental design was accomplished through a collaborative effort between BNK, NC, MU, GH, DC, OP, OP, LCSW, RC, PWK and AAP. All heterogeneous rats used in the experiments were provided by LCSW. Behavioral assays were conducted by BNK, NC, ADC, VL, ED, LS, JLH, ATR, MU, SA and AK, with analyses by BNK and NC. Both AG and DC conducted the behavioral clustering analysis. Genetic assays and analyses were conducted by ASC, KMHN, KC, DC, BP, KSZ, BL, BBJ, TMS and OP. This manuscript was written by BNK, PWK and AAP. All authors contributed edits to the manuscript.
Corresponding author
Ethics declarations
Competing interests
The authors declare no competing interests.
Ethical approval
All experimental procedures were approved by the MUSC Institutional Animal Care and Use Committee and the Italian Ministry of Health (UCAM). All testing procedures were compliant with the National Institute of Health Guide for the Care and Use of Laboratory Animals and the Assessment and Accreditation of Laboratory Animals Care, as well as the European Community Council Directive for Care and Use of Laboratory Animals. No human studies were performed.
Additional information
Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Supplementary information
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/.
About this article
Cite this article
Kuhn, B.N., Cannella, N., Chitre, A.S. et al. Genome-wide association study reveals multiple loci for nociception and opioid consumption behaviors associated with heroin vulnerability in outbred rats. Mol Psychiatry (2025). https://doi.org/10.1038/s41380-025-02922-4
Received:
Revised:
Accepted:
Published:
DOI: https://doi.org/10.1038/s41380-025-02922-4
