Multi- and Transgenerational Outcomes of an Exposure to a Mixture of Endocrine-Disrupting Chemicals (EDCs) on Puberty and Maternal Behavior in the Female Rat

Background: The effects of endocrine-disrupting chemicals (EDCs) on fertility and reproductive development represent a rising concern in modern societies. Although the neuroendocrine control of sexual maturation is a major target of EDCs, little is known about the potential role of the hypothalamus in puberty and ovulation disruption transmitted across generations. Objectives: We hypothesized that developmental exposure to an environmentally relevant dose of EDC mixture could induce multi- and/or transgenerational alterations of sexual maturation and maternal care in female rats through epigenetic reprograming of the hypothalamus. We investigated the transmission of a disrupted reproductive phenotype via the maternal germline or via nongenomic mechanisms involving maternal care. Methods: Adult female Wistar rats were exposed prior to and during gestation and until the end of lactation to a mixture of the following 13 EDCs: di-n-butyl phthalate (DnBP), di(2-ethylhexyl) phthalate (DEHP), bisphenol A (BPA), vinclozolin, prochloraz, procymidone, linuron, epoxynaxole, dichlorodiphenyldichloroethylene, octyl methoxynimmate, 4-methylbenzylidene camphor (4-MBC), butylparaben, and acetaminophen. Perinatally exposed offspring (F1) were mated with unexposed males to generate germ cell (F2) and transgenerationally exposed (F3 and F4) females. Sexual maturation, maternal behavior, and hypothalamic targets of exposure were studied across generations. Results: Germ cell (F2) and transgenerationally (F3) EDC-exposed females, but not F1, displayed delayed pubertal onset and altered folliculogenesis. We reported a transgenerational alteration of key hypothalamic genes controlling puberty and ovulation (Kiss1, Esr1, and Oxt), and we identified the hypothalamic polycomb group of epigenetic repressors as actors of this mechanism. Furthermore, we found a multigenerational reduction of maternal behavior (F1–F3) induced by a loss in hypothalamic dopaminergic signaling. Using a cross-fostering paradigm, we identified that the reduction in maternal phenotype was normalized in EDC-exposed pups raised by unexposed dams, but no reversal of the pubertal phenotype was achieved. Discussion: Rats developmentally exposed to an EDC mixture exhibited multi- and transgenerational disruption of sexual maturation and maternal care via hypothalamic epigenetic reprogramming. These results raise concerns about the impact of EDC mixtures on future generations. https://doi.org/10.1289/EHP8795


Introduction
Endocrine-disrupting chemicals (EDCs) affect populations as much as individuals, given their environmental ubiquity (Gore et al. 2015). Close to 800 compounds are known or suspected to have the ability to disrupt normal endocrine function (Bergman et al. 2013) and induce a higher risk of diseases such as reproductive failure, cancer, obesity, metabolic syndrome, and neurodevelopmental disorders (Demeneix and Slama 2019). Small-scale studies have shown that follicular fluid samples of women visiting fertility clinics contained a wide array of chemicals with endocrine-disrupting activity, such as dichlorodiphenyltrichloroethane (DDT), phthalates, bisphenol A (BPA), and perfluorinated compounds (Petro et al. 2012(Petro et al. , 2014, indicating direct exposure of maturing oocytes and their surrounding steroid-producing cells. In current human pregnancy conditions, it is estimated that a woman and her fetus are exposed to a low-dose mixture of at least 100 EDCs (Demeneix and Slama 2019). In rodents, several EDCs in combination have been shown to produce effects at doses that individually are not associated with any observable response Christiansen et al. 2009;Hass et al. 2017;Howdeshell et al. 2008;Johansson et al. 2016). Expert panels have recently urged for further studies using environmentally relevant doses (EFSA 2019; Gore et al. 2015), which are crucial from a regulatory point of view because current risk assessment is solely based on the effects of individual chemicals.
Normal puberty is essential for achievement of reproductive capacity. Several epidemiological (Aksglaede et al. 2009;Parent et al. 2015) and animal (Franssen et al. 2016;López-Rodríguez et al. 2019) studies have suggested that pubertal timing is affected by EDCs and could be an early marker of fertility issues detected later in life. The rapid secular trend in age at onset of puberty reported around the world is likely due to the effect of exogenous environmental factors, including EDCs (Biro et al. 2010;Parent et al. 2015Parent et al. , 2003. Although the causal link between EDC exposure and pubertal timing in humans remains difficult to establish, several studies (Table 1) have shown a disruption of the neuroendocrine control of puberty in rodents, sheep, and nonhuman primates after developmental exposure to single (Collet et al. 2010;Franssen et al. 2016;Gao et al. 2018;Kurian et al. 2015;Naulé et al. 2014;Ruiz-Pino et al. 2019;Zalko et al. 2016) and mixtures of EDCs Hass et al. 2010;Johansson et al. 2016) at environmentally relevant doses as summarized in Tables 1 and 2 for the compounds used in our study. In mammals, sexual maturation is driven by increased pulsatile secretion of hypothalamic gonadotropin-releasing hormone (GnRH) released into the portal vasculature, ultimately stimulating pulsatile release of luteinizing hormone (LH) from the Table 1. Literature summary showing effects of the 13 compounds used in the EDC mixture on reproductive outcomes and the hypothalamic control of reproduction. pituitary into the peripheral circulation and inducing ovarian steroidogenesis and ovulation (Ojeda and Skinner 2006). During the prepubertal period, the secretory activity of GnRH neurons is silent, under predominant trans-synaptic inhibitory control provided by GABAergic, Opiatergic, and RFamide inputs (Ojeda et al. 2003). At puberty this inhibition is lifted, and a concomitant increase in excitatory inputs to the GnRH network is provided by glutamatergic and kisspeptin neurons (Clarkson and Herbison 2006;Terasawa et al. 2018). Recent evidence suggests that this trans-synaptic activation is controlled by an epigenetic switch. This epigenetic switch coordinates the transcriptional activity of arcuate nucleus (ARC) kisspeptin neurons implicated in stimulating GnRH release (Lomniczi et al. 2013;Toro et al. 2018). Before puberty, the transcriptional activity of Kiss1, the gene encoding kisspeptins, was found to be downregulated by members of the polycomb group (PcG) of transcriptional repressors, by catalyzing the trimethylation of histone 3 at lysine 27 (H3K27me3), a histone mark associated with gene silencing (Lomniczi et al. 2013). As puberty approaches, the PcG is evicted from the Kiss1 promoter, and the Trithorax Group (TrxG) of epigenetic activators is recruited, resulting in increased histone methylation at lysine 4 (H3K4me3) and acetylation at lysine 27 (H3K27ac) to the Kiss1 promoter/enhancer regions, respectively, leading to an increase in Kiss1 mRNA transcription . Several studies have identified Kiss1 expression as a target of EDCs (Table  1, see Lopez-Rodriguez et al. 2021 for review) in the anteroventral periventricular nucleus (AVPV) (Losa et al. 2011;Naulé et al. 2014;Ruiz-Pino et al. 2019;Wang et al. 2014) and the arcuate nucleus (ARC) (Cao et al. 2012;Hu et al. 2013;Losa et al. 2011;Ruiz-Pino et al. 2019) cell populations, potentially resulting in abnormal puberty and ovulation (López- Rodríguez et al. 2019;Mahalingam et al. 2017;Patel et al. 2017;Santamaría et al. 2016), but no study so far has identified the pathways conveying epigenetic information from environmental disruptors to the hypothalamic cells controlling pubertal onset and ovulation. Although the direct long-lasting consequences of developmental exposure to EDCs have been extensively addressed (Demeneix and Slama 2019; Gore et al. 2015), the consequences across generations is only starting to be studied (Brehm and Flaws 2019;Lopez-Rodriguez et al. 2021). Exposure through maternal lineage (F0), leads to direct exposure of F1 pups in utero and/or through lactation, whereas F2 pups are directly exposed as germ cell. This phenomenon is known as mutigenerational inheritance. Transgenerational inheritance requires the study at least of the F3 generation when indirect action of the exposure is largely explained by alterations in the epigenetic reprogramming of the germline. Transgenerational consequences of EDC exposures through maternal lineage (Table  2) have been observed on pubertal timing (Rattan et al. 2018b;Zhou et al. 2017;Ziv-Gal et al. 2015), female fertility, and ovarian follicle development (Berger et al. 2016;Brehm et al. 2018;Drobna et al. 2018;Nilsson et al. 2012;Pocar et al. 2017) in rodents. These transgenerational reproductive abnormalities caused by EDCs may be mediated by alterations of the neuroendocrine control of sexual maturation and reproduction (Goldsby et al. 2017;Wolstenholme et al. 2013).
Besides transgenerational inheritance, phenotypic traits can be propagated through an experience-based nongenomic mechanism. This mechanism known as multigenerational inheritance can be transmitted by maternal care and exclusively involves somatic tissue alterations (Champagne and Curley 2016). Natural variations in rat maternal behavior, such as licking and grooming, were shown to shape the offspring hypothalamus by modifying the expression of Esr1 and Nr3c1 (Francis et al. 1999;Peña et al. 2013). This behavioral pattern of expression persisted into adulthood and was transmitted through generations (Francis et al. 1999). In rats, licking and grooming behaviors were found to be vulnerable to BPA exposure (Seta et al. 2005). In addition, natural variations in licking and grooming were found to be able to alter pubertal timing in the offspring (Cameron et al. 2008a). Whether variations in pubertal timing caused by changes in licking and grooming are produced via reprogramming of the hypothalamus is still unknown.
Here we hypothesize that developmental exposure to an EDC mixture through the maternal lineage induces transgenerational alterations of female pubertal timing and ovulation and a multigenerational reduction in maternal behavior. In the current study,  The studies used some of the EDCs found in our mixture.
Environmental Health Perspectives 087003-3 129(8) August 2021 female rats (F0) were exposed to a mixture of 13 EDCs at relevant exposure concentrations prior to and during gestation and until the end of lactation. F1 females (exposed in utero and through lactation) and subsequent generations were bred with unexposed males to generate the F2 (exposed through the germ cells), F3, and F4 generations (exposed transgenerationally). Sexual maturation and maternal behavior were evaluated across generations. We aimed to compare transcriptional and chromatin profiles after direct (F1) and transgenerational (F3) exposure. Furthermore, by using a cross-fostering paradigm, we were able to distinguish between transgenerational consequences of the EDC mixture transmitted via the germline and multigenerational consequences due to variations in maternal care.

Animals
Adult female Wistar rats purchased from the animal facility of the University of Liège, Belgium, were housed in standardized conditions (12 h inverted dark:light phase, 22.8°C, and food and water ad libitum). All animals were raised in EDC-free cages (polypropylene cages, Ref. 1291H006; Tecnilab) and fed EDC-and phytoestrogen-free chow (V135 R/Z low phytoestrogen pellets; SSNIFF Diet). Water was supplied in glass bottles. All experiments were carried out with the approval of the Belgian Ministry of Agriculture and the Ethics Committee at the University of Liège.

Experimental Design
Early adult, postnatal day 40 (P40) male and female Wistar rats were purchased from the animal facility of the University of Liège. The animals were never inbred. After habituation to the animal care facility, adult (P70) F0 females were mated for 2 wk with an adult (P70) unexposed male to generate the F1 generation. F0 females were exposed to an EDC mixture or corn oil (vehicle) starting 2 wk before mating and until the end of lactation to a mixture of 13 EDCs ( Figure 1). Depending on the timing at pregnancy during mating, adult females were exposed for a total of 57-65 d prior and during gestation and until offspring weaning. Daily exposure was done by injecting 50 lL of EDC mixture or corn oil in a wafer (Delacre N.V. S.A.) and providing one wafer to each female. Before and during mating, an EDC and one control (CTL) female were housed together. Females were monitored daily for body weight changes due to pregnancy. After confirmation of gestation, females were housed individually until weaning. When necessary, females were separated to provide the wafer containing the EDC mixture, and cages were systematically verified after 10 min to check for complete ingestion. Thereafter, females were moved back to their original cages. Females were randomly assigned to treatment. At P70, in utero EDC mixture-exposed and control females from the F1 generation were mated with an unexposed male ordered from the animal facility of the University of Liège to generate the F2 generation (exposed through the germline). Similar procedures were used to generate the F3 and F4 generations, transgenerationally exposed through the maternal lineage. Females used to generate subsequent generations were randomly selected within each group. For every generation, pups were homogenized for sex ratio (1:1) and litter size (average of 12 pups per litter). From weaning at P21 to P70, control females were paired with an EDC female and housed together. The pairing procedure was followed for all females from the F1 to the F4 generations. All behavioral and reproductive data were obtained from two independent nonrelated cohorts of animals. Samples used for all experiments originated from an original pool of 40 F0 females per group/per cohort.

Cross-Fostering
To distinguish between germ-cell vs. experience-based phenotype transmission, a cross-fostering paradigm was used. To avoid an effect of cross-fostering, a maximum of two pups per litter were cross-fostered with pups from F1 dams as follows: Control F1 dams raised either EDC-exposed pups (CE) or control pups (CC) from another dam. EDC-exposed F1 raised either EDC-exposed pups (EE) or control pups (EC) from another dam. Cross-fostering was carried out within the first 24 h after delivery and a 1-mm tail incision was done to distinguish cross-fostered pups. Samples used from all cross-fostering experiments originated from 17 EDCexposed and 16 control dams from the F1 generation.

Maternal Behavior
To assess EDC mixture effects on maternal behavior throughout generations, a set of maternal behaviors were quantified in EDC mixture and control lactating females in F0 to F3 generation. Additionally, we quantified maternal behavior in F1 dams (EDC mixture-exposed or control dam with either control or EDC mixture-exposed pups from another dam) with cross-fostered pups and in F2 cross-fostered lactating females (EDC mixture exposed or nonexposed females raised by an EDC mixture-exposed or control dam). Maternal behavior was recorded with an infrared camera (Bell and Howell DNV16HDZ-BK) for 1 h during the dark phase from lactational day 2 to 8 directly in the home cage under undisturbed conditions. Randomization was done to avoid recording every female at the same hour of the day. A set of in-nest behaviors (licking/grooming, arched-back/blanket/passive nursing, and nest building) and off-nest behaviors (eating/drinking, grooming, active or resting alone) were quantified as reviewed in Lonstein et al. (2015) by an experimenter blinded to the condition. Briefly, each video was visualized by an experimenter blinded to the condition and at the same time that the duration of each behavior was manually counted using a data sheet. Licking and grooming was determined when the female engaged in licking the anogenital regions of the pup. Nursing behavior was defined by pups attached to the dam's nipple. Arched-back nursing or blanket-posture nursing (lying down on the side) were quantified separately. Passive nursing was defined as pups attached to the dams when the dams were engaged in activities outside the nest (being active, eating/ drinking). Off-nest behaviors were quantified when the female was eating/drinking, grooming herself, or resting alone outside the nest. Total duration of every behavior was quantified within each 1-h period recording. Figure 1. Multi-and transgenerational design of an EDC mixture exposure on maternal behavior and sexual maturation via the female lineage. Adult female rats (F0) were exposed to a mixture of 13 EDCs at relevant exposure concentrations starting 2 wk before mating and until the end of lactation. Mating was carried out for period of 2 wk. Depending on the day of pregnancy during mating, adult females were exposed for a total of 49-56 d, which spanned the period before, during, and after gestation. The four subsequent generations were evaluated for sexual maturation (VO, GnRH interval interpulse, estrous cycle, and folliculogenesis) and maternal behavior (from P2 to P8). Massive parallel RNA sequencing was carried out using hypothalamic (MBH-PoA) explants from the F1 and F3 generation to decipher direct (F1) vs. transgenerational (F3) target genes of the EDC mixture exposure, followed by qPCR validation at three time points, P6, 21, and 60. Target genes were studied for histone posttranslational modifications and DNA methylation using chromatin immunoprecipitation (ChIP) and bisulfite sequencing (BS), respectively. Note: EC, estrus cycle; EDC, endocrine-disrupting chemical; VO, vaginal opening.

Pubertal Onset and Estrous Cyclicity
To assess the effect of an EDC mixture exposure on sexual maturation throughout generation, females were followed for vaginal opening (from the F1 to the F4 generation) and estrus cyclicity (from the F1 to the F3 generation) as described previously (Franssen et al. 2014;Toro et al. 2018). Briefly, from P25, females were daily inspected for vaginal opening by two experimenters. From the P40 and until P70, estrous cycle was evaluated using vaginal smears taken every morning. Regular cycles were defined as a sequence of diestrus 1, diestrus 2, proestrus, and estrus in 4 consecutive days (Goldman et al. 2007). The percentage of females having a regular cycle and the time spent in every stage of the cycle were calculated from 6-10 weeks of age.

Ovarian Histology and Uterus Weight
To assess the effect of an EDC mixture exposure on folliculogenesis throughout generations, ovaries from F1 to F3 females were removed at P70 for histological quantification of follicle development. Females were anesthetized with isoflurane (IsoFlo ® ; Zoetis) and culled by quick decapitation. After removal, ovaries were weighted together with the uterus and fixed overnight in 4% paraformaldehyde, dehydrated in 70% EtOH, and paraffin-embedded. Histological analysis was done in 8-lm coronal sections (microtome RM2245; Leica), after deparaffinization and staining with hematoxylin and eosin. For quantification, every other sections throughout the whole ovary were digitalized using an automated digital microscopy system DotSlide (Olympus BX51TF). DotSlide images taken at a magnification of 10 ×, which were in a proprietary format, were converted into a standard TIFF format, and 3-fold decimated, which made them easier to handle. Thereafter, follicles at every stage of folliculogenesis (primordial, primary, secondary, antral, and atretic), cystic follicles and corpora lutea were manually quantified by an experimenter blinded to treatment with Aperio ImageScope software (version 12.3.2.8013; SCR_014311; Leica Biosystems). Total ovarian volume was automatically calculated using an original program developed using the image analysis toolbox of the MatLab (SCR_001622; The Mathworks, Inc.) software. The algorithm developed to measure the ovarian volume used the traditional tools of signal treatment (Gonzales and Woods 2008) and mathematical morphology (Soille 1999). Briefly, images were automatically segmented (Otsu 1979), followed by visual inspection and manual correction if necessary. Next, morphological operations were carried out to fill gaps and to eliminate small noise artifacts (close, erosion, and reconstruction functions). Finally, volumetric measurements were done using binary 3D image reconstruction. The follicles were classified according to wellestablished criteria (Hirshfield and Rees Midgley 1978;Peters 1969). Briefly, primordial follicles were defined as an oocyte surrounded by a partial or complete layer of squamous granulosa cells. Primary follicles displayed a single layer of granulosa cell with cuboidal shape. Secondary follicles corresponded to follicles displaying several layers of cuboidal granulosa cells surrounded by theca cells without follicular fluid (antrum). Antral follicles were identified by their size and the presence an antrum. Atretic follicles corresponded to a shrunken oocyte surrounded by a granulosa cell layer. Cystic follicles were counted when follicles did not display an oocyte and contained a large antrum with few or no granulosa cells. Finally, corpora lutea were identified by the presence of luteal cells. Double counting of late-stage follicles was avoided by digitally marking each follicle throughout the consecutive images. Each follicle was counted once whenever the oocyte was present. As analysis were done in every other section, we applied a 2-fold correction factor for quantification of early stage follicles (primordial and primary follicles) to compensate for sections not analyzed. Measurements were expressed as number of follicles or corpora lutea per volume (cubic millimeters).

Hypothalamic Explants Incubation and GnRH Assay
To assess the effect of an EDC mixture on juvenile GnRH frequency after direct in utero and lactational (F1 generation) and transgenerational (F3 generation) exposure, GnRH interpulse interval was measured using a hypothalamic explants incubation system followed by a GnRH assay from prepubertal females (P20) of the F1 and F3 generation, as described previously (Bourguignon and Franchimont 1984;Matagne et al. 2004). Briefly, females were anesthetized with isoflurane (IsoFlo ® ; Zoetis) and culled by quick decapitation. Thereafter, brain was dissected by performing two sagittal incisions along the lateral hypothalamic sulci and two transversal incisions of 2 mm ahead from the anterior boundaries of the optic chiasm and along the caudal margin of the mammillary bodies. Once dissected, hypothalamic (MBH-PoA) explants were transferred into an individual chamber, in a static incubator, submerged in MEM. The ex vivo explant incubation chamber contained a water-saturated atmosphere of 95% O 2 , 5% CO 2 at 37.5°C. Incubation medium was then collected and renewed every 7.5 min for a period of 4 h. All hypothalamic explants were incubated at the same time during the same day for each generation. The GnRH released into the incubation medium was measured in duplicate using radioimmunoassay as described previously (Bourguignon et al. 1989a). Briefly, samples were preincubated with the highly specific CR11-B81 (AB_2687904) rabbit anti-GnRH antiserum (initial dilution 1:20,000) (provided by Dr. V.D. Ramirez, Urbana, Illinois) during 24 h at 4°C. GnRH labeled with 125I (30,000 CPM) and rabbit serum (dilution 1:100) were added for 24 h at 4°C. Finally, precipitation was induced by a solution of sheep anti-rabbit antiserum (dilution 1:200; CER Groupe), polyethylene glycol (60 g=L), tween, and cellulose. Radioactivity was counted on a gamma counter (Wallac CliniGamma). The intra-assay and interassay coefficients of variation were 7% and 10%. Undetectable values were assigned the limit of detection (5 pg=7:5-min fraction).

DNA and RNA Extraction, Reverse Transcription, and Real-Time Polymerase Chain Reaction (RT-PCR)
Expression of genes involved in the hypothalamic control of puberty, reproduction, and maternal behavior were studied by quantitative real-time polymerase chain reaction (qPCR) analysis using the half of MBH-PoA explants from females exposed and nonexposed to an EDC mixture from the F1 and F3 generations at different time points (P6, P21, and P60). Additionally, MBH-PoA explants were dissected from cross-fostered (EE, CE, EC, and CC) females at P21. All MBH-PoA explants were dissected at the same time of the day and alternating between control and EDC-exposed females. Adult P60 MBH-PoA explants were dissected on diestrus. Females were anesthetized with isoflurane (IsoFlo ® ; Zoetis) and culled by quick decapitation. Thereafter, the MBH-PoA was rapidly dissected as described in the previous section. Additionally, MBH-PoA was divided in two by sectioning along the interhemispheric fissure. Total RNA and DNA were extracted from the half MBH-PoA tissue using All Prep DNA/RNA Mini kit (Qiagen) following the manufacturer's instructions. Five hundred nanograms of bulk RNA for each sample were reverse transcribed using the Transcriptor first strand cDNA synthesis kit (Roche). For real-time quantitative PCR reactions, the cDNA of our samples were diluted 10-fold, and 4 lL were added to a mix of 5 lL FastStart Universal SYBR ® Green Master (Roche), 0:4 lL of nuclease-free water and 0:3 lL of forward and reverse primer (see primer sequences in Table  S1). Primers were synthesized by Integrated DNA Technologies, Inc. qPCR was performed using a Quant Studio 12K Flex (Applied Biosystems) under the following conditions: initialization: 2 min at at 50 C, 10 min at at 95 C, followed by 40 cycles of denaturation: 15 sec at 95 C and elongation and data collection: 60 s at 60 C. Cycle threshold (Ct) values were obtained from each individual amplification curve, and the average Ct was calculated for each target gene in each sample. Quantification of relative gene expression was performed using the DDCt method implemented with the Pfaffl equation, which takes into account reaction efficiency depending on primers (Pfaffl 2001) All assays had efficiencies between 1.9 and 2.1. b-actin was used as housekeeping gene.

Bisulfite Sequencing
Genomic bulk DNA was bisulfite converted (BC), using the EZ DNA Methylation-Gold kit (Zymo Research) according to the manufacturer's instructions. The BC DNA was used as input material for PCR amplification followed by library preparation and deep sequencing. Primers were designed to amplify a 302-base pair (bp) region of the rat Th promoter, including exon 1 (−139 to +164 bp from the Th TSS) (Table S1). Amplification was carried out on a C1000 Thermal Cycler (Bio-Rad) with 20 ng of BC DNA per reaction. The amplification conditions were: 40 cycles of 94°C for 30 sec, 55°C for 30 sec, and 72°C for 1 min; sequencing libraries were prepared using the NETFLEX ® DNA Sequencing Kit (BIOO Scientific) according to manufacturer's instructions. The libraries were evaluated using an Agilent 2100 Bioanalyzer (Agilent Technologies) and were normalized to 2 nM with 10 mM Tris-HCl. The libraries were then pooled and sequenced on a MiSeq (Illumina, Inc.) by the Molecular and Cellular Biology Core (Oregon National Primate Research Center) to generate 250 base, paired-end reads. The reads were trimmed using Trim Galore (http://www.Bioinformatics. Babraham.Ac.Uk/Projects/Trim_galore/) and aligned to the rat reference genome (Ensemble Rat Rnor_6) using Bismark (Krueger and Andrews 2011). Alignment data was converted to CpG methylation rate using the Bismark methylation extractor. Only reads with a >98% cytosine conversion and >98% alignment were used; no sample was excluded from the analysis. The CpG methylation rates were calculated as the ratio of methylated reads over the total number of reads. Methylation rates for CpGs with fewer than 10 reads were excluded from further analysis; samples had an average of 1,750 × read coverage.

Immunohistochemistry
To assess the expression of proteins related to maternal behavior, tyrosine hydroxylase (TH) immunoreactivity was measured in EDC mixture exposed and nonexposed F1 females at P21. TH-ir was measured in the substantia nigra (SN), the ventral tegmental area (VTA), and the median preoptic area (mPoA). Females were anesthetized with pentobarbital sodium (50 mg=kg, intraperitoneal) and sequentially perfused with PBS and 4% paraformaldehyde. Brains were subsequently removed and post-fixed in 4% paraformaldehyde at 4°C overnight. Thereafter, coronal sections (30 lm) were cut on a vibratome and used for immunohistofluorescence.
For immunohistochemistry, sections were incubated in phosphate-buffered saline (PBS) and blocked with 10% donkey serum for 1 h at room temperature and then were incubated overnight with the primary antibody anti-TH (1:1000, Th-mouse; 22941, ImmunoStar) for 24 h at room temperature or 4°C. Thereafter, a goat antimouse fluorophore-conjugated secondary antibody (1:400, AlexaFluor 488, ab150113) was incubated for 2 h at room temperature. A high-resolution Zeiss microscope LSM880 implemented with the fast Airyscan detector was used for visualization. Three slides per animal per region of interest (ventral midbrain sections containing the SN and VTA or mPoA sections) were used for quantification. An observer blind to condition outlined each region of interest for nuclei-specific analysis. Region area and number of immunoreactive cells were automatically quantified using Imaris 9.3. Statistics were performed on the total cell count per animal per region across slice sections.

Massively Parallel RNA Sequencing
Poly A mRNA was purified from 1 lg total bulk RNA using NETFLEX poly A beads (Perking-Elmer), followed by library preparation using the NETFLEX ® Rapid Directional RNA-seq Kit 2.0 (Perking-Elmer). In short, after fragmentation with divalent cations and heat, mRNA was used as template for reverse transcription using random hexamer primers. cDNAs were then blunted and 3 0 end A-tailed to facilitate adaptor ligation. Six base-pair Illumina adaptors were ligated, followed by 12 rounds of PCR amplification. Free dNTPs were removed using AMPure XP beads (BeckmanCoulter). Distribution of DNA sizes in the library was confirmed by Bioanalyzer analysis (Agilent). Library titer was determined by RT-PCR (Kapa Biosystems) on a Quant Studio 12K Flex Real Time PCR System (ThermoFisher). Four samples were sequenced per lane on a HiSeq 4000 (Illumina). Sequencing was done using a single-read 100-cycle protocol. The resulting base call files (.bcl) were converted to standard fastq formatted sequence files using Bcl2Fastq (Illumina). Sequencing quality was assessed using FastQC (Babraham Bioinformatics). The RNAseq procedure was carried out by the Genomics and Cell Characterization Core Facility at the University of Oregon. To determine the differential gene expression values we used the gene-level edgeR analysis package. We performed an initial trimming and adapter removal pass using Trimmomatic. After this, reads were aligned to the rn6 build of the rat genome with Bowtie2/TopHat2 and assigned to gene-level genomic features with the Rsubread featureCounts package (Liao et al. 2019) based on the Ensembl 83 annotation set. Differential expression between experimental groups was analyzed using the generalized linear modeling approaches implemented in edgeR [version 3.7; Robinson et al. (2010)]. Lists of differentially expressed genes/ transcripts were identified based on significance of pairwise comparison of experimental groups. Gene ontology and enrichment analysis were performed using the database for annotation, visualization, and integrated discovery (DAVID). RNAseq data sets are accessible through Gene Expression Omnibus (GEO) Series accession number GSE168151 (https://www.ncbi.nlm.nih.gov/ geo/query/acc.cgi?acc=GSE168151).

Chromatin Immunoprecipitation (ChIP) Assay
To assess activatory and repressive histone modifications at specific gene promoters affected by EDC mixture exposure, we performed ChIP assay using extracted chromatin from the hypothalamus of prepubertal rats at P21. As described previously (Mueller et al. 2011), ChIP procedure was carried out by crosslinking tissue in PBS containing a protease inhibitor cocktail (PI, 1 mM phenylmethylsulfonylfluoride, 7 lg=mL aprotinin, 0:7 lg=mL pepstatin A, and 0:5 lg=mL leupeptin), a phosphatase inhibitor cocktail (PhI, 1 mM b-glycerophosphate, 1 mM sodium pyrophosphate, and 1 mM sodium fluoride), and an HDAC inhibitor (20 mM sodium butyrate) at 4°C, and 1% formaldehyde for 10 min at room temperature. After two washing steps in PBS, samples were lysed with 200 lL sodium dodecyl sulfate (SDS) buffer (0.5% SDS, 50 mM Tris-HCl, and 10 mM EDTA) containing protease, phosphatase, and HDAC inhibitors and sonicated for 45 s to yield chromatin fragments of ∼ 500 bp, using the microtip of a Fisher Scientific FB 705 sonicator. Size fragmentation was confirmed by agarose gel electrophoresis. The sonicated chromatin was clarified by centrifugation at 14,000 rpm. for 10 min at 4°C, brought up to 1 mL in Chip Dilution Buffer (16:7 mM TrisHCl, pH 8.1, 150 mM NaCl, 1:2 mM EDTA, 1.1% Triton X-100, and 0.01% SDS) containing the PI and PhI cocktails and the histone deacetylase (HDAC) inhibitor described above. The samples were then stored at −80 C for subsequent immunoprecipitation. For this step, chromatin was precleared with Protein A/G beads (Dynabeads™; Invitrogen) for 1 h at 4°C. Then 25-to 50-lL aliquots of chromatin were incubated with 2-5 lg of the antibodies described in Table S2. The complexes were incubated with 25 lL of protein A or G beads solution (Dynabeads™) at 4°C overnight with mild agitation. The next day the beads were washed first with 0:5 mL low-salt wash buffer (20 mM Tris-HCl, pH 8.1, 150 mM NaCl, 2 mM EDTA, 1% Triton X-100, and 0.1% SDS), followed by high-salt wash buffer (20 mM Tris-HCl, pH 8.1, 500 mM NaCl, 2 mM EDTA, 1% Triton X-100, and 0.1% SDS), LiCl buffer (10 mM Tris-HCl, pH 8.1, 250 M LiCl, 1% Nonidet P-40, 1% sodium deoxycholate, and 1 mM EDTA), and finally with TE buffer (10 mM Tris-HCl, pH 8.0, and 1 mM EDTA). Thereafter, the immunocomplexes were eluted with 100 lL of 0:1 M NaHCO3 and 1% SDS at 65°C for 45 min. To reverse the crosslinking reaction, we added 4 lL of 5 M NaCl and incubated the samples at 95°C for 30 min. We recovered the DNA using ChIP DNA Clean & Concentrator columns (Zymo Research) and stored the resulting material at −80 C before qPCR analysis. All the chemicals mentioned above were purchased from Sigma-Aldrich.
qPCR Detection of Chromatin Immunoprecipitated DNA Genomic regions of interest were amplified by qPCR. Primer sequences and accession numbers of the genes analyzed are shown in Table S1. PCR reactions were performed using 1 lL of each immunoprecipitate (IP) or input samples, primer mix (1 lM each primer), and SYBR ® Green Power Up Master Mix™ (Thermo Fisher) in a final volume of 10 lL. Input samples consisted of 10% of the chromatin volume used for immunoprecipitation. The thermocycling conditions used were as follows: 95°C for 5 min, followed by 40 cycles of 15 sec at 95°C and 60 sec at 60°C. Data are expressed as percentage of IP signal/input signal.

Statistics
All statistical analyzes were performed using Prism (version 7.0, GraphPad). Data was subjected to a normality, and an equal variance test and parametric test were used when conditions were accomplished. Parametric test used were one-way or two-way analysis of variance (ANOVA) followed by Student-Newman-Keuls or Sidak's for multiple comparisons, respectively; or Student's t-test to compare two groups. When comparing percentages, groups were subjected to an arcsine transformation before statistical analysis to convert the values from a binomial to a normal distribution. Data that did not reach normality were analyzed using a Mann-Whitney test. When making multiple comparisons, a was adjusted by using the Bonferroni correction. The investigator was group blinded in all physiological and molecular determination. Sample sizes, reported in figure captions, were defined based on previous studies and calculation of an adequate statistical power. Report of descriptive, statistical values and effect size are available in Tables S3-S6. Effect sizes were calculated according to Cohen's delta formula. The level of statistical significance was p < 0:05.

Transgenerational Epigenetic Programing after Developmental Exposure to the EDC Mixture
We hypothesized that the transgenerational delay in pubertal onset and GnRH secretion maturation observed in F3-EDC females could be due to transcriptional and epigenetic disruption of the hypothalamic networks controlling pubertal development.
We used massive parallel RNA sequencing to identify gene regulatory pathways altered by transgenerational EDC exposure in the MBH-PoA of F3-EDC females at P21, corresponding to the juvenile phase of prepubertal development. Gene ontology analysis showed that the main down-regulated biological processes were neurohypophyseal hormone activity, neuropeptide receptor activity, and dopamine signaling. Some of the down-regulated molecular functions affected included GABA, dopamine, and glutamate signaling, whereas xenobiotic processes were upregulated ( Figure S2). mRNA expression of hypothalamic genes critical for pubertal onset was quantified at P6, P21, and P60 in transgenerationally exposed females (F3) using qPCR. We selected genes to be validated based on a functional analysis of RNAseq data and existing literature regarding the transcriptional control of female puberty in rats. Hypothalamic expression of Kiss1, Esr1, and Oxt was significantly lower in F3-EDC females at P21 (p < :05) and P60 (p < :05) ( Figure 4A; Table S3). We also identified a set of genes involved in the control of energy balance and reproduction (Cart, Pomc, Grin2d, Grid2, and Avp) and stress responsiveness (Nr3c1 and Crh) that were altered at P21 in F3-EDC females in comparison with controls ( Figure 4A; Figure S3; Table S5).
To gain insight into the potential contribution of histone modifications to the transcriptional changes caused by the transgenerational EDC exposure, we used ChIP assay to quantitate the repressive (H3K27me3 and H3K9me3) and activating (H3K4me3, K3K9ac) histone modifications at the promoter region of target genes in the MBH-PoA of F3 control and EDC transgenerationally exposed females ( Figure 4B; Figure S4; Tables S3 and S5). Transcriptional down-regulation caused by transgenerational EDC mixture exposure was consistently associated with either higher expression in repressive histone marks or lower expression of activating histone modifications at the promoter region of genes critical for the onset of puberty. We observed higher expression of the repressive marks H3K27me3 at the promoter region of the Kiss1 gene (p < :01), whereas an abundance of activating marks H3K4me3 or K3K9ac was observed to be lower at Esr1 (p < :05), Kiss1 (p < :01), and Oxt (p < :05) ( Figure 4B; Figure S4; Table S3 and S5). No difference in histone landscape was detected at the Cart and Grin2d promoter, whereas lower expression of the activatory H3K4me3 was detected at the Pomc promoter (p < :05) ( Figure 4B; Figure S4; Table S5). Furthermore, we observed a higher expression of repressive H3K9me3 and H3K27me3 together with a higher expression in activatory H3K4me3 at the Nr3c1 promoter of F3 EDC-exposed females.
We further aimed at identifying hypothalamic targets potentially involved in the alterations of maternal care caused by the EDC mixture that may explain the transmission of the maternal phenotype across generations. We compared massive parallel RNA sequencing gene ontology results obtained from MBH-PoA explants at P21 after direct in utero and lactational (F1 generation) vs. transgenerational (F3 generation) EDC mixture exposure. Gene ontology analysis of RNAseq done in MBH-PoA after direct in utero and lactational (F1) EDC exposure showed down-regulated genes related to D2 dopamine receptor binding, glutamate binding and Figure 3. Estrous cycle, ovarian follicle development, and ovarian weight across generation (F1-F3 generation) of female rats exposed to EDC mixture or vehicle (F0 generation). (A-C, left) Percentage of females showing regular cycle and average time spent in each stage of the estrous cycle in F1 (n = 20=group), F2 (n = 15=group), or F3 (n = 14-15=group) animals. Samples originate from eight different litters per group per generation. Data were analyzed using a twoway ANOVA (A-C, middle). Average number of ovarian follicles, corporal lutea, and cysts per cubic millimeter in F1, F2, and F3 females at P60 (n = 9-10=group). Samples originate from eight different litters per group per generation. Data were analyzed using a two-way ANOVA (A-C, right). Normalized ovarian weight (ovarian weight per body weight in grams) measured at P60 in F1 (n = 14=group), F2 (n = 16=group), and F3 (n = 9-10=group) females. Samples originate from eight different litters per group per generation. Data was analyzed using Student's t-test. Bars represent mean ± SEM ( * p < 0:05, ** p < 0:01, *** p < 0:001 vs. CTL). Summary data are reported in Table S3. Note: ANOVA, analysis of variance; CTL, control; D, diestrus; E, estrus; EDC, endocrine-disrupting chemical; P, proestrus; Reg. cycle, regular cycle; SEM, standard error of the mean. neurotransmitter uptake ( Figure S9). After transgenerational (F3) exposure, down-regulated genes belonged to the following enriched categories: social behavior, maternal behavior, grooming behavior, and synaptic dopaminergic transmission ( Figure S2).
Common target genes identified in F1 and F3 females and related to dopaminergic signaling were validated using RT qPCR. Critical dopaminergic signaling genes Th, Dnm1 and Darpp32, involved in maternal motivation, showed significantly lower expression in the MBH-PoA of F1-EDC animals at P21 and/or P60 ( Figure 6A; Table S3). The dopaminergic receptor 1 (Drd1) was found to be up-regulated, a possible sign of a compensatory mechanism induced by lower synaptic dopamine levels ( Figure 6A). Genes associated with stress responsiveness (Nr3c1 and Crh) were not altered ( Figure S10A).
To determine whether hypothalamic transcriptional alterations in dopaminergic signaling observed in the F1 generation involved epigenetic mechanisms, we studied Th DNA methylation and histone modifications. Although we identified no difference in DNA methylation ( Figure 6B; Table S3), we observed higher levels of the repressive histone H3K27me3 modification in the MBH-PoA of F1-EDC females (p < :05) ( Figure 6C). No differences were observed in H3K9me3, H3K4me3, or H3K9ac status. Furthermore, no alterations in histone modification landscape were detected in the 5 0 regulatory regions of Darpp32 or Drd1 ( Figure S10B and Table S5).
To further characterize whether the alterations in dopaminergic signaling are related to the observed disruption of licking and grooming behavior, we quantified TH-immunoreactivity in the SN, VTA, and mPoA of F1 females at P21. A lower number of TH-immunopositive cells was observed in the mPoA of in utero EDC-exposed females (CTL: 3:98 ± 0:61, EDC: 2:53 ± 0:36, p < :001) ( Figure 6D-E; Table S3), confirming mRNA expression results. No differences were observed in the VTA or nucleus accumbens (NAc).
To determine whether transcriptional alterations in dopaminergic signaling observed in the F1 generation are transmitted across generations and through maternal care, we studied mRNA expression levels of dopaminergic genes in the F3 generation. Th and Drd1 mRNA expression were found to be significantly lower in the MBH-PoA of transgenerationally exposed F3-EDC females (p < :05) ( Figure 6F). Th down-regulation was associated with higher DNA methylation at 3 specific CpGs of the Th promoter (p < :05), as well as enhanced repressive H3K27me3 histone mark (p < :05) ( Figure 6G-H). No differences were found in the activating histone modifications H3K4me3 and H3K9ac.

Reproductive Function and Maternal Behavior after EDC Exposure and Cross-Fostering
Our data identified a multigenerational transmission of maternal care, starting with the F1 all the way down to the F3 generation. Sexual maturation was found to be delayed in germ cell (F2) and transgenerationally exposed (F3-F4) generations and associated with epigenetic reprograming of the hypothalamus as identified in the F3 generation. Because no difference in sexual maturation, estrous cycle, or folliculogenesis were found in in utero and lactationally exposed females of the F1 generation, we hypothesized that the delay in sexual maturation in F2 and consecutive generations could be explained by hypothalamic reprograming caused by variations in maternal care, as described previously (Cameron et al. 2008a). To determine whether delayed puberty is caused via transgenerational female germ-cell inheritance or via a nongenomic multigenerational experience-based mechanism involving maternal care, a cross-fostering paradigm was carried out.
Control F2 pups raised by CTL (CC) or in utero and lactationally EDC-exposed (EC) F1 dams showed normal vaginal opening and estrous cyclicity. To the contrary, germline EDC-exposed pups (F2 generation) raised by EDC (EE) or CTL (CE) F1 dams, indistinctively showed a delay in vaginal opening (p < :01) ( Figure 7A; Table S4) and disrupted estrous cycles, characterized by a higher proportion of time spent in estrus and a lower proportion of time spent in diestrus ( Figure 7B). F2 EDC-exposed pups displayed a delay in vaginal opening, independently of being raised by an EDC dam displaying low licking and grooming or by a CTL dam displaying normal maternal care. Furthermore, adult F2 females that underwent cross-fostering and were raised by an EDC dam (EE and EC) showed less licking and grooming behavior ( Figure 7C).
To further explore the germline vs. nongenomic transmission of the phenotype, we examined the hypothalamic expression of the genes critical for pubertal onset, Kiss1, Esr1, and Oxt at P21 in females from the germline-exposed F2 generation. These animals showed a down-regulation in the hypothalamic expression of Kiss1, Esr1, and Oxt mRNA at P21, independently of being raised by an EDC (EE) or by a control (CE) dam ( Figure 7E-G, up; Table S4). Down-regulation of Kiss1, Esr1, and Oxt found in the F2 was consistently associated with a repressive chromatin state ( Figure 7E-G, down; Figure S11 and Table S6). A higher expression of the repressive mark H3K27me3 was found at the Kiss1p in EE pups. In germline EDC-exposed pups (EE and CE), Abundance of the TrxG-dependent activating marks H3K4me3 and H3K9ac and the PcG-dependent repressive mark H3K27me3 at the Kiss1, Esr1, Oxt, and Pomc promoter in the prepubertal MBH-PoA of females ancestrally exposed to a mixture of EDC (F3 generation), as measured by ChIP (n = 6=group). EDC data was normalized to control. Samples originated from six different litters per group. Dotted red lines represent repressive histone modifications (H3K27me3). Green solid lines represent activatory histone modifications (H3K27me3 and H3K9ac). Bars represent mean ± SEM ( * p < 0:05, ** p < 0:01, *** p < 0:001 vs. CTL, Student's t-test). Summary data are reported in Table S3. Note: AU, arbitrary units; ChIP, chromatin immunoprecipitation; EDC, endocrine-disrupting chemicals; Oxt, oxytocin; SEM, standard error of the mean; Th, tyrosine hydroxylase.
Esr1p was associated with lower expression of the activatory mark H3K9ac. Finally, Oxtp was associated with lower expression of the activatory H3K4me3 in EE pups.
These results show that cross-fostering could not restore normal pubertal timing in germ-cell EDC exposed animals, nor did EDC exposed dams cause delayed sexual maturation of unexposed pups. Our results strongly suggest that in these rats, delayed sexual maturation was caused by a germ-cell inheritance of the EDC mixture effects, independent of maternal care. To confirm these results, we further explored the consequences of cross-fostering on vaginal opening and hypothalamic gene expression in the transgenerationally exposed females of the F3 generation. These F3 females (EE and CE) displayed delayed vaginal opening, independent of maternal care. Furthermore, transgenerationally EDC-exposed F3 animals showed a down-regulation in the hypothalamic expression of Kiss1, Esr1, and Oxt at P21, independently of being raised by an EDC (EE) or control (CE) dam ( Figure 7I-K, up). Downregulation of Kiss1, Esr1, and Oxt found in the F3 generation was consistently associated with a repressive chromatin state ( Figure  7I-K, down). For instance, in EDC exposed pups (EE and CE) from the F3 generation, the Kiss1p was found to be associated with higher expression of the repressive H3K27me3. Furthermore, Esr1p and Oxtp were found to be associated with lower levels of H3K9ac and H3K4me3, respectively.
In our cross-fostering paradigm, maternal care was shown to be transmitted through experience of being raised by an EDC Figure 5. Maternal behavior displayed by female rats exposed to a mixture of EDC throughout four generations. (A-D, left) Time spent by dams displaying licking and grooming behavior toward pups after direct (F0 generation, n = 15=group), in utero, and through lactation (F1 generation, n = 10 − 11=group, obtained from 5 different litters per group), germ-cell (F2 generation, n = 11=group obtained from 5 different litters per group) or ancestral (F3 generation, n = 11=group, obtained from 5 different litters per group) exposure to an EDC mixture or vehicle from P2 to 8. Bar graphs show pooled time spent licking and grooming from P2-8. (A-D, right). Time spent by dams resting alone outside the nest not being involved in maternal care in F0-F3 females. Bar graphs show pooled time spent resting alone from P2-8. Plotted lines represent average time ± SEM ( * p < 0:05 vs. CTL, two-way ANOVA followed by Sidak's multiple comparisons test). Summary data are reported in Table S3. Note: ANOVA, analysis of variance; CTL, control; EDC, endocrine-disrupting chemical; SEM, standard error of the mean. EDC data from panels A, C, F, and H were normalized to control data. Samples originates from six different litters per group. Bars represent mean ± SEM ( * p < 0:05, ** p < 0:01, *** p < 0:001 vs. CTL). Data was analyzed using either two-way ANOVA (6A), followed by Sidak's multiple comparisons test or Student's t-test (6B-H). Summary data are reported in Table S3. Note: 3V, third ventricle; ANOVA, analysis of variance; AU, arbitrary units; ChIP, chromatin immunoprecipitation; CTL, control; EDC, endocrine-disrupting chemical; SN, substantia nigra; Th, tyrosine hydroxylase; VTA, ventral tegmental area. Figure 7. Sexual maturation and mRNA expression data in cross-fostered F2 offspring. (A) Average age at vaginal opening of cross-fostered germ cell (F2) EDC-exposed pups or controls raised by either in utero EDC-exposed dams or controls (n = 12-19=group). Samples originates from 6-12 F1 generation litters per group. (B) Percentage of cross-fostered females having regular cycle and time spent in each stage of the estrous cycle (n = 10=group). (C) Time spent by cross-fostered F2 dams displaying licking and grooming (n = 8=group). (D) Average age at VO in F3 females (n = 10=group). Samples originates from 6-10 litters. (E-H and M-P) Expression of Kiss1, Esr1, Oxt, and Th mRNA in the MBH-PoA of cross-fostered F2 and F3 pups at P21, as determined by qPCR (n = 6=group). Samples originates from six litters. Data from the CE, EC, and EE groups were normalized to the CC group. (I-L and Q-T) Abundance of the TrxG-dependent activating marks H3K4me3, H3K9ac, or the PcG-dependent repressive mark H3K27me3 at the Kiss1, Esr1, Oxt, and Th promoter in the prepubertal MBH-PoA of cross-fostered F2 and F3 females, as measured by ChIP (n = 6=group). Samples originates from six litters. Summary data are reported in Table S4. Bars represent mean ± SEM ( * p < 0:05, ** p < 0:01, *** p < 0:001 vs. CC; †p < 0:05, † †p < 0:01, † † †p < 0:001 vs. CE, one-way ANOVA). Note: ANOVA, analysis of variance; CC, control pup raised by control dam; CE, germ-cell EDC-exposed pup raised by control dam; ChIP, chromatin immunoprecipitation; D, diestrus; E, estrus; EC, control pup raised by in utero EDC-exposed dam; EDC, endocrine-disrupting chemical; EE, germ-cell EDC-exposed pup raised by in utero EDC-exposed dam; Oxt, oxytocin; P, proestrus; Reg. cycle, regular cycle; SEM, standard error of the mean; Th, tyrosine hydroxylase; VO, vaginal opening. dam, rather than through direct EDC mixture exposure. Because we observed lower dopaminergic signaling in the hypothalamus of transgenerationally exposed juvenile females ( Figure 6F), we further sought to determine whether Th mRNA expression was transmitted via altered maternal care or via EDC-mixture exposure-induced germline transmission. Our cross-fostering paradigm revealed that Th mRNA levels were down-regulated in F2 and F3 generations when raised by an EDC-exposed mother, independently of their ancestral (EE or EC) exposure, suggesting that the alterations in the dopaminergic system are transmitted through altered maternal care ( Figure 7H and I; Table S4). Th downregulation was associated with higher levels of the repressive H3K27me3 histone marks in the F3-exposed generation as found in the previous experiment ( Figure 7L; Figure S11).

Discussion
The potential multi-and transgenerational effects of developmental exposure to a low-dose EDC mixture on the neuroendocrine regulation of reproduction and maternal behavior via the maternal lineage in the female rat has never been addressed, to the best of our knowledge. In the present report, we provide evidence that in rats an EDC mixture exposure induced multi-and transgenerational alterations of maternal care and sexual maturation, respectively, through epigenetic reprogramming of the hypothalamus. The reproductive phenotype, found in germline-exposed F2 and transgenerationally exposed F3 females, is defined by a delay in vaginal opening, altered folliculogenesis, and estrous cyclicity. A similar transgenerational reproductive phenotype was observed after developmental exposure to pesticides, jet fuel, DEHP, and 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD) alone or in mixture. Alterations were characterized by abnormal puberty timing, lower number of primordial follicles, and accelerated follicle recruitment (Table 2; Manikkam et al. 2012aManikkam et al. , 2012bNilsson et al. 2012;Pocar et al. 2017). In our model, delayed sexual maturation in the F3 generation was accompanied by a modest but consistent delay in GnRH secretion maturation. The transgenerational consequences of the EDC mixture on sexual maturation were further confirmed by the presence of a delay in vaginal opening in the F4 generation. Furthermore, germline (F2) and transgenerationally exposed (F3) females had a significantly lower number of antral follicles and higher number of atretic follicles. This phenotype suggests a diminished ovarian efficiency, a sign of subfertility. Surprisingly, no significant reproductive effects were detected in the F1 generation that was directly exposed to the EDC mixture during embryonic and lactational development, which suggests that the reproductive phenotype may be a consequence of germline exposure. However, we also observed slightly but not significantly higher levels of atretic follicles in EDC-exposed F1 females, which might indicate a direct effect of EDCs on ovarian function as the expression of hypothalamic key players was not affected in F1 females.
The activation of the GnRH network during the juvenile period involves the transcriptional activation of KNDy neurons in the arcuate nucleus (ARC) (Han et al. 2015;Seminara et al. 2003). Developmental exposure to EDCs (Table 1) have been found to disrupt ARC kisspeptin neurons in adults, as shown for DnBP (Hu et al. 2013), DEHP Yu et al. 2020) and BPA (Losa-Ward et al. 2012;Losa et al. 2011;Ruiz-Pino et al. 2019). Here, we hypothesized that KNDy neurons could be a putative target of the EDC mixture disrupting puberty. In our study, we aimed at identifying and comparing hypothalamic transcriptional and epigenetic targets of direct in utero and lactational (F1) exposure to the EDC mixture vs. transgenerational (F3) exposure. RNAseq validation analysis using MBH-PoA hypothalamic explants from F3 females transgenerationally exposed to the EDC mixture showed that the reproductive effects were associated with lower Kiss1, Oxt, and Esr1 expression. These three genes play a crucial role in the control of the GnRH pulse generator during puberty and ovulation (Clarkson et al. 2010;Parent et al. 2008;Seminara et al. 2003;Wintermantel et al. 2006;Yeo and Herbison 2014). We found that F3 descendants of EDCtreated animals showed a disbalanced histone configuration at the Kiss1 promoter with higher expression of the repressive H3K27me3 and lower H3K4me3, a histone modification found at promoter regions of activated genes, suggesting that the delay in pubertal timing is mediated by an alteration of the epigenetic switch of PcG and TrG complexes. These results suggest that EDC mixture exposure affects the epigenetic programming of the hypothalamus transgenerationally. The EDC mixture may be specifically targeting the ARC kisspeptin population because of the quiescence of kisspeptin AVPV neurons during the juvenile period . However, because we used MBH-PoA explants in bulk, an approach that included different cell populations, futures studies including single-cell approaches should be used to characterize whether transgenerational effects of the EDC mixture target specific kisspeptin populations.
We have previously shown that metabolic cues target KNDy neurons through epigenetic modifications affecting pubertal timing (Vazquez et al. 2018). Undernutrition-induced delay in pubertal onset was found to be related to developmentally higher levels of the repression histone deacetylase sirtuin 1, responsible for lower expression of H3K9/16ac at the Kiss1 gene promoter (Vazquez et al. 2018). We showed here that the loss in Esr1 expression was associated with lower levels of H3K9ac at its promoter region, suggesting a role of the sirtuin deacetylases as putative mediators of EDC effects on reproduction. Our results indicate striking similarities between the epigenetic changes involved in the EDC-induced reprograming of the hypothalamus and those induced by undernutrition. In both cases, we identified a predominant role of the PcG and sirtuins in repressing gene expression by inducing higher levels of the repressive role of H3K27me3 and decreasing the effects of H3K9ac at the promoter region of puberty activating genes. This observation is in line with several reports suggesting that EDCs function as metabolic disruptors by sharing common intracellular pathways (Heindel et al. 2017). However, we did not find differences in body weight across generations, whereas some of the studied EDCs are known obesogens (Loganathan et al. 2021;Stahlhut et al. 2007). These EDCs, however, have never been studied in mixture at the doses used in this exposure. Surprisingly, we also found that the hypothalamic expression of Pomc and Cart, two postulated metabolic activators of the GnRH system (Parent et al. 2000;Roa and Herbison 2012), were transiently higher in exposed females. One could expect that such up-regulation of anorexigenic factors would be associated with modified weight gain. Future studies should address the expression of metabolic signals in specific nuclei.
Although delayed sexual maturation has been observed in many studies using single EDCs used in this study (Franssen et al. 2016;Naulé et al. 2014), contradicting results have also been reported following exposure to BPA at various doses (Fernandez et al. 2009;Franssen et al. 2016;Losa-Ward et al. 2012). These contradictory effects have been found to be explained by differences in cellular and molecular alterations in the hypothalamic control of reproduction (reviewed in López- Rodríguez et al. 2019). The divergence of the results may be caused by differences in dose, route of administration, windows of exposure, and, as in our case, the synergistic or additive effect of EDCs (Kortenkamp 2007).
We observed a delay in GnRH neuron maturation in the transgenerationally EDC-exposed females of the F3 generation. This hypothalamic incubation model has been shown to be a robust and reproductive tool to measure GnRH pulsatility (Bourguignon et al. 1989b;Rasier et al. 2007). Previous data have shown that GnRH secretion at this developmental stage was particularly sensitive to early exposure to EDCs (Franssen et al. 2016) and associated with consistent early or late vaginal opening. Although modest, longer GnRH interpulse interval in exposed animals were consistent with the delay in pubertal timing and the hypothalamic epigenetic reprogramming, suggesting a delayed maturation of the GnRH network.
Alterations in maternal care has been found to be transmitted through generations via a multigenerational nongenomic mechanism, affecting reproductive development by modulating hypothalamic Esr1 expression (Peña et al. 2013) and affecting the hypothalamic response to stress (Schmauss et al. 2014). In our study, in utero and lactational (F1), germ cell (F2) and transgenerationally EDC-exposed (F3) females showed diminished pup licking behavior and a longer proportion of time being active outside the nest. Only a single study has reported effects of BPA on licking and grooming behavior (Seta et al. 2005). Gene ontology analysis of hypothalamic RNAseq demonstrated that the dopaminergic system is one of the main targets of the EDC mixture in the F1 generation. In particular, EDC-exposed F1 females showed a significant loss in Th expression, suggesting a lower maternal motivation (Fang et al. 2018). This effect is limited to the mPoA because there was no TH staining loss in either the SN or the VTA. Although we found that lower levels of TH is cellspecific, it is not clear whether the low levels of TH in the mPoA is sufficient to explain the shorter time spent licking and grooming their pups. For instance, we did not observe differences in Esr1 in the F1 generations, known to be associated with low licking and grooming behavior (Peña et al. 2013). However, dopaminergic signaling alterations in the hypothalamus were found to be persistently affected transgenerationally and appeared to be mediated by epigenetic mechanisms. In exposed F3 females, Th promoter showed higher expression of the repressive H3K27me3a, and H3k9me3, as well as modestly higher levels in CpG methylation at three independent sites throughout the Th promoter. These findings show that the epigenetic reprograming of the Th locus of the F3 generation differs from F1. Moreover, although we observed lower maternal care from the F1 to the F3 generation, levels of licking and grooming were found to be variable between generations. Although we suggest that direct in utero and lactational exposure to the EDC mixture alters dopaminergic signaling and induces a multigenerational transmission of altered maternal care down to the F3 generation, further studies may address the question on whether maternal motivational brain circuits are being impaired. A maternal motivational paradigm by using a pup retrieval test could be used to determine whether alterations specifically target that component of maternal care.
Surprisingly, the reproductive phenotype appeared only after germline exposure (F2) to the EDC mixture. We did not find differences in sexual maturation in F1, as shown previously by using the same EDC mixture at high doses Johansson et al. 2016). Although we observed a trend toward higher levels of atretic follicles, we did not observe differences in estrous cyclicity which would suggest a diminished reproductive capability. Because natural variations in licking and grooming behavior have been related to sexual maturation (Cameron et al. 2008b) and our reproductive phenotype started to appear in the F2 generation, we carried out a cross-fostering experiment. The cross-fostering paradigm aimed at determining whether maternal and reproductive phenotype were related to each other. Additionally, we aimed at determining whether the reproductive phenotype was caused via germline inheritance or via nongenomic transmission. We identified that the reproductive phenotype was transmitted through the female germline and was not explained by impaired maternal behavior, because the crossfostering of the F2 generation with nonexposed mothers did not normalize pubertal timing or the epigenetically mediated downregulation of hypothalamic Kiss1, Esr1, and Oxt expression in the F2 or F3 generations. Germline-exposed (F2) pups raised by control or EDC dams displayed a delay in pubertal timing, which was also found in the transgenerationally exposed EDCs. A more modest effect in estrous cyclicity was observed in the F2 generation, suggesting that maternal care could partially restore the reproductive phenotype at adulthood, as suggested previously (Krishnan et al. 2019).
Concomitantly with the germline transmission of the reproductive phenotype, the cross-fostered EDC-exposed animals raised by control females showed a normalization of the hypothalamic dopaminergic network by higher Th expression, indicating that the reduction in maternal behavior is in fact learned and multigenerationally transmitted from the F1 generation. The transcriptional down-regulation of Th expression was found to be associated with epigenetic modifications. These results indicate that transgenerational disruption of sexual maturation and hypothalamic expression of genes critical for puberty induced by exposure to an EDC mixture in the female rat may happen through germ cell alterations. This EDC-induced epigenetic reprogramming needs to be resistant to erasure and to be transmitted across multiple generations (Jirtle and Skinner 2007). Thereafter, EDCinduced germ cell alterations flowed from germline to soma (Lim and Brunet 2013), affecting in this case the organization of the GnRH network. The mechanism by which EDCs affect germline programming and induce somatic alterations is largely unknown. It is unlikely that transgenerational alterations caused by the EDC mixture through the germline are restricted to the hypothalamus. Further studies should determine the developmental window at which the germ cells epigenome is more sensitive to environmental alterations. The combination of single-cell approaches together with transgenic models should help tracking epigenetic modifications from the germline to somatic tissues. In our model, we aimed at ensuring genetic background consistency by housing together one control and one EDC-exposed female with the same male during mating. Although this protocol allowed controlling for male genetic variance, we may have introduced a potential source of contamination through urine and feces when females were sharing housing.
Finally, we found transcriptional and epigenetic alteration of stress response genes in F3 females raised under diminished maternal care. F1-EDC females did not show differences in hypothalamic Nr3c1 or Crh expression while receiving normal maternal care from their F0 mothers. This finding confirms those of previous studies in rodents showing that maternal care induced increased responses of the HPA axis, decreasing Nr3c1 transcriptional levels via epigenetic alterations throughout generations (Liu et al. 1997;Weaver et al. 2004).
Altogether, the present study showed that developmental exposure to human-relevant doses of an EDC mixture altered the hypothalamic epigenetic programming of key puberty-activating genes in the female rat across generations, likely through the germline and maternal lineage. In addition, we showed that direct exposure to the EDC mixture down-regulated the hypothalamic dopaminergic network, affecting maternal behavior in a multigenerational and reversible manner.
conception. The authors are indebted to P. Delvenne for assistance with Papanicolaou staining and to Dr. V.D. Ramirez (Urbana, Illinois) for providing the CR11-B81 anti-GnRH antiserum. This work was supported by grants from the National institutes of Health (1R01HD084542) to A.L. and 8P51OD011092 for the operation of the Oregon National Primate Research Center as well as the FRS-FNRS (Fonds National de la Recherche Scientifique) (CDR-J013019F-33661942) and the Belgian Society for Pediatric Endocrinology and Diabetology. D.L.R. contributed to the study design and writing of the manuscript; conducted the physiological, molecular, and epigenetic experiments and statistical analysis; and evaluated and discussed the data. C.F.A. helped to perform RNA-seq and ChIP experiments; V.D., M.M., M.C., E.S., and A.G. were involved in physiological experiments and data analysis. A.L. and A.S.P. supervised the study design, analyzed the data, and wrote the manuscript, which was revised by the other authors. All the authors take full responsibility for the work.