If you don't remember your password, you can reset it by entering your email address and clicking the Reset Password button. You will then receive an email that contains a secure link for resetting your password
If the address matches a valid account an email will be sent to __email__ with instructions for resetting your password
Department of Molecular Physiology and Biophysics, Chemical and Physical Biology Program, Vanderbilt Genetics Institute, Vanderbilt Ingram Cancer Center, Vanderbilt University, Nashville, TN
Offprint requests to: Bryan J. Venters, Department of Molecular Physiology and Biophysics, Vanderbilt Genetics Institute, Vanderbilt Ingram Cancer Center, Vanderbilt University, Nashville, TN 37232
Department of Molecular Physiology and Biophysics, Chemical and Physical Biology Program, Vanderbilt Genetics Institute, Vanderbilt Ingram Cancer Center, Vanderbilt University, Nashville, TN
Chromatin immunoprecipitation-exonuclease maps for histone marks define enhancers in murine pro-erythroblasts.
•
Most pro-erythroblast enhancers are established at an earlier precursor stage.
•
Many typical, but not super-enhancer, locations are evolutionarily well conserved.
•
The hormone erythropoietin modulates chromatin dynamics of genes associated with key erythroid pathways.
The hormone erythropoietin (Epo) is required for erythropoiesis, yet its molecular mechanism of action remains poorly understood, particularly with respect to chromatin dynamics. To investigate how Epo modulates the erythroid epigenome, we performed epigenetic profiling using an ex vivo murine cell system that undergoes synchronous erythroid maturation in response to Epo stimulation. Our findings define the repertoire of Epo-modulated enhancers, illuminating a new facet of Epo signaling. First, a large number of enhancers rapidly responded to Epo stimulation, revealing a cis-regulatory network of Epo-responsive enhancers. In contrast, most of the other identified enhancers remained in an active acetylated state during Epo signaling, suggesting that most erythroid enhancers are established at an earlier precursor stage. Second, we identified several hundred super-enhancers that were linked to key erythroid genes, such as Tal1, Bcl11a, and Mir144/451. Third, experimental and computational validation revealed that many predicted enhancer regions were occupied by TAL1 and enriched with DNA-binding motifs for GATA1, KLF1, TAL1/E-box, and STAT5. Additionally, many of these cis-regulatory regions were conserved evolutionarily and displayed correlated enhancer:promoter acetylation. Together, these findings define a cis-regulatory enhancer network for Epo signaling during erythropoiesis, and provide the framework for future studies involving the interplay of epigenetics and Epo signaling.
Erythropoiesis has long served as an experimental paradigm for understanding gene regulatory mechanisms governing cell identity and development. The hormone erythropoietin (Epo) is physiologically necessary and sufficient to trigger the initiation of terminal erythroid differentiation of pro-erythroblasts. Epo initiates erythroid differentiation and gene expression patterns by binding to its cognate receptor, which activates the Jak2–Stat5 signaling axis [
], providing numerous insights into the molecular pathways that control red blood cell development. Additionally, expression profiling of Epo starvation experiments in pro-erythroblasts derived from primary fetal liver revealed that Epo signaling modulates the expression of several hundred genes involved in cell survival signaling and cell identity [
]. Enhancer segments typically are co-occupied by combinations of transcription factors (TFs) that influence gene expression programs, often integrating signals from multiple pathways. Collectively, TFs orchestrate cell type-specific transcription programs by binding to their cognate sequence motifs within specific enhancers and recruiting co-regulators and RNA polymerase II to promoters to initiate transcription [
]. Interestingly, a previous study found that broad features of chromatin states remain largely unchanged during Gata1-induced differentiation in the murine G1E-ER4 cell line despite the dramatic transcriptional changes that accompany erythropoiesis [
]. The study suggested that erythroid enhancers are established in erythroid precursor cells, but precisely when this occurs remains unclear. Although the locations of erythroid enhancers have been determined in murine and human cell systems [
], the ways in which Epo influences the enhancer landscape are currently unknown. Defining the comprehensive set of Epo-modulated enhancers will better illuminate the molecular and epigenetic pathways controlled by Epo.
To investigate the interplay of Epo signaling and the chromatin landscape during erythroid development, we used highly purified pro-erythroblasts derived from mice injected with an anemia-inducing strain of the Friend virus (FVA), as previously described [
]. Physiologically, mice injected with the FVA strain rapidly accumulate immature pro-erythroblast precursor cells in the spleen, a condition termed erythroblastosis. FVA-derived pro-erythroblast proliferation depends on the truncated form of the stem cell-derived tyrosine kinase receptor (sf-STK) but not JAK2-STAT5 activity [
]. Importantly, FVA pro-erythroblasts remain sensitive to physiological levels of Epo in an ex vivo culture system and, on Epo stimulation, will activate JAK2-STAT5 signaling and synchronously differentiate into fully mature, enucleated erythrocytes. Thus, FVA-derived pro-erythroblasts represent an ideal system for delineating the molecular action of Epo on a genomic scale and examining the epigenetic effects of Epo in a physiological context.
Here, we report the use of chromatin immunoprecipitation-exonuclease (ChIP-exo) analyses to profile the enhancer landscape during Epo-initiated erythropoiesis. Genome-wide locations for H3K4me1, H3K27ac, and H3K4me3 were examined in the FVA model system at baseline and 1 hour after Epo-initiated synchronous, terminal erythroid differentiation. This unique model cell system allowed us to uniquely define the immediate impact of Epo signaling on the epigenetic landscape during erythropoiesis.
Our work illustrates several new aspects of the role of Epo in erythroid chromatin biology. Remarkably, Epo stimulation alters the histone mark signatures across several thousand enhancer locations, highlighting Epo's underappreciated role in re-programming the epigenome. Enhancer marks were highly dynamic within the first hour of Epo stimulation, enabling classification of enhancer behavior based on the Epo-induced modulation of these histone marks. On Epo stimulation, a similar number of enhancers displayed increased (n = 1,589) and decreased (n = 1,529) acetylation. For example, loss of enhancer acetylation was linked to genes known to be downregulated during erythropoiesis, such as Cd44. Likewise, enhancers displaying increased acetylation levels in response to Epo were linked to genes involved in regulating the survival pathways, TFs, and chromatin regulators. Nevertheless, our finding that the vast majority of enhancers, including super-enhancers, were unaffected within the first hour of Epo stimulation fits with the notion that chromatin states remain largely unchanged during erythropoiesis [
Bioinformatic analyses of motif enrichment and H3K27ac enhancer:promoter correlation validated many of the enhancer locations identified in this study. Enhancer overlap with published data further validated enhancer locations, and revealed the location of evolutionarily conserved cis-regulatory modules across mouse and human erythroid cell genomes. Given that TAL1 is one of the master TF regulators of erythropoiesis, we also illustrate that TAL1 binds to thousands of enhancers, adding a layer of validation to our enhancer maps and suggesting a mechanistic link between Epo signaling, TAL1 occupancy, and chromatin dynamics. Together, these aspects of the presented study identify cis-regulatory enhancer circuits controlled by Epo signaling pathways. Importantly, the epigenetic profiles from the current study provide a framework for generating testable hypotheses about the ways in which Epo, enhancers, and TFs coordinately control erythroid expression programs.
Methods
Isolation of murine pro-erythroblasts
Highly purified pro-erythroblasts were obtained from spleens of mice infected with the Friend virus as previously described [
], with the following modifications. All animal procedures were performed in compliance with and approval from the Vanderbilt Division of Animal Care (DAC) and the Institutional Animal Care and Use Committee (IACUC). Female BALB/cJ mice (12 weeks old, Jackson Laboratories) were infected via intraperitoneal injection of ∼104 spleen focus-forming units of FVA. At 13 to 15 days postinfection, the mice were sacrificed and their spleens removed. The spleens were homogenized to a single-cell suspension by mincing and passing through a sterile 100-μm nylon mesh filter into sterile solution of 0.2% bovine serum albumin (BSA) in 1× phosphate-buffered saline (PBS). The filtrate was pipetted repeatedly to ensure a single-cell suspension. The homogenized spleen cells were size-separated by gravity sedimentation for 4 hours at 4°C in a continuous gradient of 1% to 2% de-ionized BSA. The sedimentation apparatus consisted of a 25-cm-diameter sedimentation chamber containing a 2.4-L BSA gradient, two BSA gradient chambers containing 1.2-L 1% and 2% de-ionized BSA in 1× PBS, and a cell-loading chamber (ProScience Inc.) containing the 50-mL cell suspension. After 4 hours of sedimentation, cells were collected in 50-mL fractions, with pro-erythroblasts typically enriched in fractions 5–20 of 24 total fractions. Typically, about 109 pro-erythroblasts were obtained from the separation of 1010 nucleated spleen cells (6–7 g spleen weight) across three 25-cm sedimentation chambers.
Cell culture
To study the effects of erythropoietin (Epo) on terminal erythroid differentiation, FVA-derived pro-erythroblasts were cultured at 106 cells/mL in Iscove's modified Dulbecco's medium (IMDM, Life Technologies No. 12440043), 30% heat-inactivated fetal bovine serum (Gibco, No. 26140-079), 1% penicillin–streptomycin (Gibco, No. 15140-122), 10% de-ionized BSA, and 100 μmol/L α-thioglycerol (MP Biomedicals, No. 155723). Terminal erythroid differentiation of purified pro-erythroblasts was induced by the addition of 0.4 U/mL human recombinant Epo (10kU/ml Epogen by Amgen, NDC 55513-144-10) to medium. At the desired times after the addition of Epo, cells were crosslinked by the addition of 1% formaldehyde for 10 min. Crosslinking was quenched by adding 125 mmol/L glycine. Crosslinked cells were collected by centrifugation for 5 min at 1,000g at 4°C, washed once with 1× PBS, flash frozen in liquid nitrogen, and stored at −80°C until used for ChIP analysis.
ChIP-exo and antibodies
With the following modifications, ChIP-exo was performed as previously described [
] with chromatin extracted from 50 million cells, ProteinG MagSepharose resin (GE Healthcare), and 5 μg of antibody directed against the H3K4me1, H3K4me3, H3K27ac, or Tal1 (Abcam ab8895, ab8580, ab4729, and Santa Cruz Biotech sc-12984, respectively). First, formaldehyde crosslinked cells were lysed with buffer 1 (50 mmol/L HEPES–KOH, pH 7.5; 140 mmol/L NaCl; 1 mmol/L EDTA; 10% glycerol; 0.5% Nonidet P-40; 0.25% Triton X-100) and washed once with buffer 2 (10 mmol/L Tris–HCL, pH 8; 200 mmol/L NaCl; 1 mmol/L EDTA; 0.5 mmol/L EGTA), and the nuclei were lysed with buffer 3 (10 mmol/L Tris–HCl, pH 8; 100 mmol/L NaCl; 1 mmol/L EDTA; 0.5 mmol/L EGTA; 0.1% sodium deoxycholate; 0.5% N-lauroylsarcosine). All cell lysis buffers were supplemented with fresh EDTA-free complete protease inhibitor cocktail (CPI, Roche, No. 11836153001). Purified chromatin was sonicated with a Bioruptor (Diagenode) to obtain fragments ranging from 100 bp to 500 bp. Triton X-100 was added to extract at 1% to neutralize sarcosine. Insoluble chromatin debris was removed by centrifugation, and sonication extracts were stored at −80°C until used for ChIP analysis.
Illumina sequencing and data pre-processing
Libraries were sequenced using an Illumina NextSeq500 sequencer as single-end reads 75 nucleotides in length on high output mode (Supplementary Table E1, online only, available at www.exphem.org). The sequence reads were aligned to the mouse mm10 reference genome using BWA-MEM algorithm with default parameters [
]. Because the patterns described here were evident among individual biological replicates and replicates were well correlated (Supplementary Figure E1, online only, available at www.exphem.org), we merged all tags from biological replicate data sets for final analyses.
Data visualization
Raw sequencing tags were smoothed (20bp bin, 100bp sliding window) and normalized to reads per kilobase per million (RPKM) using deepTOOLS [
] was used to display the distribution of histone marks relative to the transcription start site (TSS). RPKM normalized density plots were generated using Java TreeView [
Figure 1Experimental overview and distribution of histone marks at the β-globin locus control region. (A) Shown is the workflow for generating and isolating highly purified Epo-responsive pro-erythroblasts from a mouse injected with FVA. (B) ChIP-exo cartoon illustrates how the 5′ to 3′ lambda exonuclease trims DNA to the crosslink point, effectively footprinting the protein–DNA interaction. Mapping the 5′ ends of the sequence tags to the reference genome demarcates the exonuclease barrier and, thus, the precise site of protein-DNA crosslinking. (C) Genome browser view of H3K4me1, H3K27ac, and H3K4me3 ChIP-exo signal in murine pro-erythroblasts shown at the β-globin (Hbb) LCR. Tag distributions were smoothed and RPKM was normalized using deepTools
Figure 2ChIP-exo reveals histone modification patterns at gene promoters. (A) ChIP-exo libraries were prepared for H3K4me1, H3K27Ac, and H3K4me3 and visualized using ChAsE
. Aligned heatmaps indicate RPKM normalized number of reads across a 4-kb genomic interval in 20-bp bins relative to the TSS before and after Epo stimulation. Regions are sorted in descending order based on average row tag density for 0-hour H3K27ac. Each row represents one protein-coding gene, with 20,140 lines/genes present. Red and blue reflect high and low read densities, respectively. (B) Composite plots below each heatmap quantify the normalized tag density. The central trace denotes the average tag density for each 20-bp bin, and the orange fill reflects the standard deviation. (C) Model illustrating the nucleosome position for enriched histone marks, which are derived from the tag enrichment patterns in (A) and Supplementary Figure E2 (online only, available at www.exphem.org).
To identify genomic intervals whose patterns of histone marks were consistent with enhancers, we applied a hidden Markov modeling algorithm using ChromHMM [
]. To obtain a unified set of enhancer intervals, the data for each histone mark were merged across both time points. Briefly, the merged files were binarized using the BinarizeBed function with the signal threshold option set to 100. The minimum number of nonredundant chromatin states was heuristically determined to be six using the LearnModel function with default parameters.
Quantification and annotation of ChromHMM enhancer intervals
Chromatin states 1 and 2 (Supplementary Table E2, online only, available at www.exphem.org), whose chromatin patterns were consistent with enhancers, were annotated and quantified using the Hypergeometric Optimization of Motif EnRichment (HOMER) suite [
]. Briefly, bam files were converted to tag directories using the makeTagDirectory function with the –genome, –checkGC, and –format options. To quantify and normalize tags within enhancer regions to RPKM, the analyzeRepeats function was used with the –rpkm and –d options. The log2 fold change (from 0 to 1 hour after Epo) was then computed for each histone mark. Finally, enhancer–gene associations were inferred based on the nearest gene paradigm [
], and each enhancer interval was linked to a gene using the annotatePeaks function with the –noadj and –d options.
Figure 3Hidden Markov modeling of chromatin states identify unique enhancer signatures. (A) The cartoon illustrates distinct chromatin signatures associated with enhancers and promoters. (B) The transition matrix illustrates the probability that the state in the row (State From) will transition to the state in the column (State To) using hidden Markov modeling of chromatin states
. (C) Heatmap representation of chromatin state emissions from hidden Markov modeling that predicts enhancer intervals from H3K4me1, H3K27ac, and H3K4me3. The red box highlights states 1 and 2, which display chromatin state patterns consistent with enhancers that enrich H3K4me1 and/or H3K27ac but not H3K4me3. (D) Functional annotation inferences for each chromatin state are listed next to the number of intervals comprising the state and the average length. (E, F) Heatmaps illustrating the enrichment of each learned chromatin state ±2000 bp from the TSS and transcription termination site (TTS), respectively.
], and implemented by HOMER using the findPeaks function with the –style option set to super. For this analysis, the minimum distance to stitch peaks together (-minDist) and the local fold change thresholds (-L) were set to 6,000 and 1, respectively. Data for H3K4me1 and H3K27ac marks were combined to generate a single “enhancer” tag directory for 0 and 1 hour, individually. Enhancer intervals whose midpoints were within 1 kb of a TSS were filtered. The remaining enhancer boundaries were trimmed 3 kb distal to a TSS if their interval edges overlapped with a TSS. As with the ChromHMM enhancer analysis described previously, the super-enhancer regions were quantified and annotated using the analyzeRepeats and annotatePeaks functions.
Classification of enhancers
Similar to the criteria detailed in Ostuni et al. [
], enhancer intervals were assigned to seven enhancer classes based on H3K27ac and H3K4me1 RPKM values and H3K27ac fold-change for signals in response to Epo stimulation (Supplementary Tables E3 and E4, online only, available at www.exphem.org). Briefly, RPKM values within the lowest quartile and the upper three quartiles were designated as “low” and “high,” respectively. Occupancy changes in response to Epo were designated as “up” or “down” based on whether RPKM fold-changes increased or decreased by more than twofold, respectively. RPKM fold changes that were between these values were considered as no change. The criteria for each enhancer class are summarized below: (1) Constitutive activated enhancers displayed “high” H3K27ac and “up” in H3K27ac. (2) Constitutive not-activated enhancers displayed “high” H3K27ac and no change in H3K27ac. (3) Poised activated enhancers displayed “high” H3K4me1 levels, “low” H3K27ac, and “up” in H3K27ac. (4) Poised not-activated enhancers displayed “high” H3K4me1 levels, “low” H3K27ac, and no change in H3K27ac. (5) Latent enhancers displayed “low” H3K4me1 and H3K27ac levels and “up” in H3K4me1 and H3K27ac levels. (6) Latent to poised enhancers displayed “low” H3K4me1 and H3K27ac levels, “up” in H3K4me1, and no change or “down” in H3K27ac levels. (7) Repressed enhancers displayed “high” H3K27ac and “down” in H3K27ac. Lastly, 3,279 intervals were excluded from further analysis because they displayed biologically irrelevant characteristics, such as low H3K4me1 and/or H3K27ac levels that further decreased on Epo stimulation.
Motif discovery
De novo motif discovery for all enhancer regions was conducted using the findMotifsGenome function in the HOMER suite. To comprehensively identify locations for motifs related to erythroid cell function that were enriched in the de novo search, we conducted a directed search for the following consensus motifs [
] with zero mismatches allowed: GATA1 (WGATAR), KLF1 (YMCDCCCW), TAL1-Ebox (CANNTG), ETS (YWTCCK), and STAT5 (TTCYHDGAA). The scanMotifGenomeWide function in HOMER and intersectBED function in BEDtools [
Processed data from data analyses are available in the Supplementary Material (online only, available at www.exphem.org). Raw sequencing data are available at NCBI Sequence Read Archive (SRP082181).
Results
Experimental overview and distribution of histone marks at the β-globin locus control region
To study the molecular action of Epo, we leveraged the anemia-inducing strain of the Friend virus (FVA) murine model system that has been reported to recapitulate normal erythropoiesis, as evidenced by JAK2-STAT5 signaling, globin expression kinetics, cell morphology, cell surface marker kinetics, and cellular enucleation [
]. Indeed, the FVA-derived pro-erythroblasts were recently used as a standard for developing an improved flow cytometry sorting scheme for bone marrow-derived erythroblasts [
]. This system enables facile large-scale procurement of highly purified murine pro-erythroblast cell populations that synchronously respond to Epo (Fig. 1A) [
Figure 4Epo-induced remodeling of the erythroid enhancer landscape. (A) Classification of erythroid enhancers based on their response to Epo. (B–D) Representative genome browser views of Epo-responsive enhancers (B and D: constitutive activated and latent activated, respectively) and an Epo non-responsive enhancer (C, constitutive not activated). Red boxes denote genomic intervals identified as enhancers by hidden Markov modeling. (E) Density plots illustrating the dynamic behavior of histone marks for each enhancer interval induced by a 1-hour Epo treatment. The number of enhancer intervals shown is indicated to the left of each plot. The normalized RPKM median value for each column is represented as a bar graph below the density plots.
Friend virus infection in susceptible mice leads to a multistage disease, characterized initially by erythroblastosis and then later by erythroleukemia [
]. FVA-derived pro-erythroblasts used in this study were harvested during the initial cell expansion stage of erythroblastosis and were not erythroleukemic. The two strains of the Friend virus result in distinct phenotypic outcomes caused by amino acid differences in the Friend virus glycoprotein, GP55. The polycythemia-inducing GP55P variant requires both EpoR receptor and sf-STK, thereby constitutively activating JAK2-STAT5 signaling and Epo-independent differentiation [
]. Mechanistically, the anemia-inducing GP55A variant directly interacts with the sf-STK receptor to enable erythroid precursor proliferation in the absence of Epo (Fig. 1A). Indeed, a previous report indicated that neither EPOR nor STAT5 is required for FVA-induced erythroblastosis [
To assess the extent to which Epo influences the epigenetic landscape, we performed ChIP-exo on three key histone H3 modifications (H3K4me1, H3K27ac, and H3K4me3) before and after cells were exposed to Epo (0 and 1 hour Epo, Fig. 1B). For each antibody, ChIP-exo signals across biological replicates were highly correlated (R = 0.96–0.99), indicating high reproducibility (Supplementary Figure E1, online only, available at www.exphem.org). Importantly, these three histone marks enable complex pattern recognition algorithms, such as hidden Markov modeling, to identify candidate enhancers genome-wide and to distinguish them from promoter-specific patterns [
]. The putative functional nature of enhancers in the present study should be noted, because enhancer predictions are based on structural histone modification patterns. In Figure 1C, the histone modification patterns are illustrated at the well-characterized β-globin locus and nearby locus control region (LCR), a cluster of enhancers that control β-globin transcription [
]. Consistent with the marks found at enhancers, H3K4me1 and H3K27ac, but not H3K4me3, were enriched at the LCR. H3K4me3 was enriched at the promoters of the β-major and -minor globin genes that are expressed in adult erythroid cells. In contrast, the embryonically expressed βh1 and εγ loci lacked both H3K27ac and H3K4me3 associated with active promoters.
ChIP-exo analysis reveals histone modification patterns at gene promoters
Next, we examined the global ChIP-exo patterns for H3K4me1, H3K27ac, and H3K4me3 at promoters of protein-coding genes. The ChIP-exo signal for each mark before and after Epo treatment (0 and 1 hour, respectively) was aligned to the coding TSS and displayed as a heatmap (Fig. 2A; sorted by maximal peak position in Supplementary Figure E2, online only, available at www.exphem.org). Each histone mark was detected at nearly one half of all protein-coding genes (Supplementary Table E6, online only, available at www.exphem.org), consistent with previous reports indicating that these marks are not found at all genes, but rather, are enriched in the subset of genes actively transcribed in erythroid cells [
Figure 5Validation of candidate enhancers. (A) TF binding motifs over-represented in each enhancer class are shown as motif logos based on their respective position weight matrices (PWMs). De novo discovery of PWMs, p values, and best Jaspar motif match were generated using HOMER. (B) Bar graph illustrating percentage of each enhancer class containing at least one GATA1 consensus motif (WGATAR) with zero mismatches. Activated and not-activated subclasses are denoted by “A” and “NA,” respectively. (C) Bar graph of the percentage of enhancers, separated by class, that are nearest to genes with greater than twofold increase (green) or decrease (red) in H3K27ac at promoters (±1 kb). Activated and not-activated subclasses are denoted by “A” and “NA,” respectively. (D) Percentage of enhancer intervals occupied by ≥1 ChIP-exo Tal1 peaks and broken out by enhancer class. Activated and not-activated subclasses are denoted by “A” and “NA,” respectively. (E) Venn overlap between enhancer intervals in this study (red fill) and mouse G1E-ER4 chromatin states 2 and 3 (H3Kme1 and H3K4me1/H3K27me3, respectively) from Wu et al.
It is well established that the genome-wide distribution of histone marks at promoters is not random, but exhibits distinct preferences to specific nucleosomes within promoter regions [
]. Given the mapping precision afforded by ChIP-exo analysis, we examined the spatial distribution of histone marks relative to each other and nucleosome positions relative to the TSS. As expected, we found that H3K4me1 and H3K4me3 tend to be mutually exclusive in their preferred nucleosomes, with H3K4me3 most enriched at the −1 and +1 nucleosomes. In contrast, H3K4me1 marks were enriched at nucleosomes distal to the TSS, particularly in the +3 and +4 positions, but decayed rapidly thereafter. Notably, the H3K27ac ChIP-exo signal was enriched in a pattern similar to that of H3K4me3 at the −1 and +1 nucleosomes, reflecting the expected signature pattern at the nucleosomes flanking core promoter regions [
]. Composite plot distribution peaks (in Fig. 2B and C) further highlight the nucleosomal preferences of specific histone marks at promoters of protein-coding genes in pro-erythroblasts.
Hidden Markov modeling of chromatin states reveals unique enhancer signatures
Pioneering studies by the ENCODE consortium suggest that the mammalian genome may harbor hundreds of thousands of enhancer regions, with each of the approximately 200 different cell types containing distinct repertoires of enhancer elements [
]. Given this complexity and despite recent advances in high-throughput technologies, it remains a challenge to identify the complete set of bona fide functional enhancers in a given cell type. Nevertheless, this large-scale project and other studies have found that enhancers and promoters display distinct epigenetic signatures (Fig. 3A), which can be used to predict enhancer locations in differing cell types with a relatively high success rate [
]. This signature distinguishes enhancers from promoters, which are marked by histone 3 lysine 4 trimethylation (H3K4me3).
Thus, to systematically identify the genomic locations of chromatin signatures consistent with enhancers in FVA-derived pro-erythroblasts, we applied hidden Markov modeling (HMM) to ChIP-exo data (Fig. 3). ChromHMM employs a powerful algorithm for finding unique patterns in complex data [
]. HMM analysis of the ChIP-exo data for H3K4me1, H3K27ac, and H3K4me3 at 0- and 1-hour Epo time points revealed six non-redundant chromatin states (Fig. 3B, C). Chromatin states 1 and 2 represent nearly 60,000 candidate enhancers with an average length of ∼900 bp and signatures consistent with strong and weak enhancers, respectively—that is, having H3K4me1 and/or H3K27ac but lacking H4K4me3, respectively. As expected, these candidate enhancer regions tended to be distal to transcription start and end sites (Fig. 3E, F; quantified in Supplementary Figure E3, online only, available at www.exphem.org). In contrast, chromatin states 4–6 were strongly associated with promoter regions (Fig. 3E) and were each strongly enriched with H3K4me3 relative to states 1 and 2. Lastly, chromatin state 3 was depleted of ChIP signal for the three histone marks in this study and encompassed the majority of the genome, suggesting low representation of these marks and/or the possible presence of other histone modifications in the intervals of this state.
Epo-induced remodeling of the erythroid enhancer landscape
Tens of thousands of enhancers are scattered across mammalian genomes in a cell type-specific manner and are dynamically shaped in response to environmental stimuli [
], the ways in which Epo influences erythroid enhancers remain unclear. We hypothesized that Epo signaling reshapes the epigenetic landscape of a specific set of erythroid enhancers and drives the erythroid transcriptional program. To test this hypothesis, we systematically classified the nearly 60,000 candidate enhancer regions based on their response to Epo stimulation, in a manner similar to that previously described [
]. Poised enhancers initially lack H3K27ac but can acquire the acetyl mark with cellular stimulation from environmental cues, such as cytokines and hormones [
Overall, response to Epo stimulation was highly dynamic, enabling stratification of the ∼60,000 genomic intervals into four major classes (Fig. 4A). Candidate enhancers across each class shared several features, including commonly residing within intergenic and intronic regions of the genome, averaging about 1 kb in length, and varying quite dramatically in their distances to the nearest gene (Supplementary Figure E4, online only, available at www.exphem.org). First, the class of constitutive enhancers was associated with both H3K4me1 and H3K27ac in unstimulated pro-erythroblasts (0-hour Epo). Although most (45,297) were unaffected after 1 hour of Epo stimulation (constitutive not-activated, Fig.4E), 855 displayed increased acetylation within the first hour of Epo stimulation (constitutive activated, Fig. 4E). Figure 4B illustrates the Epo-induced dynamics of a constitutive activated enhancer at the Tnfrsf13c gene, which interestingly, encodes a receptor known to promote survival in B cells [
]. This observation supports the notion that the enhancer we identified upstream of the Tnfrsf13c locus may be mediating its transcriptional response to Epo to further promote pro-survival pathways that Epo is known to regulate [
]. As illustrated in Figure 4C, the Gata1 locus, which encodes one of the three master regulators of the erythroid expression program, exemplifies the constitutive not-activated class of enhancers. Consistent with this observation, previous studies have indicated that maximal Gata1 gene expression, and presumably increased enhancer acetylation, occurs in basophilic and orthochromatic erythroblasts [
], which are not present in the model system until 12–36 hours after Epo stimulation.
Next, the class of poised enhancers was characterized by H3K4me1, but not H3K27ac, in unstimulated pro-erythroblasts. Of the 7,841 poised enhancers we detected, 688 were Epo-responsive and displayed increased acetylation after 1 hour of Epo stimulation (poised activated, Fig. 4E). Gene ontology analysis (Supplementary Table E7, online only, available at www.exphem.org) revealed that the set of genes proximal to poised activated enhancers was most associated with positive regulation of TF activity (p = .001), including Kit, the gene encoding the receptor tyrosine kinase that activates a number of signaling pathways critical to cell survival during erythropoiesis. Also included in the group of poised activated enhancers were genes encoding Wnt signaling proteins (WNT3A and WNT10B), critical regulators that balance self-renewal, proliferation, and differentiation in hematopoietic cells [
We also identified a class of more than 100 latent enhancers, which gained H3K4me1 and/or H3K27ac on Epo stimulation (Fig. 4E). Remarkably, nearly one half of these latent enhancers (45) were activated within 1 hour of Epo stimulation. Figure 4D illustrates the Epo-induced dynamics of the class of latent activated enhancers at the Bop1 (Block of proliferation 1) gene locus. The precise timing of the Bop1 latent enhancer activation with Epo-stimulated erythropoiesis is consistent with the cellular transition to the inherently non-proliferating mode of terminal differentiation. Interestingly, as a protein involved in rRNA processing, Bop1 has been implicated in ribosomopathies [
Finally, the class of enhancers consistent with repression displayed both H3K4me1 and H3K27ac in untreated pro-erythroblasts, but lost or exhibited diminished H3K27ac within 1 hour of Epo stimulation (Fig. 4E). This class of enhancers was enriched near numerous genes encoding TF regulators (p = 10−7), including MAX. MAX is the primary partner for MYC, conferring the ability of MYC/MAX dimers to bind E-boxes. Together, MYC/MAX dimers control the opposing cell fate decisions of cell cycle progression and apoptosis [
], suggesting that Epo-induced enhancer remodeling may contribute to cell fate decisions during erythropoiesis via MAX. Consistent with the observation that Epo simulation leads to repression of an enhancer linked to the binding partner for MYC, a previous report indicated that downregulation of MYC is important during erythropoiesis [
In summary, we identified more than 3,000 Epo-responsive enhancers, of which 1,588 and 1,529 were activated and repressed, respectively, within the first hour of terminal erythroid differentiation. This suggests that Epo rapidly remodels a subset of the epigenetic landscape of pro-erythroblasts undergoing erythropoiesis. Consistent with a previous report [
] indicating that chromatin states remain largely unchanged during erythroid differentiation, the vast majority of the candidate enhancer regions identified in this study (n = 52,508) were static during the first hour of Epo exposure. This indicates that they are established at a prior precursor stage and are not Epo-responsive or that they respond to Epo at some point after 1 hour of treatment.
Validation and evolutionary conservation of candidate enhancers
The possibility of false positives is a caveat to algorithms that predict enhancer regions based on functional genomic data. To provide corroborative evidence for the validity of candidate enhancers, we applied four complementary computational and experimental approaches: (1) DNA motif analysis, (2) correlation to promoter acetylation, (3) overlap analysis with orthogonal published data sets, and (4) TF binding via ChIP-exo analysis. Transcription factors regulate cell type-specific gene expression patterns by binding to their cognate motifs within enhancer regions and recruiting the transcription machinery to target genes [
]. In erythroid cells, GATA1, KLF1, and TAL1 are master TF regulators and play a major role in coordinating the precise timing of gene expression patterns during erythropoiesis [
The LIM-only protein Lmo2 is a bridging molecule assembling an erythroid, DNA-binding complex which includes the TAL1, E47, GATA-1 and Ldb1/NLI proteins.
Therefore, we hypothesized that HMM-identified candidate enhancers that are physiologically relevant in erythroid cells would be enriched for motifs corresponding to the master regulator TFs of erythropoiesis. To test this hypothesis, we conducted de novo motif discovery analysis of all ∼60,000 enhancer regions from Figure 4. Strikingly, the Gata motif was the most commonly enriched (p = 10−13 to 10−598) across the enhancer classes (Fig. 5A). Other erythroid-related motifs found in the de novo search included the KLF1, ETS, and STAT5 consensus motifs.
Figure 6Identification of super-enhancers at erythroid genes involved in cell fate decisions. (A) Histone modification enrichment across super-enhancers in erythroid cells before and after Epo stimulation. (B) Density plots illustrating the H3K4me1, H3K27ac, and H3K4me3 ChIP-exo signal (RPKM normalized) across super-enhancer intervals before and after Epo treatment. The normalized RPKM median value for each column is represented as a bar graph below the density plots. (C–E) Genome browser views of super-enhancers before and after Epo stimulation at the Tal1, Bcl11a, and Mir144/Mir451 genes, respectively. Red boxes denote genomic intervals identified as super-enhancers by the Whyte et al. algorithm employed by HOMER. (F,G) Venn overlap of super-enhancer intervals identified in this study (red fill) compared with super-enhancer regions from Huang et al.
Because de novo motif analysis provides initial clues to over-represented motifs in the genomic regions, we next conducted a more sensitive, direct motif search. Given the connection between Epo signaling and GATA1 phosphorylation [
], we focused on identifying all occurrences of the GATA consensus motif within enhancers (Fig. 5B, Supplementary Figure E5, online only, available at www.exphem.org). Indeed, for each class of enhancers, the observed frequency of the GATA motif was significantly greater than an equivalent number of randomly sampled genomic regions (Fig. 5B, empirical p value <.001). Altogether, our motif analysis provides corroborative evidence for validating candidate enhancers.
Previous studies indicated that enhancer:promoter interactions may be inferred by correlating H3K27ac patterns [
]. Therefore, we examined the relationship between Epo-stimulated changes in H3K27ac at enhancers and those of the nearest gene promoter. In general, genes that increased in promoter acetylation were linked to nearby enhancers that also acquired H3K27ac in response to Epo (e.g., activated enhancer subclasses denoted by ‘A’) (Fig. 5C), with the repressed enhancer class least represented. Overall, the correlation between decreased promoter and enhancer acetylation was far weaker, although the repressed class of enhancers was most associated with a decrease in promoter acetylation. Not only do these results support the enhancer classification methodology (Fig. 4), they provide further corroborative evidence for the candidate enhancers in this study.
Given the previous finding that TF binding is an accurate predictor of enhancer activity [
], we performed TAL1 ChIP-exo analysis in FVA-derived pro-erythroblasts and identified 19,025 and 13,913 enriched peaks for the 0- and 1-hour Epo treatments, respectively. Strikingly, nearly 60% (20,289 of 32,938) of all TAL1 peaks resided within candidate enhancer regions identified in this study, which are separated by class in Figure 5D. In other words, TAL1 alone occupied 33% of all candidate enhancers (18,370 of 55,626 enhancers with one or more TAL1 peaks), suggesting that these enhancer regions are valid.
Lastly, we compared our set of enhancer locations to published orthogonal data sets from mouse and human studies [
] were shared with the current study (Fig. 5F), indicating substantial evolutionary conservation of erythroid cis-regulatory modules. In summary, based on evidence from the four aforementioned complementary approaches used to validate enhancers, we conclude that candidate enhancer regions identified in this study generally represent bona fide enhancers and are likely to be physiologically and functionally relevant to terminal erythroid differentiation.
Identification of super-enhancers at erythroid genes involved in cell fate decisions
Super-enhancers are an emerging class of long, clustered enhancers important for establishing cell fate and identity [
], but whether Epo stimulation affects super-enhancers remains unclear. Therefore, to investigate the extent to which Epo signaling influences super-enhancer dynamics in pro-erythroblasts, we applied the algorithm previously described for identifying super-enhancers [
] to our H3K4me1 and H3K27ac ChIP-exo data sets. We identified 395 super-enhancers dispersed throughout the mouse genome that were characterized by their length (21,192 bp on average; Supplementary Figure E4) and exceptionally high H3K4me1 and H3K27ac levels (Fig. 6A,B). We found that the super-enhancers were generally nonresponsive to Epo with respect to H3K27ac (Fig. 6B). Consistent with the notion that super-enhancers tend to be associated with genes controlling cell fate and identity, they were discovered near a number of genes known to play roles in erythroid cell physiology, such as Tal1, Bcl11a, and Mir144/451 (Fig. 6C–E, Supplementary Figure E6 for additional screen shots, and a full list of super-enhancer locations and nearby genes may be found in Supplementary Table E3, online only, available at www.exphem.org). Tal1 and Bcl11a encode TFs that are known to regulate erythroid gene expression patterns [
]. Interestingly, beyond the representative genes in Figure 6C–E, TAL1 bound to 92% of all super-enhancer regions (365 of 395 enhancers; Supplementary Figure E7, online only, available at www.exphem.org). In addition, GATA motifs were found in nearly all super-enhancer regions (99.7%, 394 of 395 super-enhancers; Supplementary Figure E5, online only, available at www.exphem.org). Indeed, the GATA motif overlapped with 95.6% of TAL1 peaks within super-enhancers (345 of 361 Tal1 peaks for 0-hour Epo; Supplementary Figure E8, online only, available at www.exphem.org). Overlap analysis between super-enhancer locations revealed that 57% (225 of 395; Fig. 6F) of super-enhancer locations from this study were shared with a previous report in mouse G1E-ER4 cells [
], thereby corroborating our super-enhancer analysis. In contrast, super-enhancer locations were poorly conserved evolutionarily when compared to super-enhancer locations from CD34+ derived pro-erythroblasts (16%, 64 of 395; Fig. 6G).
To investigate what types of genes are associated with Epo-modulated enhancers, we conducted gene ontology and pathway enrichment analysis (Fig. 7A). We found that many enhancers modulated by Epo stimulation were linked to genes involved in cell-signaling pathways, such as JAK-STAT (n = 24), PI3K (n = 45), and FOXO (n = 25) signaling. Although these pathways are known to regulate erythroid differentiation [
], it is surprising that enhancers linked to these pathways are remodeled in response to Epo. Nevertheless, this observation is consistent with the notion that Epo-driven epigenetic changes may form feed-forward loops that serve to strengthen the kinetics and robustness of Epo-responsive pathways. During terminal erythroid differentiation, the cytoskeletal structure undergoes a dramatic transformation, which involves actin and the glycosylation of a number of erythroid transmembrane proteins, such as Band 3 and Glycophorin A [
]. Indeed, within the first hour of exposure to Epo, the epigenetic landscape is altered near genes involved in regulating the actin cytoskeleton (n = 31), glycosylation (n = 153), and transport (n = 221). This suggests that cytoskeletal reorganization may be regulated in part at the epigenetic level.
Figure 7Epigenetic integration of signaling pathways by Epo-modulated enhancers. (A) Gene ontology and pathway enrichment analysis for genes linked to Epo-modulated enhancers. Top enriched pathways, biological processes, and molecular functions are represented by blue, gray, and white fill bars, respectively. (B) Model illustrating the interplay of pathways enriched in the set of genes linked to Epo-modulated enhancers.
Collectively, TFs orchestrate cell type-specific transcription programs by binding to their cognate sequence motifs within specific enhancers. TFs then recruit co-regulators and RNA polymerase II to promoters to initiate erythroid transcription programs. Strikingly, analysis of genes linked to Epo-modulated enhancers revealed an enrichment of numerous regulators of transcription (n = 284), chromatin regulators (remodeling and modification, n = 60), and phosphorylation regulators (kinase and phosphatase activity, n = 97). In particular, Epo simulation altered enhancers linked to Bach1, Bcl11a, Elf1, Klf10, Sox6, Tcf3, and Tif1γ, each of which has previously been reported to be involved in controlling erythroid expression programs [
Transcriptional regulation of Elf-1: Locus-wide analysis reveals four distinct promoters, a tissue-specific enhancer, control by PU.1 and the importance of Elf-1 downregulation for erythroid maturation.
]. Together, these observations suggest that Epo stimulation influences a network of erythroid-relevant TFs at the epigenetic level, and may serve to coordinately mobilize TF activity to orchestrate erythroid expression programs.
Discussion
How hormone signaling is manifested in chromatin dynamics remains an important, unanswered question. Here, we have illustrated that the hormone Epo reprograms a repertoire of enhancers in murine pro-erythroblasts that are linked to integrated networks of erythroid transcriptional regulators and signaling pathways. Identifying where and how these histone marks are altered in response to stimuli reveals the plasticity of erythroid enhancers and provides new insights into the molecular action of Epo.
Our experimental murine model uniquely allowed us to investigate how the enhancer landscape changes on Epo stimulation. We identified ∼60,000 enhancers in pro-erythroblasts, and in response to Epo stimulation, their histone mark dynamics enabled stratification of enhancers into four major classes in addition to the super-enhancer class. Consistent with the findings of a previous study [
], we found that the enhancer landscape remains largely unchanged before and after Epo stimulation of pro-erythroblasts. However, nearly 3,000 enhancers were remodeled within the first hour of Epo treatment, suggesting that Epo initiates terminal erythroid differentiation in part by influencing a portion of the erythroid enhancer repertoire. Whether Epo continues to influence and reshape the epigenetic landscape throughout erythropoiesis remains unclear.
Two common themes emerged among the genes linked to Epo-responsive enhancers: (1) over-representation of pathways involved in coordinating the transition from a proliferative state to a non-proliferative differentiating state (Fig. 7B), and (2) enrichment of genes encoding TFs. The transition to a non-proliferative differentiating state is exemplified by Epo-responsive enhancers linked to Tnfrsf13c and Bop1 (Fig. 4B, D) and is consistent with the physiological transitions that accompany the cell fate decision to undergo terminal erythroid differentiation. In particular, once Epo triggers the cellular commitment to undergo erythropoiesis, the cell must halt proliferative signaling pathways and commence differentiation pathways, while blocking the cellular inclination toward apoptosis. Interestingly, Epo starvation and add-back experiments resulted in strongly upregulated Tnfrsf13c expression in murine fetal pro-erythroblasts [
]. This finding is in agreement with our discovery of a constitutive enhancer upstream of Tnfrsf13c that is further activated within the first hour of Epo stimulation.
Mechanistically, this study provides insight into how Epo stimulation alters the enhancer landscape in erythroid cells, which are presumably bound by TFs to orchestrate the dynamic gene expression programs during erythropoiesis. Indeed, the top Gene Ontology category for genes linked to the Epo-responsive class of constitutive activated enhancers is the positive regulation of transcription from Pol II promoters (Supplementary Table E7, online only, available at www.exphem.org; p = .0001), which included the genes for TFs Lmo2, Klf5, and Klf13 and the Ep300 gene, which encodes an acetyltransferase. Furthermore, we found super-enhancers near a number of important genes encoding TFs that regulate erythropoiesis, such as Bcl11a and Tal1. It is interesting to note that a recent study reported a Tal1 super-enhancer specific to T-cell acute lymphoblastic leukemia cells, leading the study's authors to suggest that the Tal1 super-enhancer is oncogenic [
]. However, in the context of normal erythropoiesis, we propose that the Tal1 super-enhancer identified in the present study serves a physiological function controlling early Tal1 expression kinetics in murine pro-erythroblasts. With respect to super-enhancers, a recent report described a super-enhancer at the α-globin locus and elegantly dissected the cis-regulatory requirements controlling α-globin expression [
]. This super-enhancer also was found in the present study (Supplementary Figure E6D; online only; www.exphem.org), further supporting the validity of the enhancers presented.
Several previous studies have used epigenetic profiling to identify candidate enhancers in human and mouse erythroid cells [
]. Although these studies provided detailed insights into erythroid enhancers during erythroid differentiation, it is important to note that the inherent limitations of these model systems precluded an assessment of Epo-dependent enhancer dynamics. Furthermore, putative enhancers were identified without considering H3K4me1 [
], which would miss the dynamic behavior of activated, repressed, poised, and latent enhancer states.
Specifically, our analysis of both H3K4me1 and H3K27ac histone modifications by ChIP-exo analysis in the FVA system circumvents these limitations and provides new insights into erythroid biology by tying dynamic enhancer states directly to Epo stimulation. In addition, our finding that Epo uncovers nearly 100 latent enhancers in erythroid cells further expands the repertoire of environmental stimuli that reveal cryptic enhancer locations [
The primary strength of using chromatin state mapping to identify enhancers is that it provides insights into the plasticity of the enhancer landscape in an unbiased and genome-wide manner, particularly in response to stimuli. We acknowledge the limitations of the present study for experimental models and approaches. First, enhancer predictions based on epigenetic profiling alone do not address functional significance. Corroborative evidence presented in Figure 5 sought to mitigate this limitation, particularly with functional genomic assays, such as TAL1 ChIP-exo analysis (Fig. 5F). Still, establishing a causal link between enhancer activation or repression and a corresponding gene activity is a challenging task, particularly on a genome-wide scale. Further confounding this pairwise enhancer–promoter association is the observation made in previous studies that, together, multiple enhancers can fine-tune the expression of a single gene [
]. In agreement, the present study found that on average, a given gene was linked to approximately five discrete enhancer regions. However, whether multiple nonclustered enhancers (i.e., not super-enhancers) synergistically or antagonistically regulate a given gene locus requires a gene-by-gene analysis. Lastly, our study along with previous erythroid enhancer reports [
] inferred enhancer–gene linkages using the nearest gene paradigm. A previous report using human cancer cell lines suggested that as many as 40% of enhancers can skip the nearest gene to loop to a more distant gene [
], but the extent to which this gene skipping occurs in a normal physiological context remains unclear. Although limited by resolution, long-range chromatin interaction assays performed on a genome-wide scale would extend previous locus-specific looping studies in erythroid cells [
] and would complement the current study by furthering our understanding of the dynamics of enhancer repertoire usage during terminal erythroid differentiation.
Here, we used the FVA murine model system of erythropoiesis that uniquely allowed us to investigate how the hormone Epo shapes the enhancer landscape. This study provides novel insights into how hormone signaling pathways are able to control transcriptional programs during cellular differentiation by altering enhancer-associated marks. Together, these findings define a cis-regulatory enhancer network for Epo signaling during erythropoiesis and provide the framework for future studies involving the interplay of epigenetics and Epo signaling.
Acknowledgments
We thank the Vanderbilt Technologies for Advanced Genomics (VANTAGE) Core, Olivia Koues, and Ivo Violovich for technical and computational support with Illumina sequencing. Special thanks also to Scott Beeler for additional computational assistance.
AAP was supported by Vanderbilt Molecular Endocrinology Training Program Grant 5T32 DK07563.
The LIM-only protein Lmo2 is a bridging molecule assembling an erythroid, DNA-binding complex which includes the TAL1, E47, GATA-1 and Ldb1/NLI proteins.
Transcriptional regulation of Elf-1: Locus-wide analysis reveals four distinct promoters, a tissue-specific enhancer, control by PU.1 and the importance of Elf-1 downregulation for erythroid maturation.
30133-9&id=fx1.jpg){kind=link}
30133-9&id=gr1.jpg){kind=link}
30133-9&id=gr2.jpg){kind=link}
30133-9&id=gr3.jpg){kind=link}
30133-9&id=gr4.jpg){kind=link}
30133-9&id=gr5.jpg){kind=link}
30133-9&id=gr6.jpg){kind=link}
30133-9&id=gr7.jpg){kind=link}