Transgenerational Effects of Bisphenol A on Gene Expression and DNA Methylation of Imprinted Genes in Brain

Transgenerational Effects of Bisphenol A on Gene Expression and DNA Methylation of Imprinted... Abstract Bisphenol A (BPA) is a ubiquitous man-made endocrine disrupting compound (EDC). Developmental exposure to BPA changes behavioral and reproductive phenotypes, and these effects can last for generations. We exposed embryos to BPA, producing two lineages: controls and BPA exposed. In the third filial generation (F3), brain tissues containing the preoptic area, the bed nucleus of the stria terminalis, and the anterior hypothalamus were collected. RNA sequencing (RNA-seq) and subsequent data analyses revealed 50 differentially regulated genes in the brains of F3 juveniles from BPA vs control lineages. BPA exposure can lead to loss of imprinting, and one of the two imprinted genes in our data set, maternally expressed gene 3 (Meg3), has been associated with EDCs and neurobehavioral phenotypes. We used quantitative polymerase chain reaction to examine the two imprinted genes in our data set, Meg3 and microRNA-containing gene Mirg (residing in the same loci). Confirming the RNA-seq, Meg3 messenger RNA was higher in F3 brains from the BPA lineage than in control brains. This was true in brains from mice produced with two different BPA paradigms. Next, we used pyrosequencing to probe differentially methylated regions of Meg3. We found transgenerational effects of BPA on imprinted genes in brain. Given these results, and data on Meg3 methylation in humans, we suggest this gene may be a biomarker indicative of early life environmental perturbation. Endocrine disrupting compounds (EDCs) are man-made chemicals used to produce a large variety of daily-use materials ranging from linings of food cans to cosmetics and flame retardants. EDCs have enough structural and electrostatic similarity to steroid hormones that they are able to interfere with normal hormonal signaling. Bisphenol A (BPA) is the most widespread of the EDCs and is detected in >92% of the US population (1). Experimental studies of BPA are based largely on the fetal origins of adult disease hypothesis (2). Dams are dosed with BPA, which also exposes developing embryos, and dosing is often extended into lactation (3). The experimental subjects are the offspring, which are examined as adults. In addition to effects on the first filial generation (F1) offspring, BPA has transgenerational actions on subsequent generations (4–6). When multigenerational or transgenerational effects have been found, particularly in studies using low doses of EDCs, it has been assumed that epigenetic modifications are responsible for changes in gene expression (7). The first transgenerational study using BPA exposure to F1 embryos in utero revealed reduced fertility in the third filial generation (F3) of male rats (8, 9). More recently, in female mice, fertility decline has been recorded in F1, second filial generation (F2), and F3 exposed to BPA (6). The fact that some phenotypes persisted into the F3 makes them truly transgenerational. A mixture of phenotypes including obesity, infertility, and kidney disease were observed in F3 rats exposed ancestrally to a combination of BPA and two phthalates (10). In rats and mice, early life exposure to low doses of BPA affects many behaviors, including increased anxiety (11, 12), changes in activity (13), and reductions in some social behaviors (4, 14–17). Some of the differences in social behaviors persist in offspring for up to four generations, making them transgenerational (4, 17). Here, we asked which genes in the brain are changed transgenerationally by gestational exposure to BPA in the F3 generation. Work with the antifungal agent vinclozolin has shown changes in gene expression in F3 rats in several brain areas (5). In the current study, female mice consumed BPA incorporated into a phytoestrogen-free diet, whereas controls received the same diet without BPA. Nonsibling F1 adults were bred to produce the F2 generation, which were used to produce F3 mice. Changes found in the F3 generation are presumably transmitted via the germline and are likely to be permanent (18, 19). Using this paradigm, we have previously shown transgenerational differences between BPA and control lineages in behavior, targeted gene expression, and immunocytochemistry for estrogen receptor α (4, 17, 20). Brains from males in the F3 lineages were used for RNA sequencing (RNA-seq). We probed tissue from the combined preoptic area (POA), bed nucleus of the stria terminalis (BNST), and anterior hypothalamus (ANT HT). Because past studies have demonstrated multigenerational effects of BPA on imprinted genes (21, 22), we selected the two imprinted genes in our data set (from the Dlk1 to Dio3 domain) for subsequent analysis. One of these genes, maternally expressed gene 3 (Meg3), functions as a tumor suppressor (23, 24) and has been implicated in neurobehavioral problems (25, 26). In humans, changes in Meg3 DNA methylation have been associated with lead exposure, and in mice, methoxyclor modifies Meg3 methylation in sperm (27, 28). Precocious puberty onset in rats is correlated with Meg3 expression in brain (29). Meg3 RNA is expressed in the brain and pituitary (30, 31). Thus, changes in Meg3 expression could have transgenerational actions on behavior and reproduction via either the brain–pituitary–gonadal axis or the brain–pituitary–adrenal axis. The functions of Mirg are not known. Real-time quantitative polymerase chain reaction (qPCR) was used to measure expression of these genes in brains of F3 mice of both sexes, of several ages, and in two different mouse strains (C57BL/6J and FVB). Moreover, in addition to the C57BL/6J mice we used FVB mice bred such that females from the BPA lineage mated with unexposed males in each generation (6). This mating scheme (Fig. 1) ensures that transgenerational effects are transmitted via the dam. Finally, we examined DNA methylation for some of the 5′–C–phosphate–G–3′ (CpG) sites in the intergenic region of the Dlk1 to Dio3 domain’s intergenic differentially methylated region (IG-DMR) and the promoter of Meg3 differentially methylated region (DMR) in F3 brains. None of the CpG sites examined were significantly differentially methylated based on ancestral BPA status, and we conclude that other epigenetic mechanisms are involved in transgenerational modifications in Meg3 expression in brain. Figure 1. View largeDownload slide Cartoons illustrating the two breeding schemes used in the study. Dams are indicated by pink bows. BPA-exposed mice are indicated by gray color. (A) C57BL/6J black mice were paired, and half of the females (on left) were exposed to BPA in diet (gray mice). The subsequent generations were not exposed to additional BPA but were mated to each other to produce the BPA lineage. Black mice (right) were not exposed to BPA and served as the control lineage. (B) White FVB mice were paired, and half the pregnant females were fed BPA daily from embryonic day 11 to birth (left). The subsequent generations of females were mated with unexposed males (white). Control white mice (right) were not exposed to BPA. Figure 1. View largeDownload slide Cartoons illustrating the two breeding schemes used in the study. Dams are indicated by pink bows. BPA-exposed mice are indicated by gray color. (A) C57BL/6J black mice were paired, and half of the females (on left) were exposed to BPA in diet (gray mice). The subsequent generations were not exposed to additional BPA but were mated to each other to produce the BPA lineage. Black mice (right) were not exposed to BPA and served as the control lineage. (B) White FVB mice were paired, and half the pregnant females were fed BPA daily from embryonic day 11 to birth (left). The subsequent generations of females were mated with unexposed males (white). Control white mice (right) were not exposed to BPA. Methods and Materials Animals The C57BL/6 mice used were progeny of mice purchased from the Jackson Laboratory (stock no. 000664; Bar Harbor, ME). Adult females (between 8 and 10 weeks of age) were randomly assigned to either a phytoestrogen-free chow (Harlan Teklad, Madison, WI; catalog no. TD95092) or the identical chow supplemented with 5 mg/kg BPA diet (Harlan Teklad; catalog no.TD09386). Mice were placed on the diet 10 days before pairing with males for 2 weeks. All mice consumed food and water ad libitum; dams continued on their diets throughout gestation. During the dosing period, mice were observed daily for abnormal behavior and signs of toxicity. Lights were on a 12:12 light/dark cycle (lights off at 1200). We have previously calculated that pregnant C57BL/6J mice consume ~20 μg BPA per day at this dose and that their blood levels of BPA are within the range of those found in humans (17). As in our previous studies (4, 17), all F1 offspring were fostered within 24 hours of birth to dams on the control diet. Foster dams retained two of their own pups (these mice were not used in our studies) and received four foster pups (two of each sex). At weaning, postnatal day (PN) 21, mice were placed on standard chow (Harlan Teklad Diet; catalog no. 7912) containing phytoestrogens and group housed by litter and sex. Adult F1 males and females (nonsiblings) were mated to produce F2 mice, and these were mated to produce the F3 generation (Fig. 1A). Next, we asked whether our findings in C57BL/6J mice could generalize to other strains of inbred mice and different BPA transgenerational paradigms. We conducted an experiment with a second inbred strain, FVB (Charles River, Wilmington, MA). All mice received food (Harlan Teklad Rodent Diet containing phytoestrogens; catalog no. 8604) and water ad libitum. Using an outcrossed design (Fig. 1B), we mated female FVB mice (F0) with control male FVB mice at 12 weeks of age. When pregnancy was confirmed by the appearance of a vaginal sperm plug, females were removed from the males and individually caged. Here, BPA exposure was limited to the maternal line; sires in each generation of mating were unexposed. From embryonic day 11 to birth, female mice were orally dosed once daily in the early part of the lights-on phase of the 12:12 light/dark cycle with tocopherol-stripped corn oil containing one of three doses (0.5, 20, or 50 µg/kg body weight per day) or no BPA. Mice voluntarily consumed the corn oil administered in a pipette in the corner of their mouth (6, 32). All of these BPA doses are lower than the dose used in the experiments with C57BL/6J mice. During the dosing period, mice were observed daily for abnormal behavior and signs of toxicity. The F1 females were used to generate F2 females, and in turn, F2 females were used to generate F3 females (Fig. 1B). In these experiments, control and BPA-exposed F1 and F2 females were mated with fertility-confirmed, nonexposed males to generate the next generation. All procedures were in compliance with and approved by the University of Virginia, University of Illinois, or the North Carolina State University Animal Care and Use Committees. Brain collection and tissue punches For the RNA-seq, qPCR, and DNA methylation studies mice were euthanized, and brains were rapidly removed and frozen in powdered dry ice. From PN28 mice, we collected one piece of tissue that contained the BNST, the POA, and the ANT HT. Brains were cut (120 μM) in a cryostat (H/I Bright OTF5000; Hacker Instruments, Huntingdon, England) in the coronal plane, and tissue punches were collected as we have previously described (33). None of the mice had been used in behavioral studies. No more than one mouse of each sex per litter was used for each experiment. To collect RNA from PN0 pups, brains were removed and cut free-hand, and hypothalamus was collected and frozen. On PN4, F3 pups were euthanized, and whole brains were collected for use in the gene expression analyses described below. RNA-seq and analysis RNA from the BNST, POA, and ANT HT of each mouse (three males each from F1 and F3 controls and three each from F1 and F3 BPA lineages) was sent to Expression Analysis (Durham, NC) for sequencing. Poly-A tailed messenger RNA (mRNA) and long non-coding RNA (lncRNA) were extracted and quality tested. Complementary DNA (cDNA) was synthesized from the purified RNA and amplified and sequenced on an Illumina Hi-SEquation 2000. Between 30.7 and 46.7 million 50-base-pair paired-end reads were obtained for each of the 12 samples. The fragment size was ~175 base pairs, leaving ~75 bases between the paired-end reads unsequenced. We restricted this analysis to males to include genes from both sex chromosomes. The results of the RNA-seq were made available by Expression Analysis as 24 FASTQ files; one set of paired-end reads for each of the 12 RNA samples. The FASTQ files contained the read base calls and their associated Phred quality scores. Average Phred quality scores were >35 for all samples (34). The FASTQ files were inspected for individual base call quality and trimmed by one base. The quality-trimmed FASTQ files were then aligned to the mm9 mouse reference genome, and then those aligned reads were assigned to genes and quantified. Differential expression (DE) analysis was then done to determine which genes were differentially expressed between treatment and controls in both the F1 and F3 generations. Only the intersections of the genes from the results of two analysis approaches are presented here. Complete F3 RNA-seq data will be uploaded to National Center for Biotechnology Information GEO database. Inspecting and quality trimming the FASTQ files was done in the initial RNA-seq workflow by Expression Analysis using ea-utils, which is a proprietary tool of that company. In the initial workflow, the alignment step was done in STAR (35). STAR is an ultrafast, universal read alignment tool running on the University of Virginia 24-core Unix platform that is able to identify and map reads to splice junction sites. The assignment of the aligned sequence reads to genomic features, and the quantification of those reads was done with featureCounts (36). A principal component analysis plot of the regularized log transform of the data set indicated that one of the control F3 mice was an extreme outlier relative to all the other mice. This outlier mouse was removed from the data set, and the data set was renormalized for the subsequent downstream analyses. DE analysis of the featureCounts raw counts matrix data were conducted with DESeq2, a package in Bioconductor (37). Contrasts between BPA and control lineage mice produced a log2-fold change (lfc) measure of effect size and an adjusted P value measure of statistical significance for multiple comparisons for ~19,000 genes. Genes were ranked by both absolute value lfc and adjusted P value to determine ones that were statistically significant between two conditions. In this experiment, genes with lfc >0.5 and an adjusted P <0.05 were considered differentially expressed. The exact same workflow was executed on the FASTQ files from within the Galaxy platform (38) using software that is part of the Tuxedo Suite for RNA-seq analysis (39). The FASTQ files were first inspected for quality with FASTQC and then trimmed by one base with FASTX Trimmer (40, 41). Next, alignment to the University of California, Santa Cruz, mm9 mouse reference genome was performed in TopHat2 (42). TopHat2 alignment produced binary alignment map files, which were then input to the Tuxedo Suite DE analysis tool Cuffdiff (43). Cuffdiff output consisted of a table of DE calculations for ~19,000 genes with measurable expression levels, including FPKM effect sizes and adjusted P values, as well as other information on transcription start sites, coding sequences, promoters, and transcript-level DE. Real-time qPCR We performed follow-up studies to verify the results of RNA-seq, and here we report on the expression of imprinted genes Meg3 and Mirg. Both genes produce noncoding RNAs present in the imprinted Dlk1 to Dio3 region on mouse chromosome 12. The quantitative analyses of Meg3 and Mirg expression were conducted on samples derived from the POA, BNST, and ANT HT from PN28 (n = 4 to 9), hypothalamus from day-of-birth mice (PN0; n = 6 to 8), and PN4 whole brains (n = 5 to 6 per treatment group). As described above, all mice were in the F3 generation. For all analyses, only one male and one female per litter were used to avoid litter effects. Total RNA was extracted with QIAzol and purified with the miRNeasy Kit (Qiagen, Valencia, CA; catalog nos. 79306 and 217084) from frozen brain tissues. Concentration of RNA was measured with NanoDrop (Thermo Fisher Scientific Hanover Park, IL) and integrity checked by separation of RNA in agarose gel. Using Reverse Transcription II (Life Technology) and random primers (Life Technology, Carlsbad, CA; catalog nos. 18064-014 and 48190011), we reverse transcribed 300 ng of total RNA into cDNA according to the manufacturer’s protocol. cDNA was diluted before qPCR analysis. Each sample was analyzed in triplicate. Expression levels of two target genes (Meg3, Mirg) and endogenous control [β2 microglobulin (B2M)] were evaluated with an ABI StepOnePlus system (Thermo Fisher Scientific, Carlsbad, CA) and Fast Start SYBRGreen Master Mix (Life Technology; catalog no. 4385612). The following primers were used for amplification: Meg3 sense 5′-TGGGGATGGGTCTCTAGGTG-3′ and Meg3 antisense 5′-CCACTGACCCACAGTAACCC-3′, amplicon size 85 bp (NR_0203633.2 and NR_027651); Mirg sense 5′-TTGACTCCAGAAGATGCTCC-3′ and Mirg antisense 5′-CCTCAGGTTCCTAAGCAAGG-3′, amplicon size 170 bp (NR_028265.1); B2M sense 5′-GGCTCACACTGAATTCACCCCCAC-3′ and B2M antisense 5′-ACATGTCTCGATCCCAGTCGGT-3′, amplicon size 104 bp (NM_009735.3). Each set of primers was initially tested for efficiency (between 95% and 105%) and specificity, via melting curve analysis. For data evaluation, the comparative ΔΔ threshold cycle method was used. A calibrator sample was run on each plate to adjust for plate-to-plate variation. Samples with threshold cycle values >35 cycles, as well as outliers identified as samples with values above (or below) the 1.5-fold of the interquartile range from the third (or the first) quartile, were excluded from the analysis. DNA methylation Isolation of DNA DNA was isolated from BNST, POA, and ANT HT area of juvenile (PN28) F3 C57BL/6J mice (n = 5 or 6). Frozen brains were sectioned and collected as described previously. QIAamp DNA Micro Kit (Qiagen, Valencia, CA; catalog no. 56304) was used to isolate genomic DNA according to the manufacturer’s instructions. Bisulfite treatment and polymerase chain reaction amplification Bisulfite-treated DNA was used to evaluate the percentage of methylation in CpG sites in the promoter of Meg3 on chromosome 12 (44) (Supplemental Fig. 1). We treated 1 μg of DNA with sodium bisulfite by using the EZ DNA Methylation-Gold Kit (Zymo Research, Irvine, CA; catalog no. D5006). Briefly, 1 µg of DNA was combined with 2 µL of M-Dilution Buffer in 20 µl total volume and incubated at 37°C for 15 minutes before being combined with CT Conversion Reagent (Zymo Research). This mixture underwent bisulfite treatment at 98°C for 10 minutes and 53°C for 4 hours in a thermal cycler (Eppendorf AG, Hamburg, Germany). After desulfonation, DNA was purified with Zymo-Spin IC columns and eluted with 12 µL of M-Elution Buffer (Zymo Research; catalog no. C1004-250). As a control, commercially available Methylated Mouse DNA (Zymo Research; catalog no. D5012) and DNA isolated from matured mouse sperm were used. See Supplemental Table 1 for details. Pyrosequencing analysis of the IG-DMR and Meg3 DMR To analyze the methylation pattern of the IG-DMR, located the between Dlk1 and Gtl2 region, and the DMR in Meg3 promoter, we used a quantitative pyrosequencing assay. The regions of interest span the nucleotide positions for IG-DMR 81241 to 81540 (containing 29 CpGs) and for Meg3 DMR 94226 to 94488 (R5, containing six CpGs). In the National Center for Biotechnology Information database the gene accession number is AJ320506.1 and is referred to as the IG-DMR by Hiura et al. (45) and DMR for Meg3 promoter by Sato et al. (44). Bisulfite-treated DNA prepared as described above was amplified by polymerase chain reaction with the ZymoTaq DNA Polymerase Kit (ZymoResearch; catalog no. E2002). PCR was performed in buffer containing 5 mM of MgCl2 with 1 μL of bisulfite-treated DNA. Primers used for amplification of IG-DMR and Meg3 DMR are listed in Table 1. In both cases, reverse primers were conjugated with biotin at the 5′ end, to allow the enrichment of a single strand with streptavidin beads for the pyrosequencing reaction. The thermocycler conditions for each reaction are given in Supplemental Table 1. Single-strand amplicons were isolated with Pyrosequencing Work Station and sequenced on a Pyromark Q96 pyrosequencing instrument (Qiagen). Table 1. Sequences and Positions (AJ320506.1) of Primers Region  Forward Primer (Position)  Reverse Primer (Position)  Sequencing Primer (Position)  IG-DMR upper strand  5′-GTGGTTTGTTATGGGTAAGTTTT-3′ (81252 to 81274)  Biotin-5′-CTTCCCTCACTCCAAAAATTAAA-3′ (81545 to 81567)  5′-GGTAAGTTTTATGGTTTATTGTATA-3′ (81265 to 81289)  IG-DMR lower strand  5′-TTAGGAGTTAAGGAAAAGAAAGAAATAG-3′ (81529 to 81556)  Biotin-5′-ATCATAAACAAATCCCATAACTTACT-3′ (81259 to 81284)  5′-GTTAAGGAAAAGAAAGAAATAGT-3′ (81265 to 81289)  Meg3 DMR  5′-GTTAGTGTTGGGGATTTTTTTTTAAAG-3′ (94226 to 94252)  Biotin-5′-TCAACCACCAAAACCAAATTTTTAA-3′ (94464 to 94488)  S1: 5′-TTAGTTTGGTTTTTAGTATTTAATA-3′ (94261 to 94285)  S2: 5′-TTTATTAGGGTTTTTTTTTTTATTA-3′ (94348 to 94369)  Region  Forward Primer (Position)  Reverse Primer (Position)  Sequencing Primer (Position)  IG-DMR upper strand  5′-GTGGTTTGTTATGGGTAAGTTTT-3′ (81252 to 81274)  Biotin-5′-CTTCCCTCACTCCAAAAATTAAA-3′ (81545 to 81567)  5′-GGTAAGTTTTATGGTTTATTGTATA-3′ (81265 to 81289)  IG-DMR lower strand  5′-TTAGGAGTTAAGGAAAAGAAAGAAATAG-3′ (81529 to 81556)  Biotin-5′-ATCATAAACAAATCCCATAACTTACT-3′ (81259 to 81284)  5′-GTTAAGGAAAAGAAAGAAATAGT-3′ (81265 to 81289)  Meg3 DMR  5′-GTTAGTGTTGGGGATTTTTTTTTAAAG-3′ (94226 to 94252)  Biotin-5′-TCAACCACCAAAACCAAATTTTTAA-3′ (94464 to 94488)  S1: 5′-TTAGTTTGGTTTTTAGTATTTAATA-3′ (94261 to 94285)  S2: 5′-TTTATTAGGGTTTTTTTTTTTATTA-3′ (94348 to 94369)  Primer sequences for DNA methylation of the IG-DMR and the promoter (DMR) for Meg3. View Large Statistical analysis Statistically significant DE for the RNA-seq data set was defined in this study as a log2-fold change absolute value >0.5 and multiple-test adjusted P <0.05. The STAR/featureCounts/DESeq2 DE output and the Galaxy TopHat2/Cuffdiff DE output were examined to determine which genes fell into this category. The genes that were common to both approaches were subjected to further exploratory data analysis and visualization, including volcano plots and heat map clustering. Volcano plots and heat map clustering were performed in R Statistical Programming Language. Results from the qPCR and DNA methylation studies were evaluated with two-way analysis of variance followed by Tukey post hoc tests. Statistical significance was defined as P < 0.05. Results Identifying differentially expressed genes associated with BPA exposure We performed two separate RNA-seq analysis pipelines to identify differentially expressed genes. The STAR/featureCounts/DESeq2 analysis indicated 129 genes that were differentially expressed between BPA and control lineages in the F3 generation. The Galaxy TopHat2/Cuffdiff analysis identified 124 genes that were differentially expressed between BPA and controls in the F3 generation. Of these, 50 genes were common to both analyses, 45 were upregulated by BPA, and 5 were downregulated. The statistically significant (P < 0.05) differentially expressed genes common to both analyses are displayed in Table 2. Table 2. Differentially Expressed Genes From RNA-seqof F3 Brain Regions in Males Descended From BPA-Exposed F1 Females vs Unexposed Control Lineage F3 Gene  Gene Description  Fold Change  log2 Fold Change  P adj  Genes Upregulated in F3 due to F1 BPA Exposure         Adam8  Disintegrin and metallopeptidase domain 8  2.49  1.31  0.0001   Adamts16  Disintegrin-like and metallopeptidase (reprolysin type) with thrombospondin type 1 motif, 16  1.64  0.71  0.0175   Akap8l  Kinase (PRKA) anchor protein 8-like  1.49  0.57  0.0134   Ankrd24  Ankyrin repeat domain 24  1.50  0.59  0.0135   Atxn2l  Ataxin 2-like related protein  1.45  0.54  0.0028   Bzrap1  TSPO associated protein 1  1.46  0.55  0.0095   Celsr3  Cadherin, EGF LAG seven-pass G-type receptor 3  1.63  0.71  0.0049   Col11a1  Collagen, type XI, α 1  2.15  1.10  0.0000   Col24a1  Collagen, type XXIV, α 1  1.88  0.91  0.0153   Col4a2  Collagen, type IV, α 2  1.54  0.63  0.0134   Fam193b  Family with sequence similarity 193, member B  1.49  0.57  0.0162   Flnb  Filamin, β  1.63  0.71  0.0000   Gigyf1  GRB10 interacting GYF protein 1  1.52  0.61  0.0074   Gm14827  Predicted gene 14827  1.75  0.80  0.0095   Igsf9  Immunoglobulin superfamily, member 9  1.99  0.99  0.0029   Igsf9b  Immunoglobulin superfamily, member 9b adhesion protein  1.45  0.54  0.0287   Leng8  Leukocyte receptor cluster member 8  1.81  0.86  0.0442   Lime1  Lck interacting transmembrane adaptor 1  1.62  0.70  0.0024   Lrrc16b  Capping protein regulator and myosin 1 linker 3  1.57  0.65  0.0274   Lrrc45  Leucine rich repeat containing 45  1.64  0.71  0.0061   Malat1  Metastasis associated lung adenocarcinoma transcript 1  1.58  0.66  0.0134   Mapk11  Mitogen-activated protein kinase 11  1.53  0.62  0.0287   Meg3  Maternally expressed 3 long noncoding RNA  1.87  0.90  0.0007   Miat  Myocardial infarction associated long noncoding RNA  2.00  1.00  0.0000   Mirg  miRNA containing gene  1.49  0.58  0.0257   Nbeal2  Neurobeachin-like 2  1.65  0.72  0.0035   Neat1  Nuclear paraspeckle assembly transcript 1  1.76  0.81  0.0035   Nktr  Natural killer tumor recognition sequence  1.42  0.50  0.0061   Pan2  PAN2 poly(A) specific ribonuclease subunit  1.47  0.55  0.0100   Plekhg4  Pleckstrin homology domain containing, family G (with RhoGef domain) member 4  2.03  1.02  0.0028   Plekhn1  Pleckstrin homology domain containing, family N member 1  1.78  0.83  0.0019   Plxna3  Plexin A3  1.53  0.61  0.0017   Pnpla3  Patatin-like phospholipase domain containing 3  1.59  0.66  0.0028   Ptpru  Protein tyrosine phosphatase, receptor type, U  1.50  0.59  0.0076   Rgs11  Regulator of G-protein signaling 11  1.69  0.76  0.0098   Rps6kb2  Ribosomal protein S6 kinase, polypeptide 2  1.65  0.72  0.0122   Sema4g  Sema domain, immunoglobulin domain, transmembrane domain, and short cytoplasmic domain (semaphorin) 4G  1.53  0.61  0.0040   Shank1  SH3/ankyrin domain gene 1  1.51  0.59  0.0288   Sim1  Single-minded homolog 1  2.53  1.34  0.0050   Slc2a4rg-ps  Solute carrier family 39 (zinc transporter), member 2  1.76  0.82  0.0054   Slc39a2  SH3/ankyrin domain gene 1  2.02  1.01  0.0014   Snhg11  Small nucleolar RNA host gene 11  1.74  0.79  0.0007   Srrm2  Serine/arginine repetitive matrix 2  1.63  0.71  0.0000   Ssh3  Slingshot homolog 3  1.79  0.84  0.0028   Vwa5b2  von Willebrand factor A domain containing 5B2  1.67  0.74  0.0103  Genes Downregulated in F3 due to F1 BPA Exposure         Aldh1a1  Aldehyde dehydrogenase family 1, subfamily A1  0.66  −0.59  0.0226   Mt2  Metallothionein 2  0.60  −0.73  0.0019   Pmch  Promelanin-concentrating hormone  0.39  −1.36  0.0274   Rpl23  Ribosomal protein L23  0.66  −0.60  0.0386   Rps18  Ribosomal protein S18  0.67  −0.57  0.0050  Gene  Gene Description  Fold Change  log2 Fold Change  P adj  Genes Upregulated in F3 due to F1 BPA Exposure         Adam8  Disintegrin and metallopeptidase domain 8  2.49  1.31  0.0001   Adamts16  Disintegrin-like and metallopeptidase (reprolysin type) with thrombospondin type 1 motif, 16  1.64  0.71  0.0175   Akap8l  Kinase (PRKA) anchor protein 8-like  1.49  0.57  0.0134   Ankrd24  Ankyrin repeat domain 24  1.50  0.59  0.0135   Atxn2l  Ataxin 2-like related protein  1.45  0.54  0.0028   Bzrap1  TSPO associated protein 1  1.46  0.55  0.0095   Celsr3  Cadherin, EGF LAG seven-pass G-type receptor 3  1.63  0.71  0.0049   Col11a1  Collagen, type XI, α 1  2.15  1.10  0.0000   Col24a1  Collagen, type XXIV, α 1  1.88  0.91  0.0153   Col4a2  Collagen, type IV, α 2  1.54  0.63  0.0134   Fam193b  Family with sequence similarity 193, member B  1.49  0.57  0.0162   Flnb  Filamin, β  1.63  0.71  0.0000   Gigyf1  GRB10 interacting GYF protein 1  1.52  0.61  0.0074   Gm14827  Predicted gene 14827  1.75  0.80  0.0095   Igsf9  Immunoglobulin superfamily, member 9  1.99  0.99  0.0029   Igsf9b  Immunoglobulin superfamily, member 9b adhesion protein  1.45  0.54  0.0287   Leng8  Leukocyte receptor cluster member 8  1.81  0.86  0.0442   Lime1  Lck interacting transmembrane adaptor 1  1.62  0.70  0.0024   Lrrc16b  Capping protein regulator and myosin 1 linker 3  1.57  0.65  0.0274   Lrrc45  Leucine rich repeat containing 45  1.64  0.71  0.0061   Malat1  Metastasis associated lung adenocarcinoma transcript 1  1.58  0.66  0.0134   Mapk11  Mitogen-activated protein kinase 11  1.53  0.62  0.0287   Meg3  Maternally expressed 3 long noncoding RNA  1.87  0.90  0.0007   Miat  Myocardial infarction associated long noncoding RNA  2.00  1.00  0.0000   Mirg  miRNA containing gene  1.49  0.58  0.0257   Nbeal2  Neurobeachin-like 2  1.65  0.72  0.0035   Neat1  Nuclear paraspeckle assembly transcript 1  1.76  0.81  0.0035   Nktr  Natural killer tumor recognition sequence  1.42  0.50  0.0061   Pan2  PAN2 poly(A) specific ribonuclease subunit  1.47  0.55  0.0100   Plekhg4  Pleckstrin homology domain containing, family G (with RhoGef domain) member 4  2.03  1.02  0.0028   Plekhn1  Pleckstrin homology domain containing, family N member 1  1.78  0.83  0.0019   Plxna3  Plexin A3  1.53  0.61  0.0017   Pnpla3  Patatin-like phospholipase domain containing 3  1.59  0.66  0.0028   Ptpru  Protein tyrosine phosphatase, receptor type, U  1.50  0.59  0.0076   Rgs11  Regulator of G-protein signaling 11  1.69  0.76  0.0098   Rps6kb2  Ribosomal protein S6 kinase, polypeptide 2  1.65  0.72  0.0122   Sema4g  Sema domain, immunoglobulin domain, transmembrane domain, and short cytoplasmic domain (semaphorin) 4G  1.53  0.61  0.0040   Shank1  SH3/ankyrin domain gene 1  1.51  0.59  0.0288   Sim1  Single-minded homolog 1  2.53  1.34  0.0050   Slc2a4rg-ps  Solute carrier family 39 (zinc transporter), member 2  1.76  0.82  0.0054   Slc39a2  SH3/ankyrin domain gene 1  2.02  1.01  0.0014   Snhg11  Small nucleolar RNA host gene 11  1.74  0.79  0.0007   Srrm2  Serine/arginine repetitive matrix 2  1.63  0.71  0.0000   Ssh3  Slingshot homolog 3  1.79  0.84  0.0028   Vwa5b2  von Willebrand factor A domain containing 5B2  1.67  0.74  0.0103  Genes Downregulated in F3 due to F1 BPA Exposure         Aldh1a1  Aldehyde dehydrogenase family 1, subfamily A1  0.66  −0.59  0.0226   Mt2  Metallothionein 2  0.60  −0.73  0.0019   Pmch  Promelanin-concentrating hormone  0.39  −1.36  0.0274   Rpl23  Ribosomal protein L23  0.66  −0.60  0.0386   Rps18  Ribosomal protein S18  0.67  −0.57  0.0050  Differentially expressed genes obtained by RNA-seq of male F3 POA, BNST, and ANT HT. F3 males are descended from F1 females exposed to BPA in utero vs F3 controls from the same lineage. Fold change values and P values are from STAR/featureCounts/DESeq2 analysis and are the genes that were common to that workflow and the Galaxy TopHat2/Cuffdiff workflow having absolute value log2-fold changes >0.5 and adj P < 0.05. There were 45 upregulated and 5 downregulated common genes. The P adj is the Benjamini-Hochberg, multiple-test-adjusted P value. Abbreviations: adj, adjusted; miRNA, micro RNA. View Large Following DE analysis, the RNA-seq results were visualized to highlight especially important genes. The volcano plot in Fig. 2 shows effect size on the x-axis (log2 fold change) and a measure of statistical significance (−log10 of the adjusted P value) on the y-axis. The most statistically significant differentially expressed genes with the largest fold changes appear in the upper, outer corners of the plot. Our plot shows the 50 differentially regulated genes in the F3 generation. These genes have log2-fold change absolute values >0.5 and an adjusted P <0.05. Figure 2. View largeDownload slide Volcano plot showing the 50 differentially expressed genes in the F3 generation common to the STAR/featureCounts/DESeq2 and TopHat/Cuffdiff RNA-seq analyses. Only the section of the log-scale y-axis for which the adjusted P < 0.05 [−log10(0.05) = 1.3] is shown. Four genes (white dots on the right) were upregulated by a log2-fold change of >1 and one gene (white dot on the left) was downregulated by a log2-fold change of >1 without direct exposure to BPA. Ctl, control. Figure 2. View largeDownload slide Volcano plot showing the 50 differentially expressed genes in the F3 generation common to the STAR/featureCounts/DESeq2 and TopHat/Cuffdiff RNA-seq analyses. Only the section of the log-scale y-axis for which the adjusted P < 0.05 [−log10(0.05) = 1.3] is shown. Four genes (white dots on the right) were upregulated by a log2-fold change of >1 and one gene (white dot on the left) was downregulated by a log2-fold change of >1 without direct exposure to BPA. Ctl, control. To visualize genes whose expression levels may be correlated, genes were clustered based on a mathematical algorithm with correlation distance and average linkage, which groups genes with similar expression profiles across biological replicates. Clustered genes were depicted in a heat map (Fig. 3), where colors indicate relative expression levels of a gene across experiments. Progressively similar gene expression profiles were organized in a hierarchical tree structure dendrogram. We show the 50 genes common to both analyses that are differentially expressed between BPA and control mouse brains in the F3 generation. Finally, Ingenuity was used to determine the top canonical gene pathways (Supplemental Fig. 2), and top networks were combined to visualize unions (Supplemental Figs. 3 and 4). Because pathway analyses are more precise with larger data sets, we combined all important genes from the STAR/featureCounts/DESeq2 analysis (n = 129) and Galaxy TopHat2/Cuffdiff (n = 124). The full data set contained 203 genes because 50 genes were common to both analyses. BPA has been reported to change expression of imprinted genes, thus we selected the two imprinted genes (46, 47) in the data set, Meg3 and Mirg, for further analysis. Expression of both of these genes was higher in BPA lineage as compared with control brains. Figure 3. View largeDownload slide Heat map and clustering diagram for the 50 differentially expressed genes in the F3 generation common to the Star/featureCounts/DESeq2 and TopHat/Cuffdiff RNA-seq analyses with Benjamini-Hochberg multiple-test-adjusted P < 0.05 absolute value and lfcs > 0.5. The colors of the cells indicate the standardized number of normalized mRNA reads of the F3 BPA and control samples for that gene. In the top part of the heat map, the BPA F3 generation is more highly expressed (shades of purple) than the control F3 generation. In the bottom part of the heat map, the BPA F3 generation is less expressed (shades of light blue) than the control F3 generation. Note also that in the horizontal dendrogram the columns with three biological replicates for BPA F3 cluster together and that the columns with the two biological replicates for control F3 cluster together. Figure 3. View largeDownload slide Heat map and clustering diagram for the 50 differentially expressed genes in the F3 generation common to the Star/featureCounts/DESeq2 and TopHat/Cuffdiff RNA-seq analyses with Benjamini-Hochberg multiple-test-adjusted P < 0.05 absolute value and lfcs > 0.5. The colors of the cells indicate the standardized number of normalized mRNA reads of the F3 BPA and control samples for that gene. In the top part of the heat map, the BPA F3 generation is more highly expressed (shades of purple) than the control F3 generation. In the bottom part of the heat map, the BPA F3 generation is less expressed (shades of light blue) than the control F3 generation. Note also that in the horizontal dendrogram the columns with three biological replicates for BPA F3 cluster together and that the columns with the two biological replicates for control F3 cluster together. Real-time qPCR In F3 PN28 POA, BNST, and ANT HT the expression of Meg3 and Mirg varied by sex; in both cases, males had a significantly higher level of transcripts than females [F(1, 21) = 12.36, F(1, 16) = 22.52, respectively, P < 0.001 at least; Fig. 4A, 4D]. For Mirg, trends for an effect of lineage (P = 0.057) and an interaction between lineage and sex (P = 0.06) were noted. These trends were caused by higher gene expression in control males. In the case of Meg3, BPA lineage males had higher expression than any other group. The next two experiments extended these data to different brain regions, different developmental ages, and another mouse strain. Figure 4. View largeDownload slide Quantitative gene expression analysis of Meg3 and Mirg in (A and D) POA, BNST, and ANT HT of C57BL/6J, PN28 F3 offspring, (B and E) C57BL/6J, PN0 hypothalamus of F3 offspring, and (C and F) whole brains from PN4 F3 FVB mice. Means ± standard errors of the mean are shown. Numbers per group are listed in the individual histograms. Dark gray histograms denote males, and white denote females. Light gray histograms represent data from sexes combined, because there were no sex differences. *Significant overall sex difference, M > F, P < 0.05. **Trend for effect of lineage, BPA > control, P = 0.087. ***Significantly greater than all other male dose groups, P < 0.05. #Significantly lower than the female control group, P < 0.05. ##Significantly less than the control female group, P < 0.05. BF, BPA lineage female; BM, BPA lineage male; CF, control lineage female; CM, control lineage male. Figure 4. View largeDownload slide Quantitative gene expression analysis of Meg3 and Mirg in (A and D) POA, BNST, and ANT HT of C57BL/6J, PN28 F3 offspring, (B and E) C57BL/6J, PN0 hypothalamus of F3 offspring, and (C and F) whole brains from PN4 F3 FVB mice. Means ± standard errors of the mean are shown. Numbers per group are listed in the individual histograms. Dark gray histograms denote males, and white denote females. Light gray histograms represent data from sexes combined, because there were no sex differences. *Significant overall sex difference, M > F, P < 0.05. **Trend for effect of lineage, BPA > control, P = 0.087. ***Significantly greater than all other male dose groups, P < 0.05. #Significantly lower than the female control group, P < 0.05. ##Significantly less than the control female group, P < 0.05. BF, BPA lineage female; BM, BPA lineage male; CF, control lineage female; CM, control lineage male. The entire hypothalamus from PN0 C57BL/6J F3 mice revealed a trend for higher expression of Meg3 in both sexes from the BPA lineage [F(1, 25) = 3.17, P < 0.087] as compared with control lineage (Fig. 4B). No statistically significant effects of sex or lineage were found for Mirg (Fig. 4E). Lastly, expression of Meg3 was affected by sex [F(1, 29) = 15.06, P < 0.0006] and dose [F(3, 29) = 3.80, P < 0.02], and we noted an interaction [F(3, 29) = 14.54, P < 0.0001] between sex and dose (Fig. 4) in whole brains from FVB F3 pups. Exposure to 0.5 or 50 µg of BPA/kg per day significantly increased the expression of Meg3 in F3 FVB male brain (P < 0.03) as compared with controls. In females, Meg3 was reduced, compared with controls, only at the lowest BPA exposure level (P < 0.02) (Fig. 4C). Ancestral BPA status also affected expression of Mirg transcript [F(3, 30) = 7.36, P < 0.001], which was higher in control than in BPA lineage whole brains. The effect of lineage on Mirg expression was caused by differences between control brains, which have significantly greater expression than brains in the lowest- and the highest-dose BPA lineages (P < 0.001) (Fig. 4F). No sex differences or interactions were present. DNA methylation Because Meg3 gene expression in male brains was consistent with the RNA-seq data, we limited our DNA methylation study to this gene. We included PN28 BNST, POA, and ANT HT from both sexes to determine whether DNA methylation was correlated with sex differences in gene expression. We analyzed methylation status of 29 CpGs in the IG-DMR region (chromosome 12: 81241 to 81540), a region known to contain repeated motifs (48), and 6 CpGs in the Meg3 promoter region (chromosome 12: 94226 to 94488). We did not find any sex or exposure related differences in methylation status of IG-DMR. These data are displayed as heat maps (Fig. 5). However, in the Meg3 promoter, we found that three CpG sites were sexually dimorphic, with higher methylation in males as compared with females [F(1, 25) = 9.58, 4.35, 8.44, respectively; P < 0.05 at least] (Fig. 6). Ancestral BPA had no effect on DNA methylation in either sex. Considering that the transcript level of Meg3 was significantly lower in females, we expected CpG methylation to be higher in females than males. Thus, the observed sexually dimorphic changes in Meg3 transcript level are probably driven by mechanisms other than DNA methylation, although we cannot exclude the possibility that other CpG sites in promoter region or IG-DMR that we did not evaluate may contribute to those differences. Figure 5. View largeDownload slide Heat map and clustering diagram. (A) CpG regions 1 to 8 and 21 to 28 in the IG-DMR of the Dlk1 to Dio3 domain. (B) CpG sites in the Meg3 promoter DMR. The colors of the cells show percentage methylation, with red denoting highly methylated sites and green representing less methylated sites. Horizontal lines represent individual samples. BF, BPA lineage female; BM, BPA lineage male; CF, control female; CM, control male. Figure 5. View largeDownload slide Heat map and clustering diagram. (A) CpG regions 1 to 8 and 21 to 28 in the IG-DMR of the Dlk1 to Dio3 domain. (B) CpG sites in the Meg3 promoter DMR. The colors of the cells show percentage methylation, with red denoting highly methylated sites and green representing less methylated sites. Horizontal lines represent individual samples. BF, BPA lineage female; BM, BPA lineage male; CF, control female; CM, control male. Figure 6. View largeDownload slide (A–F) Methylation status of six CpG sites in Meg3 DMR R5 region. Means ± standard errors of the mean are graphed for CpG sites in six positions. The numbers per group are given in individual histograms. Dark gray histograms represent males, and white histograms represent females. *Significant sex difference in methylation, P < 0.05; males are hypermethylated as compared with females. BF, BPA lineage female; BM, BPA lineage male; CF, control female; CM, control male. Figure 6. View largeDownload slide (A–F) Methylation status of six CpG sites in Meg3 DMR R5 region. Means ± standard errors of the mean are graphed for CpG sites in six positions. The numbers per group are given in individual histograms. Dark gray histograms represent males, and white histograms represent females. *Significant sex difference in methylation, P < 0.05; males are hypermethylated as compared with females. BF, BPA lineage female; BM, BPA lineage male; CF, control female; CM, control male. Discussion Here we report on RNA-seq, qPCR, and DNA methylation data from mouse tissue exposed transgenerationally to an EDC. After performing two separate RNA-seq analysis pipelines we found 50 genes that were significantly different in both analyses. We selected two genes upregulated by BPA for our follow-up studies: Meg3 and Mirg. We chose these genes because BPA has documented effects on other imprinted genes in other tissues (21, 49, 50). In addition, both genes are lncRNAs that, as a class, have been implicated as important epigenetic factors involved in transgenerational inheritance (51). Meg3 and Mirg are related as they reside in the same imprinted region on mouse distal chromosome 12 (52). At least three paternally expressed genes (Dlk1, Rtl1, and Dio3) and four maternally expressed genes (Meg3, antiRlt1, Rian, and Mirg) are located in this region. Meg3 is of particular interest. Our data show that in males, Meg3, a maternally expressed lncRNA, is altered by ancestral gestational exposure to BPA, and this is via maternal inheritance. Meg3 is present in brain, pituitary, pancreas, placenta, and many other tissues (53, 54). The gene has four to five isoforms, all of which are present in brain (55). Meg3 is part of a tumor-suppressing pathway, which includes p53 (30, 56). In brain and pituitary, mutations of this gene are upregulated in several tumor types, including glioblastomas and pituitary tumors (24, 31). Expression of Meg3 is noted in 90% to 95% of pituitary cells and is not restricted by cell type (24). Meg3 is also indirectly involved in a variety of hormone-related functions including the timing of puberty in rats (29). Human epidemiology studies have shown that DNA hypermethylation of the Meg3 promoter in blood cells is correlated with poor infant temperament (25) and lead exposure–induced obesity, and hypomethylation is noted in sperm (57). The insecticide methoxychlor has partial estrogenic activity (58) and causes hypomethylation of Meg3 in sperm and liver of F1 and F2 mouse offspring from F0 treated dams (28). We speculate that this gene is sensitive to environmental variation and has the potential to be a useful biomarker in human studies. The RNA-seq data showed significantly more expression of Meg3 in male PN28 POA, BNST, and ANT HT. Real-time qPCR conducted on RNA from brains of C57BL/6J mice collected at two ages (PN28 and PN0), from two subregions (POA, BNST, and ANT HT and only the hypothalamus), and in both sexes, largely confirmed the RNA-seq data for Meg3. In all conditions, in males we noted at least a trend for more Meg3 expression in BPA vs control lineage brains. Interestingly, the qPCR data from the whole brains of FVB F3 mice confirm the data collected in C57BL/6J mice. Males from the lowest BPA dose lineage had significantly more Meg3 mRNA than controls. These significant and trending data confirm RNA-seq data collected in C57BL/6J PN28 brain, and the data generalize to another mouse strain, other BPA doses, and a different treatment and breeding paradigm. In contrast, RNA-seq results showed a significant elevation of Mirg in male BPA lineage brains, yet with qPCR an increase in expression in the BPA lineage brains was only confirmed only in the PN0 hypothalamus. In fact, the reverse effect was noted in the POA, BNST, and ANT HT of C57BL/6J mice and whole brains of PN4 FVB mice. Although we did not find BPA-related differences in Meg3 methylation, we did find small sex differences in methylation in three CpG sites. In these cases, females had lower methylation than males. This result is in contrast to our qPCR data showing higher Meg3 expression in males than females. In a microarray study conducted on cortex and hippocampal tissues from day-of-birth C57BL/6J mice, Meg3 was higher in female than in male brains (59). These data were confirmed with qPCR, and thus these results are the opposite of our qPCR data but are in agreement with the methylation results we report. A previous study from our group used microarray to compare gene expression in F3 whole embryonic brains on gestational day 18.5 exposed to the same dose of BPA used here (17). Comparison of the genes common to the two data sets revealed only two genes that overlap. Semaphorins are a class of well-characterized guidance molecules acting primarily during brain development (60). Semaphorin 4G (Scam4g) is necessary for granular cell migration during cerebellar formation in mice (61). The other gene, regulator of G protein signaling (Rgs11), has been studied in the retina and is involved in depolarization of bipolar cells and is associated with one of the glutamate receptors (62). Less is known about Rgs11 function in brain, where it is expressed throughout, and it is highest in human cerebellum (63). To date, one other study used RNA-seq to examine effects of first-generation in utero exposure to BPA (64). Two important brain areas were compared in neonatal rats: the hippocampus, important for learning and memory, and the hypothalamus, which is essential for regulation of the pituitary hormones and other reproductive functions. Interestingly, none of the important genes discovered in these areas match the data we present here. Arambula et al. (64) found differences in expression of several genes including Esr1, Esr2, and Oxt, and they validated their RNA-seq data with qPCR. These three genes were also differentially expressed in F1 mouse brains exposed to BPA in studies using other rodents and dosing paradigms (16, 17). Gene expression in F3 vinclozolin and control lineage rats has been assessed in several neural regions (65, 66). For some of these studies TaqMan qPCR arrays with 48 select genes were used. The focus of the study was on neuroendocrine-related functions. None of the genes identified by the TaqMan arrays are present in our data set (67). Microarray data from the same group did not reveal any overlap with the genes in our short list. We did not find any differences in methylation caused by ancestral BPA status in the section of the Meg3 promoter we explored in F3 brains. This finding is indirectly supported by a study limited to imprinted genes and DNA methylation (68) in which alterations in DNA methylation in F1 germ cells did not persist to F3. We speculate that other epigenetic mechanisms, such as histone modifications, regulate expression of Meg3 across generations. In fact, genomic imprinting of loci that include Meg3 is linked to increased histone acetylation (69). In human and mouse pancreatic tumors, Meg3 can be regulated by H3K4me3 and DNA methylation (70). Recently it has been shown that polycomb repressive complex 2 is required to maintain expression of maternal lncRNAs from the Meg3–Mirg locus in mouse embryonic stem cells. In its absence, the entire locus becomes silent because of de novo methylation at the IG-DMR (71). Thus, any disturbance, including environmental exposures, can lead to alterations in IG-DMR methylation, resulting in transcript-level changes. Future studies will examine other epigenetic mechanisms that may regulate gene expression transgenerationally as a result of ancestral BPA. Abbreviations: ANT HT anterior hypothalamus B2M β2 microglobulin BNST bed nucleus of the stria terminalis BPA bisphenol A cDNA complementary DNA CpG 5′–C–phosphate–G–3′ DE differential expression DMR differentially methylated region EDC endocrine disrupting compound F0 female FVB mice F1 first filial generation F2 second filial generation F3 third filial generation IG-DMR intergenic differentially methylated region lfc log-fold change lncRNA long non-coding RNA Meg3 maternally expressed gene 3 mRNA messenger RNA PN postnatal day POA preoptic area qPCR quantitative polymerase chain reaction RNA-seq RNA sequencing. Acknowledgments We thank Dr. Stephen D. Turner (University of Virginia) and Dr. David Scarr (Center for Human Health and the Environment at North Carolina State University) for technical assistance with parts of this study. Financial Support: This work was funded by National Institute of Environmental Health Sciences (NIEHS) Grants R01 ES022759 (to E.F.R.) and P01 ES022848 (to J.A.F.), and Environmental Protection Agency Grant RD-83459301 (to J.A.F.). We acknowledge help with the pyrosequencing from NIEHS Center for Human Health and the Environment Grant P30ES025128. A.D.H. gratefully acknowledges financial support from the 4-VA programs at James Madison University and the University of Virginia. Current Affiliation: J. T. Wolstenholme’s current affiliation is the Department of Pharmacology and Toxicology, Virginia Commonwealth University, Richmond, Virginia 23298. Disclosure Summary: The authors have nothing to disclose. References 1. Vandenberg LN, Chahoud I, Heindel JJ, Padmanabhan V, Paumgartten FJ, Schoenfelder G. Urinary, circulating, and tissue biomonitoring studies indicate widespread exposure to bisphenol A. Environ Health Perspect . 2010; 118( 8): 1055– 1070. Google Scholar CrossRef Search ADS PubMed  2. Barker DJ. A new model for the origins of chronic disease. Med Health Care Philos . 2001; 4( 1): 31– 35. Google Scholar CrossRef Search ADS PubMed  3. Welshons WV, Nagel SC, vom Saal FS. Large effects from small exposures. III. Endocrine mechanisms mediating effects of bisphenol A at levels of human exposure. Endocrinology . 2006; 147( 6, suppl): S56– S69. Google Scholar CrossRef Search ADS PubMed  4. Wolstenholme JT, Goldsby JA, Rissman EF. Transgenerational effects of prenatal bisphenol A on social recognition. Horm Behav . 2013; 64( 5): 833– 839. Google Scholar CrossRef Search ADS PubMed  5. Skinner MK, Anway MD, Savenkova MI, Gore AC, Crews D. Transgenerational epigenetic programming of the brain transcriptome and anxiety behavior. PLoS One . 2008; 3( 11): e3745. Google Scholar CrossRef Search ADS PubMed  6. Ziv-Gal A, Wang W, Zhou C, Flaws JA. The effects of in utero bisphenol A exposure on reproductive capacity in several generations of mice. Toxicol Appl Pharmacol . 2015; 284( 3): 354– 362. Google Scholar CrossRef Search ADS PubMed  7. Anderson OS, Kim JH, Peterson KE, Sanchez BN, Sant KE, Sartor MA, Weinhouse C, Dolinoy DC. Novel epigenetic biomarkers mediating bisphenol A exposure and metabolic phenotypes in female mice. Endocrinology . 2017; 158( 1): 31– 40. Google Scholar PubMed  8. Salian S, Doshi T, Vanage G. Perinatal exposure of rats to Bisphenol A affects the fertility of male offspring. Life Sci . 2009; 85( 21–22): 742– 752. Google Scholar CrossRef Search ADS PubMed  9. Salian S, Doshi T, Vanage G. Neonatal exposure of male rats to Bisphenol A impairs fertility and expression of sertoli cell junctional proteins in the testis. Toxicology . 2009; 265( 1–2): 56– 67. Google Scholar CrossRef Search ADS PubMed  10. Manikkam M, Tracey R, Guerrero-Bosagna C, Skinner MK. Plastics derived endocrine disruptors (BPA, DEHP and DBP) induce epigenetic transgenerational inheritance of obesity, reproductive disease and sperm epimutations. PLoS One . 2013; 8( 1): e55387. Google Scholar CrossRef Search ADS PubMed  11. Patisaul HB, Bateman HL. Neonatal exposure to endocrine active compounds or an ERbeta agonist increases adult anxiety and aggression in gonadally intact male rats. Horm Behav . 2008; 53( 4): 580– 588. Google Scholar CrossRef Search ADS PubMed  12. Rubin BS, Lenkowski JR, Schaeberle CM, Vandenberg LN, Ronsheim PM, Soto AM. Evidence of altered brain sexual differentiation in mice exposed perinatally to low, environmentally relevant levels of bisphenol A. Endocrinology . 2006; 147( 8): 3681– 3691. Google Scholar CrossRef Search ADS PubMed  13. Anderson OS, Peterson KE, Sanchez BN, Zhang Z, Mancuso P, Dolinoy DC. Perinatal bisphenol A exposure promotes hyperactivity, lean body composition, and hormonal responses across the murine life course. FASEB J . 2013; 27( 4): 1784– 1792. Google Scholar CrossRef Search ADS PubMed  14. Cox KH, Gatewood JD, Howeth C, Rissman EF. Gestational exposure to bisphenol A and cross-fostering affect behaviors in juvenile mice. Horm Behav . 2010; 58( 5): 754– 761. Google Scholar CrossRef Search ADS PubMed  15. Sullivan AW, Beach EC, Stetzik LA, Perry A, D’Addezio AS, Cushing BS, Patisaul HB. A novel model for neuroendocrine toxicology: neurobehavioral effects of BPA exposure in a prosocial species, the prairie vole (Microtus ochrogaster). Endocrinology . 2014; 155( 10): 3867– 3881. Google Scholar CrossRef Search ADS PubMed  16. Kundakovic M, Gudsnuk K, Franks B, Madrid J, Miller RL, Perera FP, Champagne FA. Sex-specific epigenetic disruption and behavioral changes following low-dose in utero bisphenol A exposure. Proc Natl Acad Sci USA . 2013; 110( 24): 9956– 9961. Google Scholar CrossRef Search ADS PubMed  17. Wolstenholme JT, Edwards M, Shetty SR, Gatewood JD, Taylor JA, Rissman EF, Connelly JJ. Gestational exposure to bisphenol a produces transgenerational changes in behaviors and gene expression. Endocrinology . 2012; 153( 8): 3828– 3838. Google Scholar CrossRef Search ADS PubMed  18. Anway MD, Skinner MK. Epigenetic transgenerational actions of endocrine disruptors. Endocrinology . 2006; 147( 6, suppl): S43– S49. Google Scholar CrossRef Search ADS PubMed  19. Jirtle RL, Skinner MK. Environmental epigenomics and disease susceptibility. Nat Rev Genet . 2007; 8( 4): 253– 262. Google Scholar CrossRef Search ADS PubMed  20. Goldsby JA, Wolstenholme JT, Rissman EF. Multi- and transgenerational consequences of bisphenol A on sexually dimorphic cell populations in mouse brain. Endocrinology . 2017; 158( 1): 21– 30. Google Scholar PubMed  21. Susiarjo M, Sasson I, Mesaros C, Bartolomei MS. Bisphenol A exposure disrupts genomic imprinting in the mouse. PLoS Genet . 2013; 9( 4): e1003401. Google Scholar CrossRef Search ADS PubMed  22. Kochmanski J, Marchlewicz EH, Savidge M, Montrose L, Faulk C, Dolinoy DC. Longitudinal effects of developmental bisphenol A and variable diet exposures on epigenetic drift in mice. Reprod Toxicol . 2017; 68: 154– 163. Google Scholar CrossRef Search ADS PubMed  23. Wang Y, Shen Y, Dai Q, Yang Q, Zhang Y, Wang X, Xie W, Luo Z, Lin C. A permissive chromatin state regulated by ZFP281-AFF3 in controlling the imprinted Meg3 polycistron. Nucleic Acids Res . 2017; 45( 3): 1177– 1185. Google Scholar CrossRef Search ADS PubMed  24. Gejman R, Batista DL, Zhong Y, Zhou Y, Zhang X, Swearingen B, Stratakis CA, Hedley-Whyte ET, Klibanski A. Selective loss of MEG3 expression and intergenic differentially methylated region hypermethylation in the MEG3/DLK1 locus in human clinically nonfunctioning pituitary adenomas. J Clin Endocrinol Metab . 2008; 93( 10): 4119– 4125. Google Scholar CrossRef Search ADS PubMed  25. Fuemmeler BF, Lee CT, Soubry A, Iversen ES, Huang Z, Murtha AP, Schildkraut JM, Jirtle RL, Murphy SK, Hoyo C. DNA methylation of regulatory regions of imprinted genes at birth and its relation to infant temperament. Genet Epigenet . 2016; 8: 59– 67. Google Scholar CrossRef Search ADS PubMed  26. Kagami M, Mizuno S, Matsubara K, Nakabayashi K, Sano S, Fuke T, Fukami M, Ogata T. Epimutations of the IG-DMR and the MEG3-DMR at the 14q32.2 imprinted region in two patients with Silver–Russell syndrome–compatible phenotype. Eur J Hum Genet . 2015; 23( 8): 1062– 1067. Google Scholar CrossRef Search ADS PubMed  27. Nye MD, King KE, Darrah TH, Maguire R, Jima DD, Huang Z, Mendez MA, Fry RC, Jirtle RL, Murphy SK, Hoyo C. Maternal blood lead concentrations, DNA methylation of MEG3 DMR regulating the DLK1/MEG3 imprinted domain and early growth in a multiethnic cohort. Environ Epigenet . 2016; 2( 1): dvv009. Google Scholar CrossRef Search ADS   28. Stouder C, Paoloni-Giacobino A. Specific transgenerational imprinting effects of the endocrine disruptor methoxychlor on male gametes. Reproduction . 2011; 141( 2): 207– 216. Google Scholar CrossRef Search ADS PubMed  29. Tao YH, Sharif N, Zeng BH, Cai YY, Guo YX. Lateral ventricle injection of orexin-A ameliorates central precocious puberty in rat via inhibiting the expression of MEG3. Int J Clin Exp Pathol . 2015; 8( 10): 12564– 12570. Google Scholar PubMed  30. Pease M, Ling C, Mack WJ, Wang K, Zada G. The role of epigenetic modification in tumorigenesis and progression of pituitary adenomas: a systematic review of the literature. PLoS One . 2013; 8( 12): e82619. Google Scholar CrossRef Search ADS PubMed  31. Zhang X, Zhou Y, Mehta KR, Danila DC, Scolavino S, Johnson SR, Klibanski A. A pituitary-derived MEG3 isoform functions as a growth suppressor in tumor cells. J Clin Endocrinol Metab . 2003; 88( 11): 5119– 5126. Google Scholar CrossRef Search ADS PubMed  32. Wang W, Hafner KS, Flaws JA. In utero bisphenol A exposure disrupts germ cell nest breakdown and reduces fertility with age in the mouse. Toxicol Appl Pharmacol . 2014; 276( 2): 157– 164. Google Scholar CrossRef Search ADS PubMed  33. Bonthuis PJ, Cox KH, Rissman EF. X-chromosome dosage affects male sexual behavior. Horm Behav . 2012; 61( 4): 565– 572. Google Scholar CrossRef Search ADS PubMed  34. Cock PJ, Fields CJ, Goto N, Heuer ML, Rice PM. The Sanger FASTQ file format for sequences with quality scores, and the Solexa/Illumina FASTQ variants. Nucleic Acids Res . 2010; 38( 6): 1767– 1771. Google Scholar CrossRef Search ADS PubMed  35. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, Batut P, Chaisson M, Gingeras TR. STAR: ultrafast universal RNA-seq aligner. Bioinformatics . 2013; 29( 1): 15– 21. Google Scholar CrossRef Search ADS PubMed  36. Liao Y, Smyth GK, Shi W. featureCounts: an efficient general purpose program for assigning sequence reads to genomic features. Bioinformatics . 2014; 30( 7): 923– 930. Google Scholar CrossRef Search ADS PubMed  37. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol . 2014; 15( 12): 550. Google Scholar CrossRef Search ADS PubMed  38. Blankenberg D, Von Kuster G, Coraor N, Ananda G, Lazarus R, Mangan M, Nekrutenko A, Taylor J. Galaxy: a web-based genome analysis tool for experimentalists. Curr Protoc Mol Biol . 2010; Chapter 19: Unit 19.10.1– 21. Google Scholar CrossRef Search ADS   39. Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, Pimentel H, Salzberg SL, Rinn JL, Pachter L. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat Protoc . 2012; 7( 3): 562– 578. Google Scholar CrossRef Search ADS PubMed  40. Kroll KW, Mokaram NE, Pelletier AR, Frankhouser DE, Westphal MS, Stump PA, Stump CL, Bundschuh R, Blachly JS, Yan P. Quality control for RNA-Seq (QuaCRS): an integrated quality control pipeline. Cancer Inform . 2014; 13( suppl 3): 7– 14. Google Scholar PubMed  41. Gordon A, Hannon G. Fastx-toolkit. FASTQ/A short-reads preprocessing tools. http://hannonlab.cshl.edu/fastx_toolkit/. 42. Trapnell C, Pachter L, Salzberg SL. TopHat: discovering splice junctions with RNA-Seq. Bioinformatics . 2009; 25( 9): 1105– 1111. Google Scholar CrossRef Search ADS PubMed  43. Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, Salzberg SL, Wold BJ, Pachter L. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol . 2010; 28( 5): 511– 515. Google Scholar CrossRef Search ADS PubMed  44. Sato S, Yoshida W, Soejima H, Nakabayashi K, Hata K. Methylation dynamics of IG-DMR and Gtl2-DMR during murine embryonic and placental development. Genomics . 2011; 98( 2): 120– 127. Google Scholar CrossRef Search ADS PubMed  45. Hiura H, Komiyama J, Shirai M, Obata Y, Ogawa H, Kono T. DNA methylation imprints on the IG-DMR of the Dlk1–Gtl2 domain in mouse male germline. FEBS Lett . 2007; 581( 7): 1255– 1260. Google Scholar CrossRef Search ADS PubMed  46. Tierling S, Dalbert S, Schoppenhorst S, Tsai CE, Oliger S, Ferguson-Smith AC, Paulsen M, Walter J. High-resolution map and imprinting analysis of the Gtl2-Dnchc1 domain on mouse chromosome 12. Genomics . 2006; 87( 2): 225– 235. Google Scholar CrossRef Search ADS PubMed  47. Miyoshi N, Wagatsuma H, Wakana S, Shiroishi T, Nomura M, Aisaka K, Kohda T, Surani MA, Kaneko-Ishino T, Ishino F. Identification of an imprinted gene, Meg3/Gtl2 and its human homologue MEG3, first mapped on mouse distal chromosome 12 and human chromosome 14q. Genes Cells . 2000; 5( 3): 211– 220. Google Scholar CrossRef Search ADS PubMed  48. Paulsen M, Takada S, Youngson NA, Benchaib M, Charlier C, Segers K, Georges M, Ferguson-Smith AC. Comparative sequence analysis of the imprinted Dlk1–Gtl2 locus in three mammalian species reveals highly conserved genomic elements and refines comparison with the Igf2–H19 region. Genome Res . 2001; 11( 12): 2085– 2094. Google Scholar CrossRef Search ADS PubMed  49. Doshi T, D’souza C, Vanage G. Aberrant DNA methylation at Igf2–H19 imprinting control region in spermatozoa upon neonatal exposure to bisphenol A and its association with post implantation loss. Mol Biol Rep . 2013; 40( 8): 4747– 4757. Google Scholar CrossRef Search ADS PubMed  50. Chao HH, Zhang XF, Chen B, Pan B, Zhang LJ, Li L, Sun XF, Shi QH, Shen W. Bisphenol A exposure modifies methylation of imprinted genes in mouse oocytes via the estrogen receptor signaling pathway. Histochem Cell Biol . 2012; 137( 2): 249– 259. Google Scholar CrossRef Search ADS PubMed  51. Rissman EF, Adli M. Minireview: transgenerational epigenetic inheritance: focus on endocrine disrupting compounds. Endocrinology . 2014; 155( 8): 2770– 2780. Google Scholar CrossRef Search ADS PubMed  52. Kobayashi S, Wagatsuma H, Ono R, Ichikawa H, Yamazaki M, Tashiro H, Aisaka K, Miyoshi N, Kohda T, Ogura A, Ohki M, Kaneko-Ishino T, Ishino F. Mouse Peg9/Dlk1 and human PEG9/DLK1 are paternally expressed imprinted genes closely located to the maternally expressed imprinted genes: mouse Meg3/Gtl2 and human MEG3. Genes Cells . 2000; 5( 12): 1029– 1037. Google Scholar CrossRef Search ADS PubMed  53. You L, Wang N, Yin D, Wang L, Jin F, Zhu Y, Yuan Q, De W. Downregulation of long noncoding RNA Meg3 affects insulin synthesis and secretion in mouse pancreatic beta cells. J Cell Physiol . 2016; 231( 4): 852– 862. Google Scholar CrossRef Search ADS PubMed  54. Murphy SK, Huang Z, Hoyo C. Differentially methylated regions of imprinted genes in prenatal, perinatal and postnatal human tissues. PLoS One . 2012; 7( 7): e40924. Google Scholar CrossRef Search ADS PubMed  55. Qu C, Jiang T, Li Y, Wang X, Cao H, Xu H, Qu J, Chen JG. Gene expression and IG-DMR hypomethylation of maternally expressed gene 3 in developing corticospinal neurons. Gene Expr Patterns . 2013; 13( 1–2): 51– 56. Google Scholar CrossRef Search ADS PubMed  56. Zhou Y, Zhong Y, Wang Y, Zhang X, Batista DL, Gejman R, Ansell PJ, Zhao J, Weng C, Klibanski A. Activation of p53 by MEG3 non-coding RNA. J Biol Chem . 2007; 282( 34): 24731– 24742. Google Scholar CrossRef Search ADS PubMed  57. Soubry A, Guo L, Huang Z, Hoyo C, Romanus S, Price T, Murphy SK. Obesity-related DNA methylation at imprinted genes in human sperm: results from the TIEGER study. Clin Epigenetics . 2016; 8( 1): 51. Google Scholar CrossRef Search ADS PubMed  58. Laws SC, Carey SA, Ferrell JM, Bodman GJ, Cooper RL. Estrogenic activity of octylphenol, nonylphenol, bisphenol A and methoxychlor in rats. Toxicol Sci . 2000; 54( 1): 154– 167. Google Scholar CrossRef Search ADS PubMed  59. Armoskus C, Moreira D, Bollinger K, Jimenez O, Taniguchi S, Tsai HW. Identification of sexually dimorphic genes in the neonatal mouse cortex and hippocampus. Brain Res . 2014; 1562: 23– 38. Google Scholar CrossRef Search ADS PubMed  60. Mann F, Chauvet S, Rougon G. Semaphorins in development and adult brain: implication for neurological diseases. Prog Neurobiol . 2007; 82( 2): 57– 79. Google Scholar CrossRef Search ADS PubMed  61. Maier V, Jolicoeur C, Rayburn H, Takegahara N, Kumanogoh A, Kikutani H, Tessier-Lavigne M, Wurst W, Friedel RH. Semaphorin 4C and 4G are ligands of Plexin-B2 required in cerebellar development. Mol Cell Neurosci . 2011; 46( 2): 419– 431. Google Scholar CrossRef Search ADS PubMed  62. Ray TA, Heath KM, Hasan N, Noel JM, Samuels IS, Martemyanov KA, Peachey NS, McCall MA, Gregg RG. GPR179 is required for high sensitivity of the mGluR6 signaling cascade in depolarizing bipolar cells. J Neurosci . 2014; 34( 18): 6334– 6343. Google Scholar CrossRef Search ADS PubMed  63. Larminie C, Murdock P, Walhin JP, Duckworth M, Blumer KJ, Scheideler MA, Garnier M. Selective expression of regulators of G-protein signaling (RGS) in the human central nervous system. Brain Res Mol Brain Res . 2004; 122( 1): 24– 34. Google Scholar CrossRef Search ADS PubMed  64. Arambula SE, Belcher SM, Planchart A, Turner SD, Patisaul HB. Impact of low dose oral exposure to bisphenol A (BPA) on the neonatal rat hypothalamic and hippocampal transcriptome: a CLARITY-BPA Consortium study. Endocrinology . 2016; 157( 10): 3856– 3872. Google Scholar CrossRef Search ADS PubMed  65. Crews D, Gillette R, Scarpino SV, Manikkam M, Savenkova MI, Skinner MK. Epigenetic transgenerational inheritance of altered stress responses. Proc Natl Acad Sci USA . 2012; 109( 23): 9143– 9148. Google Scholar CrossRef Search ADS PubMed  66. Skinner MK, Savenkova MI, Zhang B, Gore AC, Crews D. Gene bionetworks involved in the epigenetic transgenerational inheritance of altered mate preference: environmental epigenetics and evolutionary biology. BMC Genomics . 2014; 15( 1): 377. Google Scholar CrossRef Search ADS PubMed  67. Gillette R, Miller-Crews I, Skinner MK, Crews D. Distinct actions of ancestral vinclozolin and juvenile stress on neural gene expression in the male rat. Front Genet . 2015; 6: 56. Google Scholar CrossRef Search ADS PubMed  68. Iqbal K, Tran DA, Li AX, Warden C, Bai AY, Singh P, Wu X, Pfeifer GP, Szabó PE. Deleterious effects of endocrine disruptors are corrected in the mammalian germline by epigenome reprogramming. Genome Biol . 2015; 16( 1): 59. Google Scholar CrossRef Search ADS PubMed  69. Carr MS, Yevtodiyenko A, Schmidt CL, Schmidt JV. Allele-specific histone modifications regulate expression of the Dlk1–Gtl2 imprinted domain. Genomics . 2007; 89( 2): 280– 290. Google Scholar CrossRef Search ADS PubMed  70. Modali SD, Parekh VI, Kebebew E, Agarwal SK. Epigenetic regulation of the lncRNA MEG3 and its target c-MET in pancreatic neuroendocrine tumors. Mol Endocrinol . 2015; 29( 2): 224– 237. Google Scholar CrossRef Search ADS PubMed  71. Das PP, Hendrix DA, Apostolou E, Buchner AH, Canver MC, Beyaz S, Ljuboja D, Kuintzle R, Kim W, Karnik R, Shao Z, Xie H, Xu J, De Los Angeles A, Zhang Y, Choe J, Jun DL, Shen X, Gregory RI, Daley GQ, Meissner A, Kellis M, Hochedlinger K, Kim J, Orkin SH. PRC2 is required to maintain expression of the maternal Gtl2-Rian-Mirg locus by preventing de novo DNA methylation in mouse embryonic stem cells. Cell Reports . 2015; 12( 9): 1456– 1470. Google Scholar CrossRef Search ADS PubMed  Copyright © 2018 Endocrine Society http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Endocrinology Oxford University Press

Transgenerational Effects of Bisphenol A on Gene Expression and DNA Methylation of Imprinted Genes in Brain

Loading next page...
 
/lp/ou_press/transgenerational-effects-of-bisphenol-a-on-gene-expression-and-dna-wAXr4uC0H1
Publisher
Endocrine Society
Copyright
Copyright © 2018 Endocrine Society
ISSN
0013-7227
eISSN
1945-7170
D.O.I.
10.1210/en.2017-00730
Publisher site
See Article on Publisher Site

Abstract

Abstract Bisphenol A (BPA) is a ubiquitous man-made endocrine disrupting compound (EDC). Developmental exposure to BPA changes behavioral and reproductive phenotypes, and these effects can last for generations. We exposed embryos to BPA, producing two lineages: controls and BPA exposed. In the third filial generation (F3), brain tissues containing the preoptic area, the bed nucleus of the stria terminalis, and the anterior hypothalamus were collected. RNA sequencing (RNA-seq) and subsequent data analyses revealed 50 differentially regulated genes in the brains of F3 juveniles from BPA vs control lineages. BPA exposure can lead to loss of imprinting, and one of the two imprinted genes in our data set, maternally expressed gene 3 (Meg3), has been associated with EDCs and neurobehavioral phenotypes. We used quantitative polymerase chain reaction to examine the two imprinted genes in our data set, Meg3 and microRNA-containing gene Mirg (residing in the same loci). Confirming the RNA-seq, Meg3 messenger RNA was higher in F3 brains from the BPA lineage than in control brains. This was true in brains from mice produced with two different BPA paradigms. Next, we used pyrosequencing to probe differentially methylated regions of Meg3. We found transgenerational effects of BPA on imprinted genes in brain. Given these results, and data on Meg3 methylation in humans, we suggest this gene may be a biomarker indicative of early life environmental perturbation. Endocrine disrupting compounds (EDCs) are man-made chemicals used to produce a large variety of daily-use materials ranging from linings of food cans to cosmetics and flame retardants. EDCs have enough structural and electrostatic similarity to steroid hormones that they are able to interfere with normal hormonal signaling. Bisphenol A (BPA) is the most widespread of the EDCs and is detected in >92% of the US population (1). Experimental studies of BPA are based largely on the fetal origins of adult disease hypothesis (2). Dams are dosed with BPA, which also exposes developing embryos, and dosing is often extended into lactation (3). The experimental subjects are the offspring, which are examined as adults. In addition to effects on the first filial generation (F1) offspring, BPA has transgenerational actions on subsequent generations (4–6). When multigenerational or transgenerational effects have been found, particularly in studies using low doses of EDCs, it has been assumed that epigenetic modifications are responsible for changes in gene expression (7). The first transgenerational study using BPA exposure to F1 embryos in utero revealed reduced fertility in the third filial generation (F3) of male rats (8, 9). More recently, in female mice, fertility decline has been recorded in F1, second filial generation (F2), and F3 exposed to BPA (6). The fact that some phenotypes persisted into the F3 makes them truly transgenerational. A mixture of phenotypes including obesity, infertility, and kidney disease were observed in F3 rats exposed ancestrally to a combination of BPA and two phthalates (10). In rats and mice, early life exposure to low doses of BPA affects many behaviors, including increased anxiety (11, 12), changes in activity (13), and reductions in some social behaviors (4, 14–17). Some of the differences in social behaviors persist in offspring for up to four generations, making them transgenerational (4, 17). Here, we asked which genes in the brain are changed transgenerationally by gestational exposure to BPA in the F3 generation. Work with the antifungal agent vinclozolin has shown changes in gene expression in F3 rats in several brain areas (5). In the current study, female mice consumed BPA incorporated into a phytoestrogen-free diet, whereas controls received the same diet without BPA. Nonsibling F1 adults were bred to produce the F2 generation, which were used to produce F3 mice. Changes found in the F3 generation are presumably transmitted via the germline and are likely to be permanent (18, 19). Using this paradigm, we have previously shown transgenerational differences between BPA and control lineages in behavior, targeted gene expression, and immunocytochemistry for estrogen receptor α (4, 17, 20). Brains from males in the F3 lineages were used for RNA sequencing (RNA-seq). We probed tissue from the combined preoptic area (POA), bed nucleus of the stria terminalis (BNST), and anterior hypothalamus (ANT HT). Because past studies have demonstrated multigenerational effects of BPA on imprinted genes (21, 22), we selected the two imprinted genes in our data set (from the Dlk1 to Dio3 domain) for subsequent analysis. One of these genes, maternally expressed gene 3 (Meg3), functions as a tumor suppressor (23, 24) and has been implicated in neurobehavioral problems (25, 26). In humans, changes in Meg3 DNA methylation have been associated with lead exposure, and in mice, methoxyclor modifies Meg3 methylation in sperm (27, 28). Precocious puberty onset in rats is correlated with Meg3 expression in brain (29). Meg3 RNA is expressed in the brain and pituitary (30, 31). Thus, changes in Meg3 expression could have transgenerational actions on behavior and reproduction via either the brain–pituitary–gonadal axis or the brain–pituitary–adrenal axis. The functions of Mirg are not known. Real-time quantitative polymerase chain reaction (qPCR) was used to measure expression of these genes in brains of F3 mice of both sexes, of several ages, and in two different mouse strains (C57BL/6J and FVB). Moreover, in addition to the C57BL/6J mice we used FVB mice bred such that females from the BPA lineage mated with unexposed males in each generation (6). This mating scheme (Fig. 1) ensures that transgenerational effects are transmitted via the dam. Finally, we examined DNA methylation for some of the 5′–C–phosphate–G–3′ (CpG) sites in the intergenic region of the Dlk1 to Dio3 domain’s intergenic differentially methylated region (IG-DMR) and the promoter of Meg3 differentially methylated region (DMR) in F3 brains. None of the CpG sites examined were significantly differentially methylated based on ancestral BPA status, and we conclude that other epigenetic mechanisms are involved in transgenerational modifications in Meg3 expression in brain. Figure 1. View largeDownload slide Cartoons illustrating the two breeding schemes used in the study. Dams are indicated by pink bows. BPA-exposed mice are indicated by gray color. (A) C57BL/6J black mice were paired, and half of the females (on left) were exposed to BPA in diet (gray mice). The subsequent generations were not exposed to additional BPA but were mated to each other to produce the BPA lineage. Black mice (right) were not exposed to BPA and served as the control lineage. (B) White FVB mice were paired, and half the pregnant females were fed BPA daily from embryonic day 11 to birth (left). The subsequent generations of females were mated with unexposed males (white). Control white mice (right) were not exposed to BPA. Figure 1. View largeDownload slide Cartoons illustrating the two breeding schemes used in the study. Dams are indicated by pink bows. BPA-exposed mice are indicated by gray color. (A) C57BL/6J black mice were paired, and half of the females (on left) were exposed to BPA in diet (gray mice). The subsequent generations were not exposed to additional BPA but were mated to each other to produce the BPA lineage. Black mice (right) were not exposed to BPA and served as the control lineage. (B) White FVB mice were paired, and half the pregnant females were fed BPA daily from embryonic day 11 to birth (left). The subsequent generations of females were mated with unexposed males (white). Control white mice (right) were not exposed to BPA. Methods and Materials Animals The C57BL/6 mice used were progeny of mice purchased from the Jackson Laboratory (stock no. 000664; Bar Harbor, ME). Adult females (between 8 and 10 weeks of age) were randomly assigned to either a phytoestrogen-free chow (Harlan Teklad, Madison, WI; catalog no. TD95092) or the identical chow supplemented with 5 mg/kg BPA diet (Harlan Teklad; catalog no.TD09386). Mice were placed on the diet 10 days before pairing with males for 2 weeks. All mice consumed food and water ad libitum; dams continued on their diets throughout gestation. During the dosing period, mice were observed daily for abnormal behavior and signs of toxicity. Lights were on a 12:12 light/dark cycle (lights off at 1200). We have previously calculated that pregnant C57BL/6J mice consume ~20 μg BPA per day at this dose and that their blood levels of BPA are within the range of those found in humans (17). As in our previous studies (4, 17), all F1 offspring were fostered within 24 hours of birth to dams on the control diet. Foster dams retained two of their own pups (these mice were not used in our studies) and received four foster pups (two of each sex). At weaning, postnatal day (PN) 21, mice were placed on standard chow (Harlan Teklad Diet; catalog no. 7912) containing phytoestrogens and group housed by litter and sex. Adult F1 males and females (nonsiblings) were mated to produce F2 mice, and these were mated to produce the F3 generation (Fig. 1A). Next, we asked whether our findings in C57BL/6J mice could generalize to other strains of inbred mice and different BPA transgenerational paradigms. We conducted an experiment with a second inbred strain, FVB (Charles River, Wilmington, MA). All mice received food (Harlan Teklad Rodent Diet containing phytoestrogens; catalog no. 8604) and water ad libitum. Using an outcrossed design (Fig. 1B), we mated female FVB mice (F0) with control male FVB mice at 12 weeks of age. When pregnancy was confirmed by the appearance of a vaginal sperm plug, females were removed from the males and individually caged. Here, BPA exposure was limited to the maternal line; sires in each generation of mating were unexposed. From embryonic day 11 to birth, female mice were orally dosed once daily in the early part of the lights-on phase of the 12:12 light/dark cycle with tocopherol-stripped corn oil containing one of three doses (0.5, 20, or 50 µg/kg body weight per day) or no BPA. Mice voluntarily consumed the corn oil administered in a pipette in the corner of their mouth (6, 32). All of these BPA doses are lower than the dose used in the experiments with C57BL/6J mice. During the dosing period, mice were observed daily for abnormal behavior and signs of toxicity. The F1 females were used to generate F2 females, and in turn, F2 females were used to generate F3 females (Fig. 1B). In these experiments, control and BPA-exposed F1 and F2 females were mated with fertility-confirmed, nonexposed males to generate the next generation. All procedures were in compliance with and approved by the University of Virginia, University of Illinois, or the North Carolina State University Animal Care and Use Committees. Brain collection and tissue punches For the RNA-seq, qPCR, and DNA methylation studies mice were euthanized, and brains were rapidly removed and frozen in powdered dry ice. From PN28 mice, we collected one piece of tissue that contained the BNST, the POA, and the ANT HT. Brains were cut (120 μM) in a cryostat (H/I Bright OTF5000; Hacker Instruments, Huntingdon, England) in the coronal plane, and tissue punches were collected as we have previously described (33). None of the mice had been used in behavioral studies. No more than one mouse of each sex per litter was used for each experiment. To collect RNA from PN0 pups, brains were removed and cut free-hand, and hypothalamus was collected and frozen. On PN4, F3 pups were euthanized, and whole brains were collected for use in the gene expression analyses described below. RNA-seq and analysis RNA from the BNST, POA, and ANT HT of each mouse (three males each from F1 and F3 controls and three each from F1 and F3 BPA lineages) was sent to Expression Analysis (Durham, NC) for sequencing. Poly-A tailed messenger RNA (mRNA) and long non-coding RNA (lncRNA) were extracted and quality tested. Complementary DNA (cDNA) was synthesized from the purified RNA and amplified and sequenced on an Illumina Hi-SEquation 2000. Between 30.7 and 46.7 million 50-base-pair paired-end reads were obtained for each of the 12 samples. The fragment size was ~175 base pairs, leaving ~75 bases between the paired-end reads unsequenced. We restricted this analysis to males to include genes from both sex chromosomes. The results of the RNA-seq were made available by Expression Analysis as 24 FASTQ files; one set of paired-end reads for each of the 12 RNA samples. The FASTQ files contained the read base calls and their associated Phred quality scores. Average Phred quality scores were >35 for all samples (34). The FASTQ files were inspected for individual base call quality and trimmed by one base. The quality-trimmed FASTQ files were then aligned to the mm9 mouse reference genome, and then those aligned reads were assigned to genes and quantified. Differential expression (DE) analysis was then done to determine which genes were differentially expressed between treatment and controls in both the F1 and F3 generations. Only the intersections of the genes from the results of two analysis approaches are presented here. Complete F3 RNA-seq data will be uploaded to National Center for Biotechnology Information GEO database. Inspecting and quality trimming the FASTQ files was done in the initial RNA-seq workflow by Expression Analysis using ea-utils, which is a proprietary tool of that company. In the initial workflow, the alignment step was done in STAR (35). STAR is an ultrafast, universal read alignment tool running on the University of Virginia 24-core Unix platform that is able to identify and map reads to splice junction sites. The assignment of the aligned sequence reads to genomic features, and the quantification of those reads was done with featureCounts (36). A principal component analysis plot of the regularized log transform of the data set indicated that one of the control F3 mice was an extreme outlier relative to all the other mice. This outlier mouse was removed from the data set, and the data set was renormalized for the subsequent downstream analyses. DE analysis of the featureCounts raw counts matrix data were conducted with DESeq2, a package in Bioconductor (37). Contrasts between BPA and control lineage mice produced a log2-fold change (lfc) measure of effect size and an adjusted P value measure of statistical significance for multiple comparisons for ~19,000 genes. Genes were ranked by both absolute value lfc and adjusted P value to determine ones that were statistically significant between two conditions. In this experiment, genes with lfc >0.5 and an adjusted P <0.05 were considered differentially expressed. The exact same workflow was executed on the FASTQ files from within the Galaxy platform (38) using software that is part of the Tuxedo Suite for RNA-seq analysis (39). The FASTQ files were first inspected for quality with FASTQC and then trimmed by one base with FASTX Trimmer (40, 41). Next, alignment to the University of California, Santa Cruz, mm9 mouse reference genome was performed in TopHat2 (42). TopHat2 alignment produced binary alignment map files, which were then input to the Tuxedo Suite DE analysis tool Cuffdiff (43). Cuffdiff output consisted of a table of DE calculations for ~19,000 genes with measurable expression levels, including FPKM effect sizes and adjusted P values, as well as other information on transcription start sites, coding sequences, promoters, and transcript-level DE. Real-time qPCR We performed follow-up studies to verify the results of RNA-seq, and here we report on the expression of imprinted genes Meg3 and Mirg. Both genes produce noncoding RNAs present in the imprinted Dlk1 to Dio3 region on mouse chromosome 12. The quantitative analyses of Meg3 and Mirg expression were conducted on samples derived from the POA, BNST, and ANT HT from PN28 (n = 4 to 9), hypothalamus from day-of-birth mice (PN0; n = 6 to 8), and PN4 whole brains (n = 5 to 6 per treatment group). As described above, all mice were in the F3 generation. For all analyses, only one male and one female per litter were used to avoid litter effects. Total RNA was extracted with QIAzol and purified with the miRNeasy Kit (Qiagen, Valencia, CA; catalog nos. 79306 and 217084) from frozen brain tissues. Concentration of RNA was measured with NanoDrop (Thermo Fisher Scientific Hanover Park, IL) and integrity checked by separation of RNA in agarose gel. Using Reverse Transcription II (Life Technology) and random primers (Life Technology, Carlsbad, CA; catalog nos. 18064-014 and 48190011), we reverse transcribed 300 ng of total RNA into cDNA according to the manufacturer’s protocol. cDNA was diluted before qPCR analysis. Each sample was analyzed in triplicate. Expression levels of two target genes (Meg3, Mirg) and endogenous control [β2 microglobulin (B2M)] were evaluated with an ABI StepOnePlus system (Thermo Fisher Scientific, Carlsbad, CA) and Fast Start SYBRGreen Master Mix (Life Technology; catalog no. 4385612). The following primers were used for amplification: Meg3 sense 5′-TGGGGATGGGTCTCTAGGTG-3′ and Meg3 antisense 5′-CCACTGACCCACAGTAACCC-3′, amplicon size 85 bp (NR_0203633.2 and NR_027651); Mirg sense 5′-TTGACTCCAGAAGATGCTCC-3′ and Mirg antisense 5′-CCTCAGGTTCCTAAGCAAGG-3′, amplicon size 170 bp (NR_028265.1); B2M sense 5′-GGCTCACACTGAATTCACCCCCAC-3′ and B2M antisense 5′-ACATGTCTCGATCCCAGTCGGT-3′, amplicon size 104 bp (NM_009735.3). Each set of primers was initially tested for efficiency (between 95% and 105%) and specificity, via melting curve analysis. For data evaluation, the comparative ΔΔ threshold cycle method was used. A calibrator sample was run on each plate to adjust for plate-to-plate variation. Samples with threshold cycle values >35 cycles, as well as outliers identified as samples with values above (or below) the 1.5-fold of the interquartile range from the third (or the first) quartile, were excluded from the analysis. DNA methylation Isolation of DNA DNA was isolated from BNST, POA, and ANT HT area of juvenile (PN28) F3 C57BL/6J mice (n = 5 or 6). Frozen brains were sectioned and collected as described previously. QIAamp DNA Micro Kit (Qiagen, Valencia, CA; catalog no. 56304) was used to isolate genomic DNA according to the manufacturer’s instructions. Bisulfite treatment and polymerase chain reaction amplification Bisulfite-treated DNA was used to evaluate the percentage of methylation in CpG sites in the promoter of Meg3 on chromosome 12 (44) (Supplemental Fig. 1). We treated 1 μg of DNA with sodium bisulfite by using the EZ DNA Methylation-Gold Kit (Zymo Research, Irvine, CA; catalog no. D5006). Briefly, 1 µg of DNA was combined with 2 µL of M-Dilution Buffer in 20 µl total volume and incubated at 37°C for 15 minutes before being combined with CT Conversion Reagent (Zymo Research). This mixture underwent bisulfite treatment at 98°C for 10 minutes and 53°C for 4 hours in a thermal cycler (Eppendorf AG, Hamburg, Germany). After desulfonation, DNA was purified with Zymo-Spin IC columns and eluted with 12 µL of M-Elution Buffer (Zymo Research; catalog no. C1004-250). As a control, commercially available Methylated Mouse DNA (Zymo Research; catalog no. D5012) and DNA isolated from matured mouse sperm were used. See Supplemental Table 1 for details. Pyrosequencing analysis of the IG-DMR and Meg3 DMR To analyze the methylation pattern of the IG-DMR, located the between Dlk1 and Gtl2 region, and the DMR in Meg3 promoter, we used a quantitative pyrosequencing assay. The regions of interest span the nucleotide positions for IG-DMR 81241 to 81540 (containing 29 CpGs) and for Meg3 DMR 94226 to 94488 (R5, containing six CpGs). In the National Center for Biotechnology Information database the gene accession number is AJ320506.1 and is referred to as the IG-DMR by Hiura et al. (45) and DMR for Meg3 promoter by Sato et al. (44). Bisulfite-treated DNA prepared as described above was amplified by polymerase chain reaction with the ZymoTaq DNA Polymerase Kit (ZymoResearch; catalog no. E2002). PCR was performed in buffer containing 5 mM of MgCl2 with 1 μL of bisulfite-treated DNA. Primers used for amplification of IG-DMR and Meg3 DMR are listed in Table 1. In both cases, reverse primers were conjugated with biotin at the 5′ end, to allow the enrichment of a single strand with streptavidin beads for the pyrosequencing reaction. The thermocycler conditions for each reaction are given in Supplemental Table 1. Single-strand amplicons were isolated with Pyrosequencing Work Station and sequenced on a Pyromark Q96 pyrosequencing instrument (Qiagen). Table 1. Sequences and Positions (AJ320506.1) of Primers Region  Forward Primer (Position)  Reverse Primer (Position)  Sequencing Primer (Position)  IG-DMR upper strand  5′-GTGGTTTGTTATGGGTAAGTTTT-3′ (81252 to 81274)  Biotin-5′-CTTCCCTCACTCCAAAAATTAAA-3′ (81545 to 81567)  5′-GGTAAGTTTTATGGTTTATTGTATA-3′ (81265 to 81289)  IG-DMR lower strand  5′-TTAGGAGTTAAGGAAAAGAAAGAAATAG-3′ (81529 to 81556)  Biotin-5′-ATCATAAACAAATCCCATAACTTACT-3′ (81259 to 81284)  5′-GTTAAGGAAAAGAAAGAAATAGT-3′ (81265 to 81289)  Meg3 DMR  5′-GTTAGTGTTGGGGATTTTTTTTTAAAG-3′ (94226 to 94252)  Biotin-5′-TCAACCACCAAAACCAAATTTTTAA-3′ (94464 to 94488)  S1: 5′-TTAGTTTGGTTTTTAGTATTTAATA-3′ (94261 to 94285)  S2: 5′-TTTATTAGGGTTTTTTTTTTTATTA-3′ (94348 to 94369)  Region  Forward Primer (Position)  Reverse Primer (Position)  Sequencing Primer (Position)  IG-DMR upper strand  5′-GTGGTTTGTTATGGGTAAGTTTT-3′ (81252 to 81274)  Biotin-5′-CTTCCCTCACTCCAAAAATTAAA-3′ (81545 to 81567)  5′-GGTAAGTTTTATGGTTTATTGTATA-3′ (81265 to 81289)  IG-DMR lower strand  5′-TTAGGAGTTAAGGAAAAGAAAGAAATAG-3′ (81529 to 81556)  Biotin-5′-ATCATAAACAAATCCCATAACTTACT-3′ (81259 to 81284)  5′-GTTAAGGAAAAGAAAGAAATAGT-3′ (81265 to 81289)  Meg3 DMR  5′-GTTAGTGTTGGGGATTTTTTTTTAAAG-3′ (94226 to 94252)  Biotin-5′-TCAACCACCAAAACCAAATTTTTAA-3′ (94464 to 94488)  S1: 5′-TTAGTTTGGTTTTTAGTATTTAATA-3′ (94261 to 94285)  S2: 5′-TTTATTAGGGTTTTTTTTTTTATTA-3′ (94348 to 94369)  Primer sequences for DNA methylation of the IG-DMR and the promoter (DMR) for Meg3. View Large Statistical analysis Statistically significant DE for the RNA-seq data set was defined in this study as a log2-fold change absolute value >0.5 and multiple-test adjusted P <0.05. The STAR/featureCounts/DESeq2 DE output and the Galaxy TopHat2/Cuffdiff DE output were examined to determine which genes fell into this category. The genes that were common to both approaches were subjected to further exploratory data analysis and visualization, including volcano plots and heat map clustering. Volcano plots and heat map clustering were performed in R Statistical Programming Language. Results from the qPCR and DNA methylation studies were evaluated with two-way analysis of variance followed by Tukey post hoc tests. Statistical significance was defined as P < 0.05. Results Identifying differentially expressed genes associated with BPA exposure We performed two separate RNA-seq analysis pipelines to identify differentially expressed genes. The STAR/featureCounts/DESeq2 analysis indicated 129 genes that were differentially expressed between BPA and control lineages in the F3 generation. The Galaxy TopHat2/Cuffdiff analysis identified 124 genes that were differentially expressed between BPA and controls in the F3 generation. Of these, 50 genes were common to both analyses, 45 were upregulated by BPA, and 5 were downregulated. The statistically significant (P < 0.05) differentially expressed genes common to both analyses are displayed in Table 2. Table 2. Differentially Expressed Genes From RNA-seqof F3 Brain Regions in Males Descended From BPA-Exposed F1 Females vs Unexposed Control Lineage F3 Gene  Gene Description  Fold Change  log2 Fold Change  P adj  Genes Upregulated in F3 due to F1 BPA Exposure         Adam8  Disintegrin and metallopeptidase domain 8  2.49  1.31  0.0001   Adamts16  Disintegrin-like and metallopeptidase (reprolysin type) with thrombospondin type 1 motif, 16  1.64  0.71  0.0175   Akap8l  Kinase (PRKA) anchor protein 8-like  1.49  0.57  0.0134   Ankrd24  Ankyrin repeat domain 24  1.50  0.59  0.0135   Atxn2l  Ataxin 2-like related protein  1.45  0.54  0.0028   Bzrap1  TSPO associated protein 1  1.46  0.55  0.0095   Celsr3  Cadherin, EGF LAG seven-pass G-type receptor 3  1.63  0.71  0.0049   Col11a1  Collagen, type XI, α 1  2.15  1.10  0.0000   Col24a1  Collagen, type XXIV, α 1  1.88  0.91  0.0153   Col4a2  Collagen, type IV, α 2  1.54  0.63  0.0134   Fam193b  Family with sequence similarity 193, member B  1.49  0.57  0.0162   Flnb  Filamin, β  1.63  0.71  0.0000   Gigyf1  GRB10 interacting GYF protein 1  1.52  0.61  0.0074   Gm14827  Predicted gene 14827  1.75  0.80  0.0095   Igsf9  Immunoglobulin superfamily, member 9  1.99  0.99  0.0029   Igsf9b  Immunoglobulin superfamily, member 9b adhesion protein  1.45  0.54  0.0287   Leng8  Leukocyte receptor cluster member 8  1.81  0.86  0.0442   Lime1  Lck interacting transmembrane adaptor 1  1.62  0.70  0.0024   Lrrc16b  Capping protein regulator and myosin 1 linker 3  1.57  0.65  0.0274   Lrrc45  Leucine rich repeat containing 45  1.64  0.71  0.0061   Malat1  Metastasis associated lung adenocarcinoma transcript 1  1.58  0.66  0.0134   Mapk11  Mitogen-activated protein kinase 11  1.53  0.62  0.0287   Meg3  Maternally expressed 3 long noncoding RNA  1.87  0.90  0.0007   Miat  Myocardial infarction associated long noncoding RNA  2.00  1.00  0.0000   Mirg  miRNA containing gene  1.49  0.58  0.0257   Nbeal2  Neurobeachin-like 2  1.65  0.72  0.0035   Neat1  Nuclear paraspeckle assembly transcript 1  1.76  0.81  0.0035   Nktr  Natural killer tumor recognition sequence  1.42  0.50  0.0061   Pan2  PAN2 poly(A) specific ribonuclease subunit  1.47  0.55  0.0100   Plekhg4  Pleckstrin homology domain containing, family G (with RhoGef domain) member 4  2.03  1.02  0.0028   Plekhn1  Pleckstrin homology domain containing, family N member 1  1.78  0.83  0.0019   Plxna3  Plexin A3  1.53  0.61  0.0017   Pnpla3  Patatin-like phospholipase domain containing 3  1.59  0.66  0.0028   Ptpru  Protein tyrosine phosphatase, receptor type, U  1.50  0.59  0.0076   Rgs11  Regulator of G-protein signaling 11  1.69  0.76  0.0098   Rps6kb2  Ribosomal protein S6 kinase, polypeptide 2  1.65  0.72  0.0122   Sema4g  Sema domain, immunoglobulin domain, transmembrane domain, and short cytoplasmic domain (semaphorin) 4G  1.53  0.61  0.0040   Shank1  SH3/ankyrin domain gene 1  1.51  0.59  0.0288   Sim1  Single-minded homolog 1  2.53  1.34  0.0050   Slc2a4rg-ps  Solute carrier family 39 (zinc transporter), member 2  1.76  0.82  0.0054   Slc39a2  SH3/ankyrin domain gene 1  2.02  1.01  0.0014   Snhg11  Small nucleolar RNA host gene 11  1.74  0.79  0.0007   Srrm2  Serine/arginine repetitive matrix 2  1.63  0.71  0.0000   Ssh3  Slingshot homolog 3  1.79  0.84  0.0028   Vwa5b2  von Willebrand factor A domain containing 5B2  1.67  0.74  0.0103  Genes Downregulated in F3 due to F1 BPA Exposure         Aldh1a1  Aldehyde dehydrogenase family 1, subfamily A1  0.66  −0.59  0.0226   Mt2  Metallothionein 2  0.60  −0.73  0.0019   Pmch  Promelanin-concentrating hormone  0.39  −1.36  0.0274   Rpl23  Ribosomal protein L23  0.66  −0.60  0.0386   Rps18  Ribosomal protein S18  0.67  −0.57  0.0050  Gene  Gene Description  Fold Change  log2 Fold Change  P adj  Genes Upregulated in F3 due to F1 BPA Exposure         Adam8  Disintegrin and metallopeptidase domain 8  2.49  1.31  0.0001   Adamts16  Disintegrin-like and metallopeptidase (reprolysin type) with thrombospondin type 1 motif, 16  1.64  0.71  0.0175   Akap8l  Kinase (PRKA) anchor protein 8-like  1.49  0.57  0.0134   Ankrd24  Ankyrin repeat domain 24  1.50  0.59  0.0135   Atxn2l  Ataxin 2-like related protein  1.45  0.54  0.0028   Bzrap1  TSPO associated protein 1  1.46  0.55  0.0095   Celsr3  Cadherin, EGF LAG seven-pass G-type receptor 3  1.63  0.71  0.0049   Col11a1  Collagen, type XI, α 1  2.15  1.10  0.0000   Col24a1  Collagen, type XXIV, α 1  1.88  0.91  0.0153   Col4a2  Collagen, type IV, α 2  1.54  0.63  0.0134   Fam193b  Family with sequence similarity 193, member B  1.49  0.57  0.0162   Flnb  Filamin, β  1.63  0.71  0.0000   Gigyf1  GRB10 interacting GYF protein 1  1.52  0.61  0.0074   Gm14827  Predicted gene 14827  1.75  0.80  0.0095   Igsf9  Immunoglobulin superfamily, member 9  1.99  0.99  0.0029   Igsf9b  Immunoglobulin superfamily, member 9b adhesion protein  1.45  0.54  0.0287   Leng8  Leukocyte receptor cluster member 8  1.81  0.86  0.0442   Lime1  Lck interacting transmembrane adaptor 1  1.62  0.70  0.0024   Lrrc16b  Capping protein regulator and myosin 1 linker 3  1.57  0.65  0.0274   Lrrc45  Leucine rich repeat containing 45  1.64  0.71  0.0061   Malat1  Metastasis associated lung adenocarcinoma transcript 1  1.58  0.66  0.0134   Mapk11  Mitogen-activated protein kinase 11  1.53  0.62  0.0287   Meg3  Maternally expressed 3 long noncoding RNA  1.87  0.90  0.0007   Miat  Myocardial infarction associated long noncoding RNA  2.00  1.00  0.0000   Mirg  miRNA containing gene  1.49  0.58  0.0257   Nbeal2  Neurobeachin-like 2  1.65  0.72  0.0035   Neat1  Nuclear paraspeckle assembly transcript 1  1.76  0.81  0.0035   Nktr  Natural killer tumor recognition sequence  1.42  0.50  0.0061   Pan2  PAN2 poly(A) specific ribonuclease subunit  1.47  0.55  0.0100   Plekhg4  Pleckstrin homology domain containing, family G (with RhoGef domain) member 4  2.03  1.02  0.0028   Plekhn1  Pleckstrin homology domain containing, family N member 1  1.78  0.83  0.0019   Plxna3  Plexin A3  1.53  0.61  0.0017   Pnpla3  Patatin-like phospholipase domain containing 3  1.59  0.66  0.0028   Ptpru  Protein tyrosine phosphatase, receptor type, U  1.50  0.59  0.0076   Rgs11  Regulator of G-protein signaling 11  1.69  0.76  0.0098   Rps6kb2  Ribosomal protein S6 kinase, polypeptide 2  1.65  0.72  0.0122   Sema4g  Sema domain, immunoglobulin domain, transmembrane domain, and short cytoplasmic domain (semaphorin) 4G  1.53  0.61  0.0040   Shank1  SH3/ankyrin domain gene 1  1.51  0.59  0.0288   Sim1  Single-minded homolog 1  2.53  1.34  0.0050   Slc2a4rg-ps  Solute carrier family 39 (zinc transporter), member 2  1.76  0.82  0.0054   Slc39a2  SH3/ankyrin domain gene 1  2.02  1.01  0.0014   Snhg11  Small nucleolar RNA host gene 11  1.74  0.79  0.0007   Srrm2  Serine/arginine repetitive matrix 2  1.63  0.71  0.0000   Ssh3  Slingshot homolog 3  1.79  0.84  0.0028   Vwa5b2  von Willebrand factor A domain containing 5B2  1.67  0.74  0.0103  Genes Downregulated in F3 due to F1 BPA Exposure         Aldh1a1  Aldehyde dehydrogenase family 1, subfamily A1  0.66  −0.59  0.0226   Mt2  Metallothionein 2  0.60  −0.73  0.0019   Pmch  Promelanin-concentrating hormone  0.39  −1.36  0.0274   Rpl23  Ribosomal protein L23  0.66  −0.60  0.0386   Rps18  Ribosomal protein S18  0.67  −0.57  0.0050  Differentially expressed genes obtained by RNA-seq of male F3 POA, BNST, and ANT HT. F3 males are descended from F1 females exposed to BPA in utero vs F3 controls from the same lineage. Fold change values and P values are from STAR/featureCounts/DESeq2 analysis and are the genes that were common to that workflow and the Galaxy TopHat2/Cuffdiff workflow having absolute value log2-fold changes >0.5 and adj P < 0.05. There were 45 upregulated and 5 downregulated common genes. The P adj is the Benjamini-Hochberg, multiple-test-adjusted P value. Abbreviations: adj, adjusted; miRNA, micro RNA. View Large Following DE analysis, the RNA-seq results were visualized to highlight especially important genes. The volcano plot in Fig. 2 shows effect size on the x-axis (log2 fold change) and a measure of statistical significance (−log10 of the adjusted P value) on the y-axis. The most statistically significant differentially expressed genes with the largest fold changes appear in the upper, outer corners of the plot. Our plot shows the 50 differentially regulated genes in the F3 generation. These genes have log2-fold change absolute values >0.5 and an adjusted P <0.05. Figure 2. View largeDownload slide Volcano plot showing the 50 differentially expressed genes in the F3 generation common to the STAR/featureCounts/DESeq2 and TopHat/Cuffdiff RNA-seq analyses. Only the section of the log-scale y-axis for which the adjusted P < 0.05 [−log10(0.05) = 1.3] is shown. Four genes (white dots on the right) were upregulated by a log2-fold change of >1 and one gene (white dot on the left) was downregulated by a log2-fold change of >1 without direct exposure to BPA. Ctl, control. Figure 2. View largeDownload slide Volcano plot showing the 50 differentially expressed genes in the F3 generation common to the STAR/featureCounts/DESeq2 and TopHat/Cuffdiff RNA-seq analyses. Only the section of the log-scale y-axis for which the adjusted P < 0.05 [−log10(0.05) = 1.3] is shown. Four genes (white dots on the right) were upregulated by a log2-fold change of >1 and one gene (white dot on the left) was downregulated by a log2-fold change of >1 without direct exposure to BPA. Ctl, control. To visualize genes whose expression levels may be correlated, genes were clustered based on a mathematical algorithm with correlation distance and average linkage, which groups genes with similar expression profiles across biological replicates. Clustered genes were depicted in a heat map (Fig. 3), where colors indicate relative expression levels of a gene across experiments. Progressively similar gene expression profiles were organized in a hierarchical tree structure dendrogram. We show the 50 genes common to both analyses that are differentially expressed between BPA and control mouse brains in the F3 generation. Finally, Ingenuity was used to determine the top canonical gene pathways (Supplemental Fig. 2), and top networks were combined to visualize unions (Supplemental Figs. 3 and 4). Because pathway analyses are more precise with larger data sets, we combined all important genes from the STAR/featureCounts/DESeq2 analysis (n = 129) and Galaxy TopHat2/Cuffdiff (n = 124). The full data set contained 203 genes because 50 genes were common to both analyses. BPA has been reported to change expression of imprinted genes, thus we selected the two imprinted genes (46, 47) in the data set, Meg3 and Mirg, for further analysis. Expression of both of these genes was higher in BPA lineage as compared with control brains. Figure 3. View largeDownload slide Heat map and clustering diagram for the 50 differentially expressed genes in the F3 generation common to the Star/featureCounts/DESeq2 and TopHat/Cuffdiff RNA-seq analyses with Benjamini-Hochberg multiple-test-adjusted P < 0.05 absolute value and lfcs > 0.5. The colors of the cells indicate the standardized number of normalized mRNA reads of the F3 BPA and control samples for that gene. In the top part of the heat map, the BPA F3 generation is more highly expressed (shades of purple) than the control F3 generation. In the bottom part of the heat map, the BPA F3 generation is less expressed (shades of light blue) than the control F3 generation. Note also that in the horizontal dendrogram the columns with three biological replicates for BPA F3 cluster together and that the columns with the two biological replicates for control F3 cluster together. Figure 3. View largeDownload slide Heat map and clustering diagram for the 50 differentially expressed genes in the F3 generation common to the Star/featureCounts/DESeq2 and TopHat/Cuffdiff RNA-seq analyses with Benjamini-Hochberg multiple-test-adjusted P < 0.05 absolute value and lfcs > 0.5. The colors of the cells indicate the standardized number of normalized mRNA reads of the F3 BPA and control samples for that gene. In the top part of the heat map, the BPA F3 generation is more highly expressed (shades of purple) than the control F3 generation. In the bottom part of the heat map, the BPA F3 generation is less expressed (shades of light blue) than the control F3 generation. Note also that in the horizontal dendrogram the columns with three biological replicates for BPA F3 cluster together and that the columns with the two biological replicates for control F3 cluster together. Real-time qPCR In F3 PN28 POA, BNST, and ANT HT the expression of Meg3 and Mirg varied by sex; in both cases, males had a significantly higher level of transcripts than females [F(1, 21) = 12.36, F(1, 16) = 22.52, respectively, P < 0.001 at least; Fig. 4A, 4D]. For Mirg, trends for an effect of lineage (P = 0.057) and an interaction between lineage and sex (P = 0.06) were noted. These trends were caused by higher gene expression in control males. In the case of Meg3, BPA lineage males had higher expression than any other group. The next two experiments extended these data to different brain regions, different developmental ages, and another mouse strain. Figure 4. View largeDownload slide Quantitative gene expression analysis of Meg3 and Mirg in (A and D) POA, BNST, and ANT HT of C57BL/6J, PN28 F3 offspring, (B and E) C57BL/6J, PN0 hypothalamus of F3 offspring, and (C and F) whole brains from PN4 F3 FVB mice. Means ± standard errors of the mean are shown. Numbers per group are listed in the individual histograms. Dark gray histograms denote males, and white denote females. Light gray histograms represent data from sexes combined, because there were no sex differences. *Significant overall sex difference, M > F, P < 0.05. **Trend for effect of lineage, BPA > control, P = 0.087. ***Significantly greater than all other male dose groups, P < 0.05. #Significantly lower than the female control group, P < 0.05. ##Significantly less than the control female group, P < 0.05. BF, BPA lineage female; BM, BPA lineage male; CF, control lineage female; CM, control lineage male. Figure 4. View largeDownload slide Quantitative gene expression analysis of Meg3 and Mirg in (A and D) POA, BNST, and ANT HT of C57BL/6J, PN28 F3 offspring, (B and E) C57BL/6J, PN0 hypothalamus of F3 offspring, and (C and F) whole brains from PN4 F3 FVB mice. Means ± standard errors of the mean are shown. Numbers per group are listed in the individual histograms. Dark gray histograms denote males, and white denote females. Light gray histograms represent data from sexes combined, because there were no sex differences. *Significant overall sex difference, M > F, P < 0.05. **Trend for effect of lineage, BPA > control, P = 0.087. ***Significantly greater than all other male dose groups, P < 0.05. #Significantly lower than the female control group, P < 0.05. ##Significantly less than the control female group, P < 0.05. BF, BPA lineage female; BM, BPA lineage male; CF, control lineage female; CM, control lineage male. The entire hypothalamus from PN0 C57BL/6J F3 mice revealed a trend for higher expression of Meg3 in both sexes from the BPA lineage [F(1, 25) = 3.17, P < 0.087] as compared with control lineage (Fig. 4B). No statistically significant effects of sex or lineage were found for Mirg (Fig. 4E). Lastly, expression of Meg3 was affected by sex [F(1, 29) = 15.06, P < 0.0006] and dose [F(3, 29) = 3.80, P < 0.02], and we noted an interaction [F(3, 29) = 14.54, P < 0.0001] between sex and dose (Fig. 4) in whole brains from FVB F3 pups. Exposure to 0.5 or 50 µg of BPA/kg per day significantly increased the expression of Meg3 in F3 FVB male brain (P < 0.03) as compared with controls. In females, Meg3 was reduced, compared with controls, only at the lowest BPA exposure level (P < 0.02) (Fig. 4C). Ancestral BPA status also affected expression of Mirg transcript [F(3, 30) = 7.36, P < 0.001], which was higher in control than in BPA lineage whole brains. The effect of lineage on Mirg expression was caused by differences between control brains, which have significantly greater expression than brains in the lowest- and the highest-dose BPA lineages (P < 0.001) (Fig. 4F). No sex differences or interactions were present. DNA methylation Because Meg3 gene expression in male brains was consistent with the RNA-seq data, we limited our DNA methylation study to this gene. We included PN28 BNST, POA, and ANT HT from both sexes to determine whether DNA methylation was correlated with sex differences in gene expression. We analyzed methylation status of 29 CpGs in the IG-DMR region (chromosome 12: 81241 to 81540), a region known to contain repeated motifs (48), and 6 CpGs in the Meg3 promoter region (chromosome 12: 94226 to 94488). We did not find any sex or exposure related differences in methylation status of IG-DMR. These data are displayed as heat maps (Fig. 5). However, in the Meg3 promoter, we found that three CpG sites were sexually dimorphic, with higher methylation in males as compared with females [F(1, 25) = 9.58, 4.35, 8.44, respectively; P < 0.05 at least] (Fig. 6). Ancestral BPA had no effect on DNA methylation in either sex. Considering that the transcript level of Meg3 was significantly lower in females, we expected CpG methylation to be higher in females than males. Thus, the observed sexually dimorphic changes in Meg3 transcript level are probably driven by mechanisms other than DNA methylation, although we cannot exclude the possibility that other CpG sites in promoter region or IG-DMR that we did not evaluate may contribute to those differences. Figure 5. View largeDownload slide Heat map and clustering diagram. (A) CpG regions 1 to 8 and 21 to 28 in the IG-DMR of the Dlk1 to Dio3 domain. (B) CpG sites in the Meg3 promoter DMR. The colors of the cells show percentage methylation, with red denoting highly methylated sites and green representing less methylated sites. Horizontal lines represent individual samples. BF, BPA lineage female; BM, BPA lineage male; CF, control female; CM, control male. Figure 5. View largeDownload slide Heat map and clustering diagram. (A) CpG regions 1 to 8 and 21 to 28 in the IG-DMR of the Dlk1 to Dio3 domain. (B) CpG sites in the Meg3 promoter DMR. The colors of the cells show percentage methylation, with red denoting highly methylated sites and green representing less methylated sites. Horizontal lines represent individual samples. BF, BPA lineage female; BM, BPA lineage male; CF, control female; CM, control male. Figure 6. View largeDownload slide (A–F) Methylation status of six CpG sites in Meg3 DMR R5 region. Means ± standard errors of the mean are graphed for CpG sites in six positions. The numbers per group are given in individual histograms. Dark gray histograms represent males, and white histograms represent females. *Significant sex difference in methylation, P < 0.05; males are hypermethylated as compared with females. BF, BPA lineage female; BM, BPA lineage male; CF, control female; CM, control male. Figure 6. View largeDownload slide (A–F) Methylation status of six CpG sites in Meg3 DMR R5 region. Means ± standard errors of the mean are graphed for CpG sites in six positions. The numbers per group are given in individual histograms. Dark gray histograms represent males, and white histograms represent females. *Significant sex difference in methylation, P < 0.05; males are hypermethylated as compared with females. BF, BPA lineage female; BM, BPA lineage male; CF, control female; CM, control male. Discussion Here we report on RNA-seq, qPCR, and DNA methylation data from mouse tissue exposed transgenerationally to an EDC. After performing two separate RNA-seq analysis pipelines we found 50 genes that were significantly different in both analyses. We selected two genes upregulated by BPA for our follow-up studies: Meg3 and Mirg. We chose these genes because BPA has documented effects on other imprinted genes in other tissues (21, 49, 50). In addition, both genes are lncRNAs that, as a class, have been implicated as important epigenetic factors involved in transgenerational inheritance (51). Meg3 and Mirg are related as they reside in the same imprinted region on mouse distal chromosome 12 (52). At least three paternally expressed genes (Dlk1, Rtl1, and Dio3) and four maternally expressed genes (Meg3, antiRlt1, Rian, and Mirg) are located in this region. Meg3 is of particular interest. Our data show that in males, Meg3, a maternally expressed lncRNA, is altered by ancestral gestational exposure to BPA, and this is via maternal inheritance. Meg3 is present in brain, pituitary, pancreas, placenta, and many other tissues (53, 54). The gene has four to five isoforms, all of which are present in brain (55). Meg3 is part of a tumor-suppressing pathway, which includes p53 (30, 56). In brain and pituitary, mutations of this gene are upregulated in several tumor types, including glioblastomas and pituitary tumors (24, 31). Expression of Meg3 is noted in 90% to 95% of pituitary cells and is not restricted by cell type (24). Meg3 is also indirectly involved in a variety of hormone-related functions including the timing of puberty in rats (29). Human epidemiology studies have shown that DNA hypermethylation of the Meg3 promoter in blood cells is correlated with poor infant temperament (25) and lead exposure–induced obesity, and hypomethylation is noted in sperm (57). The insecticide methoxychlor has partial estrogenic activity (58) and causes hypomethylation of Meg3 in sperm and liver of F1 and F2 mouse offspring from F0 treated dams (28). We speculate that this gene is sensitive to environmental variation and has the potential to be a useful biomarker in human studies. The RNA-seq data showed significantly more expression of Meg3 in male PN28 POA, BNST, and ANT HT. Real-time qPCR conducted on RNA from brains of C57BL/6J mice collected at two ages (PN28 and PN0), from two subregions (POA, BNST, and ANT HT and only the hypothalamus), and in both sexes, largely confirmed the RNA-seq data for Meg3. In all conditions, in males we noted at least a trend for more Meg3 expression in BPA vs control lineage brains. Interestingly, the qPCR data from the whole brains of FVB F3 mice confirm the data collected in C57BL/6J mice. Males from the lowest BPA dose lineage had significantly more Meg3 mRNA than controls. These significant and trending data confirm RNA-seq data collected in C57BL/6J PN28 brain, and the data generalize to another mouse strain, other BPA doses, and a different treatment and breeding paradigm. In contrast, RNA-seq results showed a significant elevation of Mirg in male BPA lineage brains, yet with qPCR an increase in expression in the BPA lineage brains was only confirmed only in the PN0 hypothalamus. In fact, the reverse effect was noted in the POA, BNST, and ANT HT of C57BL/6J mice and whole brains of PN4 FVB mice. Although we did not find BPA-related differences in Meg3 methylation, we did find small sex differences in methylation in three CpG sites. In these cases, females had lower methylation than males. This result is in contrast to our qPCR data showing higher Meg3 expression in males than females. In a microarray study conducted on cortex and hippocampal tissues from day-of-birth C57BL/6J mice, Meg3 was higher in female than in male brains (59). These data were confirmed with qPCR, and thus these results are the opposite of our qPCR data but are in agreement with the methylation results we report. A previous study from our group used microarray to compare gene expression in F3 whole embryonic brains on gestational day 18.5 exposed to the same dose of BPA used here (17). Comparison of the genes common to the two data sets revealed only two genes that overlap. Semaphorins are a class of well-characterized guidance molecules acting primarily during brain development (60). Semaphorin 4G (Scam4g) is necessary for granular cell migration during cerebellar formation in mice (61). The other gene, regulator of G protein signaling (Rgs11), has been studied in the retina and is involved in depolarization of bipolar cells and is associated with one of the glutamate receptors (62). Less is known about Rgs11 function in brain, where it is expressed throughout, and it is highest in human cerebellum (63). To date, one other study used RNA-seq to examine effects of first-generation in utero exposure to BPA (64). Two important brain areas were compared in neonatal rats: the hippocampus, important for learning and memory, and the hypothalamus, which is essential for regulation of the pituitary hormones and other reproductive functions. Interestingly, none of the important genes discovered in these areas match the data we present here. Arambula et al. (64) found differences in expression of several genes including Esr1, Esr2, and Oxt, and they validated their RNA-seq data with qPCR. These three genes were also differentially expressed in F1 mouse brains exposed to BPA in studies using other rodents and dosing paradigms (16, 17). Gene expression in F3 vinclozolin and control lineage rats has been assessed in several neural regions (65, 66). For some of these studies TaqMan qPCR arrays with 48 select genes were used. The focus of the study was on neuroendocrine-related functions. None of the genes identified by the TaqMan arrays are present in our data set (67). Microarray data from the same group did not reveal any overlap with the genes in our short list. We did not find any differences in methylation caused by ancestral BPA status in the section of the Meg3 promoter we explored in F3 brains. This finding is indirectly supported by a study limited to imprinted genes and DNA methylation (68) in which alterations in DNA methylation in F1 germ cells did not persist to F3. We speculate that other epigenetic mechanisms, such as histone modifications, regulate expression of Meg3 across generations. In fact, genomic imprinting of loci that include Meg3 is linked to increased histone acetylation (69). In human and mouse pancreatic tumors, Meg3 can be regulated by H3K4me3 and DNA methylation (70). Recently it has been shown that polycomb repressive complex 2 is required to maintain expression of maternal lncRNAs from the Meg3–Mirg locus in mouse embryonic stem cells. In its absence, the entire locus becomes silent because of de novo methylation at the IG-DMR (71). Thus, any disturbance, including environmental exposures, can lead to alterations in IG-DMR methylation, resulting in transcript-level changes. Future studies will examine other epigenetic mechanisms that may regulate gene expression transgenerationally as a result of ancestral BPA. Abbreviations: ANT HT anterior hypothalamus B2M β2 microglobulin BNST bed nucleus of the stria terminalis BPA bisphenol A cDNA complementary DNA CpG 5′–C–phosphate–G–3′ DE differential expression DMR differentially methylated region EDC endocrine disrupting compound F0 female FVB mice F1 first filial generation F2 second filial generation F3 third filial generation IG-DMR intergenic differentially methylated region lfc log-fold change lncRNA long non-coding RNA Meg3 maternally expressed gene 3 mRNA messenger RNA PN postnatal day POA preoptic area qPCR quantitative polymerase chain reaction RNA-seq RNA sequencing. Acknowledgments We thank Dr. Stephen D. Turner (University of Virginia) and Dr. David Scarr (Center for Human Health and the Environment at North Carolina State University) for technical assistance with parts of this study. Financial Support: This work was funded by National Institute of Environmental Health Sciences (NIEHS) Grants R01 ES022759 (to E.F.R.) and P01 ES022848 (to J.A.F.), and Environmental Protection Agency Grant RD-83459301 (to J.A.F.). We acknowledge help with the pyrosequencing from NIEHS Center for Human Health and the Environment Grant P30ES025128. A.D.H. gratefully acknowledges financial support from the 4-VA programs at James Madison University and the University of Virginia. Current Affiliation: J. T. Wolstenholme’s current affiliation is the Department of Pharmacology and Toxicology, Virginia Commonwealth University, Richmond, Virginia 23298. Disclosure Summary: The authors have nothing to disclose. References 1. Vandenberg LN, Chahoud I, Heindel JJ, Padmanabhan V, Paumgartten FJ, Schoenfelder G. Urinary, circulating, and tissue biomonitoring studies indicate widespread exposure to bisphenol A. Environ Health Perspect . 2010; 118( 8): 1055– 1070. Google Scholar CrossRef Search ADS PubMed  2. Barker DJ. A new model for the origins of chronic disease. Med Health Care Philos . 2001; 4( 1): 31– 35. Google Scholar CrossRef Search ADS PubMed  3. Welshons WV, Nagel SC, vom Saal FS. Large effects from small exposures. III. Endocrine mechanisms mediating effects of bisphenol A at levels of human exposure. Endocrinology . 2006; 147( 6, suppl): S56– S69. Google Scholar CrossRef Search ADS PubMed  4. Wolstenholme JT, Goldsby JA, Rissman EF. Transgenerational effects of prenatal bisphenol A on social recognition. Horm Behav . 2013; 64( 5): 833– 839. Google Scholar CrossRef Search ADS PubMed  5. Skinner MK, Anway MD, Savenkova MI, Gore AC, Crews D. Transgenerational epigenetic programming of the brain transcriptome and anxiety behavior. PLoS One . 2008; 3( 11): e3745. Google Scholar CrossRef Search ADS PubMed  6. Ziv-Gal A, Wang W, Zhou C, Flaws JA. The effects of in utero bisphenol A exposure on reproductive capacity in several generations of mice. Toxicol Appl Pharmacol . 2015; 284( 3): 354– 362. Google Scholar CrossRef Search ADS PubMed  7. Anderson OS, Kim JH, Peterson KE, Sanchez BN, Sant KE, Sartor MA, Weinhouse C, Dolinoy DC. Novel epigenetic biomarkers mediating bisphenol A exposure and metabolic phenotypes in female mice. Endocrinology . 2017; 158( 1): 31– 40. Google Scholar PubMed  8. Salian S, Doshi T, Vanage G. Perinatal exposure of rats to Bisphenol A affects the fertility of male offspring. Life Sci . 2009; 85( 21–22): 742– 752. Google Scholar CrossRef Search ADS PubMed  9. Salian S, Doshi T, Vanage G. Neonatal exposure of male rats to Bisphenol A impairs fertility and expression of sertoli cell junctional proteins in the testis. Toxicology . 2009; 265( 1–2): 56– 67. Google Scholar CrossRef Search ADS PubMed  10. Manikkam M, Tracey R, Guerrero-Bosagna C, Skinner MK. Plastics derived endocrine disruptors (BPA, DEHP and DBP) induce epigenetic transgenerational inheritance of obesity, reproductive disease and sperm epimutations. PLoS One . 2013; 8( 1): e55387. Google Scholar CrossRef Search ADS PubMed  11. Patisaul HB, Bateman HL. Neonatal exposure to endocrine active compounds or an ERbeta agonist increases adult anxiety and aggression in gonadally intact male rats. Horm Behav . 2008; 53( 4): 580– 588. Google Scholar CrossRef Search ADS PubMed  12. Rubin BS, Lenkowski JR, Schaeberle CM, Vandenberg LN, Ronsheim PM, Soto AM. Evidence of altered brain sexual differentiation in mice exposed perinatally to low, environmentally relevant levels of bisphenol A. Endocrinology . 2006; 147( 8): 3681– 3691. Google Scholar CrossRef Search ADS PubMed  13. Anderson OS, Peterson KE, Sanchez BN, Zhang Z, Mancuso P, Dolinoy DC. Perinatal bisphenol A exposure promotes hyperactivity, lean body composition, and hormonal responses across the murine life course. FASEB J . 2013; 27( 4): 1784– 1792. Google Scholar CrossRef Search ADS PubMed  14. Cox KH, Gatewood JD, Howeth C, Rissman EF. Gestational exposure to bisphenol A and cross-fostering affect behaviors in juvenile mice. Horm Behav . 2010; 58( 5): 754– 761. Google Scholar CrossRef Search ADS PubMed  15. Sullivan AW, Beach EC, Stetzik LA, Perry A, D’Addezio AS, Cushing BS, Patisaul HB. A novel model for neuroendocrine toxicology: neurobehavioral effects of BPA exposure in a prosocial species, the prairie vole (Microtus ochrogaster). Endocrinology . 2014; 155( 10): 3867– 3881. Google Scholar CrossRef Search ADS PubMed  16. Kundakovic M, Gudsnuk K, Franks B, Madrid J, Miller RL, Perera FP, Champagne FA. Sex-specific epigenetic disruption and behavioral changes following low-dose in utero bisphenol A exposure. Proc Natl Acad Sci USA . 2013; 110( 24): 9956– 9961. Google Scholar CrossRef Search ADS PubMed  17. Wolstenholme JT, Edwards M, Shetty SR, Gatewood JD, Taylor JA, Rissman EF, Connelly JJ. Gestational exposure to bisphenol a produces transgenerational changes in behaviors and gene expression. Endocrinology . 2012; 153( 8): 3828– 3838. Google Scholar CrossRef Search ADS PubMed  18. Anway MD, Skinner MK. Epigenetic transgenerational actions of endocrine disruptors. Endocrinology . 2006; 147( 6, suppl): S43– S49. Google Scholar CrossRef Search ADS PubMed  19. Jirtle RL, Skinner MK. Environmental epigenomics and disease susceptibility. Nat Rev Genet . 2007; 8( 4): 253– 262. Google Scholar CrossRef Search ADS PubMed  20. Goldsby JA, Wolstenholme JT, Rissman EF. Multi- and transgenerational consequences of bisphenol A on sexually dimorphic cell populations in mouse brain. Endocrinology . 2017; 158( 1): 21– 30. Google Scholar PubMed  21. Susiarjo M, Sasson I, Mesaros C, Bartolomei MS. Bisphenol A exposure disrupts genomic imprinting in the mouse. PLoS Genet . 2013; 9( 4): e1003401. Google Scholar CrossRef Search ADS PubMed  22. Kochmanski J, Marchlewicz EH, Savidge M, Montrose L, Faulk C, Dolinoy DC. Longitudinal effects of developmental bisphenol A and variable diet exposures on epigenetic drift in mice. Reprod Toxicol . 2017; 68: 154– 163. Google Scholar CrossRef Search ADS PubMed  23. Wang Y, Shen Y, Dai Q, Yang Q, Zhang Y, Wang X, Xie W, Luo Z, Lin C. A permissive chromatin state regulated by ZFP281-AFF3 in controlling the imprinted Meg3 polycistron. Nucleic Acids Res . 2017; 45( 3): 1177– 1185. Google Scholar CrossRef Search ADS PubMed  24. Gejman R, Batista DL, Zhong Y, Zhou Y, Zhang X, Swearingen B, Stratakis CA, Hedley-Whyte ET, Klibanski A. Selective loss of MEG3 expression and intergenic differentially methylated region hypermethylation in the MEG3/DLK1 locus in human clinically nonfunctioning pituitary adenomas. J Clin Endocrinol Metab . 2008; 93( 10): 4119– 4125. Google Scholar CrossRef Search ADS PubMed  25. Fuemmeler BF, Lee CT, Soubry A, Iversen ES, Huang Z, Murtha AP, Schildkraut JM, Jirtle RL, Murphy SK, Hoyo C. DNA methylation of regulatory regions of imprinted genes at birth and its relation to infant temperament. Genet Epigenet . 2016; 8: 59– 67. Google Scholar CrossRef Search ADS PubMed  26. Kagami M, Mizuno S, Matsubara K, Nakabayashi K, Sano S, Fuke T, Fukami M, Ogata T. Epimutations of the IG-DMR and the MEG3-DMR at the 14q32.2 imprinted region in two patients with Silver–Russell syndrome–compatible phenotype. Eur J Hum Genet . 2015; 23( 8): 1062– 1067. Google Scholar CrossRef Search ADS PubMed  27. Nye MD, King KE, Darrah TH, Maguire R, Jima DD, Huang Z, Mendez MA, Fry RC, Jirtle RL, Murphy SK, Hoyo C. Maternal blood lead concentrations, DNA methylation of MEG3 DMR regulating the DLK1/MEG3 imprinted domain and early growth in a multiethnic cohort. Environ Epigenet . 2016; 2( 1): dvv009. Google Scholar CrossRef Search ADS   28. Stouder C, Paoloni-Giacobino A. Specific transgenerational imprinting effects of the endocrine disruptor methoxychlor on male gametes. Reproduction . 2011; 141( 2): 207– 216. Google Scholar CrossRef Search ADS PubMed  29. Tao YH, Sharif N, Zeng BH, Cai YY, Guo YX. Lateral ventricle injection of orexin-A ameliorates central precocious puberty in rat via inhibiting the expression of MEG3. Int J Clin Exp Pathol . 2015; 8( 10): 12564– 12570. Google Scholar PubMed  30. Pease M, Ling C, Mack WJ, Wang K, Zada G. The role of epigenetic modification in tumorigenesis and progression of pituitary adenomas: a systematic review of the literature. PLoS One . 2013; 8( 12): e82619. Google Scholar CrossRef Search ADS PubMed  31. Zhang X, Zhou Y, Mehta KR, Danila DC, Scolavino S, Johnson SR, Klibanski A. A pituitary-derived MEG3 isoform functions as a growth suppressor in tumor cells. J Clin Endocrinol Metab . 2003; 88( 11): 5119– 5126. Google Scholar CrossRef Search ADS PubMed  32. Wang W, Hafner KS, Flaws JA. In utero bisphenol A exposure disrupts germ cell nest breakdown and reduces fertility with age in the mouse. Toxicol Appl Pharmacol . 2014; 276( 2): 157– 164. Google Scholar CrossRef Search ADS PubMed  33. Bonthuis PJ, Cox KH, Rissman EF. X-chromosome dosage affects male sexual behavior. Horm Behav . 2012; 61( 4): 565– 572. Google Scholar CrossRef Search ADS PubMed  34. Cock PJ, Fields CJ, Goto N, Heuer ML, Rice PM. The Sanger FASTQ file format for sequences with quality scores, and the Solexa/Illumina FASTQ variants. Nucleic Acids Res . 2010; 38( 6): 1767– 1771. Google Scholar CrossRef Search ADS PubMed  35. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, Batut P, Chaisson M, Gingeras TR. STAR: ultrafast universal RNA-seq aligner. Bioinformatics . 2013; 29( 1): 15– 21. Google Scholar CrossRef Search ADS PubMed  36. Liao Y, Smyth GK, Shi W. featureCounts: an efficient general purpose program for assigning sequence reads to genomic features. Bioinformatics . 2014; 30( 7): 923– 930. Google Scholar CrossRef Search ADS PubMed  37. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol . 2014; 15( 12): 550. Google Scholar CrossRef Search ADS PubMed  38. Blankenberg D, Von Kuster G, Coraor N, Ananda G, Lazarus R, Mangan M, Nekrutenko A, Taylor J. Galaxy: a web-based genome analysis tool for experimentalists. Curr Protoc Mol Biol . 2010; Chapter 19: Unit 19.10.1– 21. Google Scholar CrossRef Search ADS   39. Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, Pimentel H, Salzberg SL, Rinn JL, Pachter L. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat Protoc . 2012; 7( 3): 562– 578. Google Scholar CrossRef Search ADS PubMed  40. Kroll KW, Mokaram NE, Pelletier AR, Frankhouser DE, Westphal MS, Stump PA, Stump CL, Bundschuh R, Blachly JS, Yan P. Quality control for RNA-Seq (QuaCRS): an integrated quality control pipeline. Cancer Inform . 2014; 13( suppl 3): 7– 14. Google Scholar PubMed  41. Gordon A, Hannon G. Fastx-toolkit. FASTQ/A short-reads preprocessing tools. http://hannonlab.cshl.edu/fastx_toolkit/. 42. Trapnell C, Pachter L, Salzberg SL. TopHat: discovering splice junctions with RNA-Seq. Bioinformatics . 2009; 25( 9): 1105– 1111. Google Scholar CrossRef Search ADS PubMed  43. Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, Salzberg SL, Wold BJ, Pachter L. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol . 2010; 28( 5): 511– 515. Google Scholar CrossRef Search ADS PubMed  44. Sato S, Yoshida W, Soejima H, Nakabayashi K, Hata K. Methylation dynamics of IG-DMR and Gtl2-DMR during murine embryonic and placental development. Genomics . 2011; 98( 2): 120– 127. Google Scholar CrossRef Search ADS PubMed  45. Hiura H, Komiyama J, Shirai M, Obata Y, Ogawa H, Kono T. DNA methylation imprints on the IG-DMR of the Dlk1–Gtl2 domain in mouse male germline. FEBS Lett . 2007; 581( 7): 1255– 1260. Google Scholar CrossRef Search ADS PubMed  46. Tierling S, Dalbert S, Schoppenhorst S, Tsai CE, Oliger S, Ferguson-Smith AC, Paulsen M, Walter J. High-resolution map and imprinting analysis of the Gtl2-Dnchc1 domain on mouse chromosome 12. Genomics . 2006; 87( 2): 225– 235. Google Scholar CrossRef Search ADS PubMed  47. Miyoshi N, Wagatsuma H, Wakana S, Shiroishi T, Nomura M, Aisaka K, Kohda T, Surani MA, Kaneko-Ishino T, Ishino F. Identification of an imprinted gene, Meg3/Gtl2 and its human homologue MEG3, first mapped on mouse distal chromosome 12 and human chromosome 14q. Genes Cells . 2000; 5( 3): 211– 220. Google Scholar CrossRef Search ADS PubMed  48. Paulsen M, Takada S, Youngson NA, Benchaib M, Charlier C, Segers K, Georges M, Ferguson-Smith AC. Comparative sequence analysis of the imprinted Dlk1–Gtl2 locus in three mammalian species reveals highly conserved genomic elements and refines comparison with the Igf2–H19 region. Genome Res . 2001; 11( 12): 2085– 2094. Google Scholar CrossRef Search ADS PubMed  49. Doshi T, D’souza C, Vanage G. Aberrant DNA methylation at Igf2–H19 imprinting control region in spermatozoa upon neonatal exposure to bisphenol A and its association with post implantation loss. Mol Biol Rep . 2013; 40( 8): 4747– 4757. Google Scholar CrossRef Search ADS PubMed  50. Chao HH, Zhang XF, Chen B, Pan B, Zhang LJ, Li L, Sun XF, Shi QH, Shen W. Bisphenol A exposure modifies methylation of imprinted genes in mouse oocytes via the estrogen receptor signaling pathway. Histochem Cell Biol . 2012; 137( 2): 249– 259. Google Scholar CrossRef Search ADS PubMed  51. Rissman EF, Adli M. Minireview: transgenerational epigenetic inheritance: focus on endocrine disrupting compounds. Endocrinology . 2014; 155( 8): 2770– 2780. Google Scholar CrossRef Search ADS PubMed  52. Kobayashi S, Wagatsuma H, Ono R, Ichikawa H, Yamazaki M, Tashiro H, Aisaka K, Miyoshi N, Kohda T, Ogura A, Ohki M, Kaneko-Ishino T, Ishino F. Mouse Peg9/Dlk1 and human PEG9/DLK1 are paternally expressed imprinted genes closely located to the maternally expressed imprinted genes: mouse Meg3/Gtl2 and human MEG3. Genes Cells . 2000; 5( 12): 1029– 1037. Google Scholar CrossRef Search ADS PubMed  53. You L, Wang N, Yin D, Wang L, Jin F, Zhu Y, Yuan Q, De W. Downregulation of long noncoding RNA Meg3 affects insulin synthesis and secretion in mouse pancreatic beta cells. J Cell Physiol . 2016; 231( 4): 852– 862. Google Scholar CrossRef Search ADS PubMed  54. Murphy SK, Huang Z, Hoyo C. Differentially methylated regions of imprinted genes in prenatal, perinatal and postnatal human tissues. PLoS One . 2012; 7( 7): e40924. Google Scholar CrossRef Search ADS PubMed  55. Qu C, Jiang T, Li Y, Wang X, Cao H, Xu H, Qu J, Chen JG. Gene expression and IG-DMR hypomethylation of maternally expressed gene 3 in developing corticospinal neurons. Gene Expr Patterns . 2013; 13( 1–2): 51– 56. Google Scholar CrossRef Search ADS PubMed  56. Zhou Y, Zhong Y, Wang Y, Zhang X, Batista DL, Gejman R, Ansell PJ, Zhao J, Weng C, Klibanski A. Activation of p53 by MEG3 non-coding RNA. J Biol Chem . 2007; 282( 34): 24731– 24742. Google Scholar CrossRef Search ADS PubMed  57. Soubry A, Guo L, Huang Z, Hoyo C, Romanus S, Price T, Murphy SK. Obesity-related DNA methylation at imprinted genes in human sperm: results from the TIEGER study. Clin Epigenetics . 2016; 8( 1): 51. Google Scholar CrossRef Search ADS PubMed  58. Laws SC, Carey SA, Ferrell JM, Bodman GJ, Cooper RL. Estrogenic activity of octylphenol, nonylphenol, bisphenol A and methoxychlor in rats. Toxicol Sci . 2000; 54( 1): 154– 167. Google Scholar CrossRef Search ADS PubMed  59. Armoskus C, Moreira D, Bollinger K, Jimenez O, Taniguchi S, Tsai HW. Identification of sexually dimorphic genes in the neonatal mouse cortex and hippocampus. Brain Res . 2014; 1562: 23– 38. Google Scholar CrossRef Search ADS PubMed  60. Mann F, Chauvet S, Rougon G. Semaphorins in development and adult brain: implication for neurological diseases. Prog Neurobiol . 2007; 82( 2): 57– 79. Google Scholar CrossRef Search ADS PubMed  61. Maier V, Jolicoeur C, Rayburn H, Takegahara N, Kumanogoh A, Kikutani H, Tessier-Lavigne M, Wurst W, Friedel RH. Semaphorin 4C and 4G are ligands of Plexin-B2 required in cerebellar development. Mol Cell Neurosci . 2011; 46( 2): 419– 431. Google Scholar CrossRef Search ADS PubMed  62. Ray TA, Heath KM, Hasan N, Noel JM, Samuels IS, Martemyanov KA, Peachey NS, McCall MA, Gregg RG. GPR179 is required for high sensitivity of the mGluR6 signaling cascade in depolarizing bipolar cells. J Neurosci . 2014; 34( 18): 6334– 6343. Google Scholar CrossRef Search ADS PubMed  63. Larminie C, Murdock P, Walhin JP, Duckworth M, Blumer KJ, Scheideler MA, Garnier M. Selective expression of regulators of G-protein signaling (RGS) in the human central nervous system. Brain Res Mol Brain Res . 2004; 122( 1): 24– 34. Google Scholar CrossRef Search ADS PubMed  64. Arambula SE, Belcher SM, Planchart A, Turner SD, Patisaul HB. Impact of low dose oral exposure to bisphenol A (BPA) on the neonatal rat hypothalamic and hippocampal transcriptome: a CLARITY-BPA Consortium study. Endocrinology . 2016; 157( 10): 3856– 3872. Google Scholar CrossRef Search ADS PubMed  65. Crews D, Gillette R, Scarpino SV, Manikkam M, Savenkova MI, Skinner MK. Epigenetic transgenerational inheritance of altered stress responses. Proc Natl Acad Sci USA . 2012; 109( 23): 9143– 9148. Google Scholar CrossRef Search ADS PubMed  66. Skinner MK, Savenkova MI, Zhang B, Gore AC, Crews D. Gene bionetworks involved in the epigenetic transgenerational inheritance of altered mate preference: environmental epigenetics and evolutionary biology. BMC Genomics . 2014; 15( 1): 377. Google Scholar CrossRef Search ADS PubMed  67. Gillette R, Miller-Crews I, Skinner MK, Crews D. Distinct actions of ancestral vinclozolin and juvenile stress on neural gene expression in the male rat. Front Genet . 2015; 6: 56. Google Scholar CrossRef Search ADS PubMed  68. Iqbal K, Tran DA, Li AX, Warden C, Bai AY, Singh P, Wu X, Pfeifer GP, Szabó PE. Deleterious effects of endocrine disruptors are corrected in the mammalian germline by epigenome reprogramming. Genome Biol . 2015; 16( 1): 59. Google Scholar CrossRef Search ADS PubMed  69. Carr MS, Yevtodiyenko A, Schmidt CL, Schmidt JV. Allele-specific histone modifications regulate expression of the Dlk1–Gtl2 imprinted domain. Genomics . 2007; 89( 2): 280– 290. Google Scholar CrossRef Search ADS PubMed  70. Modali SD, Parekh VI, Kebebew E, Agarwal SK. Epigenetic regulation of the lncRNA MEG3 and its target c-MET in pancreatic neuroendocrine tumors. Mol Endocrinol . 2015; 29( 2): 224– 237. Google Scholar CrossRef Search ADS PubMed  71. Das PP, Hendrix DA, Apostolou E, Buchner AH, Canver MC, Beyaz S, Ljuboja D, Kuintzle R, Kim W, Karnik R, Shao Z, Xie H, Xu J, De Los Angeles A, Zhang Y, Choe J, Jun DL, Shen X, Gregory RI, Daley GQ, Meissner A, Kellis M, Hochedlinger K, Kim J, Orkin SH. PRC2 is required to maintain expression of the maternal Gtl2-Rian-Mirg locus by preventing de novo DNA methylation in mouse embryonic stem cells. Cell Reports . 2015; 12( 9): 1456– 1470. Google Scholar CrossRef Search ADS PubMed  Copyright © 2018 Endocrine Society

Journal

EndocrinologyOxford University Press

Published: Jan 1, 2018

There are no references for this article.

You’re reading a free preview. Subscribe to read the entire article.


DeepDyve is your
personal research library

It’s your single place to instantly
discover and read the research
that matters to you.

Enjoy affordable access to
over 18 million articles from more than
15,000 peer-reviewed journals.

All for just $49/month

Explore the DeepDyve Library

Search

Query the DeepDyve database, plus search all of PubMed and Google Scholar seamlessly

Organize

Save any article or search result from DeepDyve, PubMed, and Google Scholar... all in one place.

Access

Get unlimited, online access to over 18 million full-text articles from more than 15,000 scientific journals.

Your journals are on DeepDyve

Read from thousands of the leading scholarly journals from SpringerNature, Elsevier, Wiley-Blackwell, Oxford University Press and more.

All the latest content is available, no embargo periods.

See the journals in your area

DeepDyve

Freelancer

DeepDyve

Pro

Price

FREE

$49/month
$360/year

Save searches from
Google Scholar,
PubMed

Create lists to
organize your research

Export lists, citations

Read DeepDyve articles

Abstract access only

Unlimited access to over
18 million full-text articles

Print

20 pages / month

PDF Discount

20% off