Weighted Gene Correlation Network Analysis (WGCNA) Reveals Novel Transcription Factors Associated With Bisphenol A Dose-Response

Despite Bisphenol-A (BPA) being subject to extensive study, a thorough understanding of molecular mechanism remains elusive. Here we show that using weighted gene correlation network analysis (WGCNA), which takes advantage of a graph theoretical approach to understanding correlations amongst genes and grouping genes into modules that typically have co-ordinated biological functions and regulatory mechanisms, that despite some commonality in altered genes, there is minimal overlap between BPA and estrogen in terms of network topology. We confirmed previous findings that ZNF217 and TFAP2C are involved in the estrogen pathway, and are implicated in BPA as well, although for BPA they appear to be active in the absence of canonical estrogen-receptor driven gene expression. Furthermore, our study suggested that PADI4 and RACK7/ZMYNDB8 may be involved in the overlap in gene expression between estradiol and BPA. Lastly, we demonstrated that even at low doses there are unique transcription factors that appear to be driving the biology of BPA, such as SREBF1. Overall, our data is consistent with other reports that BPA leads to subtle gene changes rather than profound aberrations of a conserved estrogen signaling (or other) pathways.

Introduction

Bisphenol A (BPA) is an industrial chemical used in the manufacture of polycarbonate plastic found in a number of consumer products such as thermal paper, canned foods and epoxy resins (Rubin, 2011) – although many of these uses are being phased out (Zimmerman and Anastas, 2015). Amongst the general population, exposure to BPA is widespread, with very low levels of BPA present in the majority of urinary samples taken in the general population (Calafat et al., 2008), although serum levels are estimated to be lower (Teeguarden et al., 2013). Release of BPA to the environment exceeds one million pounds per year (Rubin, 2011).

Bisphenol-A has been subjected to a high-level of scrutiny – “bisphenol A” returns over 11,500 abstracts in PubMed, with over 700 articles per year being published every year since 2013 (“Pubmed Bisphenol A, n.d.). Within HSDB (the Hazardous Substance Database), there are over 79 peer-reviewed animal studies (TOXNET, n.d.). The CLARITY study, a three generation chronic study with low-levels of BPA used to mimic population exposures, involved 3,500 rats: while it resulted in no revision of safety standards by the FDA, it still failed to bring about a consensus as to a safe level (Academics Urge Caution in Interpreting Clarity-Bpa Results, n.d.). Despite the overwhelming amount of data, the mechanism(s) by which BPA may exert adverse effects remains unclear, nor is there a widely agreed upon endpoint on which to base a safe dose.

Bisphenol-A was presumed to have potentially estrogenic activity, as well as potential carcinogenicity, based on its structural similarity to DES (Diethylstilbestrol) and other synthetic estrogens, as well as appearing to trigger gene expression similar to estrogen receptor agonists, despite its relatively low binding affinity for estrogen receptors (LaPensee et al., 2009). Two different hypotheses have been put forward to explain this discrepancy: one, BPA may bind to different domains of ESR1 or ESR2 and recruit different co-regulators (Safe et al., 2002), or alternatively, BPA may exert its effects through non-classical estrogen receptors, such as membrane-bound ER (GPR30) (LaPensee et al., 2009) or ERRγ (ERRG) (Okada et al., 2008), which is one of several “orphan” receptors that are classified as estrogen-related receptors (Horard and Vanacker, 2003). The question of BPA’s ultimate molecular initiating event is not academic – on the presumption that BPA’s effects are mediated via estrogen receptors, several alternatives were proposed, such as Bisphenol F and Bisphenol S, but both compounds have proven equally problematic (Rochester and Bolden, 2015).

In our previous work for the Mapping the Human Toxome project (Kleensang et al., 2014; Bouhifd et al., 2015), we demonstrated that using non-inferential statistical methods that did not depend on existing annotations such as IDEA (Pendse et al., 2017) and WGCNA (Maertens et al., 2015) offered a powerful method to untangle possible regulatory mechanisms and providing insight into possible Pathways of Toxicity compared to either inferential-based methods or approaches such as pathway enrichment analysis that depend exclusively on annotations. Building upon our previous work using WGCNA applied to in vitro transcriptomic data to more fully understand the transcription factors that are driving the biology of estradiol (Pendse et al., 2017), we used WGCNA to examine a previously published transcriptomic dataset (Shioda et al., 2013). Briefly, Shioda et al. (2013) aimed to study the sensitivities of estrogen responsive genes to various endocrine disrupting chemicals (EDCs) based on the transcriptomic profile of MCF-7 cells exposed to either estrogen or several xenoestrogens (including BPA) over a dose-response curve ranging from picomolar to micromolar concentrations for a 48 h time period. Based on their analysis, they found that a gene signature of “estrogen-responsive genes” allowed the estrogenic substances to be ranked in terms of potency. Additionally, the heat map of BPA-inducible genes demonstrated a weak transcriptional activation at very low BPA concentration as well as a strong peak at high concentration. However, BPA has differences as well as similarities to estrogen in terms of gene signatures: therefore, we sought to explore specifically the differences between estrogen and BPA as well as the differences between low-dose BPA and high-dose BPA for possible regulatory mechanisms.

Our analysis shows that while there is substantial overlap between genes altered by BPA and estrogen, which might imply that BPA is indeed “estrogenic,” there are important differences in network topology as well as biological function, and that the overlap appears to be driven by transcription factors such as ZNF217, TFAP2C, PADI4, and RACK7/ZMYND8 rather than the estrogen receptor per se. Furthermore, BPA (even at the lower end of the dose response curve - defined here as less than 12.5 μM) has pathways that are likely not mediated by estrogen receptors, but instead by other transcription factors, such as SREBF1. Moreover, our data is consistent with other reports that BPA leads to subtle, diffuse gene changes that are comparatively difficult to capture with inferential methods, and that low-dose BPA has distinct effects compared to higher doses.

Materials and Methods

Data

Dataset GSE50705, a comprehensive analysis of estrogen and xenoestrogen dose-response curves on MCF-7 cells after 48 h of exposure, was downloaded from GEO via GEOQuery (Davis and Meltzer, 2007) as normalized data and all analyses were performed with R/Bioconductor (Gentleman et al., 2004).

Weighted Gene Correlation Network Analysis

A WGCNA network (Langfelder and Horvath, 2007) was generated for several subsets of the data: Estrogen (n = 36), BPA (n = 44), and low dose BPA (n = 32), as well as a consensus network for estrogen and BPA together (n = 80) using the 10,000 most highly expressed genes for each subset of the data as determined by rank means expression, the approach in Maertens et al. (2015); consensus networks and module statistics followed overall the approach in Langfelder et al. (2008). Briefly, the network was derived based on a signed Spearman correlation using a β of 10 as a weight function. The topological overlap metric (TOM) (Yip and Horvath, 2007) was derived from the resulting adjacency matrix, and was used to cluster the modules using the blockwiseModules function (blockwiseConsensusModules, for the consensus modules) and the dynamic tree cut algorithm (Langfelder et al., 2008) with a height of 0.25 and a deep split level of 2, a reassign threshold of 0.2 and a minimum module size of 30 (100 for the consensus network). The eigenmodules— essentially the first principal component of the modules, which can be used as a “signature” of the modules gene expression —were then correlated with dose, and each module that was correlated with the dose-response curve with a p-value < 0.01 (p-value < 0.05 for the consensus network) was considered statistically significant.

Transcription Factor Analysis

All statistically significant modules were analyzed in EnrichR (Chen et al., 2013) using the CHEA dataset (Lachmann et al., 2010) restricted to MCF-7/10 cells –as well as the ARCHS4 TF-Coexpression dataset with an adjusted p-value less than 0.01 based on Fisher’s exact test.

Functional Annotation Analysis

Module “hubs” were defined as having high-ranking kME (which ranks the connectivity of genes) within the module and predicted as high-degree within the STRING database (Szklarczyk et al., 2017) of protein–protein interactions. The three highest ranking modules were analyzed in STRING for the enrichment of predicted protein interactions as well as functional annotation via GO Biological Process and Molecular Function. All analysis with STRING was done with medium stringency settings, and included all possible interactions (text-mining, database, experiments, co-expression, neighborhood, gene fusion, and co-occurrence).

TCGA Data

Expression and methylation data for FIZ1 was correlated with clinical attributes (estrogen receptor, progesterone receptor, and solid tumor vs. normal vs. metastatic tumor) using MEExpress (Koch et al., 2015) based on the TCGA BRCA dataset.

Results

Consensus Network Analysis Indicates Minimal Overlap Between Estrogen and BPA

We began by analyzing the dose-response curve of the BPA and estrogen dataset combined using WGCNA (which takes advantage of correlations amongst genes and groups genes into modules using network topology) to look for a “consensus network”-a common pattern of genes that are correlated in all conditions. The consensus network identified (Figure 1) had clearly delineated modules, and the modules identified were significantly correlated with both estrogen (Figure 2) and BPA (Figure 3). However, estrogen clearly had a stronger signal in comparison to BPA and quite possibly overwhelmed the signal from BPA. More strikingly, however, the majority of modules in the consensus analysis when analyzed for correlation with both BPA and estrogen showed virtually no similarity – most modules had opposite directions of correlation, and of the few modules with similar correlations, the coefficient of correlation was very weak – only one module (“yellow”) was significant with a p-value < 0.05 (Table 1), indicating that while there may be some overlap in genetic signatures, from a network topology perspective there is minimal conservation. When analyzed for transcription factors against the CHEA dataset – a collection of ChIP-chip, ChIP-seq, ChIP-PET, and DamID studies collected into a database to infer transcriptional regulation (Lachmann et al., 2010) – the common module was enriched for E2F1, ZNF217, and RACK7, but not ESR1 or ESR2 (Table 2).

FIGURE 1

FIGURE 1. Consensus network from BPA and Estrogen dose-response curve. Gene expression similarity is determined using a pair-wise weighted correlation metric, and clustered according to a topological overlap metric into modules; assigned modules are colored on bottom, gray genes are unassigned to a module.

FIGURE 2

FIGURE 2. Consensus network modules correlated with estrogen dose using the eigenmodule (the first principal component of the module). Correlation coefficient along with p-value in parenthesis underneath; color-coded according to correlation coefficient (legend at right).

FIGURE 3

FIGURE 3. Consensus network modules correlated with BPA dose using the eigenmodule (the first principal component of the module). Correlation coefficient along with p-value in parenthesis underneath; color-coded according to correlation coefficient (legend at right).

TABLE 1

TABLE 1. Consensus network modules associated with BPA and estrogen.

TABLE 2

TABLE 2. Enriched transcription factors in conserved module in consensus network.

Estrogen and BPA Network Overlap With Transcription Factors, Including ESR1 and ESR2, but With Different Network Topologies and Different Biological Processes

Next, we derived the de novo networks individually for the entire dose-response curve of estrogen and BPA to examine the network topology, common hubs, and biological role of the modules in each network separately. Within the estrogen network (Supplementary Figure 1A), there was one large module highly correlated with estrogen dose (“turquoise”) (Table 3), which was also enriched for ESR1 and ESR2 in addition to E2F1, ZNF217, TFAP2C amongst others (Table 4); furthermore, ESR1 was a hub within the module, and the module was predominately enriched with terms related to cell-cycle as well as poly(A) RNA binding (Supplementary Table 1).

TABLE 3

TABLE 3. Estrogen modules correlated with dose.

TABLE 4

TABLE 4. Enriched transcription factors in estrogen modules.

In comparison, the top module correlated with BPA dose (“lightcyan1”) (Supplementary Figure 1B) (Table 5), was a relatively small module not enriched for any transcription factors, however, SSR1 (Signal Sequence Receptor Subunit 1) was the top hub; annotation analysis revealed that the module was enriched for genes in the GO category of “response to endoplasmic reticulum stress” and “endoplasmic reticulum unfolded protein response” (p-value of 2.93E-17 and 2.02E-12, respectively) (Supplementary Table 1). The second module correlated with dose (“royal blue”) was enriched for ESR1 and ESR2 genes (Table 6), however, neither ESR1 nor ESR2 was present in the module (or any module correlated positively or negatively with dose) and instead the main hub was TOP2A (Topoisomerase IIA). The module had a weak over-representation of genes involved in development (p-value 0.00416) and cytoskeleton organization (p-value 0.0155) (Supplementary Table 1). The other module enriched for ESR1 and ESR2 genes (“dark gray”) was also annotated to “response to unfolded protein” (p-value of 9.52E-05) (Supplementary Table 1).

TABLE 5

TABLE 5. BPA modules associated with dose.

TABLE 6

TABLE 6. Enriched transcription factors in BPA modules.

Low-Dose BPA Network Shows No Enrichment of ESR1 or ESR2 Genes

It has been speculated that BPA at low doses has fundamentally different effects than at high doses; in the original study of the dataset, the authors detected a weak, but distinct, transcriptional activity peak at low doses. Therefore, we restricted the BPA network to doses below 12.5 μM (leaving a highest dose of 6.25 μM, and most of the dose-response curve in the nanomolar/picomolar range) and calculated a network specific for this lower dose range. Despite the smaller sample size, the network still produced several modules that were significantly correlated with dose (Supplementary Figure 1C and Table 7). This low-dose BPA network shows consistent transcription factors (ZNF217, TFAP2C, RACK7/ZMYND8, and PADI4) with the larger BPA network as well as the estrogen network, but no modules were enriched for genes with ESR1 or ESR2 with a p-value cut-off of < 0.01 (Table 8).

TABLE 7

TABLE 7. Low-dose BPA modules associated with dose.

TABLE 8

TABLE 8. Enriched transcription factors in low-dose BPA modules.

The module with the highest correlation with dose (“turquoise”) was comparatively dense for predicted protein-protein interactions (p-value of < 1.0e-16, average node degree 10.3) as well as genes related to cellular macromolecule metabolic process and poly(A) RNA binding (p-value of 1.58E-16 and 2.17eE-18, respectively) (Supplementary Table 1) in contradistinction to the estrogen network, where the dominant module was enriched overwhelmingly with cell-cycle genes. Moreover, the module included both ZNF217 and TFAP2C, but neither ESR1 nor ESR2 were in the module, much less hubs. The second module correlated with dose (“Dark Green”) showed no enrichment for transcription factors, although it was enriched for protein–protein interactions and the molecular function “enzyme binding”; the third module (“dark red”), also strongly correlated with dose, showed no enrichment for transcription factors or protein–protein interactions, though it was weakly enriched for the KEGG pathway Insulin Signaling (p-value 0.0028); this module may simply be an artifact, reflect diffuse alterations that are difficult to detect, or an unknown regulatory mechanism. An additional fairly large module correlated with dose (“brown”) was enriched for both E2F1 and PADI4, strongly enriched for protein-protein interactions (p-value < 1.0E-16, average node degree 6.56) and as well as cellular metabolic process (p-value 1.32E-10) and poly(A) RNA binding (p-value 2.5E-16) (Supplementary Table 1); but, also in contrast to the estrogen network, was not enriched for cell-cycle genes.

ZNF217 has previously been shown by our work (Pendse et al., 2017) and others (Frietze et al., 2014) to be a critical component of estrogen signaling and an important prognostic factor for breast cancer (Vendrell et al., 2012). Similarly, TFAP2C is known to modulate ESR1 and GPR30 expression, and attenuate the expression of several estrogen-targeted genes (Woodfield et al., 2007). Given the presence of both ZNF217 and TFAP2C in the network as well as the strong enrichment genes targeted by these transcription factors, this suggests that these genes are indeed central to mediating BPAs phenotypic effects; however, our study shows little evidence that they are in exerting their effect in tandem with ESR1 or ESR2.

Moreover, both ZNF217 and TFAP2C were shown independently to be altered by bisphenol A in a rat seminiferous tubule culture model (Ali et al., 2014). The same study also showed alterations (albeit subtle) in PADI4 and RACK7 (official gene symbol: ZYMND8) mRNA as well as RACK7/ZMYND8 methylation; neither of these genes were in our dataset so their role in observed changes remains speculative. However, RACK7/ZYMND8 binds a large set of active enhancers, including almost all “super-enhancers,” and is therefore expected to have sweeping transcriptional effects (Shen et al., 2016). Although little is known about its role in breast cancer, it is thought to inhibit HIF-dependent breast-cancer progression (Chen et al., 2018). PADI4 is known to be implicated in cancer and is thought to respond to estrogen-simulation in MCF-7 cells through both genomic and non-genomic mechanisms (Dong et al., 2007). In breast cancer specifically, it is implicated in the ELK1/C-Fos pathway (Zhang et al., 2011). Moreover, BPA was shown to increase protein levels of PADI4 via a reactive oxygen species mechanism in neuroblastoma cells (Park et al., 2012).

Low-Dose BPA Network Had Unique Transcription Factors Not Present in the Estrogen Dataset

To further delineate possible transcription factors unique to BPA signaling compared to estrogen, we examined the list of genes in all modules statistically significantly associated with the low-dose BPA network that were not present in the estrogen network, a total of 1,901 genes. Analyzed against the CHEA dataset, the genes were again enriched for RACK7/ZMYND8, in addition to ELK1 and HIF1A (Supplementary Table 2). In order to expand our search for transcription factors that may not have been studied in MCF-7 cells in the CHEA data set, we also analyzed the list of genes for enrichment against the ARCHS4 database (Lachmann et al., 2018), which correlates transcription factor expression against gene expression in a combined database of over 20,000 RNASeq samples. Of the top 50 transcription factors identified as significantly correlated with the gene list, 18 were also present in the low-dose BPA network (Table 9). The highest-ranking transcription factor, FIZ1, is zinc-finger protein with a largely unknown biological role (Wolf and Rohrschneider, 1999) - it has a relatively poor literature base, with only 8 citations in PubMed. However, FIZ1 expression in breast cancer is statistically associated with progesterone receptor status, estrogen receptor status, and sample subtype, and it undergoes extensive CpG-island methylation (Supplementary Figure 2), and it is therefore an intriguing candidate for further study. The second highest-ranking transcription factor, SREBF1 is comparatively better characterized: it is known to be central to lipid homeostasis, regulating the LDL receptor gene as well as related fatty acid and cholesterol synthesis genes. Furthermore, SREBF1 mRNA was identified as upregulated in adipocytes by BPA (Boucher et al., 2014). Neither SREBF1 nor FIZ1 were present in the estrogen dataset, and SREBF1- and FIZ1-correlated genes were not enriched in the subset of estrogen-only genes. It is therefore plausible that these two transcription factors are more central to BPAs effects than estrogen, however, because enrichment for transcription factors motifs/regulated genes in any gene list often produce false-positives, understanding their role would require further study.

TABLE 9

TABLE 9. Transcription factors unique to low-dose BPA network.

Discussion

Estrogen signaling is unique amongst nuclear receptors in that substantial number of the genes altered by estrogen do not have canonical estrogen response elements (Miller et al., 2017) – estrogen signaling takes place within a transcriptomic and epigenomic context that markedly influences receptor activation. Our examination of the estrogen dose response curve network both confirmed several of the transcription factors identified previously, such as E2F1, ZNF217 and TFAP2C, as well as suggested other transcriptional factors such as PADI4 and RACK7/ZYMND8 that may impact estrogen signaling.

Our study is consistent with other findings that the assumption that BPA works exclusively or even predominantly on canonical ESR1 or ESR2 gene regulation may be misleading or an oversimplification (Delfosse et al., 2012; MacKay and Abizaid, 2018). To be sure, one can find gene patterns similar to those found in estrogen-induced cells, but the leap from that observation to the presumption that such changes are estrogen-mediated may not be warranted. While this study cannot determine conclusively the ultimate chain of events that leads from the molecular initiating event to the phenotypic consequences, it does suggest some hypotheses that are more probable. The lack of overlap in the consensus network indicates that despite similarity of genes, there is minimal conservation of network topology, and the one conserved module was not enriched for ESR1 or ESR2 genes. In networks drawn separately from dose-response curves for estrogen and BPA, the substantial differences in network topology, the absence of ESR1 as a hub gene in the BPA network, and the differences in biological function of the modules suggest that even at high-doses, BPAs effects are fundamentally different than estradiol. The lack of estrogen receptor target genes in the low dose BPA network in the presence of a clear signature of other transcription factors suggests that at low doses BPA’s effects are driven by mechanisms other than direct estrogen receptor activation. Additionally, regardless of molecular initiating event, assessing BPAs dose-response by looking at estrogen gene-signatures may miss interesting and important biology, such as the likely role of SREBF1. Furthermore, our study is consistent with other findings that BPA’s effects are subtle and phenotypic changes likely reflect modest effects at multiple different points (Porreca et al., 2016) and that analyzing the effects of low-dose BPA can reveal effects that are obscured at higher doses (Shioda et al., 2013). This does not necessarily lead to a “non-monotonic” dose response curve– this could be due to technical reasons, or higher-doses could cause non-specific changes that are the result of cellular stress, as evidenced by the identification of modules associated with unfolded protein response. It does, however, point to a need to consider the doses chosen for an in vitro study carefully and to not presume linear effects.

This study is certainly not a definitive study of BPA molecular mechanisms: our conclusions cannot confidently be extrapolated to other tissue types, as BPA may have tissue specific effects; MCF-7 cells are prone to artifacts (Kleensang et al., 2016); and our study did not focus on epigenetic mechanisms which are speculated as significantly underpinning much of the observed adverse events seen with BPA exposure, especially at a low dose (Singh and Li, 2012). While using the CHEA dataset and restricting candidate transcription factors to those observed in MCF-7 cells eliminates many of the false-positives intrinsic to such approaches, it also limits findings to those transcription factors that have been studied, and this may miss some important biology. Extending our analysis with the ARCHS4 database added interesting candidates, but all correlation-based approaches must be treated with caution and viewed as “hypothesis-generating,” and all exploratory data analysis techniques such as WGCNA require further targeted studies to confirm suggested molecular networks.

Nonetheless, our study does indicate that transcriptomics, especially given a high-dimensional dataset and the use of non-inferential methods, can likely aid toxicologists in having a better understanding of probable molecular targets as well as the complexity of perturbed networks - clearly, understanding BPAs effects will require a systems level approach (Hartung et al., 2017) as well as better characterization of genes that are not as yet confidently mapped as to biological function. More generally speaking, this points to the pitfall of trying to design “greener” substitutes (Maertens et al., 2014; Maertens and Hartung, 2018) in the absence of a clear, comprehensive understanding of molecular mechanism.

Author Contributions

AM: main author of the paper. VT: support programming and data analysis. AK: planned the work and revised the manuscript. TH: mentoring, revision of the manuscript, and PI Human Toxome project laying the conceptual ground.

Funding

This work was supported by an NIH Transformational Research Grant, “Mapping the Human Toxome by Systems Toxicology” (RO1 ES 020750). VT was supported by NIEHS training grant (T32 ES007141). This work was also supported by the EU-ToxRisk project (An Integrated European “Flagship” Program Driving Mechanism-Based Toxicity Testing and Risk Assessment for the 21st Century) funded by the European Commission under the Horizon 2020 program (Grant Agreement No. 681002).

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2018.00508/full#supplementary-material

FIGURE S1 | (A) Estrogen and (B) BPA Network dendrogram.

FIGURE S2 | Methylation pattern for FIZ1 in BRCA TCGA data set.

TABLE S1 | Protein-protein interaction and enrichment from STRING database.

TABLE S2 | Transcription factors for genes present only in low-dose BPA network.

References

Academics Urge Caution in Interpreting Clarity-Bpa Results (n.d.). Chemical Watch. Available at: https://chemicalwatch.com/64449/academics-urge-caution-in-interpreting-clarity-bpa-results [accessed July 13, 2018].

Google Scholar

Ali, S., Steinmetz, G., Montillet, G., Perrard, M., Loundou, A., Durand, P., et al. (2014). Exposure to low-dose bisphenol A impairs meiosis in the rat seminiferous tubule culture model: a physiotoxicogenomic approach. PLoS One 9:e106245. doi: 10.1371/journal.pone.0106245

PubMed Abstract | CrossRef Full Text | Google Scholar

Boucher, J. G., Husain, M., Rowan-Carroll, A., Williams, A., Yauk, C. L., and Atlas, E. (2014). Identification of mechanisms of action of bisphenol a-induced human preadipocyte differentiation by transcriptional profiling. Obesity 22, 2333–2343. doi: 10.1002/oby.20848

PubMed Abstract | CrossRef Full Text | Google Scholar

Bouhifd, M., Andersen, M. E., Baghdikian, C., Boekelheide, K., Crofton, K. M., Fornace, A. J., et al. (2015). The human toxome project. ALTEX 32:112. doi: 10.14573/altex.1502091

PubMed Abstract | CrossRef Full Text | Google Scholar

Calafat, A. M., Ye, X., Wong, L., Reidy, J. A., and Needham, L. L. (2008). Exposure of the U.S. population to bisphenol A and 4-tertiary-octylphenol: 2003-2004. Environ. Health Perspect. 116, 39–44. doi: 10.1289/ehp.10753

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, E. Y., Tan, C. M., Kou, Y., Duan, Q., Wang, Z., Meirelles, G. V., et al. (2013). Enrichr: interactive and collaborative HTML5 gene list enrichment analysis tool. BMC Bioinformatics 14:128. doi: 10.1186/1471-2105-14-128

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, Y., Zhang, B., Bao, L., Jin, L., Yang, M., Peng, Y., et al. (2018). ZMYND8 acetylation mediates HIF-dependent breast cancer progression and metastasis. J. Clin. Invest. 128, 1937–1955. doi: 10.1172/JCI95089

PubMed Abstract | CrossRef Full Text | Google Scholar

Davis, S., and Meltzer, P. S. (2007). GEOquery: a bridge between the Gene Expression Omnibus (GEO) and BioConductor. Bioinformatics 23, 1846–1847. doi: 10.1093/bioinformatics/btm254

PubMed Abstract | CrossRef Full Text | Google Scholar

Delfosse, V., Grimaldi, M., Pons, J., Boulahtouf, A., le Maire, A., Cavailles, V., et al. (2012). Structural and mechanistic insights into bisphenols action provide guidelines for risk assessment and discovery of bisphenol A substitutes. Proc. Natl. Acad. Sci. U.S.A. 109, 14930–14935. doi: 10.1073/pnas.1203574109

PubMed Abstract | CrossRef Full Text | Google Scholar

Dong, S., Zhang, Z., and Takahara, H. (2007). Estrogen-Enhanced Peptidylarginine Deiminase Type IV Gene (PADI4) Expression in MCF-7 Cells Is Mediated by Estrogen Receptor-$\alpha$-Promoted Transfactors Activator Protein-1. Nuclear Factor-Y and Sp1. Mol. Endocrinol. 21, 1617–1629. doi: 10.1210/me.2006-0550

PubMed Abstract | CrossRef Full Text | Google Scholar

Frietze, S., O’Geen, H., Littlepage, L. E., Simion, C., Sweeney, C. A., Farnham, P. J., et al. (2014). Global analysis of ZNF217 chromatin occupancy in the breast cancer cell genome reveals an association with ERalpha. BMC Genomics 15:520. doi: 10.1186/1471-2164-15-520

PubMed Abstract | CrossRef Full Text | Google Scholar

Gentleman, R. C., Carey, V. J., Bates, D. M., Bolstad, B., Dettling, M., Dudoit, S., et al. (2004). Bioconductor: open software development for computational biology and bioinformatics. Genome Biol. 5:R80. doi: 10.1186/gb-2004-5-10-r80

PubMed Abstract | CrossRef Full Text | Google Scholar

Hartung, T., FitzGerald, R. E., Jennings, P., Mirams, G. R., Peitsch, M. C., Rostami-Hodjegan, A., et al. (2017). Systems toxicology: real world applications and opportunities. Chem. Res. Toxicol. 30, 870–882. doi: 10.1021/acs.chemrestox.7b00003

PubMed Abstract | CrossRef Full Text | Google Scholar

Horard, B., and Vanacker, J. (2003). Estrogen receptor-related receptors: orphan receptors desperately seeking a ligand. J. Mol. Endocrinol. 31, 349–357. doi: 10.1677/jme.0.0310349

PubMed Abstract | CrossRef Full Text | Google Scholar

Kleensang, A., Maertens, A., Rosenberg, M., Fitzpatrick, S., Lamb, J., Auerbach, S., et al. (2014). t4 workshop report: Pathways of Toxicity. ALTEX 31, 53–61. doi: 10.14573/altex.1309261

PubMed Abstract | CrossRef Full Text | Google Scholar

Kleensang, A., Vantangoli, M. M., Odwin-DaCosta, S., Andersen, M. E., Boekelheide, K., Bouhifd, M., et al. (2016). Genetic variability in a frozen batch of MCF-7 cells invisible in routine authentication affecting cell function. Sci. Rep. 6:28994. doi: 10.1038/srep28994

PubMed Abstract | CrossRef Full Text | Google Scholar

Koch, A., De Meyer, T., Jeschke, J., and Van Criekinge, W. (2015). MEXPRESS: visualizing expression, DNA methylation and clinical TCGA data. BMC Genomics 16:636. doi: 10.1186/s12864-015-1847-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Lachmann, A., Torre, D., Keenan, A. B., Jagodnik, K. M., Lee, H. J., Wang, L., et al. (2018). Massive mining of publicly available RNA-seq data from human and mouse. Nat.Commun. 9, 1366. doi: 10.1038/s41467-018-03751-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Lachmann, A., Xu, H., Krishnan, J., Berger, S. I., Mazloom, A. R., and Ma’ayan, A. (2010). ChEA: transcription factor regulation inferred from integrating genome-wide ChIP-X experiments. Bioinformatics 26, 2438–2444. doi: 10.1093/bioinformatics/btq466

PubMed Abstract | CrossRef Full Text | Google Scholar

Langfelder, P., and Horvath, S. (2007). Eigengene networks for studying the relationships between co-expression modules. BMC Syst. Biol. 1:54. doi: 10.1186/1752-0509-1-54

PubMed Abstract | CrossRef Full Text | Google Scholar

Langfelder, P., Zhang, B., and Horvath, S. (2008). Defining clusters from a hierarchical cluster tree: the Dynamic Tree Cut package for R. Bioinformatics 24, 719–720. doi: 10.1093/bioinformatics/btm563

PubMed Abstract | CrossRef Full Text | Google Scholar

LaPensee, E. W., Tuttle, T. R., Fox, S. R., and Ben-Jonathan, N. (2009). } % The entry below contains non-ASCII chars that could not be converted % to a LaTeX equivalent.). Bisphenol A at Low Nanomolar Doses Confers Chemoresistance in Estrogen Receptor-$\alpha$–Positive and –Negative Breast Cancer Cells. Environ. Health Perspect. 117, 175–180. doi: 10.1289/ehp.11788

PubMed Abstract | CrossRef Full Text | Google Scholar

MacKay, H., and Abizaid, A. (2018). A plurality of molecular targets: The receptor ecosystem for bisphenol-A (BPA). Horm. Behav. 101, 59–67. doi: 10.1016/j.yhbeh.2017.11.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Maertens, A., Anastas, N., Spencer, P. J., Stephens, M., Goldberg, A., and Hartung, T. (2014). Green toxicology. ALTEX 31, 243–249. doi: 10.14573/altex.1406181

PubMed Abstract | CrossRef Full Text | Google Scholar

Maertens, A., and Hartung, T. (2018). Green Toxicology-know early about and avoid toxic product liabilities. Toxicol. Sci. 161, 285–289. doi: 10.1093/toxsci/kfx243

PubMed Abstract | CrossRef Full Text

Maertens, A., Luechtefeld, T., Kleensang, A., and Hartung, T. (2015). MPTP’s pathway of toxicity indicates central role of transcription factor SP1. Arch. Toxicol. 89, 743–755. doi: 10.1007/s00204-015-1509-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Miller, M. M., McMullen, P. D., Andersen, M. E., and Clewell, R. A. (2017). Multiple receptors shape the estrogen response pathway and are critical considerations for the future of in vitro-based risk assessment efforts. Crit. Rev. Toxicol. 47, 564–580. doi: 10.1080/10408444.2017.1289150

PubMed Abstract | CrossRef Full Text | Google Scholar

Okada, H., Tokunaga, T., Liu, X., Takayanagi, S., Matsushima, A., and Shimohigashi, Y. (2008). Direct evidence revealing structural elements essential for the high binding ability of bisphenol A to human estrogen-related receptor-gamma. Environ. Health Perspect. 116, 32–38. doi: 10.1289/ehp.10587

PubMed Abstract | CrossRef Full Text | Google Scholar

Park, B., Rhee, D., and Pyo, S. (2012). Apoptotic mechanism of bisphenol a in human neuroblastoma. FASEB J. 26.

Google Scholar

Pendse, S. N., Maertens, A., Rosenberg, M., Roy, D., Fasani, R. A., Vantangoli, M. M., et al. (2017). Information-dependent enrichment analysis reveals time-dependent transcriptional regulation of the estrogen pathway of toxicity. Arch. Toxicol. 91, 1749–1762. doi: 10.1007/s00204-016-1824-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Porreca, I., Ulloa Severino, L., D’Angelo, F., Cuomo, D., Ceccarelli, M., Altucci, L., et al. (2016). “Stockpile” of Slight Transcriptomic Changes Determines the Indirect Genotoxicity of Low-Dose BPA in Thyroid Cells. PLoS One 11:e0151618. doi: 10.1371/journal.pone.0151618

PubMed Abstract | CrossRef Full Text | Google Scholar

Pubmed Bisphenol A (n.d.). PubMED. Available at: https://www.ncbi.nlm.nih.gov/pubmed/?term = Bisphenol + A [accessed January 7, 2018].

Google Scholar

Rochester, J. R., and Bolden, A. L. (2015). Bisphenol and F: a systematic review and comparison of the hormonal activity of bisphenol a substitutes. Environ. Health Perspect. 123, 643–650. doi: 10.1289/ehp.1408989

PubMed Abstract | CrossRef Full Text | Google Scholar

Rubin, B. S. (2011). Bisphenol A: an endocrine disruptor with widespread exposure and multiple effects. J. Steroid Biochem. Mol. Biol. 127, 27–34. doi: 10.1016/j.jsbmb.2011.05.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Safe, S. H., Pallaroni, L., Yoon, K., Gaido, K., Ross, S., and McDonnell, D. (2002). Problems for risk assessment of endocrine-active estrogenic compounds. Environ. Health Perspect. 110(Suppl. 6), 925–929. doi: 10.1289/ehp.02110s6925

PubMed Abstract | CrossRef Full Text | Google Scholar

Shen, H., Xu, W., Guo, R., Rong, B., Gu, L., Wang, Z., et al. (2016). Suppression of enhancer overactivation by a RACK7-histone demethylase complex. Cell 165, 331–342. doi: 10.1016/j.cell.2016.02.064

PubMed Abstract | CrossRef Full Text | Google Scholar

Shioda, T., Rosenthal, N. F., Coser, K. R., Suto, M., Phatak, M., Medvedovic, M., et al. (2013). Expressomal approach for comprehensive analysis and visualization of ligand sensitivities of xenoestrogen responsive genes. Proc. Natl. Acad. Sci. U.S.A. 110, 16508–16513. doi: 10.1073/pnas.1315929110

PubMed Abstract | CrossRef Full Text | Google Scholar

Singh, S., and Li, S. S. (2012). Epigenetic effects of environmental chemicals bisphenol A and phthalates. Int. J. Mol. Sci. 13, 10143–10153. doi: 10.3390/ijms130810143

PubMed Abstract | CrossRef Full Text | Google Scholar

Szklarczyk, D., Morris, J. H., Cook, H., Kuhn, M., Wyder, S., Simonovic, M., et al. (2017). The STRING database in 2017: quality-controlled protein–protein association networks, made broadly accessible. Nucleic Acids Res. 45, D362–D368. doi: 10.1093/nar/gkw937

PubMed Abstract | CrossRef Full Text | Google Scholar

Teeguarden, J., Hanson-Drury, S., Fisher, J. W., and Doerge, D. R. (2013). Are typical human serum BPA concentrations measurable and sufficient to be estrogenic in the general population? Food Chem. Toxicol. 62, 949–963. doi: 10.1016/j.fct.2013.08.001

PubMed Abstract | CrossRef Full Text | Google Scholar

TOXNET (n.d.). TOXNET Databases. Available at: https://toxnet.nlm.nih.gov/cgi-bin/sis/search2/r?dbs + hsdb:@term + @DOCNO + 513[accessed July 11, 2018].

Google Scholar

Vendrell, J. A., Thollet, A., Nguyen, N. T., Ghayad, S. E., Vinot, S., Bi\’eche, I., et al. (2012). ZNF217 is a marker of poor prognosis in breast cancer that drives epithelial-mesenchyme transition and invasion. Cancer Res. 72, 3593–3606. doi: 10.1158/0008-5472.CAN-11-3095

PubMed Abstract | CrossRef Full Text | Google Scholar

Wolf, I., and Rohrschneider, L. R. (1999). Fiz1, a novel zinc finger protein interacting with the receptor tyrosine kinase Flt3. J. Biol. Chem. 274, 21478–21484. doi: 10.1074/jbc.274.30.21478

PubMed Abstract | CrossRef Full Text | Google Scholar

Woodfield, G. W., Horan, A. D., Chen, Y., and Weigel, R. J. (2007). TFAP2C controls hormone response in breast cancer cells through multiple pathways of estrogen signaling. Cancer Res. 67, 8439–8443. doi: 10.1158/0008-5472.CAN-07-2293

PubMed Abstract | CrossRef Full Text | Google Scholar

Yip, A. M., and Horvath, S. (2007). Gene network interconnectedness and the generalized topological overlap measure. BMC Bioinformatics 8:22. doi: 10.1186/1471-2105-8-22

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, X., Gamble, M. J., Stadler, S., Cherrington, B. D., Causey, C. P., Thompson, P. R., et al. (2011). Genome-wide analysis reveals PADI4 cooperates with Elk-1 to activate c-Fos expression in breast cancer cells. PLoS Genet. 7:e1002112. doi: 10.1371/journal.pgen.1002112

PubMed Abstract | CrossRef Full Text | Google Scholar

Zimmerman, J. B., and Anastas, P. T. (2015). Chemistry. Toward substitution with no regrets. Science 347, 1198–1199. doi: 10.1126/science.aaa0812

PubMed Abstract | CrossRef Full Text | Google Scholar