EDCs Reorganize Brain-Behavior Phenotypic Relationships in Rats.

All species, including humans, are exposed to endocrine-disrupting chemicals (EDCs). Previous experiments have shown behavioral deficits caused by EDCs that have implications for social competence and sexual selection. The neuromolecular mechanisms for these behavioral changes induced by EDCs have not been thoroughly explored. Here, we tested the hypothesis that EDCs administered to rats during a critical period of embryonic brain development would lead to the disruption of normal social preference behavior, and that this involves a network of underlying gene pathways in brain regions that regulate these behaviors. Rats were exposed prenatally to human-relevant concentrations of EDCs (polychlorinated biphenyls [PCBs], vinclozolin [VIN]), or vehicle. In adulthood, a sociosexual preference test was administered. We profiled gene expression of in preoptic area, medial amygdala, and ventromedial nucleus. Prenatal PCBs impaired sociosexual preference in both sexes, and VIN disrupted this behavior in males. Each brain region had unique sets of genes altered in a sex- and EDC-specific manner. The effects of EDCs on individual traits were typically small, but robust; EDC exposure changed the relationships between gene expression and behavior, a pattern we refer to as dis-integration and reconstitution. These findings underscore the effects that developmental exposure to EDCs can have on adult social behavior, highlight sex-specific and individual variation in responses, and provide a foundation for further work on the disruption of genes and behavior after prenatal exposure to EDCs.

Environmental contamination with endocrine-disrupting chemicals (EDCs) perturbs hormones and their actions in virtually all species and ecosystems [1]. Prenatal EDC exposures pose a particular risk due to the extreme sensitivity of the developing brain to gonadal hormones, which are required for sex-typical differentiation and development of neural circuits, and the manifestation of behaviors. In the hypothalamus of male rodents and other mammals, prenatal and early postnatal testicular hormones masculinize and defeminize neural circuits. In females, the ovary is relatively quiescent with low gonadal hormone secretion; together with alpha-fetoprotein preventing estrogen crossing the blood–brain barrier, the brain undergoes feminization and demasculinization under these female-typical conditions [2-5].

The effects of developmental EDC exposure on sexually dimorphic social behaviors and gene expression patterns in different brain regions have been described for several classes of chemicals. Although individual EDCs are not pure hormone agonists or antagonists, some (such as certain polychlorinated biphenyls [PCBs] and bisphenol A [BPA]) mimic or disrupt estrogen signaling [6], and others (vinclozolin [VIN] and phthalates) are antiandrogenic [7, 8]. Polychlorinated biphenyls, widespread industrial chemical contaminants, alter gene expression in the hypothalamus [9-11] and change interactions of adult rats with conspecifics [12-16]. Bisphenol A from plastic, and the fungicide VIN also change brain gene expression (BPA: [10, 17–19]) and sociosexual behavior (BPA: [20-22], VIN: [23, 24]). Prenatal exposure to phthalates causes long-lasting changes to gene expression in the hypothalamus and beyond [25, 26]. Phthalate exposure early in life also cause deficits in cognitive and social behaviors [25, 27]. In most cases, outcomes are dependent on the dose, timing, and length of exposure, as well as the sex of the animal. This is not surprising considering the dynamic nature of endogenous hormone signaling as the brain develops, and the vulnerability of estrogenic and androgenic pathways to EDCs.

Reproductive success is contingent upon sex-appropriate differentiation of the brain during early life. For an individual to reproduce successfully, appropriate dyadic interactions with another sexually mature potential mate of the opposite sex are required. This process involves assessment of an opposite-sex animal’s fitness through a variety of physical and behavioral cues, including hormonal status, as well-documented in rats [28-30]. There is plasticity in this behavior, with the decision-making process affected by prior sexual experience of both individuals, estrous cycle stage, hormone levels, and other factors. Within the brain, a complex social decision–making network [31] comprising hypothalamic [eg, ventromedial nucleus (VMN), preoptic area (POA)] and extra-hypothalamic (eg, medial amygdala [MeA]) regions expresses specific genes and proteins that modulate these behaviors [32].

Here, we tested the hypothesis that prenatal EDC exposures would cause disruptions to the pattern of expression of a suite of genes in 3 brain regions in the social decision–making network (VMN, POA, MeA) and that this underlies functional deficits in an ethologically-relevant sociosexual behavioral task. Previous work has shown that EDCs have small but significant effects on individual phenotypic traits, but the complex inter-relationships of these phenotypes have not been examined in detail. Therefore, the goal of the current study was to determine whether these relationships would break down (become “dis-organized”) and/or become reconstituted into novel patterns, in a sexually dimorphic manner.

Methods and Materials

Experimental design

All rat procedures were conducted in compliance with protocols approved by IACUC at The University of Texas at Austin. Sprague-Dawley rats purchased from Envigo (Houston) were housed in colony rooms with a consistent temperature (22°C) and light cycle (14:10 dark:light, lights off at 11:00 am). All rats had ad libitum access to water and were fed a low phytoestrogen rat chow (Teklad 2019; Envigo, Indianapolis, IN).

To generate experimental rats, virgin females were mated with sexually experienced males. Successful mating was indicated by the presence of sperm in a vaginal smear. The day after mating overnight was termed embryonic day 1 (E1). Pregnant rats received intraperitoneal (i.p.) injections of 1 of 3 treatments daily from E8-E18: (1) Vehicle (6% dimethylsulfoxide (DMSO) in sesame oil), (2) the PCB mixture Aroclor 1221 (A1221, 1 mg/kg), or (3) VIN (1 mg/kg). Each dam was exposed to the same treatment daily and received a total of 11 injections. The route, timing of treatment, and the dosages were selected to match prior work, based on ecological relevance, and to span the period of hypothalamic neurogenesis, fetal gonadal development, and the early stages of brain sexual differentiation [24, 33–36]. The choice of i.p. injection was made to enable us to replicate prior work [24, 33, 34, 37, 38]; since the time that this work was initiated we have switched to feeding the EDCs to pregnant rats on a cookie [39, 40] and have had similar results.

Aroclor 1221, while best studied as a weakly estrogenic PCB mixture due to its lightly chlorinated structures [6], can also act to disrupt androgen, thyroid, and other hormonal and nonhormonal (eg, neurotransmitter) pathways [41-43]. Vinclozolin is considered antiandrogenic, as some but not all of its effects are mimicked by the antiandrogenic pharmaceutical, flutamide [7, 8]. Because neither A1221 nor VIN is a pure hormone receptor agonist or antagonist, we did not include a “positive” control for estrogenic (A1221) or antiandrogenic (VIN) effects.

A subset of offspring (30 male and 29 female) from 9 DMSO, 10 A1221 (PCB), and 10 VIN dams were included in this study. No more than 2 same-sex rats per litter were used. We measured body weight and anogenital distance (AGD) on days P7 and P14 to calculate the Anogenital Index (). The 5 males and 5 females with the median intrasex AGI measurements were used for the subsequent experiments. The pups were weaned at P21 and re-housed in same-sex groups of 2 to 3. Beginning on the day of vaginal opening, daily vaginal smears were collected from females and cell cytology was examined as an indication of estrous cyclicity. Timing of pubertal development did not vary across treatment groups (analysis of variance [ANOVA]; male age at preputial separation: DMSO 43.5 ± 0.5, PCB 44.2 ± 0.7, VIN 44.1 ± 0.7. Female age at vaginal opening: DMSO 36.4 ± 0.7, PCB 35.3 ± 0.9, VIN 34.6 ± 0.6). Rats were euthanized at ~P120 between 8:00 and 10:00 am (lights out at 11:00 am), with females in proestrus, by rapid decapitation and brains removed and processed as described below.

Stimulus Sprague-Dawley rats for the mate preference test were purchased as virgin adults. Males were gonadectomized (GDX) and females ovariectomized (OVX) under isoflurane anesthesia in aseptic conditions [44, 45]. During the surgery, stimulus animals assigned to the hormone-replaced group also had a 1.5 cm silastic capsule containing testosterone (males: 100% testosterone (T); GDX + T) or a 1.0 cm silastic capsule containing 17β-estradiol (females: 5% E2/95% cholesterol; OVX + E2) implanted subcutaneously into the nape of the neck [44]. All rats recovered from surgery for at least 1 week prior to use in behavioral tests. Thirty-two GDX males (no hormone replacement), 32 GDX + T males, 32 OVX females (no hormone replacement), and 32 OVX + E2 females were used as stimuli throughout the study. For the latter group, on the day of use, these E2-treated females were primed for sexual receptivity by a subcutaneous injection of progesterone (P4, 0.6 mg) in sesame oil 4 hours before the experiments started.

Sociosexual preference behavior

A 1 x 1 m 3-chambered apparatus (Stoelting, Wood Dale, IL, USA) was used as the testing arena [15, 46]. Testing was conducted under dim red light approximately 2 hours into the dark phase of the light-dark cycle. Each test utilized an experimental (EDC or vehicle exposed) rat at ~3 months of age. This paradigm was selected as a standard model of sexual preference wherein a rat is given a choice between opposite-sex partners of varying attractiveness [11]; we previously found this behavior to be affected by prenatal EDCs [16]. Two opposite-sex stimulus rats, one with and one without hormone replacement, were used. Each stimulus rat was placed inside a 7 x 15 cm cylindrical cage positioned in the 2 far opposite corners of the apparatus. These cages have spaced vertical bars, allowing for limited tactile interactions between rats. The position of stimulus rats was randomized between trials and with respect to hormone status. The bars of the stimulus cage allowed for visual, olfactory, auditory, and minimal tactile interaction between the confined stimulus rat and the freely-moving experimental rat. Each trial began with the 2 stimulus rats already in position in their cylindrical cages. An experimental rat was placed in the center chamber of the apparatus with closed doors preventing entry into either side chamber for a 5-minute habituation period. After habituation, the doors were removed and the experimental rat was allowed to freely explore the entire arena for 10 minutes. Each test was recorded by overhead video. ANY-Maze (Stoelting, Wood Dale, IL, USA) was used to track the position, speed, and distance traveled of the experimental rat in each compartment of the chamber [16, 44, 46]. Recordings of the tests were scored by a trained investigator blinded to treatment for the following behaviors: nose touching (direct nose-to-nose contact between the experimental rat and a stimulus rat) and stimulus investigation (all other investigation by the experimental rat of a stimulus rat or stimulus cage). The time the experimental animal spent within 1 body length of either stimulus cage without engaging with the stimulus animal or cage (time within one body length [time nose touching + time investigating]) was defined as “time near.” To avoid testing fatigue, stimulus rats were used for no more than 3 rounds of testing per day and had 10 minutes of rest with access to food and water between each round. Stimulus rats had 2 days of rest between each day of testing. The entire apparatus was cleaned using 70% ethanol between each test subject.

Hormone radioimmunoassay

Serum concentrations of testosterone and corticosterone (CORT) were measured in duplicate samples, and estradiol (E2) in single samples (due to larger serum volume needed for this assay) using radioimmunoassays (testosterone: MP Biomedicals #07189102; CORT: MP Biomedicals #07120102; E2: Beckman Coulter #DSL-4800). Assay parameters were: CORT, limit of detection 7.7 ng/ml, intra-assay coefficient of variance (CV) 2.5%; testosterone, limit of detection 30 pg/ml, intra-assay CV 3.7%; E2: limit of detection 2.2 pg/ml, intra-assay CV 16.8%.

TaqMan Low Density qPCR Array

TaqMan low density array (TLDA) cards were designed and purchased (ThermoFisher Scientific, Waltham, MA) to include genes that fell into functional categories such as sex steroid hormone signaling, glucocorticoid stress axis, nonapeptides, gonadotropin-releasing hormone (GnRH)-related genes, neurotrophins, neurotransmission, epigenetics, and clock genes and other transcription factors (Supplemental Table 1 [47]).

Brains from experimental rats were rapidly removed, chilled on ice, and then coronally sliced at 1 mm using a chilled brain matrix. These slices were placed on slides and stored at -80 until all samples were collected. Bilateral punches were taken of the POA, MeA, and VMN using a 1-mm Palkovits punch [48]. Ribonucleic acid (RNA) from frozen POA, MeA, and VMN punches was extracted using AllPrep RNA/DNA Mini Kit (80204; Qiagen, Germantown, MD) according to the manufacturer’s protocol. To determine the integrity and purity, a subset of samples was run on a Bioanalyzer 2100 (Agilent, RNA Pico Kit 5067-1513). All samples had an RIN of 8.4 or above. RNA (200 ng) was then converted to single stranded complementary DNA using high-capacity complementary DNA reverse transcriptase kit (4374966; Life Technologies, Carlsbad, CA) according to the manufacturer’s protocol. Run parameters for the qPCR TLDA cards were: 95°C for 10 minutes, 50 cycles of 95°C for 15 seconds, and 60°C for 1 minute [11]. Gene expression cycle threshold (Ct) values were normalized using the ΔΔCt method. First, each target gene value was normalized to the expression level of the reference gene Gapdh within each subject to generate ΔCt. To standardize between subjects, the ΔCt of each gene was normalized to the median value of a control group (DMSO females) to generate ΔΔCt. Data are reported as . Two genes (Cyp11a1 & Hsd3b1) did not amplify and were excluded, leaving 44 target genes and 2 housekeeping genes (Gapdh, 18s). In all cases, significance was set at P < 0.05 after appropriate corrections for multiple comparisons.

Statistics

Gene expression data were analyzed by ANOVA, with corrections by Sidak’s multiple comparisons test for multiple comparisons. For behaviors, analyses were performed separately for each sex. Those behaviors involving a choice based on the hormone status of the stimulus rat were analyzed by a two-way ANOVA (treatment x stimulus hormone status). Other behaviors (eg, center time, distance traveled of experimental rat) were analyzed by one-way ANOVA. To explore sex differences, a two-way ANOVA for treatment x sex was used for stimulus-independent behaviors between the sexes. Hormone concentrations, body weight, and puberty timing within each sex were analyzed by one-way ANOVA. Reported p-values of multiple comparisons were adjusted using Sidak’s multiple comparisons test.

A hormone preference score was calculated as . Social preference was calculated as a ratio of time within one body length of both stimulus animals out of the total test time:. Linear regressions were used to determine correlations with significantly nonzero slopes.

Principal components analysis

To characterize coordinated phenotypic response to EDC exposure, we performed a Principal Components Analysis (PCA) on morphological, physiological, and behavioral measures, including body weight, CORT, E2, T (the latter males only), activity, social preference, hormone preference, and social activity (time spent investigating and interacting with stimulus rats) using the prcomp function in R. All variables were centered and scaled prior to PCA.

Weighted gene co-expression network analysis

To capture coordinated gene expression changes associated with EDC exposure we performed a Weighted Gene Co-expression Network Analysis (WGCNA) on the 44 target genes measured [49, 50]. WGCNA was performed independently for the 2 sexes and 3 brain regions, with a minimum module size of 5 genes. Expression values of each module were summarized as module eigengenes (ie, the first principal component of each gene co-expression module). Thus, each eigengene is the linear combination of gene expression values that explains the most variation in the expression levels of the genes contained in the module. We assessed coordinated changes in phenotypes and gene expression across treatments using general linear models.

Co-variance patterns between neural gene expression, physiology, behavior, and morphology

To assess systems-level response to EDC treatment and the potential dis-integration and/or reconstitution, defined as loss or change of correlations by EDCs, respectively, we examined co-variance patterns among all neural gene expression and physiological, behavioral, and morphological measures for each control and treatment conditions separately and visualized changes in the correlation structure across treatments. We calculated Spearman’s rank correlations between all pairwise variables. Variables were clustered using 1-correlation scores as distance variables. To visually assess the extent of integration or re-organization (or lack thereof) for each EDC treatment across levels of biological organization, neural gene expression and phenotype clustering of the control condition was maintained for treatment animals for each sex and brain region.

Results

Embryonic exposure to EDCs affected sociosexual behavior in a sex-dependent manner

The mate preference task was performed on 29 females (9 DMSO, 10 PCB, and 10 VIN) and 30 males (10 DMSO, 10 PCB, 10 VIN). Twelve male rats (4 DMSO, 4 PCB, and 4 VIN) failed to investigate both of the stimulus rat options during the allotted 10 minutes. We refer to these males as “nonresponders” in all subsequent analyses. All males were included regardless of responder status in analyses of sex differences, PCA analysis, and gene expression analysis. However, Fig. 1 shows analyses of only the responder males (6 DMSO, 6 PCB, and 6 VIN), as that test required rats to interact with both opposite-sex stimulus animals to calculate a score. Other figures are inclusive of the entire cohort of males, regardless of responder status.

Figure 1.

Time spent with each stimulus rat (social time) during mate preference is shown for females (A) and males (B) as median (bar within the box), quartiles (upper and lower limits of the box), and range (whisker) for the 10-minute mate preference test. The same graphing conventions are used for other box-and-whisker graphs. Data for males include responders only. Asterisks indicate a significant difference between the social time spent with the two stimulus options within a group. Two-way ANOVA was followed by Sidak’s multiple comparisons test. ** P < 0.01. Abbreviations: E2, estradiol; GDX, gonadectomized male; OVX, ovariectomized female; T, testosterone; P4, progesterone.

There was a main effect of hormone status of the stimulus rat on the time experimental females spent associating (within 1 body length) with the stimulus rats (F(1, 52)  = 18.77; P < 0.0001). Females prenatally exposed to DMSO or VIN preferred the hormone-replaced male. Prenatal exposure to PCB, however, abolished this preference in females (Fig. 1A). There was no effect of treatment on the total time that the experimental rats spent investigating both stimulus rats (Social Time).

In males, there was a significant interaction between treatment and hormone status of the stimulus rat on the time spent associating with the stimulus rats (F(2, 30)  = 7.113; P < 0.01). Males exposed prenatally to DMSO spent more time investigating the stimulus female with hormone replacement over the one without. However, the time males exposed prenatally to PCB or VIN spent near the 2 stimulus rat options did not differ significantly within each treatment group (Fig. 1B). These findings replicated those in our recent publication [16].

There were significant sex differences in several behavioral measures (Fig. 2). Throughout the test duration, females traveled significantly farther than males (Fig. 2A; F(1, 52)  = 55.94; P < 0.0001). Females spent more time in close proximity to but not interacting with both stimulus rats (Fig. 2B; F(1, 52)  = 10.49; P < 0.01). Females spent more time directly investigating the stimulus rats than males (Fig. 2C; F(1, 52)  = 10.43; P < 0.01). For hormone preference score, females preferred the hormone-replaced stimulus animals more strongly than did males (Fig. 2D; F(1, 52)  = 13.21; P < 0.001). VIN males spent less time in close proximity to the stimulus rats without interacting than DMSO males, shown in Fig. 2E as the social preference score (F(1, 52)  = 20.66; P < 0.0001).

Figure 2.

Sexually dimorphic behaviors in the 10-minute mate preference test are shown for all rats. A: Total distance traveled. B: Total time the experimental rat spent near (within one body length) both stimulus cages without interacting with the stimulus rats. C: Total time the experimental rat spent exploring (sniffing, touching, etc.) the stimulus cages and stimulus rats (but not nose-touching). D: Hormone preference score (social time with hormone-replaced stimulus rat/social time with both stimulus rats). E: Social preference score (social time with both stimulus rats/time in the remote portions of the side chambers [further than one body length away from the stimulus cage]). For males, nonresponders (rats who failed to venture near 1 or both stimulus options) are indicated with solid black circles, here and in subsequent figures. Two-way ANOVA was followed by Sidak’s multiple comparisons test. Main effects of sex are indicated as: ** P < 0.01, **** P < 0.0001. Effects of treatment within sex are indicated as + P < 0.05, VIN < DMSO.

EDCs did not affect circulating steroid hormone levels, but PCBs resulted in reduced body weight in males and females

We measured body weight and serum hormone concentrations (CORT, E2, and T) of behaviorally characterized rats on the day of euthanasia (Fig. 3). Prenatal PCB exposure significantly reduced the body weight of both female (F(2, 26)  = 3.506; P < 0.05; Fig. 3A) and male (F(2, 27)  = 6.080; P < 0.01; Fig. 3D) rats. There were no significant effects of treatment on hormone concentrations within males or females (Fig. 3B and 3C [females], Fig. 3E–3G [males]). In females, there was a nonsignificant trend for PCB exposure to increase E2 levels (F(2, 26)  = 2.916; P  = 0.07; Fig. 3C).

Figure 3.

Body weight and serum hormones are shown for female (A–C) and male (D–G) experimental rats. Note differences in y-axes for body weight and serum E2 for the two sexes. Asterisk indicates significant difference from DMSO. One-way ANOVA, main effect of treatment followed by Sidak’s multiple comparisons test. * P < 0.05, ** P < 0.01. Abbreviations: DMSO, dimethylsulfoxide; CORT, corticosterone; E2, estradiol; T, testosterone.

Integration across levels of biological organization revealed high levels of individual variation with treatment groups

A principal components analysis of behavioral and physiological measures in females (12 measures) and males (13 measures) revealed variation between the sexes in coordinated phenotypic response to EDC treatment, as the principal components (PCs) characterized different aspects of the phenotype (Fig. 4). Eigenvectors describing the loadings, or contribution, of phenotypes to PC variation revealed patterns of coordinated phenotypic variation that differed in strength (absolute value of the eigenvector) and directionality (positive or negative).

Figure 4.

Principal components analysis (PCA) of behavioral and physiological measures in females (A, B) and males (C, D) revealed variation between the sexes in coordinated phenotypic response to EDC treatment. The contribution of each behavioral, morphological, and physiological measure to each principal component (PC) is indicated by the eigenvectors, and the percentage each PC contributed to the variance is indicated on the x- and y-axes. Definitions of phenotypes are: StimExpl, stimulus explore time (seconds); Dist, total distance traveled (m); CentTime, center time in the open field (seconds); SocTime, total social time (within one body length distance to the stimulus rats in the chambers); SocPref, social preference; NonSocTime, total nonsocial time (time spent greater than one body length away from the stimulus rats in the side chambers); NearTime, total time near stimulus rats (seconds); HormPref, hormone preference; NoseTouch, total nose touch time (seconds); Cort, corticosterone.

In females, the first 4 PCs described 78% of the total variation in morphology, behavior, and physiology within and across treatments (Fig. 4A and 4B). Social Time, Social Preference, and Stimulus Explore loaded strongly and concordantly on PC1 (35%), indicating that PC1 primarily represents variation in time spent engaging in sociosexual behavior. With strong and concordant loadings of Hormone Preference and Near Time, PC2 (20%) represents Social Preference and Social Interaction. PC3 (12%) is strongly loaded by sex steroid (E2) levels and time spent near a stimulus rat with opposing effects. Finally, with strong and opposing loadings of CORT and body weight, PC4 (11% of the variation) may be an indicator of condition. Within treatment, females varied in their integrated response to EDC treatment; however, there were no significant differences across treatments for the first 4 PCs (Supplemental Figures 1 and 2 [47]).

In males, the first 4 PCs described 75% of the total variation in morphology, behavior, and physiology within and across treatments (Fig. 4C and 4D). PC1 (38%) primarily described the variation among males in the time they spent in the center of the apparatus, and thus represents differences between responder and nonresponder males. PC2 (15% of the variation) primarily represented, in opposing fashion, CORT levels and time spent engaging in nonsocial activity. PC3 (12% of the variation) characterized opposing variation in body weight and T among males. With strong loadings of nose touch, stimulus explore, and hormone preference, PC4 (10% of the variation) characterized variation in preference and social interaction. Because nonresponders failed to approach 1 or both of the stimulus rat options during the allotted 10 minutes biasing their hormone preference and other sociosexual behavioral scores, responder and nonresponder males differed across PC4.

Genes and their relationships in VMN, POA, and MeA

The effect of EDC treatment and sex on the expression of the 44 detectable candidate genes in 3 brain regions was examined using the TLDA qPCR array. A subset of genes was significantly affected by treatment after correction for multiple comparisons (Table 1). Raw data for the entire dataset are provided in Supplemental Table 1 [47].

Table 1.

Genes significantly affected by treatment

Females	Gene	DMSO			PCB			VIN
Region		Mean	SEM	N	Mean	SEM	N	Mean	SEM	N
VMN	Cyp19a1	2.14	0.80	9	11.93	6.99	10	2.96	0.97	10
	Hsd17b1	2.94	1.94	9	1.63	0.50	10	7.74	5.87	10
	Oxt	0.93	0.18	9	5.74	2.23	10	4.50	1.75	10
	Avp	0.93	0.16	9	4.40	2.01	10	2.87	1.04	10
	Kiss1	2.75	1.10	9	6.38	1.78	10	5.00	1.35	10
POA	Grin2b	1.04	0.06	9	4.26	2.18	10	2.60	1.53	10
MeA	Oxt	3.15	1.37	9	1.83	0.70	10	4.78	1.72	10
	Kiss1	1.06	0.32	9	2.85	0.45	10	1.68	1.12	10
Males	Gene	DMSO			PCB			VIN
Region		Mean	SEM	N	Mean	SEM	N	Mean	SEM	N
VMN	Cyp19a1	23.68	14.16	10	3.89	1.13	9	4.81	1.47	10
POA	Grin2b	2.19	1.19	10	1.83	0.78	10	5.22	4.11	10
MeA	Kiss1	6.08	1.57	10	5.84	1.22	10	4.94	0.97	10

Bolded numbers indicate significantly different from DMSO after adjusting for multiple comparisons. Abbreviations: DMSO, dimethylsulfoxide; MeA, medial amygdala; PCB, polychlorinated biphenyl; POA, preoptic area; VIN, vinclozolin; VMN, ventromedial nucleus.

In the VMN of females, PCB-exposed rats had higher expression of Cyp19a1, Oxt, Avp, and Kiss1 than DMSO females. In VIN females, Hsd17b1 and Oxt were higher than levels in DMSO females. In the female POA, only Grin2b was changed significantly by EDCs (higher in both PCB and VIN than DMSO females). Two genes were affected in the female MeA: Kiss1 (PCB > DMSO) and Oxt (PCB < DMSO; VIN > DMSO).

Three genes were significantly affected in males, 1 in each region (Table 1). These genes were: in the male VMN, Cyp19a1 (PCB < DMSO, VIN < DMSO); in the POA, Grin2b (VIN > DMSO); and in the MeA, Kiss1 (VIN < DMSO).

Four significant sex differences were identified in control DMSO rats after correcting for multiple comparisons (Table 2); the entire dataset is in Supplemental Table 2 [47]. In the VMN, males had significantly higher expression of Cyp19a1 compared with females. In the POA, females had significantly higher Kiss1 expression, while males had higher Grin2b expression. In the MeA, males had significantly higher Kiss1 expression than females.

Table 2.

Sex differences in gene expression in vehicle (DMSO) rats

Region	Gene	Females		Males		Directionality
		Mean	SEM	Mean	SEM
VMN	Cyp19a1	2.14	0.79	23.68	14.16	Male > Female
POA	Kiss1	2.14	0.87	0.42	0.22	Female > Male
	Grin2b	1.04	0.06	2.19	1.19	Male > Female
MeA	Kiss1	1.06	0.32	6.08	1.57	Male > Female

Sex differences in gene expression in DMSO rats. Abbreviations: DMSO, dimethylsulfoxide; MeA, medial amygdala; PCB, polychlorinated biphenyl; POA, preoptic area; VIN, vinclozolin; VMN, ventromedial nucleus.

To identify potential relationships among correlated genes, we conducted a weighted gene co-expression network analysis (WGCNA), which allowed us to identify modules of correlated genes and calculated eigengene values [50]. This analysis showed that co-expression modules spanned functional gene groups and varied across brain regions and sexes (Supplemental Figure 3 [47]). However, treatments did not differ in eigengene expression for either sex or brain region after adjusting for multiple hypothesis testing (Supplemental Figure 4 [47]).

Embryonic exposure to EDCs caused dis-integration and reconstitution across levels of organization

The individual variation in behavior and gene expression allowed us to examine co-variance patterns between phenotypic measures (behavior, body weight, hormones) and gene expression, in order to identify any systems-level effects of EDC treatment and ask if behavioral, physiological, and neuromolecular correlations are maintained across treatments. Correlation strengths and phenotypic clustering are illustrated as heatmaps. An example of such a heatmap is shown Fig. 5 (in this case for gene expression in the VMN of females). Heatmaps for the POA and MeA in females, and for all 3 regions in males, are shown in Supplemental Figure 5 [47]. To interpret these heatmaps, the ordering of phenotypic traits along the axes was determined by hierarchical clustering based on correlation relationships among traits in the female DMSO group as shown in Fig. 5A. This ordering of traits in 5A was preserved to illustrate the extent of “dis-integration,” as shown in the heatmaps in panels 5B (PCB) and 5D (VIN), in which the traits in the DMSO control females have lost those correlations in EDC females. To generate “reconstitution” heatmaps shown in panels 5C (PCB) and 5D (VIN), hierarchical clustering was done independently for each treatment, resulting in a novel pattern of clusters that is in a different order from the female DMSO control group. These same analyses were applied to the 3 regions and both sexes, with the remaining heatmaps provided in the Supplement [47]. In all cases, both dis-integration and reconstitution were evident.

Figure 5.

Representative correlation heatmap from the VMN of females, illustrating the dis-integration and reconsolidation of behavioral, hormonal, and neuromolecular phenotypes caused by EDC treatment. Gene expression, behavioral, and physiological measures were clustered using Spearman’s rank correlations for the DMSO control samples (A). The same ordering of phenotypes in the control group was applied to both the PCB (B) and VIN (D) treatment groups to illustrate dis-integration associated with EDC treatment. To illustrate reorganization of phenotypes, behavioral, hormonal, and neuromolecular phenotypes were clustered for using Spearman’s rank correlations for PCB (C) and VIN (E) females, resulting in a novel ordering of phenotypic traits. Correlation strength is indicated by intensity of color with yellow indicating positive correlations and indigo indicating negative correlations. Labels for the reconstitution heat mapes (C and E) are in larger font in Supplemental Figure 5 [47]. Correlations within the remaining brain regions and in both sexes are provided in Supplemental Figure 5 [47]. Abbreviation: VMN, ventromedial nucleus.

For each set of heatmaps, we observed that correlations between gene expression levels are stronger (females: median r = 0.3–0.43; males: median r = 0.28–0.43) than between genes and behavior (females: median r = 0.2–0.28; males: median r = 0.2–0.32); this result may be partially influenced by the larger number of variables in the gene expression (44 traits) than the behavioral/hormonal (15 traits) analysis. Nevertheless, for both sexes and all brain regions the co-variance structure was strongly integrated in control animals (in terms of number and size of robust clusters), whereas with EDC treatment these patterns appeared to “dis-integrate” and/or reorganize. While the present study is underpowered to perform a quantitative analysis of this “dis-integration hypothesis,” the pattern of correlation loss and reorganization in both sexes and in all three brain regions reveals striking qualitative differences between the DMSO and EDC-treated rats.

To further illustrate this dis-integration effect we conducted correlation analyses between individual traits, with an example shown here for total distance traveled during the mate preference test and VMN Grin2b expression in females (Fig. 6). A robust correlation between these factors was found for the DMSO group (F(1,7)  = 19.60; P < 0.01; R2  = 0.74; Fig. 6A) but not the PCB (R2  = 0.005; Fig. 6B) or VIN (R2  = 0.0008; Fig. 6C) females. Conversely, we found strong evidence of reconstitution in males, where the direction of correlations present in the DMSO group were reversed by EDC treatment. Figure 6D shows one such example: the negative correlation between hormone preference score and Hsd17b1 expression in the POA of DMSO males (F(1,8)  = 13.56; P < 0.01; R2  = 0.63). The direction of this correlation was reversed by PCB treatment (F(1,7)  = 17.68; P < 0.01; R2  = 0.72; Fig. 6E) and abolished by VIN treatment (R2  = 0.03; Fig. 6F). Furthermore, there was considerable individual variation in both gene expression and behavior, and the correlation heatmaps demonstrate that there are many stronger correlations between genes than between genes and behavior.

Figure 6.

Example of treatment-induced dis-integration and reconstitution of behavioral measures and gene expression. Correlation between total distance traveled during the mate preference test and Grin2b expression in the VMN of female DMSO rats (A). B: No correlation between measures after PCB exposure. C: No correlation between measures after VIN exposure. D: Negative correlation between hormone preference score and Hsd17b1 expression in the POA of male DMSO rats. E: Positive correlation between hormone preference score and Hsd17b1 expression in the POA of male PCB rats. F: No correlation between measures after VIN exposure. Abbreviations: DMSO, dimethylsulfoxide; PCB, polychlorinated biphenyl; POA, preoptic area; VIN, vinclozolin; VMN, ventromedial nucleus.

Discussion

The current study is an innovative analysis of the effects of prenatal EDCs on phenotypic relationships among neural gene expression and behavior. Our finding that there is wholesale dis-integration and reconstitution among related sociosexual behavioral measures and neuromolecular networks indicate that while individual traits may be modestly affected, the gestalt of those traits is fundamentally different in treated vs control rats. Furthermore, the results were both sex-, EDC- (VIN and PCB) and region- (POA, VMN, MeA) specific. Our work underscores the importance of evaluating multiple phenotypic traits across levels of organization in assessing the long-term outcomes of EDC exposures. What determines the pattern of dis-organization or reorganization is not known, but these responses may reflect compensatory allostatic mechanisms in the face of external perturbations.

Prenatal EDCs changed mate preference in a treatment- and sex-specific manner

Due to their high investment in reproduction, female mammals seek out mates with the best likelihood of producing fertile offspring. For example, males with the typical adult range of concentrations of testosterone [32], and odors from such males, are preferred by sexually active female rats over low- or no-testosterone male counterparts [51, 52]. In the current study, sociosexual preference behaviors were impaired by PCBs in both sexes, and by VIN in males. This difference between the EDCs is interesting and may reflect the different modes of action by which the EDCs act [2–5, 7, 8].

More specifically, prenatal exposure to PCBs abolish the females’ preference for a stimulus male with testosterone replacement over a male without testosterone replacement, thereby replicating results from another study [16]. This outcome could translate into compromised reproductive success if a female were unable to discriminate between optimal and suboptimal males in more naturalistic conditions. Interestingly, VIN treatment had no effect in the females.

While males tend to be less choosy about mates, the process of mating still involves mutual interactions and coordination between both members of the dyad. Females exhibit proceptive behaviors to solicit the sexual attention of males, and males are also able to discriminate the odor of urine from receptive females [30, 53–55]. In the current study, and unlike in females, exposure of experimental males to both classes of EDCs (PCB or VIN) abolished the preference for the hormone-primed female over the female without hormone-replacement. The male rat brain develops under the influence of relatively high concentrations of both androgens and estrogens [2-5], perhaps conferring greater sensitivity to disruption of these pathways by both VIN and PCBs, respectively.

Previous work on prenatal A1221 exposure shows disrupted sex behavior in female rats, and decreased sexual motivation in male rats [11, 12]. Exposure to other PCBs also cause reduced sexual motivation and receptivity in females and altered sexual behavior in males [10, 13, 14, 56, 57]. In male rabbits [58, 59] exposure to VIN during the prenatal period and during postnatal life (E14 to adulthood) results in a lack of sexual motivation and deficits in sexual performance (reduced erections and ejaculations). Our current finding adds to this literature on sex-specific effects of EDCs on sociosexual behaviors. As a whole, the perturbations by EDCs of conspecific interactions have implications for social preference and sexual selection [60].

Prenatal EDCs altered the neuromolecular phenotype in the hypothalamus and amygdala

Prenatal EDCs affected the expression of a small number of genes in the VMN, POA, and MeA. It is notable that the identified genes, kisspeptin (Kiss1), nonapeptides (Avp, Oxt), steroidogenic enzymes (Cyp19a1, Hsd17b1), and the glutamatergic NMDA receptor subunit 2b (Grin2b) have previously been shown to be disrupted by EDCs. These genes also play roles in social behaviors and are hormone-sensitive. For example, the hypothalamic kisspeptin system is highly responsive to the estrogenic milieu [61, 62], making it an obvious target for estrogenic EDCs. Consistent with the current results, other studies have shown that PCBs and other EDCs affect kisspeptin protein and gene expression [63, 64]. Our finding that VIN affects Kiss1 in males is also consistent with this neuropeptide’s regulation by androgens [65].

Oxytocin gene expression was decreased in the MeA of PCB females; this was the same group of rats that showed deficits in the mate preference test. The MeA has been characterized as an important target of oxytocin in mate and odor preference behaviors [66]; although it has relatively sparse oxytocin fibers, there is evidence for a role of oxytocin expression in the amygdala in social behavior [67]. Oxytocin knockout mice have deficits in social recognition associated with reduced activity in the MeA and its projection targets [68]. By contrast, in the current study Oxt is increased by both PCBs and VIN in the female VMN and vasopressin by PCBs in the female VMN. Other labs have reported effects of EDCs on the nonapeptides vasopressin, oxytocin, and their receptors in several brain regions [69-71], implicating these as targets for perinatal endocrine disruption.

Both PCB and VIN males had lower Cyp19a1 (aromatase) expression in the VMN compared to DMSO control males. This region exhibits some of the highest levels of aromatase in the brains of rats along with the POA and the bed nucleus of the stria terminalis (BNST) [72]. In the hypothalamus, aromatase expression and activity is sexually dimorphic, with males having denser expression and higher activity [72]. Prenatal exposure to a similar PCB, Aroclor 1254, reduces aromatase activity in the hypothalamus of neonatal male rats [73]. Prenatal exposure to another EDC, the phthalate di(2-ethylhexyl)phthalate (DEHP), also reduces Cyp19a1 expression in the hypothalamus of neonatal rats [73].The N-methyl-D-aspartate (NMDA) glutamate receptor subunit 2b (Grin2b) is expressed widely throughout the rat hypothalamus [26], and its presence and abundance affects functional properties of NMDA receptors. In the POA, PCBs resulted in the overxpression of Grin2b in both sexes; VIN also increased Grin2b in the female POA. Hypothalamic Grin2b expression is sensitive to circulating estradiol levels and naturally decreases during reproductive senescence [74]. The activation of GnRH neurons in the POA by glutamate is necessary for reproductive function, and administration of a specific antagonist of the NMDAR2b subunit alter GnRH and downstream luteinizing hormone (LH) release in rats [75]. Limited work also suggests that EDCs may change Grin2b expression [63, 76]. Our finding provides further support that glutamatergic neurotransmission may be altered by prenatal EDC exposure.

EDC treatment can dis-integrate and/or reconstitute the relationships of behavior and gene expression phenotypes

The concept of an “essential phenotype” of an individual was introduced in the context of transgenerational effects of EDCs as a functional readout of the combination of underlying traits (genetic, hormonal, physiological, etc.) acquired from multiple measures [77]. This idea is equally applicable to the current study on prenatally exposed rats. Here, it is not each component part of the phenotype that defines an individual, but, rather, it is the unique way that each part interacts with one another across the developmental trajectory. Our evaluation of the essential phenotype of exposed individuals in adulthood revealed that each treatment had specific effects on females as well as males. This was illustrated as the disorganization of traits normally correlated with one another in control animals; and the emergence of novel relationships in the EDC-exposed animals, ie, the emergence of a novel essential phenotype in each sex, and specific to each EDC.

Previous work has demonstrated the emergence of unexpected relationships caused by EDCs interacting with other factors within an animal’s life cycle or across generations. In a transgenerational model of EDCs, Crews et al reported a history of EDC exposure in a transgenerational model, together with stress in the F3 descendants’ own adolescence, elicited a new phenotype referred to as the “synchronicity” of the 2 insults [48, 78]. Similarly, Bell et al using a two-hit model showed that prenatal and postnatal (adolescent) PCB exposure led to outcomes on behavior, physiology, and gene expression that could not be predicted by either hit alone [15, 79]. Our dis-integration/reconstitution analysis in the current study adds to this body of work by providing a new way of ascertaining complex relationships caused by EDCs.