The closely related transcription factors (TFs), estrogen receptors ERα and ERβ, regulate divergent gene expression programs and proliferative outcomes in breast cancer. Utilizing breast cancer cells with ERα, ERβ, or both receptors as a model system to define the basis for differing response specification by related TFs, we show that these TFs and their key coregulators, SRC3 and RIP140, generate overlapping as well as unique chromatin‐binding and transcription‐regulating modules. Cistrome and transcriptome analyses and the use of clustering algorithms delineated 11 clusters representing different chromatin‐bound receptor and coregulator assemblies that could be functionally associated through enrichment analysis with distinct patterns of gene regulation and preferential coregulator usage, RIP140 with ERβ and SRC3 with ERα. The receptors modified each other's transcriptional effect, and ERβ countered the proliferative drive of ERα through several novel mechanisms associated with specific binding‐site clusters. Our findings delineate distinct TF‐coregulator assemblies that function as control nodes, specifying precise patterns of gene regulation, proliferation, and metabolism, as exemplified by two of the most important nuclear hormone receptors in human breast cancer.
To define how the estrogen receptors α and β control specific responses in breast cancer cells, genome‐wide patterns of chromatin binding of the ERα and ERβ receptors and their coregulators, SRC3 and RIP140, were determined and integrated with gene expression data and functional analyses.
The closely related transcription factors, estrogen receptors ERα and ERβ, can elicit differential cellular responses.
To understand the basis of this specificity, chromatin binding of ERs and key coregulators, and gene expression, were analyzed genome wide in human breast cancer cells containing ERα only, ERα+ERβ, and ERβ only.
A clustering‐based combinatorial analysis of ChIP‐Seq and gene expression data was used to parse genes into groups, specifying their mode of functional regulation in a particular cell background.
Through this analysis, RIP140 was identified as an ERβ‐preferential cofactor regulating cell proliferation, apoptosis, and adipogenesis programs.
A 20‐gene ERβ and RIP140 signature was developed, which predicted outcome and disease‐free survival in breast cancer patients.
Members of the superfamily of nuclear hormone receptors (NHRs) have pivotal roles in almost all aspects of physiology, including reproduction, metabolism, and immune function (Bookout et al, 2006). The estrogen receptor‐α (ERα) and ‐β (ERβ) belong to this superfamily of transcription factors (TFs) and mediate the diverse actions of estrogens, including their involvement in several disease states, such as breast cancer (Katzenellenbogen and Katzenellenbogen, 2000; Nilsson et al, 2001; Deroo and Korach, 2006). Although encoded by different genes, ERα and ERβ share the same general modular protein structure common to the NHR superfamily, and have nearly 100% amino acid sequence homology in their DNA‐binding domain. As such, it has been previously shown that ERα and ERβ are able to bind to similar DNA elements. They also show 59% amino acid homology in their ligand‐binding domains, and can heterodimerize as well as act as homodimers (Cowley et al, 1997; Pace et al, 1997).
It has been well documented that ERα controls the transcription of hundreds of genes (Frasor et al, 2003; Hah et al, 2011), and that its presence in breast cancer cells drives enhanced proliferation in response to estrogens (Frasor et al, 2003; Chang et al, 2008). Therefore, ERα is the main therapeutic target in ER‐positive breast tumors, which represent two‐thirds of all breast cancers. The precise role of ERβ and the manner in which it can impact estrogen actions in breast cancer, however, are less clear. ERβ is expressed in ca. 70% of human breast tumors, often along with ERα, with some human breast tumors expressing only ERβ (Kurebayashi et al, 2000; Speirs et al, 2004; Saji et al, 2005; Skliris et al, 2006). Although several reports have implicated ERβ as having net antiproliferative effects in breast cancer cells (Lazennec et al, 2001; Paruthiyil et al, 2004; Strom et al, 2004; Chang et al, 2006; Lin et al, 2007a; Williams et al, 2008), elucidation of the mechanistic basis for the seemingly contrasting actions of ERα and ERβ in breast cancer cells, including delineating the manner in which the genes involved are differentially selected for regulation by ERα and ERβ, and mapping of the signaling pathways utilized, remain critical issues.
When ERα and ERβ bind their ligand, 17β‐estradiol (E2), they undergo conformational changes that release heat shock proteins, enhancing receptor dimerization, interactions with coregulators (Skliris et al, 2006; Xu et al, 2009), and binding to the regulatory regions of target genes. ERs can be targeted to chromatin by direct recognition of estrogen response elements (EREs) through the agency of pioneer factors (e.g., FOXA1, GATA3, and PBX1) that modify the chromatin environment to a more permissive state, or via tethering to other TFs (e.g., Sp1 and AP1; Ali and Coombes, 2000; Glass and Rosenfeld, 2000; McKenna and O'Malley, 2002; Fullwood et al, 2009; Stender et al, 2010; Rosell et al, 2011; Jozwik and Carroll, 2012).
Given the fact that both ERs can potentially recognize similar chromatin‐binding sites, interact with a largely overlapping set of coregulators, and form both homo‐ and heterodimers in order to regulate gene expression and cell phenotypic properties, we explored how estradiol can elicit contrasting phenotypic outcomes—proproliferative versus antiproliferative—through these two closely related TFs. In this report, we have undertaken an integrative genomic approach to map in a comprehensive manner the chromatin‐binding interactions of ERα and ERβ, and their key coregulators, SRC3 and RIP140 (Cavailles et al, 1995; Glass and Rosenfeld, 2000; Xu et al, 2000; Rosell et al, 2011), in the same cell background when the receptors are present alone or together. The use of novel clustering algorithms enabled us to associate the distinct chromatin‐binding landscapes of these receptor and coregulator modules with ER‐regulated gene sets that delineate the specific cellular pathways and regulatory programs underlying the distinct phenotypic outcomes induced by hormone working through these two important NHRs in breast cancer cells.
These integrative and clustering approaches, delineating distinct genome‐wide patterns of chromatin binding of receptors and coregulators with gene expression behavior and functional outcomes, can be applied broadly to elucidate the molecular underpinnings for the transcriptional regulation and physiological effects of any TF in response to extrinsic or temporally modulated stimuli.
Genome‐wide analysis of ERα, ERβ, SRC3 and RIP140 chromatin binding by ChIP‐seq
Although ERα and ERβ have high structural and sequence homology, especially in their DNA‐binding domains, it is not known whether these closely related receptors, in the same cell background, would substitute for one another when present alone, whether they would synergize or antagonize each other at different regulatory gene sites when present together, and how their utilization of coregulators might contribute to their specification of activities at the many gene regulatory sites to which these ERs bind. To compare genome‐wide cartographies of ERα and ERβ, and their modulation of gene expression in these contexts, we utilized MCF‐7 breast cancer cells that endogenously express only ERα, or cells expressing only ERβ (adenovirally expressed ERβ with knockdown of ERα via RNAi), or both ERα and ERβ at similar levels (termed ERα/ERβ cells). The presence of ERα and/or ERβ protein in the three cell types used for genome‐wide analyses was confirmed by western blot (Supplementary Figure 1). ChIP samples for ChIP‐seq analysis and total RNA for cDNA microarray analysis under vehicle control or 10 nM estradiol treatments were prepared in parallel in these three cell contexts (see Figure 1A for the experimental scheme).
For cistrome analysis cells were treated with hormone for 45 min, whereas for gene expression studies total RNA was collected after an early (4 h) and a late (24 h) time of hormone treatment to gain insight into both rapid and more delayed responses. Whole‐genome binding sites profiled by the ChIP‐Seq approach were analyzed using MACS (Model‐based Analysis of ChIP‐Seq; Zhang et al, 2008), to identify binding sites with high confidence. Multiple rounds of ChIP‐seq data analysis were then carried out using the integrated ChIP‐seq data interpretation program, seqMINER (Ye et al, 2011), to identify, by clustering algorithms, subgroups of genomic loci for ERα, ERβ, and the coregulators SRC3 and RIP140 having similar compositional features. Further analyses were then done to relate subgroups of genomic binding sites (clusters) with binding sites for other TFs, and with genes regulated by estradiol in the three receptor cell backgrounds.
We first compared ERα and ERβ cistromes when the two receptors were present alone (ERα versus ERβ cells). As shown by the Venn diagram in Figure 1B, out of ca. 38 000 high‐confidence binding sites for each receptor in ERα or ERβ cells, only about 40% (15 000) overlapped, indicating that despite homologous DNA‐binding domains the two ERs bind to largely distinct chromatin regions. Notably, when both ERs were present together in cells, the number of binding sites for each ER was markedly reduced (to ca. 20 000), suggesting that under this condition a more complex chromatin interplay causes a mutual restriction in the total number of sites accessible to ERα and ERβ. Approximately half (7600) of the DNA regions bound by both ERα and ERβ in the ERα/ERβ cells overlapped with binding sites occupied by the receptors in the ERα‐ or ERβ‐only cells. Moreover, we identified ca. 2500 binding sites for ERα and ca. 3000 sites for ERβ, which were unique to the ERα/ERβ cells, further emphasizing that cells expressing both receptors present a different cistromic landscape.
To examine whether the concurrent presence of key coregulators for ERs would provide better stratification of ER‐binding sites, we performed ChIP‐seq for the key coregulators SRC3 and RIP140, known to be involved in estrogen‐controlled mechanisms of transcriptional regulation (Cavailles et al, 1995; Augereau et al, 2006; Lanz et al, 2010; Rosell et al, 2011; Madak‐Erdogan and Katzenellenbogen, 2012). In case of RIP140 (Figure 1C), we observed >2000 binding sites in ERα and ERα/ERβ cells, with an overlap of ca.1400. Interestingly, in the ERβ cells the number of RIP140‐binding sites was more than fourfold increased, to ca. 9000, suggestive of a major shift in the role of RIP140 when ERβ is present alone. On the other hand, this major change with ERβ was not seen for SRC3, where only a modest (ca.1.5‐fold) increase in the number of SRC3‐binding sites was observed in the ERβ cells (Figure 1D).
We then overlaid the coregulator cistromes with the ER cistromes (Figures 1E and F). As apparent from the Venn diagrams, we observed that the great majority of coregulator binding sites (>80%) were associated with genomic regions bound by ERα and/or ERβ, clearly showing that ERs are the major hubs for chromatin localization of these coregulators after hormone treatment.
Categorization of binding‐site clusters based on ERα, ERβ, SRC3, and RIP140 cistromes
To parse the wealth of data obtained from ChIP‐seq analyses, we examined whether the combinatorial usage of binding sites for ERs and coregulators might contain more physiologically and/or mechanistically relevant information. We therefore utilized a clustering approach using the seqMINER software (Ye et al, 2011), which provided a more comprehensive view of genome occupancy of ERs and coregulators across the different cell backgrounds (Figure 2A). This software compared the presence of the multiple factors at a given chromosomal location within a 300‐bp window in both directions and clustered together those binding sites that shared a similar pattern of factor localization.
In some of the clusters, we observed a consistent set ofca. 100‐bp steps in the location of the binding sites for the factors with respect to ER. As the 300‐bp window in both directions encompasses nearly three nucleosomes these factors are presumed to be recruited into the same complexes, but might have been crosslinked by the ChIP process to neighboring nucleosomes with respect to ER. These displacements are highly consistent in magnitude among different data sets, even when our data are compared with that of different factors from different laboratories (such as FOXA1 and GATA3; Supplementary Figure 2). This suggests that we are observing a biologically significant phenomenon that is inherent to these complexes and the recruitment process.
From these analyses, we identified 11 distinct clusters representing groups of loci having similar receptor and coregulator composition (Figure 2A). Cluster 1 represents ERα‐binding sites observed in ERα cells that were lost upon introduction of ERβ (ERα/ERβ cells) and were not bound by ERβ in ERβ‐only cells. Clusters 2 and 3 contained ERα‐ and/or ERβ‐binding sites in three (C2) or two (C3) of the cell contexts. Clusters 4–7 represent ER‐binding sites that were occupied by either receptor in all of the cell backgrounds. Clusters 8 and 9 contained binding sites predominantly in cells containing ERβ (ERα/ERβ or ERβ‐only), and Clusters 10 and 11 contained binding sites that were occupied by ERβ in ERβ‐only cells. A very interesting observation was that the presence of either SRC3 and/or RIP140 indicated that a binding site would be occupied by either ERα and/or ERβ in all the cell backgrounds (Clusters 4 and 5), whereas the two coregulators were virtually absent at binding sites in Clusters 1–3 and Clusters 10 and 11, which represent ERα and ERβ unique sites, respectively, and these coregulators were also reduced in Clusters 6–9. ER bound at sites in these clusters might use alternative coregulators. As discussed below, Clusters 4–7 are associated with many of the classic estrogen‐regulated genes, such as TFF1/pS2 (trefoil factor 1) and GREB1.
Associations of clusters with FOXA1‐ and GATA3‐binding sites
With the public availability of other ChIP‐seq data sets from MCF‐7 ERα cells, we combined our cistrome analysis with genome‐wide data on FOXA1 (GSE data set 26 831 (Joseph et al, 2010)) and GATA3 (GSE 29 073 (Kong et al, 2011)), which are reported to be important pioneering TFs for ER‐binding site determination and epithelial lineage differentiation (Carroll et al, 2005, 2006; Eeckhoute et al, 2007; Lin et al, 2007b; Lupien et al, 2008; Hurtado et al, 2011; Kong et al, 2011). By doing so, we were able to further stratify Cluster 1 into three smaller clusters (Clusters 1a, 1b, and 1c; Figure 2A) that are distinguished by having (Clusters 1b and 1c) or lacking (Cluster 1a) FOXA1 and/or GATA3 binding, overlapping with ERα‐binding sites (Supplementary Figure 2). FOXA1‐ and GATA3‐binding sites were also associated with Clusters 4 and 6 but, interestingly, not with Clusters 2, 3, 5, 7, or 8, despite ERα binding in these clusters. We also observed no association of GATA3 and FOXA1 with Clusters 9, 10, and 11 specific for ERβ DNA binding, implying that FOXA1‐ and GATA3‐binding sites in ERα cells do not predict ERβ binding to chromatin in ERβ cells.
Enriched TF‐binding motifs and binding‐site mapping in the clusters
To further examine these clusters, we performed TF motif analysis prediction for each cluster in Figure 2A, using Seqpos and CEAS (cis‐regulatory element annotation system; Shin et al, 2009) softwares ( http://liulab.dfci.harvard.edu/CEAS/usermanual.html). Table I shows the most highly enriched binding‐site motifs for each cluster. As expected, the ER motif was the most highly enriched motif, but other motifs were characteristic of specific clusters (Table I, and see Supplementary Tables 1 and 2 for more detailed analysis of enriched TF motifs in the ER‐binding sites in Clusters 1–11, based on SeqPos and CEAS analysis), suggesting a diversity of pathways likely to be represented after Gene Ontology (GO) analysis (described below). Also consistent with our cistromes overlay analysis, FOXA1‐ and GATA3‐binding motifs were enriched only at ERα‐binding sites (Table I). Interestingly, motif enrichment analysis failed to show any major differences in the composition of TF‐binding motifs or strength of EREs between different clusters, suggesting that the selectivity for binding might be due to certain epigenetic marks or other factors (such as FOXA1 and GATA3) that render these sites accessible for ER binding with or without SRC3 and/or RIP140 in different cell backgrounds (Supplementary Figure 2).
Regarding the genomic location of the receptor binding sites, ERα and ERβ were found to bind preferentially to distal enhancer regions, followed by intronic regions, in almost all clusters (Figure 2B). Notably, however, almost 50% of ER‐binding sites in C8, 16% in C9, and 8% in C10 mapped to proximal promoter regions. These clusters showed enrichment for cell cycle‐related genes, and E2F predicted binding site motifs (Table I). Another notable feature of cluster 8 was its high sequence conservation among different species, which was the highest of all clusters (Figure 2C).
Patterns of E2‐regulated gene expression in ERα, ERα/ERβ, and ERβ cells, and their association with ER and coregulator cistromes and cellular pathways
To understand the relative importance and functional outcomes of each of the identified clusters of ER and coregulator binding sites in each cell background, we compared gene expression cDNA microarrays in the three (ERα, ERβ, and ERα/ERβ) cell backgrounds at two time points (4 h and 24 h). After data normalization and clustering analysis of the three cell transcriptomes (Figure 3A and B), we found that at the 4‐h time point estradiol regulated, both up or down, a similar number of genes in ERα cells (657) or ERβ cells (674). In cells with both receptors, there was a greater overlap with ERα cells (∼30%), and more importantly there was a significant number of genes that were uniquely regulated in ERα cells (253), in ERβ cells (555), or in ERα/ERβ cells (243). Hence, ERα and ERβ, alone and together, regulate some common but also a large number of unique transcripts.
We observed that the number of genes regulated in all cell backgrounds was greater at 24 h versus 4 h (1296 versus 657 in ERα cells, 1133 versus 658 in ERα/ERβ cells, and 988 versus 674 in ERβ cells, respectively). As an alternative way of describing the data, we generated scatter plots of the regulated genes in each cell background (Figure 3C), which highlight two main trends in this data. One is that the direction of regulation by ERα (up (shown in red) or down (shown in blue)) does not change with the copresence of ERβ or in the ERβ‐only cells (i.e., if a gene is regulated in all the cell backgrounds, it will always change in the same direction); and two, the magnitude of the E2 response varies among cell backgrounds, generally being highest in ERα cells, followed by ERα/ERβ cells and then ERβ cells.
To investigate how changes in hormone‐mediated gene expression are related to receptor and coregulator binding events, we mapped the closest E2‐regulated gene that resides within a 20‐kb window of an ER‐binding site. We found that two‐thirds of all of the E2‐regulated genes harbored at least one receptor binding site in its proximity (2861/4190 genes; Figure 4A), and often harbored multiple sites for different clusters (Figure 4B and C). This allowed us to use powerful statistical means to evaluate the predictive power of the ER‐binding site clusters, described in Figure 2A, in determining the direction and magnitude of E2‐dependent gene regulation with different cell background and duration of ligand treatment (4 h and 24 h; Figure 4B). For this, we calculated enrichment of E2‐regulated genes relative to all the genes with a binding site in the respective cluster, and then performed clustering.
From this analysis (Figure 4B) we found that, overall, upregulated genes in the three cell contexts at both times were most highly enriched in binding sites from Clusters 2, 4, and 6. By contrast, downregulated genes overall were enriched in sites from Clusters 1b and 1c, and in most cases showed deficits in binding sites from several other clusters (e.g., C4, C3, and C7). Other patterns of binding‐site enrichment and deficits could be associated with particular upregulation or downregulation of gene expression in different cell contexts or times. Thus, this clustering approach not only is able to stratify ER‐binding sites based upon co‐occupancy of these TFs and coregulators, but also has predictive value for the direction and time of E2‐mediated gene regulation for any given cell background.
To further examine the relationship between hormone‐regulated genes and the binding sites from different clusters, we clustered the BED files for the genes that map to binding sites in each cluster with those for hormone‐regulated genes at different times of ligand treatment in different cell backgrounds that contained at least one ER‐binding site (Figure 4C). This schematically illustrates that hormone‐regulated genes are highly likely to have more than one type of ER‐binding site, often having representatives from several different clusters.
Cell receptor background‐specific effects on profiles of gene regulation and ERα‐, ERβ‐, SRC3‐, and RIP140‐binding sites
We investigated the cell background‐specific effects on selected E2‐regulated genes and ER and coregulator binding sites. For each example provided, we show a Genome Browser view of the genomic region of interest with the highlightedER‐binding site patterns in each cell background, a dose–response curve of mRNA expression, and ChIP assays for ERα, ERβ, SRC3 and RIP140 in the three cell backgrounds. For example, the well‐known estrogen‐regulated gene TFF1/pS2 has binding sites from Clusters 1, 2, and 4 (Figure 4D); its mRNA (Figure 4E) is upregulated by its natural ligand, E2, in all cell backgrounds, albeit to a different magnitude, with the highest stimulation being in the ERα/ERβ cells. Also, the E2 dose–response analysis reveals a change in sensitivity to E2, evident in the ERα/ERβ cells, where the EC50 for E2 is close to 10−12 M, compared with 10−11 to 10−10 M in ERα or ERβ cells. By ChIP assays (Figure 4F), it is also evident that there is an increase in ERα recruitment when ERβ is present, and a slight increase in ERβ recruitment when present with ERα, and this is mirrored also in elevated SRC3 and RIP140 recruitment in ERα/ERβ cells. For the progesterone receptor (PgR) gene, which harbors binding sites from Clusters 1, 2, and 7 (Figure 4G), there were no major differences in magnitude of response or recruitment of ERs in the three cell backgrounds (Figure 4H and I); however, a change in sensitivity to E2 was observed but opposite in direction to that seen for TFF1 (EC50 10−10 M for ERα cells, 10‐9 M for ERα/ERβ or ERβ cells). Finally, for the otubain2 (OTUB2) gene, which we previously identified as an ERβ target gene (Chang et al, 2006, 2008), as seen in Figure 4J–L, we confirmed little regulation of OTUB2 mRNA in ERα cells, and little ERα recruitment and no stimulation of SRC3 or RIP140 recruitment in these cells; by contrast, maximal stimulation of the OTUB2 gene expression was seen in ERβ and ERα/ERβ cells (Figure 4K and Supplementary Figure 3). Interestingly, ChIP assays showed that ERα could be recruited to the gene regulatory site when ERβ was also present (ERα/ERβ cells), and this correlated with increased recruitment of SRC3 and RIP140 in this cell background (Figure 4L). ERβ showed increased recruitment in response to E2 in ERα/ERβ and ERβ cells that was accompanied by increased recruitment of both coregulators in these two cell backgrounds (Figure 4L). At the Ki‐67 gene, RIP140 recruitment was strongly associated with E2 stimulation in ERβ‐only cells, consistent with our ChIP‐Seq data (Supplementary Figure 3).
Pathway analysis of the E2‐regulated genes and effects on hormone‐regulated proliferation and their association with ERα‐ and ERβ‐binding site clusters
The GO analysis of the binding‐site clusters was performed using web‐based Metacore, DAVID, Cytoscape ClueGO plugin and Genespring software programs (Dennis et al, 2003; Huang da et al, 2009), and the main results are shown in Figure 5 and Table I. Of note, the major group of genes found to be differentially regulated between ERα and ERβ cells was associated with regulation of the M‐phase of the cell cycle (Figure 5A–C) (a complete description of GO terms associated with all of the binding‐site clusters is given in Table I).
The GO terms ‘cell cycle related’ were associated withE2‐stimulated genes selectively in ERα cells, and this was abrogated in ERα/ERβ cells. Interestingly, genes involved in the G1–S transition of the cell cycle were upregulated by E2 in both ERα‐only and ERβ‐only cells, but not in ERα/ERβ cells, and genes involved in the G2–M transition were upregulated only in ERα‐containing cells but not in the other two cell backgrounds, highlighting profound differences in regulation of cell cycle stages by each receptor, either alone or in combination (Figure 5, Table I, and Supplementary Figure 4). We therefore examined in more detail the impact of ERα and ERβ on cell proliferation in the three cell backgrounds, with a specific interest in understanding the role of ERβ when present alone or in combination with ERα.
As seen in Figure 6A, ERα cells showed a stimulation of proliferation in response to E2. The copresence of ERβ greatly reduced proliferation, and cells expressing only ERβ showed little proliferative stimulation by E2. Because of the interesting differences in regulation of cell cycle progression by each ER that were highlighted by the combination of gene expression and cistrome analyses, we went on to validate this also by flow cytometry (Figure 6B). As predicted, we observed anE2‐dependent increase in S phase in both ERα and ERβ cells, but not in ERα/ERβ cells. Notably, even though there was a basal increase in G2–M phase in ERβ‐containing cells, ligand‐dependent changes in the proportion of G2–M phase cells were observed exclusively in ERα cells, and not in the other two cell types, indicating that the entire cell cycle machinery is activated only when ERα is present alone. The changes in cell proliferation and phases of the cell cycle were dependent on the RIP140 coregulator function, as knockdown of RIP140 abrogated the E2‐mediated increase in cell number (Supplementary Figure 4A) and percent of cells in G2–M phase in all cell backgrounds (Supplementary Figure 4B). Interestingly, in ERβ‐containing cells RIP140 knockdown also increased the percent of S‐phase cells, consistent with a report suggesting a corepressor function for RIP140 action with E2F1 (Docquier et al, 2010). There was also an increase in the subG1 population when we introduced ERβ, and this was further enhanced with RIP140 knockdown, implying a corepressor function for RIP140. However, RIP140 overall acted as a coactivator for E2‐stimulated cell proliferation mediated by ERα.
Our gene expression microarray data suggest that the effect of ERβ might be due to reduced expression of important G2–M phase activators (Figure 6C and Supplementary Figure 5), such as FOXM1, as shown in ERα/ERβ cells previously (Chang et al, 2006), and Ki‐67 and BIRC5, which all showed a RIP140‐dependent increase in expression level only in ERα cells (Supplementary Figure 4C).
Analysis of gene networks and pathways associated with the influence of ERβ: impact on cell proliferation, MAPK activation, adipogenesis, and the involvement of RIP140
Our examination of the role of ERβ in cell proliferation revealed that ERβ affects cell proliferation in at least three important ways. As demonstrated in Figure 6, ERβ changed the cell cycle profile and was unable to drive the cells into G2/M. This first mechanism could be due to the loss of central factors, such as FOXM1(Supplementary Figure 4C; Chang et al, 2006), and a direct effect on cell cycle genes via the loss of binding of ERα when ERβ is present, as shown by Cluster 1, and gain of genomic binding of ERβ, as shown by Cluster 8, which is enriched for cell cycle and apoptosis‐related genes.
Second, we observed that ERα binds to sites near important pathway modulator genes (such as MAP2K1‐MEK1 and MAP3K1‐MEKK1). In our previous studies (Madak‐Erdogan et al, 2011) we showed that MAPK signaling was essential for E2‐mediated cell proliferation. To test the impact of ERβ on MAPK signaling, we monitored ERK1/2 activation after E2 treatment of cells containing both receptors. Our findings (Figure 7A) show that E2 increased pERK1/2 in cells containing ERα only (C, control cells with ERα but no ERβ) and that increasing concentration of ERβ abrogated pERK1/2 activation.
A third mechanism we uncovered was the major influence of ERβ on metabolic pathways (Figure 7 B–D) due to its preferential coupling with RIP140. We observed in our cistromes analysis that the number of genomic binding sites for RIP140, a well‐known factor in metabolic regulation (Rosell et al, 2011; Siersbaek et al, 2012), was dramatically increased in E2‐treated cells containing ERβ only. GO analysis of these binding sites (Table I, and see Supplementary Table 3 for functional enrichment of pathway and process networks for ERβ‐modulated genes) revealed enrichment for genes associated with fatty acid metabolism, regulation of NF‐κB and inflammation, IGF‐1 signaling, and mammary gland epithelial development.
RIP140 is an estrogen‐induced protein (Cavailles et al, 1995; Augereau et al, 2006) and is known to be a repressor of metabolic pathways (Christian et al, 2005; Rosell et al, 2011), especially lipolysis in adipocytes. RIP140 was 3‐fold upregulated by E2 in ERα and ERα/ERβ cells but, most remarkably, RIP140 was more than 100‐fold upregulated by hormone in ERβ cells (Figure 7B). Of note, we observed greatly elevated expression of key adipogenic factors, such as DDIT3, FOXO1A, and MEF2A, in ERβ cells treated with E2, consistent with the binding data for ERs and coregulators (Figure 7C and Supplementary Figure 6). This gene regulation was dependent on RIP140 and SRC3, as knockdown of either coregulator blocked the mRNA induction (Supplementary Figure 7). Interestingly, recruitment of each coregulator was reduced when the other coregulator was depleted, suggesting a cooperative role for each coregulator in the recruitment of a functional complex. When we monitored mRNA expression and coregulator recruitment for other genes not associated with adipogenesis, such as pS2 or OTUB2, we saw that SRC3 was required for mRNA induction and RIP140 recruitment in all cell backgrounds, while RIP140 showed a gene and cell background‐specific function. For TFF1, RIP140 knockdown increased mRNA induction and interestingly increased recruitment of SRC3 in all receptor cell backgrounds. For OTUB2, RIP140 knockdown in ERα/ERβ cells increased mRNA induction and SRC3 recruitment, whereas both of these were blocked by the knockdown of RIP140 in ERβ cells. These results suggest that although RIP140 acts as a coactivator or corepressor in a gene and cell background‐dependent manner, adipogenesis requires high levels of RIP140 that we observe in ERβ cells.
To further assess the effect of ERβ and RIP140 on adipogenic properties of the breast cancer cells, we monitored the cellular level of triglycerides (Figure 7D), which revealed a threefold increase in triglyceride levels in ERβ cells. This effect was dependent on RIP140.
These findings in MCF‐7 cells were also verified in another ERα‐positive cell line, T47D, where we introduced ERβ and generated cells with ERα and/or ERβ. Findings from cell proliferation (Supplementary Figure 8A), adipogenesis (Supplementary Figure 8B), and adipogenesis‐related gene expression (Supplementary Figure 8C) studies in these T47D cells were consistent with our observations in MCF‐7 cells. The results rule out phenotypic differences and transcriptional regulatory effects being unique to the MCF‐7 cells. The observations in these MCF‐7 and T47D breast cancer models also highlight several novel gene networks involved in delineating the phenotypic properties of the ERβ cells and the suppression of proliferation by this receptor.
A gene signature associated with RIP140 predicts better prognosis in breast cancer patients
To examine possible importance of our gene targets in breast cancer patient prognosis, we used the Oncomine database to obtain a 20‐gene signature associated with ERβ‐ and RIP140‐binding sites, and found that this signature predicted better prognosis for patients in several large clinical data sets. The initial gene list that we used to derive the signature included genes whose expression was modulated by ERβ in our gene expression microarrays and harbored binding sites for ERβ and RIP140 within 20 kb of the transcription start site of the gene. The gene signature predicted clinical outcome in 23 data sets from Oncomine (Figure 7E and Supplementary Figure 9). Moreover, this signature was most predictive for outcome in breast cancer as compared with other cancers, as we observed an enrichment of concepts associated with breast cancer but not other cancer types (Supplementary Figure 9). The signature genes included RIP140; oxidoreductases ALDH3A2, GLUD1, PRDX3, and GPX4; ribosomal proteins RPL22, RPLP1, and RPS15; negative regulators of gene transcription BNIP3L, DICER1, and TCEAL1; cell cycle regulators RAD17, CREBL2, and CGRRF1; and ubiquitin modifiers TTC3 and UBE3A. These genes were coexpressed in breast tumors (Figure 7F), and their overexpression marked tumors that had a better prognosis in terms of overall survival and reduced recurrence (Figure 7G and Supplementary Figure 10). Thus, our integrative genomics approach revealed genes that are associated with better prognosis in breast cancer patients.
The regulation of transcription in eukaryotes involves the assembly of multiprotein complexes at specific sites on chromosomal DNA that occurs in a temporally coordinated and cell‐specific manner in response to specific stimuli (Katzenellenbogen et al, 2000; Katzenellenbogen and Katzenellenbogen, 2002; Welboren et al, 2009; Madak‐Erdogan et al, 2011; Ross‐Innes et al, 2012). Herein we have shown that multiple rounds of data analysis using clustering algorithms have enabled us to delineate a hierarchy of ER‐ and coregulator‐binding site subgroups that demarcate and govern gene expression programs and lead to the distinct physiological outcomes controlled by the closely related TFs, ERα and ERβ. The two ERs have very specific cistromes and differentially engage SRC3 and RIP140 at different chromatin sites to regulate distinct sets of target genes that can explain the ERα proliferative and ERβ antiproliferative activities of these receptors. In fact, we find that ERα is capable of regulating genes involved with all phases of the cell cycle when expressed alone, whereas ERβ when coexpressed with ERα dampens this response. When ERβ is expressed alone, it appears capable of driving the cells into S phase, but is then incapable of successfully completing the cell cycle. Moreover, this action of ERβ on cell cycle genes is complemented by ERβ cell‐selective effects on lipid metabolism genes that can be linked to an ERβ preferential usage of the coregulator RIP140.
Given the clinical importance of ERα and ERβ in breast cancer development and in the response of patients to breast cancer treatments (Ali and Coombes, 2000; Katzenellenbogen and Katzenellenbogen, 2000; Speirs et al, 2004; Saji et al, 2005; Deroo and Korach, 2006; Skliris et al, 2006; Chang et al, 2008; Nilsson and Gustafsson, 2011), these findings illuminate the molecular mechanisms that underlie the distinct differences in the biological phenotypes and in the endocrine therapy responses of ERα‐ and/or ERβ‐positive breast cancers. Our analyses of gene expression data on human breast tumors where we could divide tumors into those expressing ERα or ERβ only, or both receptors, showed distinct gene expression profiles for these tumors that were consistent with our observations in MCF‐7 cells having these three different complements of ERs. In this human tumor gene expression microarray study, major differences were observed in cell proliferation genes that were correlated with clinical outcome and patient survival (Gruvberger‐Saal et al, 2007). Patients in this cohort, whose tumors were ERα‐negative and ERβ‐positive, had a better outcome on tamoxifen, consistent with a reduced proliferation status. Moreover, when we combined our gene expression data and binding site analysis to obtain an ERβ‐specific gene signature, this included RIP140 as well as other genes whose overexpression predicted better prognosis for breast cancer patients. Although speculative, utilization of ERβ‐selective agonists (Meyers et al, 2001; Harris, 2006; Charn et al, 2010) and activation of ERβ gene programs in ERα plus ERβ‐positive or ERβ‐only‐containing breast cancers, or increasing levels of RIP140, which could have a growth suppressive activity through ERβ, might provide a benefit for patients with ERβ‐positive breast tumors.
Although there have been reports on the binding of ERα and SRC3 in cells containing only ERα (Lanz et al, 2010), and we and others have shown effects of introduction of ERβ on the expression of some ERα‐regulated genes (Paruthiyil et al, 2004; Chang et al, 2006; Lin et al, 2007a; Chang et al, 2008; Charn et al, 2010; Grober et al, 2011; Nassa et al, 2011), this is the first study to comprehensively compare the binding of ERs and coregulators, and the effects of hormone on gene expression in breast cancer cells containing ERα and ERβ alone and together. To the best of our knowledge, this is also the first report of genome‐wide ChIP‐seq binding of RIP140 in any MCF‐7 ER context, and the first for the SRC3 cistrome in ERα and ERβ, and ERβ‐only cells. We selected these two central ER coregulators for the study, as SRC3, also known as amplified in breast cancer‐1 (Anzick et al, 1997), is a coactivator that promotes estrogen‐stimulated transcriptional activity by ERα (Anzick et al, 1997; Suen et al, 1998; Azorsa et al, 2001; Louie et al, 2004); RIP140 is an especially interesting coregulator, as it can act either as a coactivator or a corepressor (Cavailles et al, 1995; Rosell et al, 2011).
ER‐positive breast cancers differ in their content of ERα and ERβ, and as the cancer progresses to increasingly malignant character there is usually a decline in the level of ERβ relative to ERα (Speirs et al, 2004; Saji et al, 2005). To examine how estrogen action in breast cancer proceeds through the two ERs, in terms of chromatin‐binding site selection and gene regulation, to result in distinct cellular phenotypes, we constructed a model that was based on MCF‐7 cells. We chose this as the cell background because it is the most well‐studied breast cancer cell line in which much is known about ERα‐binding sites as well as cistromes for certain other cofactors, such as FOXA1 and GATA3, providing a broad context for our studies (Carroll et al, 2005, 2006; Eeckhoute et al, 2007; Lupien et al, 2008; Hurtado et al, 2011; Liu et al, 2011; Ross‐Innes et al, 2012). In addition, the introduction of ERβ results in a decrease in proliferation and a change in other cellular properties that mirror the effects of ERβ in human breast cancers, based on gene microarray and patient outcome data (Gruvberger‐Saal et al, 2007).
A multilayer model of progressive specification that underlies the distinctive activity of estradiol (E2) through ERα versus ERβ in breast cancer cells
In this study, as schematized in the model in Figure 8, we have examined at multiple levels the ways by which estradiol (E2) acts through its closely related receptors in breast cancer cells to generate distinct phenotypic and proliferative outcomes. At the chromatin‐binding level, we show that ERα and ERβ select distinct, only partially overlapping sets of chromatin‐binding sites when present alone; when present together, they compete for these sites, and restrict and shift each other's binding. There are also characteristic differences in the ERα and ERβ sites in terms of TF‐binding motif enrichments and genomic location.
At a second level, the coregulators SRC3 and RIP140, known to have important roles in the mammary gland and in breast cancer (Cavailles et al, 1995; Font de Mora and Brown, 2000; Xu et al, 2000; Augereau et al, 2006), become recruited in a distinctive manner to some sites at which ERα and/or ERβ were bound to chromatin, and some hormone‐induced genes were strongly associated with ER‐binding sites that were also co‐occupied by these coregulators. Of note was the increased (ca. fourfold) number of RIP140‐binding sites in the ERβ cells, and also the greatly increased expression of RIP140 in ERβ cells, highlighting a potentially preferential action of this coregulator with ERβ. Both coregulators were highly recruited to chromatin together with ER after E2 treatment, as shown by a >80% overlap between the binding sites of the coregulators and the binding sites of the ERs. The majority of ER‐binding sites, however, did not show any colocalization with SRC3 or RIP140; at these sites, ER might be using other coregulators (McKenna and O'Malley, 2002; Green and Carroll, 2007). Recent studies have elegantly demonstrated the promoter‐specific involvement of distinct coregulators in the regulation of estrogen‐responsive genes in MCF‐7 cells. Thus, SRC3 was present at nine of the ten genes examined in MCF‐7 ERα cells, whereas fewer coregulators were involved in estrogen stimulation of other genes, such as the PgR gene (Won Jeong et al, 2012).
The overlay of our ERα and ERβ cistromes with our SRC3 and RIP140 coregulator cistromes, and with genome‐wide cistrome analysis of the FOXA1 and GATA3 pioneering factors by others (Carroll et al, 2005, 2006; Eeckhoute et al, 2007; Lupien et al, 2008), followed by clustering, enabled a stratification of binding‐site patterns that predict time and direction of gene expression and phenotypic changes in each cell background. FOXA1 and GATA3 are important DNA‐bound cofactors of ERα (Carroll et al, 2005, 2006; Eeckhoute et al, 2007; Lupien et al, 2008), and they have been shown to be necessary for ERα chromatin binding and subsequent activation of many ERα target genes. This additional layer of analysis allowed us to organize the ER‐binding sites into 11 clusters having distinct characteristics in terms of predicted TF‐binding site motifs, genomic location, conservation, associated genes and enriched pathways. Finally, the direction of E2 regulation (stimulation or repression) of the target genes could be defined based upon clustering analysis of ER and coregulator cistromes. Of note, each cluster of binding sites could be associated with distinctive functional gene ontologies that can explain avenues through which ERα and ERβ modulate cell proliferation and metabolism.
Specification of estrogen action at the chromatin layer
Despite sharing great sequence and structural homologyin the DNA‐binding domain (Katzenellenbogen and Katzenellenbogen, 2000; Nilsson et al, 2001; Deroo and Korach, 2006), ERα and ERβ select distinct, only partially overlapping chromatin‐binding sites, with more sites being accessed by each ER when present by itself. When both are together, mutual competition between them restricts the number of sites that are occupied by each ER subtype and contracts the number of overlapping sites; moreover, we also observe the appearance of many new binding sites for ERα and ERβ when both are present, which indicates a profound redistribution of cistromes when the receptors are copresent.
Although the two ERs can both homodimerize and heterodimerize, and both homo and heterodimers can be observed in vitro and in cells (Powell and Xu, 2008; Paulmurugan et al, 2011; Powell et al, 2012), the potential number of chromatin‐binding sites through which ERα‐ERβ heterodimers might function (12 000 shared sites) is only a fraction (43%) of the total number of ER‐binding sites that can be occupied by either ER in the ERα/ERβ cells (28 000). In addition, most of the new sites present only in ERα/ERβ cells are associated uniquely with either ERα or ERβ, and thus can be occupied only by homodimers of ERα and ERβ. Clearly then, the context of a chromatin‐binding site imposes a significant specification of whether ER homo or heterodimers bind, regardless of their cellular abundance.
Specification of ER action at the gene function–phenotypic layer
Among the nuclear receptor family, it has proved challenging to associate ER chromatin‐binding sites with estrogen‐regulated genes, because many of the ER regulatory sites are in far‐distal enhancer regions rather than in proximal promoters (Carroll et al, 2005, 2006; Lin et al, 2007b; Madak‐Erdogan et al, 2011). Nevertheless, through our genome‐wide location and expression analyses, we were able to define several novel mechanisms underlying the effects of ERβ on cell proliferation. By analyzing the effect of E2 on the cell cycle, we found that E2 directly regulated genes involved in all the phases of the cell cycle in ERα cells; this was dramatically altered and reduced, however, with ERβ copresence. Very interestingly, ERβ, when present alone, was able to increase the percent of cells in S phase, but cells were then unable to complete the cell cycle because they became arrested at the G2/M transition.
GO analysis of ER‐binding site clusters revealed an intriguing model for ERβ antiproliferative action. First, ERβ impinges on cell cycle genes in multiple ways: by reducing ERα binding when ERβ is present (cluster 1), by binding directly to promoters of cell cycle‐regulated genes (cluster 8), and by reducing the levels of factors implicated in G2/M transition (e.g., FOXM1). Cluster 8 was particularly interesting, as it was enriched in promoter proximal binding sites that contained overrepresented motifs for E2F TFs. The E2F family members consist of activators and repressors, and the fact that our findings show that ERβ‐binding sites in this specific cluster are more enriched with E2F motifs (Table I) raises the possibility that ERβ might be preferentially binding to the chromatin near E2F‐binding sites and utilizing some of the E2F members to affect the regulation of its target genes, consistent with their suggested involvement from prior studies (Strom et al, 2004; Stender et al, 2007).
ERβ also exerted an extensive effect on MAPK signaling, which we previously showed to be essential for stimulation of cell proliferation by E2 in ERα‐containing breast cancer cells (Madak‐Erdogan et al, 2011). Moreover, ERβ altered cell metabolism that is central to cell growth, where it was found to act in several reinforcing ways. ERβ increased the genomic binding and expression level of RIP140, a critical metabolic regulator, in ERβ‐only cells. These unique binding sites were related to genes involved in IGF‐1 signaling, inflammation, and fatty acid metabolism. Functional outcomes associated with increased binding of RIP140 included the ability of ERβ to directly bind and increase the expression of adipogenesis genes, elevating triglyceride levels in ERβ cells.
An integrative, whole‐genome and clustering analysis for elucidating the progressive specification of the biological activities of NHR complexes
Through integration of cistromes and transcriptomes, and by application of clustering analyses, we have been able to define a combinatorial pattern of regulation by these factors, in response to hormone stimulation, that specifies gene regulation across binding site clusters. Importantly, gene expression thereby becomes predictable and allows association of multiple regulatory sequences with a pattern of control defined by the cell's receptor context and the time of hormonal exposure.
Work in many laboratories, in both mammalian systems (Farnham, 2009) and in yeast (Beer and Tavazoie, 2004), has highlighted the difficulty in associating TF cistrome binding with gene expression outcomes. Our approach using clustering analyses, which we have used to elucidate patterns of binding of multiple regulatory factors and their association with gene regulation by two of the most important NHRs in human breast cancer, should be of help in characterizing the layers of specification that determine the regulation of transcription and distinct functional outcomes by any TF of interest.
Materials and methods
Cell culture, adenovirus, lentivirus, and siRNA, and ligand treatments
MCF‐7 cells were grown in minimal essential medium (MEM) (Sigma, St Louis, MO), supplemented with 5% calf serum (HyClone, Logan, UT), and 100 μg/ml penicillin/streptomycin (Invitrogen, Carlsbad, CA) (Chang et al, 2006, 2008). For experiments, the cells were maintained in phenol red‐free MEM plus 5% charcoal‐dextran‐treated calf serum for at least 3 days, and were then seeded at a density of 3 × 105 cells per 10‐cm tissue culture dish (Corning, Corning, NY) for 2 days before adenovirus infection. Recombinant adenoviruses were constructed and prepared as described (Chang et al, 2006). Adenoviral infection conditions were those described previously (Chang et al, 2006; Frasor et al, 2006; Charn et al, 2010) and were used to generate MCF‐7 cells expressing levels of ERβ equal to that of the endogenously expressed ERα. ERβ‐only cells were generated from these cells by knockdown of ERα in parental cells using the following siERα sequences from Dharmacon: forward, 5′‐ UCAUCGCAUUCCUUGCAAAdTdT‐3′, and reverse, 5′‐ UUUGCAAGGAAUGCGAUGAdTdT‐3′. siRNA experiments were performed as previously described, and resulted in knocknown of ERα mRNA and protein by greater than 95% (Chang et al, 2008). Briefly, cells were transfected with 20 nM siCtrl or siERα for 48 h after infection. Then, cells were treated with 0.1% EtOH (Veh) or 10 nM E2 (Sigma‐Aldrich) for the indicated times. For coregulator knockdown experiments, 24 h after control or ERβ adenovirus infection, cells were infected with control, SRC3 or RIP140 shRNA lentiviral particles according to the manufacturer's suggestions (Santa Cruz Biotechnology). All experiments were conducted with three or more replicates.
ChIP assay was carried out as described (Barnett et al, 2008) for ERα, ERβ, SRC3, and RIP140. Briefly, the different cells were treated with 0.1% EtOH (Veh) or 10 nM E2 for 45 min. Chromatin was crosslinked using 1% formaldehyde for 15 min at room temperature. Cells were washed with PBS and were collected. After three times 10 s sonication in ChIP lysis buffer, the samples were centrifuged for 10 min at 4 °C. The antibodies used for immunoprecipitation of DNA/protein complexes were ERα antibody HC‐20 (Santa Cruz Biotechnology); ERβ antibodies were a combination with equal parts of CWK‐F12 (produced in our lab (Choi et al, 2001)), GTX70182 (GeneTex), GR40 (Calbiochem), and PA1‐311 (Affinity Bioreagents); SRC3 antibody SC‐9119 (Santa Cruz Biotechnology); and RIP140 antibody SC‐8997 (Santa Cruz Biotechnology). Precipitated complexes were washed with RIPA buffer (three times) and TE buffer (2 times). After overnight incubation at 68 °C, ChIP DNA was isolated using QIAGEN PCR purification kit as per the manufacturer's suggestions. The ChIP DNA was used forChIP‐seq analysis and quantitative real‐time PCR.
ChIP‐seq analysis and clustering
The ChIP DNA was prepared into libraries according to Illumina Solexa ChIP‐seq sample processing methods, and single read sequencing was performed using the Illumina Solexa Genomic Analyzer. Sequences generated were mapped uniquely onto the human genome (hg18) by ELAND. MACS algorithm (Zhang et al, 2008) was used to identify enriched peak regions (default settings) with a P‐value cutoff of 6.0e−7 and FDR of 0.01. ChIP‐seq data will be available from the Gene Expression Omnibus database.
The seqMINER density array method with a 300‐bp window in both directions was used for the generation of clusters, i.e., groups of loci having similar compositional features (Ye et al, 2011). This ChIP‐seq data interpretation platform allows the comparison and integration of multiple ChIP‐seq data sets, and their extraction and visualization of specific patterns within data sets. A meta binding BED file was generated using ChIP‐seq tool set, which included ER‐binding sites from all three cell backgrounds (Blahnik et al, 2010). As there were only a relatively small number of coregulator‐specific binding sites without any receptor binding, these sites were not included in the clustering step. This reference file was then used to map the binding sites for each factor from different cell backgrounds. BED files for each cluster were used for further analysis with Galaxy Cistrome integrative analysis tools (Venn diagram, conservation, CEAS; Liu et al, 2011).
Gene Chip mRNA transcriptional profiling microarrays
Total RNA was used to generate cRNA, which was labeled with biotin according to techniques recommended by Affymetrix. All analyses were done for three or more samples for each treatment. The biotin‐labeled cRNA was then hybridized to Affymetrix U133 plus 2.0 GeneChips, which contain oligonucleotide probe sets for over 47 000 transcripts. After washing, the chips were scanned and analyzed using Affymetrix processing software. CEL files were processed using GeneSpring GX 11.0 software (Agilent) to obtain fold‐change and P‐value with Benjamini and Hochberg multiple test correction (Hochberg and Benjamini, 1990) for each gene for each treatment relative to the vehicle control. We considered genes with fold‐change >1.8 and P‐value <0.05 as statistically significant, differentially expressed.
Motif and GO category analysis
Overrepresented GO biological processes were determined by the web‐based DAVID Bioinformatics Resources database (Dennis et al, 2003; Huang da et al, 2009) and the GeneSpring and web‐based GREAT (Genomic Regions Enrichment of Annotations Tool) software (McLean et al, 2010). Motif‐enrichment analysis was done using Seqpos, which uses TRANSFAC PWMs, JASPAR, and de novo motif discovery tools (Liu et al, 2011).
Cell proliferation and flow cytometry
MCF‐7 cells were seeded at 30 000 cells/well in six‐well plates. On the second day, the cells were infected with AdGal or AdERβ at multiplicity of infection 20 for 4 h and then the medium was changed. After 24 h, the cells were transfected with 20 nM siCtrl or siERα. At 24 h, cells were treated with Veh (0.01% EtOH) or 10 nM E2 (day 0) and again at day 2. Cell number was monitored using the MTS assay (Bergamaschi et al, 2011). Flow cytometry and analysis of cell cycle stage were conducted as described using BD‐FACS Canto. Cells were fixed in 70% ethanol, stained for 30 min with 20 μg/ml propidium iodide (PI; Molecular Probe) in Triton‐X (Sigma) in the presence of DNAse‐free RNAse A (Sigma), and PI staining was measured (Bergamaschi et al, 2011; Bhatt et al, 2012).
Adipogenesis was monitored in cells grown in full media and used the adipogenesis assay kit from BioVision (Milpitas, CA) that measures total triglycerides. Absorbance was measured at 570 nm using a BioRad 680 Microplate Reader, and all assays were performed in triplicate.
Tumor data sets and data analysis
Gene lists for generation of the ERβ and RIP140 signature were selected based on the modulation of gene expression in ERβ‐containing cells, and the presence of ERβ‐ and RIP140‐binding sites within 20 kb of the transcription start site of the genes. This gave us a list of 345 genes. We further examined this list using the Oncomine database to find genes that predicted molecular outcome in breast cancer data sets with an odds ratio of at least 2 and P‐value of 0.0001. The log2 median‐centered intensity expression values for signature genes were obtained from the Oncomine database. This analysis yielded a 20‐gene ERβ and RIP140 signature. Hierarchical clustering of data was performed and displayed using Cluster 3.0 and Java TreeView software for analysis and visualization.
For breast cancer patient survival analysis, patients were stratified according to average expression value of the genes in the signature, and the top 30% and bottom 30% of patients were used for computation of Kaplan–Meier curves by the Cox–Mantel log‐rank test and Gehan–Breslow–Wilcoxon test in Graphpad Prism.
Gene expression and ChIP‐Seq data will be available from the Gene Expression Omnibus database with the following accession numbers: meta study, GSE42349; gene expression, GSE42347; and ChIP‐Seq, GSE42348.
This study was supported by NIH grants P50 AT006268 from the National Center for Complementary and Alternative Medicine (NCCAM), the Office of Dietary Supplements (ODS), and the National Cancer Institute (NCI) (BSK); also by NIH R37DK015556 (JAK) and a grant from The Breast Cancer Research Foundation (BSK). ZME received partial support from NIH T32 ES07326. THC was supported by an A*STAR graduate fellowship from The Singapore Agency for Science, Technology and Research. ETL was supported by A*STAR of Singapore and EU‐grant CRESCENDO (FP6‐018652). We thank Luke Petry for technical assistance.
Author contributions: ZME, THC, YJ, ETL, JAK, and BSK conceived and designed the experiments. ZME, THC, and YJ performed the experiments. ZME, THC, YJ, ETL, JAK, and BSK analyzed the data. ZME, THC, YJ, ETL, JAK, and BSK wrote the paper. ZME, JAK, and BSK proofread the manuscript.
Conflict of Interest
The authors declare that they have no conflict of interest.
Supplementary Figures 1‐10 and Supplementary Tables 1‐3 [msb201328-sup-0001.pdf]
Supplementary Table 1
CEAS analysis of transcription factor binding motif enrichment for binding site clusters (C1‐11) identified through SeqMINER analysis. [msb201328-sup-0002.xls]
Supplementary Table 2
SeqPos analysis of transcription factor binding motif enrichment for binding site clusters (C1‐11) identified through SeqMINER analysis. [msb201328-sup-0003.xls]
Supplementary Table 3
Functional enrichment of categories for ERβ‐modulated genes. [msb201328-sup-0004.xls]
This is an open‐access article distributed under the terms of the Creative Commons Attribution License, which permits distribution, and reproduction in any medium, provided the original author and source are credited. This license does not permit commercial exploitation without specific permission.
- Copyright © 2013 EMBO and Macmillan Publishers Limited