Epigenome Mapping Reveals Distinct Modes of Gene Regulation and Widespread Enhancer Reprogramming by the Oncogenic Fusion Protein EWS-FLI1

Reference Epigenome Mapping in a Cellular Model of EWS-FLI1 Dependence

Epigenome mapping is a powerful method for cataloging functional elements throughout the genome (Bernstein et al., 2012), and it can provide insights into the regulatory mechanisms that underlie changes of cell state (Rivera and Ren, 2013). To investigate the effect of EWS-FLI1 expression on epigenetic cell states, we mapped the epigenome of an Ewing sarcoma cell line (A673) that has emerged as a standard model for systems biology in Ewing sarcoma (http://www.ucd.ie/sbi/asset). This cell line is EWS-FLI1 dependent and was previously engineered to harbor a doxycycline-inducible small hairpin RNA against EWS-FLI1 (Carrillo et al., 2007). These cells retain a low level of EWS-FLI1 expression after knockdown, thus allowing us to compare, in an isogenic setting, the epigenomes of EWS-FLI1-high and EWS-FLI1-low cell states without causing extensive cell death (Figure 1A).

We established reference epigenome maps comprising DNA methylation, seven histone marks (H3K4me3, H3K4me1, H3K27ac, H3K56ac, H3K36me, H3K27me3, and H3K9me3), open chromatin states, and RNA levels. DNA methylation was mapped both by whole-genome bisulfite sequencing (WGBS) and reduced representation bisulfite sequencing (RRBS), histone marks by chromatin immunoprecipitation sequencing (ChIP-seq), open chromatin states by assay for transposase-accessible chromatin with high-throughput sequencing (ATAC-seq), and RNA levels by strand-specific RNA sequencing (RNA-seq). Each experiment was performed in two biological replicates each for both the EWS-FLI1-high state (cells in normal growth conditions) and the EWS-FLI1-low state (53 hr after adding doxycycline to induce knockdown of EWS-FLI1). We also performed RNA-seq in cells that were treated with each of three histone deacetylase inhibitors to further dissect the relevance of histone acetylation in our model. In total, we generated 61 datasets and sequenced 2.9 billion reads comprising 221 billion base pairs (Table S1).

To facilitate data access and reuse by other researchers, we developed a website for online exploration and data download (http://tomazou2015.computational-epigenetics.org). Our dataset can be used to view and investigate the epigenetic changes in EWS-FLI1-dependent cells one gene at a time (Figure 1B) or to study patterns of epigenetic and transcriptional deregulation on a more global level, for example by annotating known and novel non-coding transcripts and alternative promoters that are EWS-FLI1 dependent (Figure S1).

Hierarchical clustering of the ChIP-seq and ATAC-seq data showed high consistency between experiments of the same type (Figure 1C). We also observed the expected separation into repressive histone marks (H3K9me3 and H3K27me3), transcription-associated histone marks (H3K36me3), and open-chromatin-associated marks (H3K4me3, H3K4me1, H3K27ac, H3K56ac, and ATAC-seq signal). For the DNA methylation data, which were combined from WGBS and RRBS experiments, we observed high correlation between replicates (Figure 1D) and broad coverage throughout the genome (Figure 1E). The correlation among biological replicates was equally high for the RNA-seq data (Figure 1F), and we were able to identify a small number of transcripts and splicing isoforms that appear to be specific to this cell type (Figures 1G and S1). In aggregate, these observations confirm the high technical data quality of the reference epigenome maps that we generated.

Significant Association of Gene Expression, EWS-FLI1 Binding, and Open Chromatin

EWS-FLI1 is known to act as a transcriptional activator (Bailly et al., 1994; Ohno et al., 1993), but it also represses genes—indirectly by activating transcriptional repressor genes and potentially also directly by recruiting repressive protein complexes (Hahm et al., 1999; Prieur et al., 2004; Sankar et al., 2013). To investigate the genome-wide association among gene expression, EWS-FLI1 binding, and chromatin states, we stratified all genes by their RNA levels (Figure 2A, top) and by the distance from the transcription start site (TSS) to the nearest EWS-FLI1 binding peaks (Figure 2A, bottom), which were mapped previously (Bilke et al., 2013), and we plotted average ChIP-seq intensities for each histone mark around the TSS (Figures 2A and S2). Promoters of highly expressed genes and those bound by EWS-FLI1 had high levels of the open-chromatin-associated marks H3K4me3, H3K27ac, and H3K56ac. These marks also showed the characteristic dip at the TSS that is indicative of nucleosome-free regions (ENCODE Project Consortium, 2012), and they overlapped with a peak of open chromatin in the ATAC-seq data. H3K4me1 was similarly enriched in the wider promoter region but depleted in immediate vicinity of the TSS, reflecting its mutual exclusivity with H3K4me3. Furthermore, H3K36me3 levels were high in genic regions directly downstream of the TSS of genes that were highly expressed and bound by EWS-FLI1, whereas the repressive histone marks H3K9me3 and H3K27me3 were depleted in the promoters of these genes (Figures 2A and S2). Promoters with distal EWS-FLI binding (5–40 kb from the TSS) were no more enriched for the open-chromatin marks H3K4me3, H3K27ac, and H3K56ac than those showing no EWS-FLI1 binding within 40 kb from the TSS, but they were more strongly depleted for the repressive-histone marks H3K27me3 and H3K9me3.

Combining gene expression data with EWS-FLI1 binding, we observed that proximal and distal binding of EWS-FLI1 were both associated with higher RNA levels (Figure 2B), suggesting that EWS-FLI1 acts as a transcriptional activator not only when it binds in direct vicinity to the TSS but also at enhancers outside of the promoter region. Proximal binding sites of EWS-FLI1 carried the characteristic histone patterns of promoter regions (Figure 2C), whereas distal EWS-FLI1 binding sites showed the chromatin signature of active enhancers with high levels of H3K27ac and H3K4me1, low levels of H3K4me3, and a dip in H3K27ac levels at the binding site (Figure 2D). H3K27ac was consistently associated with EWS-FLI1 binding, both for promoter regions and for distal enhancer elements. To identify DNA methylation patterns associated with EWS-FLI1 binding, we thus compared DNA methylation levels of H3K27ac peaks with and without EWS-FLI1 binding, and we observed significantly lower DNA methylation levels among the EWS-FLI1-bound regions (Figure 2E). This result is consistent with recent evidence that active transcription factors can reduce DNA methylation at their binding sites (Stadler et al., 2011).

Collectively, our results support a strong association of EWS-FLI1 binding with high gene expression, open chromatin, nucleosome depletion, and reduced DNA methylation. These patterns hold true for both promoters and enhancers and for both proximal and distal binding. We did not observe widespread co-localization of EWS-FLI1 with repressive chromatin, suggesting that the repressive effect of EWS-FLI1 observed for a sizable fraction of its target genes (Figures 3A and 3B) is either indirect or caused by the depletion of active histone marks.

Chromatin-Based Clusters of EWS-FLI1-Dependent Transcripts

Several studies mapped EWS-FLI1 target genes using expression microarrays and/or ChIP-seq analysis of EWS-FLI1 binding (Bilke et al., 2013; Kauer et al., 2009; Prieur et al., 2004; Riggi et al., 2008; Smith et al., 2006; Tirode et al., 2007). Both approaches are useful but limited in their insight into the underlying regulatory mechanisms. We hypothesized that epigenome maps could help identify distinct modes of transcriptional regulation, and we devised a bioinformatic approach to classify genes in a way that takes biological knowledge about chromatin regulation into account (Figure S3; Supplemental Experimental Procedures).

First, we identified all genes that were significantly downregulated or upregulated upon EWS-FLI1 knockdown, according to our RNA-seq data (Figure 3A). To emphasize that the regulatory relationship with EWS-FLI1 may be direct or indirect, we refer to these genes as “EWS-FLI correlated” and “EWS-FLI1 anticorrelated” rather than “EWS-FLI1 activated” and “EWS-FLI1 repressed.” In total, 1,287 transcripts were significantly correlated with EWS-FLI1 (i.e., more highly expressed in the EWS-FL1-high state than in the EWS-FLI1-low state, which is equivalent to downregulation upon EWS-FLI1 knockdown), and 1,446 transcripts were significantly anticorrelated (i.e., more highly expressed in the EWS-FLI1-low state, which is equivalent to upregulation upon EWS-FLI1 knockdown) (Figure 3B).

In a second step, we annotated all correlated and anticorrelated transcripts with ChIP-seq intensity levels in their promoter regions, focusing on the five histone marks that are generally associated with gene promoters (Figure 3C). We then clustered the transcripts based on their histone marks in the state where they are more lowly expressed (correlated transcripts: EWS-FLI1-low state after EWS-FLI1 knockdown; anticorrelated transcripts: EWS-FLI1-high state before EWS-FLI1 knockdown). To that end, an initial step of unsupervised clustering was followed by expert classification of selected transcripts and supervised prediction of cluster membership for all transcripts (Figure S3). Cluster 1 comprises transcripts carrying all four active-promoter marks (H3K4me3, H3K4me1, H3K27ac, and H3K56ac), but not the repressive H3K27me3 mark; cluster 2 transcripts carry both H3K4me3 and H3K4me1 but lack the histone-acetylation marks; cluster 3 transcripts are marked by H3K4me1 only; and cluster 4 is characterized by the repressive mark H3K27me3 in the presence or absence of H3K4 methylation. Using a logistic regression model for predicting cluster membership, most correlated and anticorrelated transcripts (82% and 80%) were placed unambiguously in one of the four clusters (Figures 3D and 3E). Examples of EWS-FLI1-anticorrelated transcripts for each of the four clusters are shown in Figures 3F–3I.

Functional Characteristics and Chromatin Dynamics of EWS-FLI1-Regulated Transcripts

To validate the biological relevance of our chromatin-based gene clusters, we analyzed gene set enrichment and overlap of their promoters with catalogs of regulatory elements. Each cluster showed specific enrichment patterns, which were similar among EWS-FLI1-correlated and anticorrelated transcripts. Cluster 1 promoters were highly enriched for polymerase binding and transcription initiation activity in a broad range of cell types (Figure 4A), and they were also enriched for CpG islands (Figures 3D and 3E), which tend to associate with strong and widely active promoters. Cluster 2 promoters had lower overlap with CpG islands and a strong association with tissue-specific patterns of open versus repressed chromatin and with binding of tissue-specific transcription factors. Cluster 3 comprised tissue-specific genes that lack CpG islands. Finally, cluster 4 transcripts were strongly enriched for Polycomb repressive complex 2 binding across many cell types.

In addition to the enrichment patterns that were shared between EWS-FLI1-correlated and anticorrelated clusters, we also observed interesting differences. For example, EWS-FLI1 promoter binding was more common for EWS-FLI1-correlated transcripts (i.e., those that are downregulated after knockdown of EWS-FLI1) than for anticorrelated transcripts (Figures 3D, 3E, and 4A). Furthermore, when comparing each EWS-FLI1-correlated cluster to the corresponding EWS-FLI1-anticorrelated cluster (Figure 4B), we observed a striking enrichment for binding of E2F transcription factors among the EWS-FLI1-correlated transcripts, reinforcing a previously reported link between EWS-FLI1 and E2F (Bilke et al., 2013). The EWS-FLI1-correlated transcripts were also enriched for proliferation gene sets, consistent with the ability of EWS-FLI1 to accelerate cell growth. In contrast, anticorrelated transcripts were enriched for gene sets that are characteristic of adult stem cells and for binding of both components of the AP-1 transcription factor, FOS and JUN.

The chromatin-based gene clusters not only showed distinct patters of functional enrichment but also responded differently to EWS-FLI1 knockdown (Figure 4C). Among the EWS-FLI1-correlated transcripts, H3K27ac was most strongly reduced at cluster 1 promoters after EWS-FLI1 knockdown, whereas H3K4me3 and H3K4me1 signals were largely retained even in the EWS-FLI1-low state. In contrast, clusters 2 and 3 showed reduced levels for all three marks (H3K4me3, H3K4me1, and H3K27ac), and cluster 4 promoters underwent relatively minor changes upon EWS-FLI1 knockdown. The EWS-FLI1-anticorrelated clusters showed a more uniform response across clusters with increases in H3K4me3 and H3K27ac (Figure 4D). Together, these observations indicate that EWS-FLI1 expression promotes transcription of its target genes in at least two distinct ways: first, by further increasing H3K27ac levels and transcription of genes that are already widely expressed in proliferating cells (correlated cluster 1); and second, by establishing or maintaining the chromatin signature of active promoters (H3K4me3, H3K4me1, and H3K27ac) at cell type specific genes that would otherwise be silent in these cells (correlated clusters 2 and 3). Furthermore, the relatively uniform response among the anticorrelated clusters indicates that EWS-FLI1 mediated repression involves mechanisms that decrease H3K27ac in similar ways across all four EWS-FLI1 anticorrelated clusters.

Based on these promoter-centric analyses, H3K27ac and, to a lesser degree, H3K4me3 and H3K4me1 emerged as the histone marks that were most strongly affected by EWS-FLI1 knockdown. This result was corroborated by western blots for all seven histone marks in additional Ewing sarcoma cell lines (TC252, SK-N-MC, and ST-ET-7.2), where we observed a global reduction in H3K4me3 and H3K27ac levels upon knockdown of EWS-FLI1 (Figure 4E). These data indicate that the global and EWS-FLI1-dependent increase of active chromatin may represent a key mechanism of transcriptional rewiring in Ewing sarcoma cells.

Cluster-Specific Transcriptional Response to Treatment with Histone Deacetylase Inhibitors

Given that H3K27ac was most variable between the chromatin-based gene clusters and also most dynamic upon EWS-FLI1 knockdown (Figures 3 and 4), we hypothesized that the clustered genes should respond differently to induced global changes of acetylation, which would support the clusters’ usefulness for identifying distinct modes of gene regulation among EWS-FLI1 target genes. To test this hypothesis, we treated the cells with histone acetyltransferase (HAT) inhibitors and histone deacetylase (HDAC) inhibitors. These experiments were done for both the EWS-FLI1-high and the EWS-FLI1-low cell state, and they used drug concentrations that were low enough to avoid broad toxicity. None of the three HAT inhibitors that we tested (C646, CPTH2, and anacardic acid) had an effect on cell survival or gene expression at the chosen concentrations (data not shown); hence, we focused our analysis primarily on the HDAC inhibitors (Figures 5 and S4).

HDAC inhibitors are being evaluated as drug candidates for Ewing sarcoma (Arnaldez and Helman, 2014), but here we used them solely as a chemical biology tool for validating our gene clustering, and at much lower concentrations than if our goal were to selectively kill EWS-FLI1-expressing cells. For two out of three HDAC inhibitors tested—namely, romidepsin (FK-228) and entinostat (MS-275), but not Vorinostat (SAHA)—we were able to identify concentrations that fulfilled all three criteria that were crucial for our experiment: (1) no effect on cell viability (Figures 5A and S4A), (2) global increase in histone acetylation levels (Figure 5B), and (3) no effect on EWS-FLI1 expression levels as a potential confounder (Figure 5B).

We performed RNA-seq to measure the transcriptome response to each HDAC inhibitor at the chosen concentrations (Figures 5C–5E, S4B, and S4C). Both romidepsin and entinostat specifically upregulated the expression of EWS-FLI1-anticorrelated transcripts in cells with high EWS-FLI1 levels, suggesting that the (direct or indirect) repressive effect of EWS-FLI1 is partially reversed by treating EWS-FLI1-high cells with HDAC inhibitors. Among the anticorrelated transcripts, cluster 1 was least affected by treatment with HDAC inhibitors, most likely because cluster 1 promoters already carried high levels of H3K27ac and may therefore be less sensitive to HDAC inhibition (Figures 5D and S4C). In contrast, clusters 2, 3, and 4 carried low levels of H3K27ac in the EWS-FLI1-high state and responded more strongly to treatment with HDAC inhibitors. Examples for each of the four clusters are shown in Figure 5F.

Among the EWS-FLI1-correlated genes, there was no global trend in either direction (Figure 5C), but we observed strong cluster-specific differences in the EWS-FLI1-low state (Figures 5E and S4B). Clusters 2 and 4 were consistently upregulated in response to HDAC inhibition, whereas no such trend was observed for cluster 3, and cluster 1 genes were even slightly downregulated. Specific examples are shown in Figure 5G. These results indicate that EWS-FLI1 activates certain genes in cluster 2 and 4 by interfering with HDAC function, and treatment with HDAC inhibitors appears to recapitulate this effect in EWS-FLI1-low cells.

Widespread Reprogramming of Enhancers and Super-Enhancers by EWS-FLI1

Our promoter-centric analysis (Figures 3 and 4) established H3K4me3 and H3K27ac as the two histone marks that were most strongly associated with EWS-FLI1 expression. H3K4me3 is a prototypic promoter mark, and the observed changes in H3K4me3 were indeed located almost exclusively in promoter regions. By contrast, not only does H3K27ac occur in promoter regions (where it tends to overlap with H3K4me3), but it is also a defining mark of active enhancers when occurring in the absence of H3K4me3 (Creyghton et al., 2010; Rada-Iglesias et al., 2011). We therefore investigated the effect of EWS-FLI1 expression on H3K27ac beyond gene promoters, comparing H3K27ac patterns before and after knockdown of EWS-FLI1 across the genome. In total, 15,300 H3K27ac peaks showed significantly lower intensity in EWS-FLI1-low cells (“EWS-FLI1 correlated peaks”), and 18,727 H3K27ac peaks showed significantly higher intensity (“EWS-FLI1 anticorrelated peaks”). 27% of correlated peaks and only 6% of anticorrelated peaks overlapped with promoter regions (Figures 6A and S5A), indicating that the majority of EWS-FLI1-associated H3K27ac peaks are located at enhancer elements outside of promoter regions.

We observed striking differences between EWS-FLI1-correlated and anticorrelated H3K27ac peaks, which became apparent when we stratified all peaks by the strength of H3K27ac signal difference upon EWS-FLI1 knockdown. Among EWS-FLI1-correlated enhancers, only 10% of the top 20% most differential H3K27ac peaks were located in promoters, whereas 62% of the bottom 20% peaks (which are still significantly different for H3K27ac between the EWS-FLI1-high and low states) were located in promoters (Figure 6B). No such differences were observed for EWS-FLI1-anticorrelated peaks, of which 5% and 7% were located in promoters. We also identified strong differences in EWS-FLI1 binding, with 77% of the top hundred EWS-FLI1-correlated H3K27ac peaks being EWS-FLI bound, while the percentage for EWS-FLI1 binding fell to 20% when calculated across all EWS-FLI1-correlated peaks and to only 1% for EWS-FLI1-anticorrelated peaks (Figure 6C). Finally, for the other histone marks tested we observed much lower overlap of EWS-FLI1 (Figure S5B), indicating that EWS-FLI1 primarily drives H3K27ac at selected promoters and enhancers throughout the genome.

Our observation that the strongest EWS-FLI1-correlated enhancers were qualitatively different from other H3K27ac peaks is reminiscent of super-enhancers (Hnisz et al., 2013; Lovén et al., 2013; Whyte et al., 2013). We therefore investigated whether super-enhancers occur in EWS-FLI1 expressing cells and how they might respond to EWS-FLI1 knockdown. Using a published method for H3K27ac-based annotation of super-enhancers (Lovén et al., 2013; Whyte et al., 2013), we identified 697 super-enhancers in the EWS-FLI1-high state and annotated them with their neighboring genes (http://tomazou2015.computational-epigenetics.org). Overall, more than 10% of H3K27ac peaks in EWS-FLI1-high cells were located in super-enhancers (Figure S5C), and super-enhancers had much higher cumulative H3K27ac enrichment than typical enhancers (Figure 6D). This strong enrichment of H3K27ac was largely due to the exceeding length of super-enhancers (median length of 27,818 base pairs as compared to 1,196 base pairs for typical enhancers; Figure S5D) and not due to higher H3K27ac signal intensity per base pair (Figure S5E). As illustrated in Figure 6E, super-enhancers were essentially agglomerates of co-located H3K27ac peaks over a range of several dozen kilobases, while typical enhancers rarely comprised more than two H3K27ac peaks and were typically shorter than 2 kb. Several of the top super-enhancers were located near genes with known roles in cell proliferation, differentiation, and apoptosis. We also observed widespread loss of H3K27ac upon knockdown of EWS-FLI1 for some super-enhancers (Figure 6E), while others displayed more localized changes affecting only those H3K27ac peaks that overlapped with EWS-FLI1 binding sites (Figure S5F). Interestingly, EWS-FLI1 binding sites and EWS-FLI1-correlated H3K27ac peaks were both more likely than other H3K27ac peaks to be located in super-enhancers (Figure S5G); in addition, more than half of the super-enhancers in EWS-FLI1-expressing cells contained at least one EWS-FLI1 binding event, and a similar fraction contained at least one EWS-FLI1-correlated H3K27ac peak (Figure S5H).

To assess the cell-type specificity of EWS-FLI1-correlated H3K27ac peaks, we studied their association with active regulatory elements in 72 human cell types, based on DNase hypersensitivity data (Sheffield et al., 2013). Our analysis uncovered a striking trend toward high cell-type specificity among the strongest EWS-FLI1-correlated, EWS-FLI1-bound, and super-enhancer-associated H3K27ac peaks (Figure 6F). This trend was already visible for H3K27ac peaks in promoters and in typical enhancers, but the association with super-enhancers and overlap with EWS-FLI1 binding sites further increased the cell-type specificity (Figure 6G). Most notably, 37% of the 2,322 strongest EWS-FLI1-correlated H3K27ac enhancer peaks overlapped with DNase-hypersensitive sites present only in the one EWS-FLI1 expressing cell line that was part of the DNase dataset (SK-N-MC) and were thus unique to Ewing sarcoma. In summary, these results provide strong evidence that EWS-FLI1 not only regulates genes and promoters but also establishes highly cell-type-specific enhancers and super-enhancers in EWS-FLI1-expressing cells.

Global Epigenome Dynamics in EWS-FLI1-Dependent Cells

In this study, we used a Ewing sarcoma cell line with tunable EWS-FLI1 expression to uncover connections between this oncogenic fusion protein and the epigenome. Based on both ChIP-seq and western blots, H3K27ac emerged as the mark that was most strongly affected by EWS-FLI1 knockdown. 33,170 H3K27ac peaks were significantly altered throughout the genome (Figure 7A), and 80% of the most strongly EWS-FLI1-correlated H3K27ac peaks were bound by EWS-FLI1 (Figure 6C). In contrast, H3K27ac levels in regions that were not bound by EWS-FLI1 were either correlated or anticorrelated with EWS-FLI1 expression and likely driven by indirect mechanisms of transcription regulation. Similar but weaker patterns of EWS-FLI1 association were also observed for H3K4me3 and H3K4me1 (Figure 7A), while the other studied histone marks did not show widespread changes upon EWS-FLI1 knockdown. Differences in response to EWS-FLI1 knockdown were also notably absent for DNA methylation (Figure 7B). Although this mark was anticorrelated with EWS-FLI1 binding (Figure 2E), no changes in DNA methylation patterns were observed 53 hr after knockdown—at a time point that should be sufficient for de novo methylation to have at least started. Hence, it seems possible that Ewing sarcoma cells retain an epigenetic memory of prior EWS-FLI1 binding as part of their DNA methylation patterns.

Promoter Regulation by EWS-FLI1

Chromatin-based promoter analysis identified four clusters of EWS-FLI1-correlated transcripts and four clusters of EWS-FLI1-anticorrelated transcripts. These clusters define three different regulatory modes (Figure 4). First, EWS-FLI1-correlated cluster 1 genes, such as the cell-cycle regulator and proto-oncogene CCND1 (Sanchez et al., 2008), were expressed even in the EWS-FLI1-low state and carried all histone marks of active promoters, but their H3K27ac and transcription levels were higher in the EWS-FLI1-high state (Figure 7C). Many of the correlated cluster 1 genes were associated with cell proliferation, suggesting that increased H3K27ac at widely expressed genes may drive the rapid proliferation of EWS-FLI1-expressing cells. Second, the genes in EWS-FLI1-correlated clusters 2 and 3 (Figure 4A) were more tissue specific than genes in correlated cluster 1, and they lacked active promoter marks such as H3K27ac and/or H3K4me3 in the EWS-FLI1-low state. EWS-FLI1-dependent establishment of the missing histone marks may aberrantly activate these genes, as illustrated by the hematopoietic proto-oncogene MYB that was specifically expressed in the EWS-FLI1-high state (Figure 7D). Third, EWS-FLI1-anticorrelated genes had higher H3K27ac levels in the EWS-FLI1-low state independent of their cluster association and H3K27ac level in the EWS-FLI-high state (Figure 4D), which is illustrated by the AP-1 component FOSL1 (Figure 7E). Intriguingly, we observed broad enrichment for AP-1 binding in other cell types among EWS-FLI1-anticorrelated gene promoters (Figure 4B), suggesting that EWS-FLI1 may repress a sizable fraction of its target genes by interfering with the activating role of AP-1. HDACs seem to play a role in the process, given the upregulation of EWS-FLI1-anticorrelated genes (clusters 2, 3, and 4) upon treatment with HDAC inhibitors specifically in the EWS-FLI1-high state (Figure 5D).

Enhancer and Super-Enhancer Regulation by EWS-FLI1

Our analysis of H3K27ac resulted in a genome-wide catalog of EWS-FLI1-associated enhancers and super-enhancers. Strikingly, the most dynamic enhancers upon EWS-FLI1 knockdown were also the most cell-type specific (Figure 6), which is illustrated by the absence of any ENCODE signal for H3K27ac at the EWS-FLI1-bound enhancers upstream of the CCND1 promoter (Figure 7C). A large percentage of the observed EWS-FLI1-correlated H3K27ac peaks overlapped exclusively with regulatory elements in SK-N-MC cells, which is the one EWS-FLI1-expressing cell line in a published catalog with DNase hypersensitivity data for 72 human cell types (p < 10⁻¹⁰⁰, Fisher’s exact test). EWS-FLI1 also appears to contribute to cell-type-specific enhancer signatures by converting broadly active enhancers to a poised state, as illustrated by the enhancer element upstream of the FOSL1 locus (Figure 7E). While our study was in revision, Riggi et al. published conclusive evidence that EWS-FLI1 binds cis-regulating elements shared by Ewing sarcoma cell lines and primary tumors (Riggi et al., 2014), suggesting that our observations are not only limited to cell lines but also relevant for Ewing sarcoma biology in patients.

In summary, the establishment and analysis of reference epigenomes for EWS-FLI1-high and EWS-FLI1-low cells not only constitutes a first comprehensive epigenome dataset for a solid cancer that is driven by an oncogenic fusion protein, but it also provides relevant insights into the mechanisms by which a single oncogenic event has shaped a conductive epigenome. The dataset is available for interactive exploration and data download (http://tomazou2015.computational-epigenetics.org), thus providing a genome-wide reference that can guide future work on regulatory mechanisms and epigenetic drug effects in Ewing sarcoma.

Epigenome Mapping

A673 cells were cultured as previously described (Carrillo et al., 2007). DNA, chromatin, RNA, and protein were collected following standard procedures from untreated cells (EWS-FLI1-high state) and 53 hr after adding doxycycline to the media (EWS-FLI1-low state). To assess DNA methylation, WGBS and RRBS experiments were performed using custom protocols (see the Supplemental Experimental Procedures for details), and the data of both assays were merged. ChIP-seq experiments were performed using the iDeal ChIP-seq kit (Diagenode) according to manufacturer instructions. ATAC-seq was performed as previously described (Buenrostro et al., 2013) with minor adaptations for A673 cells. RNA-seq used Illumina kits for strand-specific library preparation. All sequencing was performed by the Biomedical Sequencing Facility at CeMM using Illumina HiSeq 2000 sequencers.

The epigenome maps as well as the raw and processed data underlying the presented analyses are available online at http://tomazou2015.computational-epigenetics.org.

Epigenetic Gene Clustering

We clustered all EWS-FLI1-associated genes based on the ChIP-seq signal intensity of five promoter-associated histone marks in the vicinity of the TSS (excluding H3K36me3 and H3K9me3, which were rarely present in promoter regions), using a two-step semi-supervised clustering method. Multinomial logistic regression models were trained based on expert-curated examples and used to classify EWS-FLI1-regulated genes derived from the RNA-seq data. We trained separate classifiers for identifying clusters among the EWS-FLI1-correlated genes and among the EWS-FLI1-anticorrelated genes.

Enrichment Analysis

To functionally annotate the gene clusters, we identified significant overlap between gene sets derived from our dataset and gene sets obtained from several public databases (see the Supplemental Experimental Procedures for details). In addition to gene-based enrichments, we also developed a method for region-based enrichment analysis called location overlap analysis (LOLA; http://lola.computational-epigenetics.org), and we identified sets of functional elements that significantly co-localized with the promoters of a given gene cluster. We used Fisher’s exact test to establish significance, and odds ratios as well as p values were used to rank the results. Comparisons were done (1) against the background of all genes identified by our RNA-seq analysis to identify enrichment patterns among EWS-FLI1-regulated genes, (2) against the background of genes in other clusters to identify cluster-specific enrichment patterns separately for EWS-FLI1-correlated and EWS-FLI1-anticorrelated genes, and (3) against the background of genes in the same cluster but with opposite direction of EWS-FLI1 regulation (correlated or anticorrelated) to identify upregulation- versus downregulation-specific enrichment patterns separately for each cluster.

Document S1. Supplemental Experimental Procedures, Figures S1–S5

Table S1. Epigenome Sequencing Statistics, Related to Figure 1

Document S2. Article plus Supplemental Information

1. Arnaldez F.I., Helman L.J. New strategies in ewing sarcoma: lost in translation? Clin. Cancer Res. 2014;20:3050–3056. doi: 10.1158/1078-0432.CCR-13-0633. [DOI] [PMC free article] [PubMed] [Google Scholar]
2. Bailly R.A., Bosselut R., Zucman J., Cormier F., Delattre O., Roussel M., Thomas G., Ghysdael J. DNA-binding and transcriptional activation properties of the EWS-FLI-1 fusion protein resulting from the t(11;22) translocation in Ewing sarcoma. Mol. Cell. Biol. 1994;14:3230–3241. doi: 10.1128/mcb.14.5.3230. [DOI] [PMC free article] [PubMed] [Google Scholar]
3. Bernstein M., Kovar H., Paulussen M., Randall R.L., Schuck A., Teot L.A., Juergens H. Ewing’s sarcoma family of tumors: current management. Oncologist. 2006;11:503–519. doi: 10.1634/theoncologist.11-5-503. [DOI] [PubMed] [Google Scholar]
4. Bernstein B.E., Birney E., Dunham I., Green E.D., Gunter C., Snyder M., ENCODE Project Consortium An integrated encyclopedia of DNA elements in the human genome. Nature. 2012;489:57–74. doi: 10.1038/nature11247. [DOI] [PMC free article] [PubMed] [Google Scholar]
5. Bilke S., Schwentner R., Yang F., Kauer M., Jug G., Walker R.L., Davis S., Zhu Y.J., Pineda M., Meltzer P.S., Kovar H. Oncogenic ETS fusions deregulate E2F3 target genes in Ewing sarcoma and prostate cancer. Genome Res. 2013;23:1797–1809. doi: 10.1101/gr.151340.112. [DOI] [PMC free article] [PubMed] [Google Scholar]
6. Brohl A.S., Solomon D.A., Chang W., Wang J., Song Y., Sindiri S., Patidar R., Hurd L., Chen L., Shern J.F. The genomic landscape of the Ewing Sarcoma family of tumors reveals recurrent STAG2 mutation. PLoS Genet. 2014;10:e1004475. doi: 10.1371/journal.pgen.1004475. [DOI] [PMC free article] [PubMed] [Google Scholar]
7. Buenrostro J.D., Giresi P.G., Zaba L.C., Chang H.Y., Greenleaf W.J. Transposition of native chromatin for fast and sensitive epigenomic profiling of open chromatin, DNA-binding proteins and nucleosome position. Nat. Methods. 2013;10:1213–1218. doi: 10.1038/nmeth.2688. [DOI] [PMC free article] [PubMed] [Google Scholar]
8. Carrillo J., García-Aragoncillo E., Azorín D., Agra N., Sastre A., González-Mediero I., García-Miguel P., Pestaña A., Gallego S., Segura D., Alonso J. Cholecystokinin down-regulation by RNA interference impairs Ewing tumor growth. Clin. Cancer Res. 2007;13:2429–2440. doi: 10.1158/1078-0432.CCR-06-1762. [DOI] [PubMed] [Google Scholar]
9. Chen J., Odenike O., Rowley J.D. Leukaemogenesis: more than mutant genes. Nat. Rev. Cancer. 2010;10:23–36. doi: 10.1038/nrc2765. [DOI] [PMC free article] [PubMed] [Google Scholar]
10. Creyghton M.P., Cheng A.W., Welstead G.G., Kooistra T., Carey B.W., Steine E.J., Hanna J., Lodato M.A., Frampton G.M., Sharp P.A. Histone H3K27ac separates active from poised enhancers and predicts developmental state. Proc. Natl. Acad. Sci. USA. 2010;107:21931–21936. doi: 10.1073/pnas.1016071107. [DOI] [PMC free article] [PubMed] [Google Scholar]
11. Crompton B.D., Stewart C., Taylor-Weiner A., Alexe G., Kurek K.C., Calicchio M.L., Kiezun A., Carter S.L., Shukla S.A., Mehta S.S. The genomic landscape of pediatric Ewing sarcoma. Cancer Discov. 2014;4:1326–1341. doi: 10.1158/2159-8290.CD-13-1037. [DOI] [PubMed] [Google Scholar]
12. Delattre O., Zucman J., Melot T., Garau X.S., Zucker J.M., Lenoir G.M., Ambros P.F., Sheer D., Turc-Carel C., Triche T.J. The Ewing family of tumors—a subgroup of small-round-cell tumors defined by specific chimeric transcripts. N. Engl. J. Med. 1994;331:294–299. doi: 10.1056/NEJM199408043310503. [DOI] [PubMed] [Google Scholar]
13. ENCODE Project Consortium An integrated encyclopedia of DNA elements in the human genome. Nature. 2012;489:57–74. doi: 10.1038/nature11247. [DOI] [PMC free article] [PubMed] [Google Scholar]
14. Hahm K.B., Cho K., Lee C., Im Y.H., Chang J., Choi S.G., Sorensen P.H., Thiele C.J., Kim S.J. Repression of the gene encoding the TGF-beta type II receptor is a major target of the EWS-FLI1 oncoprotein. Nat. Genet. 1999;23:222–227. doi: 10.1038/13854. [DOI] [PubMed] [Google Scholar]
15. Hnisz D., Abraham B.J., Lee T.I., Lau A., Saint-André V., Sigova A.A., Hoke H.A., Young R.A. Super-enhancers in the control of cell identity and disease. Cell. 2013;155:934–947. doi: 10.1016/j.cell.2013.09.053. [DOI] [PMC free article] [PubMed] [Google Scholar]
16. Huether R., Dong L., Chen X., Wu G., Parker M., Wei L., Ma J., Edmonson M.N., Hedlund E.K., Rusch M.C. The landscape of somatic mutations in epigenetic regulators across 1,000 paediatric cancer genomes. Nat. Commun. 2014;5:3630. doi: 10.1038/ncomms4630. [DOI] [PMC free article] [PubMed] [Google Scholar]
17. Kauer M., Ban J., Kofler R., Walker B., Davis S., Meltzer P., Kovar H. A molecular function map of Ewing’s sarcoma. PLoS ONE. 2009;4:e5415. doi: 10.1371/journal.pone.0005415. [DOI] [PMC free article] [PubMed] [Google Scholar]
18. Kovar H. Downstream EWS/FLI1 - upstream Ewing’s sarcoma. Genome Med. 2010;2:8. doi: 10.1186/gm129. [DOI] [PMC free article] [PubMed] [Google Scholar]
19. Lawrence M.S., Stojanov P., Polak P., Kryukov G.V., Cibulskis K., Sivachenko A., Carter S.L., Stewart C., Mermel C.H., Roberts S.A. Mutational heterogeneity in cancer and the search for new cancer-associated genes. Nature. 2013;499:214–218. doi: 10.1038/nature12213. [DOI] [PMC free article] [PubMed] [Google Scholar]
20. Lovén J., Hoke H.A., Lin C.Y., Lau A., Orlando D.A., Vakoc C.R., Bradner J.E., Lee T.I., Young R.A. Selective inhibition of tumor oncogenes by disruption of super-enhancers. Cell. 2013;153:320–334. doi: 10.1016/j.cell.2013.03.036. [DOI] [PMC free article] [PubMed] [Google Scholar]
21. Macconaill L.E., Garraway L.A. Clinical implications of the cancer genome. J. Clin. Oncol. 2010;28:5219–5228. doi: 10.1200/JCO.2009.27.4944. [DOI] [PMC free article] [PubMed] [Google Scholar]
22. Mitelman F., Johansson B., Mertens F. Fusion genes and rearranged genes as a linear function of chromosome aberrations in cancer. Nat. Genet. 2004;36:331–334. doi: 10.1038/ng1335. [DOI] [PubMed] [Google Scholar]
23. Ohno T., Rao V.N., Reddy E.S. EWS/Fli-1 chimeric protein is a transcriptional activator. Cancer Res. 1993;53:5859–5863. [PubMed] [Google Scholar]
24. Prieur A., Tirode F., Cohen P., Delattre O. EWS/FLI-1 silencing and gene profiling of Ewing cells reveal downstream oncogenic pathways and a crucial role for repression of insulin-like growth factor binding protein 3. Mol. Cell. Biol. 2004;24:7275–7283. doi: 10.1128/MCB.24.16.7275-7283.2004. [DOI] [PMC free article] [PubMed] [Google Scholar]
25. Rada-Iglesias A., Bajpai R., Swigut T., Brugmann S.A., Flynn R.A., Wysocka J. A unique chromatin signature uncovers early developmental enhancers in humans. Nature. 2011;470:279–283. doi: 10.1038/nature09692. [DOI] [PMC free article] [PubMed] [Google Scholar]
26. Riggi N., Suvà M.L., Suvà D., Cironi L., Provero P., Tercier S., Joseph J.M., Stehle J.C., Baumer K., Kindler V., Stamenkovic I. EWS-FLI-1 expression triggers a Ewing’s sarcoma initiation program in primary human mesenchymal stem cells. Cancer Res. 2008;68:2176–2185. doi: 10.1158/0008-5472.CAN-07-1761. [DOI] [PubMed] [Google Scholar]
27. Riggi N., Knoechel B., Gillespie S.M., Rheinbay E., Boulay G., Suvà M.L., Rossetti N.E., Boonseng W.E., Oksuz O., Cook E.B. EWS-FLI1 utilizes divergent chromatin remodeling mechanisms to directly activate or repress enhancer elements in Ewing sarcoma. Cancer Cell. 2014;26:668–681. doi: 10.1016/j.ccell.2014.10.004. [DOI] [PMC free article] [PubMed] [Google Scholar]
28. Rivera C.M., Ren B. Mapping human epigenomes. Cell. 2013;155:39–55. doi: 10.1016/j.cell.2013.09.011. [DOI] [PMC free article] [PubMed] [Google Scholar]
29. Sanchez G., Bittencourt D., Laud K., Barbier J., Delattre O., Auboeuf D., Dutertre M. Alteration of cyclin D1 transcript elongation by a mutated transcription factor up-regulates the oncogenic D1b splice isoform in cancer. Proc. Natl. Acad. Sci. USA. 2008;105:6004–6009. doi: 10.1073/pnas.0710748105. [DOI] [PMC free article] [PubMed] [Google Scholar]
30. Sankar S., Bell R., Stephens B., Zhuo R., Sharma S., Bearss D.J., Lessnick S.L. Mechanism and relevance of EWS/FLI-mediated transcriptional repression in Ewing sarcoma. Oncogene. 2013;32:5089–5100. doi: 10.1038/onc.2012.525. [DOI] [PMC free article] [PubMed] [Google Scholar]
31. Sheffield N.C., Thurman R.E., Song L., Safi A., Stamatoyannopoulos J.A., Lenhard B., Crawford G.E., Furey T.S. Patterns of regulatory activity across diverse human cell types predict tissue identity, transcription factor binding, and long-range interactions. Genome Res. 2013;23:777–788. doi: 10.1101/gr.152140.112. [DOI] [PMC free article] [PubMed] [Google Scholar]
32. Smith R., Owen L.A., Trem D.J., Wong J.S., Whangbo J.S., Golub T.R., Lessnick S.L. Expression profiling of EWS/FLI identifies NKX2.2 as a critical target gene in Ewing’s sarcoma. Cancer Cell. 2006;9:405–416. doi: 10.1016/j.ccr.2006.04.004. [DOI] [PubMed] [Google Scholar]
33. Stadler M.B., Murr R., Burger L., Ivanek R., Lienert F., Schöler A., van Nimwegen E., Wirbelauer C., Oakeley E.J., Gaidatzis D. DNA-binding factors shape the mouse methylome at distal regulatory regions. Nature. 2011;480:490–495. doi: 10.1038/nature10716. [DOI] [PubMed] [Google Scholar]
34. Tirode F., Laud-Duval K., Prieur A., Delorme B., Charbord P., Delattre O. Mesenchymal stem cell features of Ewing tumors. Cancer Cell. 2007;11:421–429. doi: 10.1016/j.ccr.2007.02.027. [DOI] [PubMed] [Google Scholar]
35. Tirode F., Surdez D., Ma X., Parker M., Le Deley M.C., Bahrami A., Zhang Z., Lapouble E., Grossetête-Lalami S., Rusch M., St. Jude Children’s Research Hospital–Washington University Pediatric Cancer Genome Project and the International Cancer Genome Consortium Genomic landscape of Ewing sarcoma defines an aggressive subtype with co-association of STAG2 and TP53 mutations. Cancer Discov. 2014;4:1342–1353. doi: 10.1158/2159-8290.CD-14-0622. [DOI] [PMC free article] [PubMed] [Google Scholar]
36. Whyte W.A., Orlando D.A., Hnisz D., Abraham B.J., Lin C.Y., Kagey M.H., Rahl P.B., Lee T.I., Young R.A. Master transcription factors and mediator establish super-enhancers at key cell identity genes. Cell. 2013;153:307–319. doi: 10.1016/j.cell.2013.03.035. [DOI] [PMC free article] [PubMed] [Google Scholar]

Supplementary Materials