Sex-Specific Gene Expression in the Mouse Nucleus Accumbens Before and After Cocaine Exposure.

The nucleus accumbens plays a major role in the response of mammals to cocaine. In animal models and human studies, the addictive effects of cocaine and relapse probability have been shown to be greater in females. Sex-specific differential expression of key transcripts at baseline and after prolonged withdrawal could underlie these differences. To distinguish between these possibilities, gene expression was analyzed in four groups of mice (cycling females, ovariectomized females treated with estradiol or placebo, and males) 28 days after they had received seven daily injections of saline or cocaine. As expected, sensitization to the locomotor effects of cocaine was most pronounced in the ovariectomized mice receiving estradiol, was greater in cycling females than in males, and failed to occur in ovariectomized/placebo mice. After the 28-day withdrawal period, RNA prepared from the nucleus accumbens of the individual cocaine- or saline-injected mice was subjected to RNA sequencing analysis. Baseline expression of 3% of the nucleus accumbens transcripts differed in the cycling female mice compared with the male mice. Expression of a similar number of transcripts was altered by ovariectomy or was responsive to estradiol treatment. Nucleus accumbens transcripts differentially expressed in cycling female mice withdrawn from cocaine exhibited substantial overlap with those differentially expressed in cocaine-withdrawn male mice. A small set of transcripts were similarly affected by cocaine in the placebo- or estradiol-treated ovariectomized mice. Sex and hormonal status have profound effects on RNA expression in the nucleus accumbens of naive mice. Prolonged withdrawal from cocaine alters the expression of a much smaller number of common and sex hormone-specific transcripts.

Human females have been shown to be more susceptible to the effects of cocaine and other addictive drugs than males [1-4], and these sex-specific differences can be mimicked in rodents. Female rodents are more likely than males to choose cocaine administration over food [5]. Female rodents also have exaggerated locomotor responses to psychostimulants, accelerate voluntary drug administration more rapidly, and show a greater tendency to reinstatement than males [5-13]. Sex-specific differences in gene expression during development and in the adult, along with the effects of circulating hormones, presumably lead to these differences.

Across somatic tissues, >10,000—nearly one-half of the genome—show sex-biased expression (expression biased by genetic sex and/or hormonal status) [14]. Important aspects of the basic wiring pattern in the brain are sex-specific, with girls and young women showing stronger connections between hemispheres, and boys and young men showing stronger connections within hemispheres [8, 9, 11]. Striking sex-specific differences such as these have also been observed in those with Parkinson disease and a number of other pathologic entities centered in the striatum [15]. One region of the ventral striatum, the nucleus accumbens or “reward center,” plays a central role in the response to addictive drugs such as cocaine [9]. Sex-specific differences in the nucleus accumbens transcriptome of drug-naive males and cycling females would be expected to affect their initial response to cocaine administration and their response to withdrawal from cocaine. To the best of our knowledge, sex hormonal status-dependent differences in the nucleus accumbens transcriptome have not yet been addressed.

In addition to its critical developmental role, estradiol exerts tissue-specific effects on gene expression in the adult [15, 16]. Many of the behavioral deficits observed in ovariectomized mice will be rapidly reversed by the administration of estradiol [12]. In female mice, the spine density on the dendrites of hippocampal pyramidal neurons varies in synchrony with the estrous cycle [17, 18]. To identify sex hormonal status-dependent differences in the nucleus accumbens transcriptome of naive mice and to explore sex-specific responses to cocaine, we examined the nucleus accumbens of naive and cocaine-withdrawn ovariectomized mice implanted with estradiol or placebo pellets. Our goal was to test the hypothesis that the expression of a subset of transcripts would be altered in a sex-specific, and perhaps hormone-specific, and, certainly, cocaine withdrawal-specific pattern in the nucleus accumbens. Our further expectation was that the differentially expressed transcripts would control and explain the sex-specific differences in the behavioral responses to cocaine [13, 19]. Our identification of sex-specific and cocaine-specific nucleus accumbens transcripts validated our fundamental hypothesis. In addition, our study enumerated sex-specific nucleus transcripts not identified in RNA prepared from larger brain regions (available at: https://www.ncbi.nlm.nih.gov/gene/) and revealed profound sex hormonal status-dependent differences in the nucleus accumbens transcriptome of naive mice.

1. Methods

A. Mice

All experiments were performed using male and female C57BL/6J wild-type mice (Jackson Laboratories, Bar Harbor, ME). The mice were group housed (four per cage and one sex per cage) in the University of Connecticut Health Center animal facility with a 12-hour light/dark cycle (lights on at 7:00 am). Food and water were available ad libitum. Ovariectomized mice implanted with estradiol or placebo pellets were housed together. The saline- and cocaine-treated mice were group housed but not mixed together. All experiments were executed in accordance with the University of Connecticut Health Center institutional animal care and use committee guidelines.

B. Ovariectomy

Female mice (postnatal day, 75 to 90) were ovariectomized (day 0) as previously described [20, 21]. Additionally, a small incision superficial to the muscle layer was made at the nape of the neck. A 60-day continuous-release pellet (Innovative Research of America, Inc., Sarasota, FL) was placed ~1 cm away from the incision site [20]. The pellets contained either a placebo or 0.01 mg of estradiol with a nominal release rate of 165 ng/d. The incision was closed with a suture. The mice were allowed to recover for 10 days before any experimental testing; uterine weights were obtained at sacrifice. Previous work had established that estradiol levels were undetectable in ovariectomized mice and that the 0.01-mg estradiol pellets restored estradiol levels and uterine weights to normal physiological levels [20].

C. Cocaine Sensitization

Sets of mice were given intraperitoneal injections of either 0.9% saline (10 mg/kg/d) or cocaine (10 or 20 mg/kg/d). The mice were injected on the top of their home cage and immediately placed in a 16-in. × 16-in. Plexiglas chamber (PAS Open Field System; San Diego Instruments, San Diego, CA) to monitor their locomotor activity for 45 minutes. On days 1 to 3, all the mice received one intraperitoneal injection of saline (S1 to S3). Control animals received a daily injection of saline for the next 7 days (S4 to S10). Cocaine (provided by the National Institute on Drug Abuse, Bethesda, MD) was dissolved in 0.9% saline. On the first and seventh days of cocaine administration, the mice received 10 mg/kg of cocaine (C1 and C7); 20 mg/kg of cocaine was administered on the intervening days (C2 to C6). This cocaine sensitization paradigm was adopted from the study by Pierce et al. [22] and Mazzone et al. [20]. Stereotypy was not evident using this dosing paradigm. After completion of the saline and cocaine injections, the mice remained in their home cage without treatment for 28 days. Previous studies from our laboratory, and others, have consistently demonstrated that cocaine sensitization will be maintained for >4 weeks [20, 23–26]. Because we wanted to analyze the transcriptome after withdrawal, continued sensitization to cocaine was confirmed by challenging a single mouse from each treatment group with 10 mg/kg of cocaine after 27 days of withdrawal. After validation of sensitization by recording the open field activity of the challenged mice, the saline controls and mice withdrawn from cocaine were euthanized, and coronal sections (1 mm) containing the nucleus accumbens were prepared. RNA samples from the challenged mice were sequenced but were excluded from the pooled analyses of RNA from the cocaine-withdrawn mice. Three sets of ovariectomized mice were studied for locomotion (n = 15) plus two sets of cycling females and males (n = 10).

D. RNA Sequencing, Differential Expression, and Pathway Analyses

RNA was extracted from tissue punches containing the nucleus accumbens using TRIzol; (Invitrogen, Carlsbad, CA) barcoded libraries were constructed as described previously [23]. Sequencing was accomplished using a HiSeq sequence (Illumina San Diego, CA) with eight samples per lane (11.3 ± 5 × 106 reads per sample), using 76-nt paired end reads. Alignment to the Mus musculus genome (December 2011; GRCm38/mm10) was performed as previously described [23]; 70% to 78% of the reads mapped to mm10 transcripts. Differential expression was evaluated using three methods: Cuffdiff (Center for Genome Innovation Storrs, Storrs, CT; available at: http://bioinformatics.uconn.edu/), Limma (available at: http://genexplain.cam.uchc.edu:8080/bioumlweb/), and DESeq2 (Center for Genome Innovation Storrs, Storrs, CT; available at: http://bioinformatics.uconn.edu/) [27]. Quantitative PCR (qPCR) was performed as described previously [23, 27]. All primers had a calculated melt temperature of 60°C to 61.5°C, and all products were 120 ± 4 nt [28]. RNA sequencing (RNASeq) data were analyzed using ingenuity pathway analysis (IPA; Qiagen, Germantown, MD) [27], focusing on pathways depicting direct effects of the upstream regulatory protein on downstream targets. The primary sequencing data have been deposited in the Gene Expression Omnibus under the series record number GSE111600.

E. Statistical Analysis

Simple comparisons were evaluated using paired t tests (in Excel [Microsoft, Redmond, WA] or Apache Open Office); locomotion was evaluated during the 7-day saline-or-cocaine test with repeated measures ANOVA (Prism 7; GraphPad, San Diego, CA). One-way ANOVA was performed for the effects of sex/estradiol using Excel. Two-way ANOVA was performed for the combined effects of sex (placebo or estradiol) and drug (saline or cocaine) and any interactions, using SPSS, version 25 (IBM, Armonk, NY). The cutoff for significance was P < 0.05.

2. Results

A. Locomotor Response to Cocaine Is Sensitive to Sex and Estradiol

Locomotor sensitization has been widely accepted as a reliable measure of the behavioral plasticity that occurs in response to repeated cocaine exposure [29]. In an earlier study of male mice receiving daily injections of cocaine, we observed sensitization and used high throughput sequencing to identify transcripts whose expression remained altered after 1 week of withdrawal [23]. In the present study, we extended this approach to freely cycling female mice and ovariectomized mice receiving a placebo implant or an estradiol implant [20]. The withdrawal period was extended to 28 days to focus on the long-lasting consequences of cocaine administration [19, 26], and male mice were again studied. Each of the four groups was divided equally into mice receiving 10 days of saline injections (S1 to S10) and mice receiving 3 days of saline injections (S1 to S3), followed by 7 days of cocaine injections (C1 to C7; Fig. 1A).

Figure 1.

Experimental paradigm, uterine weights, and locomotor response to cocaine. (A) Experimental timeline. On day 1, mice were ovariectomized or handled. A 10-day recovery period followed before daily injections of saline began (S1 to S3). Daily saline (S4 to S10) or cocaine (C1 to C7) injections occurred on each of the next 7 days. Locomotor activity was assessed immediately after each of the 10 injections. Nucleus accumbens were harvested 28 days after the last saline or cocaine injection. (B) Uterine weights were determined for 26 female mice per group [one-way ANOVA, F(2,75) = 240; P < 0.0001). (C) Locomotor responses were documented for all mice. Averaged data for ovariectomized females (n = 15) with a placebo implant and ovariectomized females (n = 15) with a slow-release estradiol implant shown. Repeated measures ANOVA for C1 to C7 vs S4 to S10 (both ovariectomized plus placebo and ovariectomized plus estradiol), P = 0.0005; and for C1 to C7, placebo vs estradiol, P = 0.005.

With the inclusion of a 10-day recovery period after ovariectomy and a 28-day period of withdrawal, the experimental protocol extended over a 7-week period. To determine whether the replacement dose of estradiol was adequate, uterine weights were measured at the time of sacrifice (Fig. 1B). This measurement verified the success of the ovariectomies and the efficacy of the estradiol replacement paradigm during the study period [20].

Daily locomotor data were collected for individual mice in the saline- and cocaine-treated male, cycling female, ovariectomized plus placebo, and ovariectomized plus estradiol groups. Data for the ovariectomy plus placebo and ovariectomy plus estradiol groups are shown in Fig. 1C. For all four groups, the number of ambulations declined slightly after the first saline injection, reaching a steady level that was maintained through S10 (Fig. 1C). A comparison of the daily data for these two groups revealed a substantially increased number of ambulations in the cocaine-injected mice with an estradiol implant compared with those with the placebo implant (Fig. 1C; P = 0.005, repeated measures ANOVA).

At the end of the 28-day withdrawal period, when the cocaine-injected mice were still sensitized, the nucleus accumbens was harvested from the saline- and cocaine-injected mice in all four groups. Bar-coded samples (four or five per group) were subjected to high throughput RNA sequencing. We first examined differences in the transcriptomes of the saline-treated control mice and then assessed their transcriptomes after withdrawal from cocaine treatment. The complete data set has been provided in an online repository [28].

B. Sex- and Estradiol-Dependent Expression of Transcripts in the Nucleus Accumbens

To assess sex-specific differences in gene expression in the nucleus accumbens, we first examined the expression of X- and Y-chromosome–encoded genes in all four groups of saline-treated mice (Fig. 2A and 2B). A total of 36 Y-chromosome transcripts have been identified in mm10 (GRCm38; December 2011). Transcripts encoding only seven of these genes were detected in the male nucleus accumbens. Because fragments per kilobase of transcript per million mapped reads (FPKM) values for these transcripts all were far <1 FPKM in females, we set FPKM of 1 as our detection limit. Only four Y chromosome transcripts were expressed in the male nucleus accumbens at a level >1 FPKM (Fig. 2A). Transcripts encoding Eif2s3y (eukaryotic translation initiation factor 2 subunit 3 y-linked), Ddx3y (DEAD box helicase 3 y-linked), Kdm5d (lysine-specific demethylase 5D), and Uty (ubiquitously transcribed tetratricopeptide repeat containing y-linked) were each readily detected in the male nucleus accumbens. The Y-chromosome transcripts expressed at FPKM levels of <0.2 in the male nucleus accumbens were Uba1y, Usp9y, and Zfy1. Consistent with the observation that Sry expression in the adult rodent is limited to catecholaminergic neurons [15, 16], Sry transcripts were not detected in any nucleus accumbens samples. The other 28 transcripts encoded on the Y-chromosome were undetectable in any nucleus accumbens samples [28].

Figure 2.

Expression of Y- and X-chromosome–encoded genes and sex-specific transcripts in the nucleus accumbens. (A) FPKM data for four Y-chromosome–encoded transcripts in the nucleus accumbens of all four groups of saline-injected mice. These data were used to set FPKM = 1 as the lower limit of detection for these samples. Data shown as mean ± SEM. (B) FPKM data for four X-chromosome–encoded transcripts in the nucleus accumbens of all four groups of saline-injected animals. Xist values were not significantly different statistically across the three female samples, and male Xist levels were less than the FPKM detection limit. (C) FPKM data from the nucleus accumbens for six autosome-encoded transcripts known to be sex-regulated in other tissues; note the different y-axis scales. One-way ANOVA with n = 4 in each group. *P < 0.05; **P < 0.01; ***P < 0.001; ****P < 0.0001. F(3,12) values: Eif2s3y, 43; Ddx3y, 28; Kdm5d, 29; Uty, 18; Xist, 18; Eif3s3x, 18; Gpr101, 4.9; Gabrq, 3.9; Cck, 2.1; Slc17a7, 1.0; Oxt, 2.2; Ttr, 9.6; Glp1r, 15; Esr1, 5.9. Statistically significant pairwise comparisons shown (two-tailed Student t test). NS, not significant.

A total of 995 X-chromosome transcripts were identified in mm10, and transcripts encoding 90 of these genes were detected (FPKM >1) in the naive female nucleus accumbens. As expected, Xist (X inactive specific transcript) expression in males was essentially undetectable [30] (Fig. 2B). The Xist transcript levels in female mice did not vary substantially with ovariectomy or estradiol replacement. Eif2s3x (eukaryotic translation initiation factor 2 subunit 3 x-linked) was evaluated so it could be compared with Eif2s3y in males. Eif2s3x transcript levels in the three groups of female mice were approximately double the level observed in the males and were not altered by ovariectomy (Fig. 2B). Escape from X-inactivation is thought to contribute to the roughly twofold greater levels of Eif2s3x and several other X-chromosome–encoded transcripts [30]. Expression of Gpr101 (G-protein coupled receptor 101), 1 of the 10 most highly expressed X-linked transcripts, and Gabrq, a subunit of the major inhibitory neurotransmitter receptor in the vertebrate brain, were also assessed (Fig. 2B). For both of these genes, the transcript levels in the cycling females were substantially greater than those in the ovariectomized females and males, and estradiol treatment did not reverse the effects of ovariectomy.

The expression of six autosomally encoded transcripts known to be expressed in a sex-specific manner in other tissues was next examined in the nucleus accumbens (Fig. 2C). Cck (cholecystokinin) and Slc17a7 (glutamate vesicular transporter; Vglut1) expression was low in cycling females and much greater in ovariectomized females and males. Oxt (oxytocin) and Ttr (transthyretin) transcript levels were substantially greater in the cycling females than in the males and were reduced by ovariectomy [31]. Glp1r [glucagon-like peptide 1 receptor] and Esr1 [estrogen receptor 1] transcript levels were greater than background levels in the nucleus accumbens of the cycling females and decreased in the ovariectomized and placebo mice. Expression of Glp1r and Esr1 was <1 FPKM in the male mice. Based on this limited sampling, genes identified as sex-specific in other tissues often exhibit a sex-specific pattern of expression in the nucleus accumbens (available at: https://www.ncbi.nlm.nih.gov/gene/; a comparison of these data with the data in the National Center for Biotechnology Information has been provided in an online repository [28]).

C. qPCR Validation of RNASeq Data Set

qPCR, which provides a method of verifying and extending observations made using RNASeq, requires the use of a normalizer gene appropriate for the tissue of interest [27, 32, 33]. Analysis of the FPKM data for the 40 mice included in the present study identified Gapdh as the most reliable normalizer tested (Fig. 3A). Using an independent set of RNAs prepared from the nucleus accumbens of similarly treated groups of mice, the levels of several of the X-and Y-chromosome–encoded transcripts were evaluated using qPCR (Fig. 3B and 3C). Analysis of Y-chromosome–encoded transcripts and Xist demonstrated that the signal/noise ratio for RNASeq data was better than for that for the qPCR data. The qPCR analysis of several autosomally encoded sex-responsive transcripts in this new set of RNAs (Fig. 3D) yielded data very similar to those reported in the RNASeq analysis of the first set of RNAs (Fig. 2).

Figure 3.

qPCR validation with separate samples. (A) Mean FPKM data for all 40 mice (pooled by sex and drug treatment) for seven transcripts commonly used as normalizers in qPCR experiments [27, 32, 33] and their relative variance (SD as percentage of the mean). Based on this information, Gapdh was selected as the best normalizer for the nucleus accumbens samples. Using RNA prepared from a separate set of five mice, qPCR data were obtained for the (B) Y-chromosome–encoded transcripts, (C) X-chromosome-encoded transcripts, and (D) autosome-encoded transcripts shown in Fig. 2. Based on Glp1r and Esr1, the qPCR limit of detection was ~0.0008 with respect to Gapdh. (B–D) One-way ANOVAs with n = 4 in each group. *P < 0.05; **P < 0.01; ***P < 0.001; ****P < 0.0001. (B,C) F(3,12): Eif2s3y, 56; Ddx3y, 34; Uty, 6.8; Xist, 18; Eif3s3x, 15; Gpr101, not significant (NS); Gabrq, 3.2. (D) F(2,9): Cck, 15; Slc17a7, NS; Oxt, NS; Ttr, 6.0; Glp1r, NS; Esr1, 4.8. Statistically significant pairwise comparisons shown (two-tailed Student t test). (E) Pairwise correlation of qPCR and FPKM data for the transcripts analyzed in Fig. 2B and 2C in comparisons of different groups: female saline (FS)/male saline (MS); FS/ovariectomized, placebo-treated, saline (OPS); MS/OPS. Y-chromosome–encoded transcripts and Xist were not included in this comparison. The range of values observed using FPKM was larger than the range observed using qPCR, as expected [23]. The slopes of the regression lines were 0.75 (FS/MS), 0.43 (FS/OPS), and 0.57 (MS/OPS).

When relative changes in gene expression determined using RNASeq for the first set of samples and qPCR for the second set were compared, an excellent correlation was found for the data for female to male mice, female to ovariectomized/placebo mice, and male to ovariectomized/placebo mice (Fig. 3E). Reflecting the greater dynamic range of the sequencing approach, the magnitudes of the relative changes in transcript levels seen with RNASeq were larger than those seen with qPCR; a similar difference was observed in a previous study [23].

D. Regulation of Nucleus Accumbens Transcripts by Biological Sex (Female Saline/Male Saline Comparison)

Three programs were used to identify transcripts expressed at substantially different levels in the nucleus accumbens of cycling females and males (Fig. 4A). A total of 340 transcripts were identified; 79 transcripts were identified as differentially expressed by all three programs, and 148 were identified by at least two of the programs (Fig. 4A). Just as in other studies comparing methods of analyzing differential expression [27, 34], DESeq identified the smallest number (n = 102) of transcripts. We focused our analysis on the 148 transcripts identified as differentially regulated by at least two of the three analytical programs used; this group included one-half (75 transcripts) more highly expressed in cycling females and one-half (73 transcripts) more highly expressed in males. Only 15 of the transcripts more highly expressed in females were encoded by the X-chromosome, and only 4 of the transcripts more highly expressed in males were encoded by the Y-chromosome. The complete list of 148 differentially expressed transcripts has been provided in an online repository [28].

Figure 4.

Comparison of FPKM data from cycling female and male nucleus accumbens. (A) Venn diagram of transcripts deemed differentially expressed using Cuffdiff, DESeq2, and Limma. A core of 148 transcripts were identified using two or more methods. (B) Heatmap comparing FPKM data for these 148 transcripts. To combine data over a wide range of FPKM values, data for each transcript were normalized to the average across the eight samples for that transcript and plotted as a logarithm. The 148 columns were ordered from the highest to lowest expression in males. The data cover almost a 10,000-fold range of expression. (C) IPA parameters. Transcripts affected by each upstream regulator are identified; red and green symbols indicate increased and decreased expression, respectively, in our cycling female vs male data set. Colored lines indicate the predicted effect of each upstream regulator on expression of these transcripts (red, activation; blue, inactivation; gray, no clear prediction; dashed, indirect effect). (D–F) Relationships and log2 FPKM ratios shown for the three pathways with the highest z-scores.

A comparison of the four cycling female and four male nucleus accumbens samples from the saline-injected mice (Fig. 4B) revealed good agreement within each set of samples, especially when comparing the left (high in males) and right (high in females) 20% of the samples. For these 148 transcripts, the average relative variance (SD divided by the mean) within each group was 4.2%.

IPA of the group of 148 differentially expressed transcripts identified regulatory networks of potential interest. The three networks with the highest z-scores were driven by transcription factors (Myc, Hif1a, Hic1) [35], and the IPA predictions and actual experimental data were in full agreement (Fig. 4C). Of the 14 transcripts regulated in these three pathways, three (Eno1, Fhl2, Nov) occurred in more than one pathway. The Myc (myelocytomastosis oncogene) network includes Cyth2, an inorganic tripolyphosphatase [36], Dctpp1 (dCTP pyrophosphatase1), a general deoxynucleotide triphosphate pyrophosphatase, which maintains the balance in the cytoplasmic deoxynucleotide triphosphate pool [37], and two common metabolic enzymes, Eno1 (enolase 1), a glycolytic enzyme, and Ldha (lactate dehydrogenase) (Fig. 4D). The higher level of Myc transcript observed in cycling females has been predicted to decrease expression of each of these target genes, as observed.

The Hif1a (hypoxia inducible factor 1a) network includes Fhl2 (four and a half LIM domains 2), a negative regulator of adrenergic signaling [38], Nov (nephroblastoma overexpressed gene), an insulin-like-growth factor-binding-protein [39], Il15 (interleukin 15), and Ctgf (connective tissue growth factor), in addition to the two metabolic enzymes in the Myc network (Fig. 4E). As for Myc, the higher level of Hif1a transcript observed in the cycling females could have contributed to the observed decrease in the levels of all six transcripts (Fig. 4E).

The Hic1 (hypermethylated in cancer 1) network includes Fhl2 and Nov, along with Fst (follistatin), which binds and inactivates activin and other TGFβ family proteins, preventing uncontrolled cellular proliferation and follicle-stimulating hormone secretion, and Tmem47 (also termed TM4SF10 or transmembrane tetraspanin family 10), a claudin-like cell junction protein that is transiently expressed during podocyte development [40] (Fig. 4F). Although Hic1 transcript levels did not differ in cycling females and males, a decrease in its activity would be consistent with the observed increase in the expression of Fst and Tmem47 and the observed decrease in expression of Nov and Fhl2 (Fig. 4F).

E. Effect of Ovariectomy on the Nucleus Accumbens Transcriptome (Ovariectomized, Placebo-Treated, Saline/Female Saline Comparison)

Circulating estradiol and progesterone levels are low in ovariectomized females, and the uterine weight will decrease to <15% of normal (Fig. 1B) [20, 41, 42]. DESeq2, Limma, and Cuffdiff were used to identify transcripts expressed at substantially different levels in ovariectomized/placebo vs cycling female mice. A total of 343 transcripts were recognized by at least two of these programs (Fig. 5A); 236 of these transcripts were more highly expressed in cycling females, and 107 were more highly expressed in ovariectomized/placebo mice. The complete list of 343 differentially expressed transcripts has been provided in an online repository [28]. Transcripts expressed in a sex-specific manner in the preoptic area were also found in these nucleus accumbens samples: Eif2s3x, Esr1, Kdm5c, and Syt4 [30].

Figure 5.

Comparison of FPKM data from ovariectomized/placebo and cycling female nucleus accumbens. (A) Venn diagram of transcripts deemed differentially expressed by Cuffdiff, DESeq2, and Limma; the 343 core transcripts were identified using at least two methods. (B) Cross-correlation of 21 differentially expressed transcripts identified both in the cycling female to male comparison (Fig. 4) and in the ovariectomized/placebo to cycling female comparison (slope of regression line, 0.98). For these 21 transcripts, the direction and magnitude of the changes in expression were very similar in both comparisons. (C–E) IPA of the 343 core differentially expressed transcripts identified in ovariectomized/placebo vs cycling female nucleus accumbens. The three pathways with the highest z-score are depicted. For each, FPKM data are shown.

A comparison of the differentially expressed genes identified when comparing the cycling females and males (148 transcripts) (Fig. 4A) and cycling females and ovariectomized/placebo mice (343 transcripts) (Fig. 5A) identified 21 common transcripts (Fig. 5B; and the online repository [28]). A plot of the transcript ratios for the 21 genes common to these two data sets yielded a straight line with a slope of 0.98 (Fig. 5B). Expression of Slc17a7 (Vglut1) was 16-fold greater in the ovariectomized/placebo females and the males than in the cycling females. Expression of Esr1 (estrogen receptor α) was almost fourfold lower in ovariectomized/placebo females and males than in the cycling females.

Using the 343 transcripts identified as differentially expressed in the ovariectomized and placebo vs cycling females, IPA identified six pathways with activation scores of 1.9 or greater. The three pathways involving the direct action of an upstream regulator on target genes are depicted in Fig. 5C–5E. The Huntingtin (Htt) network includes transcripts encoding Agrn (agrin), an extracellular proteoglycan involved in synaptic development [43, 44] and three collagen genes (Col16a1, Col4a2, Col6a1) expressed together in various parts of the nervous system [45]. The observed decrease in Htt expression in ovariectomized/placebo mice could have caused the observed decrease in the expression of these genes. Other potential Htt targets include Grin2d (NMDA receptor subunit 2D), Hap1 (Huntington associated protein 1), which binds Kalirin [46, 47], Fhl1 (four and a half LIM domains 1), Pdxk (pyridoxal kinase), a gene whose loss is associated with an increased risk of Parkinson disease [48], and Pnoc (prepronociceptin), a peptide precursor strongly implicated in neonatal abstinence syndrome [49]. Expression of each of these transcripts was reduced in ovariectomized/placebo vs cycling females (Fig. 5C). The observed increase in the expression of Hmgb2 (high mobility group box 2), Ndufs3 (nicotinamide adenine dinucleotide hydrate/ubiquinone oxidoreductase subunit 3), and Pcp4 (Purkinje cell protein 4), which modulates calcium-calmodulin signaling and calcium-mediated toxicity [50] could also be mediated by Htt.

Nupr1 (nuclear protein transcription regulator 1) is known to increase the pancreatic β-cell mass and promote the growth of certain tumors [51]. The Nupr1 network includes Agrn, Abcc5 (ATP-binding cassette subfamily c member 5), Bcl9 (B-cell lymphoma 9), Dhtkd1 (mitochondrial function), and Zmym3 (zinc finger MYM-type protein 3), a DNA binding protein with a role in cytoskeletal organization. The levels of each had decreased, as predicted by the observed increase in Nupr1 (Fig. 5D). Fcf1 (RNA processing protein) and Sat1 (spermidine N1-acetyltransferase 1), targeted by anti-Zika virus drugs [52], increased with ovariectomy, as expected.

Ahr (aryl hydrocarbon receptor), a transcriptional activator, has roles in xenobiotic metabolism, cell cycle, and circadian regulation. This network includes the collagen genes found in the Htt network, along with Timp2 (metallopeptidase inhibitor), another component of the extracellular matrix. The levels of each of these transcripts, along with Esr1, declined, as predicted by the IPA (Fig. 5E).

F. Effect of Estradiol on the Nucleus Accumbens Transcriptome of Ovariectomized Mice (Ovariectomized, Estradiol-Treated, Saline/Ovariectomized, Placebo-Treated, Saline Comparison)

The timed-release estradiol pellets used in the present study maintained the uterine weight in the normal range for nearly 2 months (Fig. 1B) [20]. At least two of the three differential expression programs identified 162 differentially expressed transcripts in the ovariectomized/estradiol vs ovariectomized/placebo nucleus accumbens (Fig. 6A); 116 of these transcripts were more highly expressed in the estradiol-treated mice, and 46 were more highly expressed in the placebo-treated mice. The complete list has been provided in an online repository [28].

Figure 6.

Comparison of FPKM data from ovariectomized/estradiol and ovariectomized/placebo nucleus accumbens. (A) Venn diagram of transcripts deemed differentially expressed by Cuffdiff, DESeq2, and Limma. The 162 core transcripts were identified using at least two methods. IPA of these core transcripts identified three pathways with high z-scores. (B) Nupr1 was not differentially expressed. Expression of transcripts encoding (C) Myc and (D) Erg (ETS transcription factor) was reduced in ovariectomized/estradiol mice compared with ovariectomized/placebo mice. (E) Cross-correlation of the 15 differentially expressed transcripts identified both in the ovariectomized to cycling female comparison (Fig. 5) and in the ovariectomized to estradiol-replaced ovariectomized female (slope of regression line, 0.59).

The three pathways identified by the IPA program that consist entirely of direct effects are shown in Fig. 6B–6D. Two of the upstream network controllers were identified in previous pairings. A Nupr1 network was identified in the cycling female to ovariectomized/placebo comparison (Fig. 5D), and a Myc network was identified in the cycling female to male comparison (Fig. 4D). Although the Nupr1 transcript levels were not changed, an alteration in its activity could have contributed to the increased levels of the four transcripts listed. Nupr1 is known to form functional complexes with several transcriptional regulators, including estrogen receptors [51, 53]. Hcn1 (hyperpolarization activated cyclic nucleotide gated potassium channel 1) is being explored as a candidate target for pain treatment [54].

In this Myc network, Acvr2a (activin A receptor 2A) is increased by estradiol treatment, as is a heat shock protein (Hsp90aa1), along with Ipo7 (importin 7), Pdgfra (platelet-derived growth factor receptor-α), and Tiam1 (T-cell invasion and metastasis 1), which is important in dendritic spine function [55] (Fig. 6C).

Erg (ETS transcription factor) was expressed at low levels (1.4–2 FPKM) in the nucleus accumbens samples and declined in ovariectomized mice treated with estradiol. Magi1 (membrane-associated guanylate kinase 1), which was recently found to have copy number variants associated with schizophrenia [56], is a downstream target in two of the estradiol-regulated pathways. Dctpp1 was also identified as a downstream target of Myc in the male-female pairing (Fig. 4D). Phactr2 (phosphatase and actin regulator 2) expression in the hypothalamus is increased in females given ethinyl estrogen [57].

The cross-correlation of the 15 differentially expressed transcripts identified both in the ovariectomized to cycling female comparison (Fig. 5) and in the ovariectomized to estradiol-replaced ovariectomized female (Fig. 6E; slope of the regression line, 0.59) showed that these 15 transcripts were altered by ovariectomy and partially (60%) recovered with estradiol replacement to the cycling female levels. Ten times as many transcripts (n = 147) did not appear in both the ovariectomy and the estradiol replacement transcript lists.

G. Nucleus Accumbens Transcriptome of Cycling Females Is Altered After 4 Weeks of Withdrawal From Chronic Cocaine (Female Cocaine/Female Saline Comparison)

At least two of the differential expression programs used identified 50 transcripts in the nucleus accumbens of cycling females withdrawn from experimenter-administered cocaine 28 days earlier compared with saline-injected control females (Fig. 7A); 40% of the transcripts were more highly expressed in the cocaine-withdrawn females, and 60% were more highly expressed in the saline-injected control mice.

Figure 7.

For cycling female mice, the nucleus accumbens transcripts were altered by withdrawal from chronic cocaine. (A) Venn diagram of transcripts deemed differentially expressed by Cuffdiff, DESeq2, and Limma in the cocaine-injected [female cocaine (FC)] vs saline-injected [female saline (FS)] cycling female comparison. The 50 core transcripts were identified using at least two methods. (B) For each transcript, the average saline value was used to normalize all samples for that transcript; the log2 values were plotted as a heat map using Excel. **Transcripts previously identified as cocaine-regulated in a study of male mice examined 1 wk after withdrawal from a similar treatment with cocaine [23]. (C) IPA of the single statistically significant transcription factor pathway,

To assess any additional variability between the mice caused by cocaine, the data for these 50 transcripts in the individual mice are shown as a heat map (Fig. 7B). For these 50 transcripts, the relative variance (SD divided by the mean) within each group was 24%, much greater than the relative variance in the sets of untreated mice (4%; Fig. 4B). This might reflect the much smaller range of levels of transcripts regulated by cocaine (23-fold in Fig. 7B) compared with the differences between the males and females (>10,000-fold in Fig. 4B). Twenty-five of these transcripts were identified in an earlier study of differential gene expression in the nucleus accumbens of male mice receiving the same chronic cocaine treatment but withdrawn from cocaine for only 7 days [23] (Fig. 7B).

IPA of these 50 differentially expressed transcripts identified Trp53 (cellular tumor antigen, p53), a tumor suppressor, as the only high probability upstream regulator (Fig. 7C). Two of its downstream targets, Igfbp4 and Wnt2, are well-studied hormonal regulators. Igfbp4 binding to insulin-like growth factor alters its activity, and Wnt2 binding to its frizzled-family receptor affects Wnt/β-catenin signaling. Other secretory products included in this group are Fst, Olfm2 (noelin-2), Fgf3, and Fibin (Fin bud initiation factor homolog). In addition, multiple G protein coupled receptors (Htr1a, Adra2a) and ion channels [Hvcn1 (voltage gated hydrogen channel 1), Gabrd (GABA receptor subunit-δ), and Kcnc4 (voltage-gated potassium channel)] appear in this group.

H. A Common Set of Transcripts Remains Altered After Withdrawal From Cocaine in Estradiol- and Placebo-Treated Ovariectomized Mice (Ovariectomized, Placebo-Treated, Cocaine/Ovariectomized, Placebo-Treated, Saline and Ovariectomized, Estradiol-Treated, Cocaine/Ovariectomized, Estradiol-Treated, Saline Comparisons)

Differentially expressed transcripts were also identified in the nucleus accumbens of ovariectomized/placebo and ovariectomized/estradiol mice withdrawn from cocaine vs their saline controls. In the nucleus accumbens of the ovariectomized/placebo group, withdrawal from cocaine affected the expression of 43 transcripts (Fig. 8A). Cocaine withdrawal affected the expression of a larger number of transcripts (n = 156) in the nucleus accumbens of the ovariectomized/estradiol group (Fig. 8B). None of the transcripts regulated by cocaine in ovariectomized females was included in the list of 50 transcripts regulated by cocaine in cycling females (Fig. 7B).

Figure 8.

Effect of withdrawal from chronic cocaine on ovariectomized/placebo and ovariectomized/estradiol mice. Venn diagrams of transcripts deemed differentially expressed after withdrawal from cocaine by Cuffdiff, DESeq2, and Limma in (A) ovariectomized/placebo mice [ovariectomized, placebo-treated, cocaine/ovariectomized, placebo-treated, saline (OPC/OPS)] and (B) ovariectomized/estradiol mice [ovariectomized, estradiol-treated, cocaine/ovariectomized, estradiol-treated saline (OEC/OES)]. When the core cocaine-regulated transcripts (43 in placebo and 156 in estradiol) were compared, only 5 transcripts were found in both data sets. (C) FPKM values for these five transcripts are plotted; based on FPKM and qPCR, withdrawal from cocaine reduced levels of all five transcripts to ∼60% of the corresponding saline value. Two-way ANOVA (placebo or estradiol; saline or cocaine) showed no statistically significant effect of sex hormone and no significant interaction of hormone with cocaine. The main effect of cocaine was significant in every case, F(3,12): Acvr1c, P = 0.049, F = 4.5; Hhip, P = 0.037, F = 5.2; Kcnq3, P = 0.047, F = 4.86; Mme, P = 0.023, F = 6.3; Prkcb, P = 0.048, F = 4.6.

Only five transcripts appeared in both the ovariectomized, placebo-treated, cocaine/ovariectomized, placebo-treated, saline and ovariectomized, estradiol-treated, cocaine/ovariectomized, estradiol-treated saline groups of cocaine-regulated transcripts. Expression of each of these five transcripts was decreased in the ovariectomized mice withdrawn from chronic cocaine. As shown in Fig. 8C, the response to cocaine was remarkably similar for all five transcripts. Two-way ANOVA revealed a main effect of cocaine treatment (P = 0.023-0.048) but no significant effect of sex or hormone status (P = 0.486-0.939) or sex or hormone and cocaine interaction (P = 0.275-0.989; four per group). The average qPCR ratio (cocaine/saline) was 0.61 on the separate set of samples not used for the RNASeq analyses (Fig. 8C; qPCR data not shown), in excellent agreement with the average fold change from the FPKM data (cocaine/saline = 0.57). The levels of each of these five transcripts were all increased in the male mice withdrawn from chronic cocaine for 1 week [23], and they had decreased after cocaine withdrawal in these females.

The transcripts in this group are of special interest. Acvr1c (activin receptor type-1C) encodes a member of the TGFβ receptor subfamily known to be expressed in the brain [58]. Hhip (hedgehog-interacting protein), a secreted glycoprotein that interacts with all three hedgehog family members modulates hedgehog signaling and is widely expressed in fetal and adult tissues (available at: https://www.omim.org/entry/606178). Kcnq3, a voltage-gated potassium channel (also known as Kv7.3), forms heterotetramers thought to be critical for establishing the subthreshold excitability of neurons and their synaptic inputs [59]. Mme (neutral endopeptidase 24.11; also called neprilysin or enkephalinase, EC:3.4.24.11), a type II membrane protease, degrades opioid peptides such as Met- and Leu-enkephalin, along with angiotensins and atrial natriuretic peptide [60, 61]. Prkcb (protein kinase C-β, EC:2.7.11.13], a calcium-activated, phospholipid- and diacylglycerol-dependent serine/threonine-protein kinase, is involved in insulin and many other signaling pathways.

3. Discussion

Activation of the dopaminergic pathway from the ventral tegmental area to the nucleus accumbens is an essential first step in the process of addiction [62]. Although not a perfect model for addiction in humans, animal models have provided important insights into the complex effects of cocaine on the nervous system, exhibiting behavioral and biochemical changes after lengthy withdrawal from daily cocaine exposure, such as seen with human patients [8, 26, 63, 64]. When rodents are exposed to cocaine multiple times in a novel environment, they develop locomotor sensitization. However, even among inbred rodents, the responses have varied, and only those rodents that show locomotor sensitization exhibit the biochemical readouts indicative of drug addiction [22]. Locomotor sensitization to cocaine does not occur in mice lacking the dopamine plasma membrane transporter (DAT) or in mice lacking the D1 dopamine receptor (Drd1a), and these mice exhibit weak or no acquisition of cocaine self-administration [65-69]. Tests of cocaine self-administration in rodents, which might provide the best insight into human addiction, indicate an essential role for the Drd1a receptor [70]. Elimination of Drd2 yields increased self-administration and increased locomotor sensitization to chronic cocaine injections [69]. Transcripts encoding the dopamine receptors and DAT were not among those identified as differentially regulated at baseline or after withdrawal from cocaine (cycling female vs male; ovariectomized; placebo or estradiol pellets).

Sex-specific baseline differences in gene expression are key to understanding the effects of drugs of abuse. With only one other analysis of female nucleus accumbens RNASeq [71], it is clear that more data on male/female differences in specific brain regions are needed. RNASeq studies can provide an overwhelming amount of transcriptomic data; however, integrating the results from different studies has been difficult. For example, RNASeq analyses performed 2 hours, 24 hours, and 7 days after cocaine exposure yielded profoundly different results [72]. Differences in the depth of sequencing analysis and the completeness of the reporting of the sequencing data have made comparisons even more challenging. A major limitation of the present study was the limited number of mice in each of the eight groups (five mice for treatment, four mice for RNASeq), less than one-half the number included in studies with fewer groups of experimental animals [23, 73].

Nine studies report RNASeq data for the nucleus accumbens, striatum, or prefrontal cortex of mice [13, 19, 23, 26, 72–75] or monkeys [76]; however, not all had been treated with psychostimulants. The number of differentially expressed transcripts varied widely, with 4513 identified in cocaine-treated mice [74] and 328 identified in monkeys [76]. Although self-administration [19, 26, 76] is more like the human situation, the more common experimenter-administered approach [13, 23, 72, 74, 75] allows for analysis of more animals, a key to the reliable detection of reproducible differences. The drug dosages have also varied, as have the times of administration (7 to 100 days) and withdrawal (2 hours to 30 days).

A. Differences in the Nucleus Accumbens Transcriptome of Naive Female and Male Mice

Our understanding of the reasons the responses of men and women to many drugs of abuse differ is quite limited [8–10, 19]. Women and female experimental animals have been shown to be more sensitive to mood-altering drugs than men and male experimental animals by several criteria, including a more marked initial physiological response, quicker progression to repeated self-administration, and faster reacquisition of drug-seeking behavior after a prolonged withdrawal period [1–4, 13, 71, 77]. The focus of these studies has been the identification of differences in gene expression between naive female and male rodents and after prolonged cocaine withdrawal. Because the neurochemical mechanisms underlying these responses can differ in females and males, a treatment found to be beneficial for males might not be beneficial for females [9, 10, 13]. For example, different cell death pathways are triggered by stroke in females and males [21].

To understand why males and females respond to cocaine in such a different way, we first tested the hypothesis that baseline gene expression differed substantially in the nucleus accumbens of normally cycling female and male mice, because it is commonly accepted that differential gene expression underlies most sexually dimorphic traits [77]. A total of 148 genes were differentially expressed in the cycling female vs the male nucleus accumbens, with one-half greater in the female and one-half greater in the male. The differentially expressed genes encode metabolic enzymes such as Eno1 and Ldha, along with secreted factors such as Il15, connective tissue growth factor, and follistatin. The baseline differences identified in multiple signaling networks suggest that the initial response of the female brain to cocaine might differ substantially from that of the male brain.

At baseline, only 10% of the genes encoded on the X- and Y-chromosomes were expressed in the adult nucleus accumbens. The X-chromosome encodes 20% of the genes that are more highly expressed in females than in males. A comparison of gene expression in the nucleus accumbens of ovariectomized mice to gene expression in cycling females identified 343 differentially regulated transcripts. Two-thirds were more highly expressed in cycling females. Y chromosome transcripts (Ddx3y, Eif2s3y, Kdm5d, Uty) were abundant in the male NAc, as in the few previous studies that reported these transcripts [23, 71, 73], as were Xist transcripts in females [71]. Similar levels of the male-specific transcripts (Ddx3y, Eif2s3y, Kdm5d, Uty), common X-chromosome transcripts (Eif2s3x, Gabrq, Gpr101), and sex-regulated autosomal transcripts (Cck, Esr1, Glp1r, Slc17a7, Ttr) from Fig. 2 were all identified by Bottomly et al. [73] in the striatum, which includes the NAc. It is unclear why many of these transcripts have been largely absent from a large fraction of the reported data on the NAc transcriptome [19, 74, 75].

The set of 21 differentially expressed transcripts common to cycling females vs males and cycling females vs ovariectomized females exhibited an extremely high cross-correlation (R = 0.96), suggesting regulation by similar factors in ovariectomized females and males. This group of transcripts included Esr1 and Vglut1 (Slc17a7). With 16-fold lower levels of Vglut1 transcript in the nucleus accumbens of cycling females than in the nucleus accumbens of ovariectomized or male mice, one would expect to find functionally important differences in glutamatergic signaling by medium spiny neurons. Similarly robust changes in the expression of Cck and Nov would also be expected to affect signaling.

We next evaluated the ability of estradiol to counteract the changes in gene expression observed after ovariectomy. The estradiol implants used maintained the uterine weight at the level observed in cycling females [20]. A total of 162 differentially expressed transcripts were identified in the nucleus accumbens of estradiol-treated vs placebo-treated ovariectomized mice. Pathway analysis of these transcripts highlighted several of the same pathways identified in the comparison of cycling females and ovariectomized/placebo and male mice. Only 15 of the 162 transcripts were also identified in a comparison of ovariectomized/placebo mice vs cycling females; for this subset of transcripts, the ratios were highly correlated (R = 0.85), perhaps suggesting common regulatory pathways.

B. Sex-Dependent Effects of Prolonged Cocaine Withdrawal on the Nucleus Accumbens Transcriptome

A total of 50 differentially expressed transcripts were identified in the nucleus accumbens of cycling females withdrawn from cocaine for 28 days. One-half of these genes were identified in an earlier study of differential gene expression in the nucleus accumbens of male mice withdrawn from cocaine for 7 days [23]. Many of the genes whose long-term expression was altered in response to cocaine were clearly affected in both males and females. Pathway analysis identified Trp53, a tumor suppressor, as a potential regulatory factor; several of its downstream targets (Igfbp4 and Wnt2) and other members of this group are secretory products [Fst, Olfm2 (noelin-2), Fgf3, and Fibin].

In a previous study [23], cocaine withdrawal was shown to decrease Eif2s3y expression and increase Uty, Kdm5d, and Ddx3y expression in male mice. However, the smaller number of animals in our study precluded substantial detection of these changes. In the present study, pathway analysis highlighted changes in the extracellular matrix components such as agrin, collagen, and Timp2, along with a few molecules relevant to signaling such as Pnoc (prepronociceptin) and receptor subunits such as Grin2D that would be expected to play a role in signaling in the nucleus accumbens.

None of the transcripts affected by cocaine withdrawal of ovariectomized/progesterone (43 transcripts) or ovariectomized/estradiol (156 transcripts) mice appeared in the list of 50 transcripts affected by cocaine withdrawal of cycling female mice. The expression of five transcripts was decreased by cocaine withdrawal in both groups of ovariectomized mice. All five transcripts exhibited an increase in expression in male mice withdrawn from cocaine for 1 week [23]. Included in this intriguing group of genes is a member of the TGFβ family, a hedgehog interacting protein and neprilysin (enkephalinase), which degrades opioid peptides. Secretory proteins whose expression remains altered after withdrawal from cocaine might prove useful as biomarkers for monitoring changes that occur during withdrawal.

4. Conclusions

Our studies have revealed the existence of profound sex-specific differences in gene expression in the nucleus accumbens of naive cycling female and male mice. These baseline differences will play a critical role in sex-specific differences to an initial exposure to cocaine and to the long-lasting effects of drug withdrawal. Many of the genes whose expression was altered in male mice withdrawn from cocaine for 1 week [23] were affected in a similar manner in female mice withdrawn from cocaine for 28 days, providing a list of target genes that merit further study. Our data have demonstrated the existence of sex-specific, cocaine-sensitive metabolic pathways, secreted factors, and extracellular matrix components in the mouse nucleus accumbens that can be used to guide the development of sex-specific approaches to drug addiction.