Combinatorial assembly of developmental stage-specific enhancers controls gene expression programs during human erythropoiesis

Abstract

Gene-distal enhancers are critical for tissue-specific gene expression, but their genomic determinants within a specific lineage at different stages of development are unknown. Here we profile chromatin state maps, transcription factor occupancy, and gene expression profiles during human erythroid development at fetal and adult stages. Comparative analyses of human erythropoiesis identify developmental stage-specific enhancers as primary determinants of stage-specific gene expression programs. We find that erythroid master regulators GATA1 and TAL1 act cooperatively within active enhancers but confer little predictive value for stage specificity. Instead, a set of stage-specific coregulators collaborates with master regulators and contributes to differential gene expression. We further identify and validate IRF2, IRF6, and MYB as effectors of an adult-stage expression program. Thus, the combinatorial assembly of lineage-specific master regulators and transcriptional coregulators within developmental stage-specific enhancers determines gene expression programs and temporal regulation of transcriptional networks in a mammalian genome.

Figure 1. Comparative Genomic Analyses of Human Erythropoiesis

(A) Fetal and adult CD34+ HSCPs were differentiated into ProEs ex vivo. Cells at matched stages of differentiation were collected for gene expression profiling and ChIP-seq analyses. (B) Expression of human embryonic (ε-), fetal (γ-), and adult (β-) globin mRNAs. Results are means ± SD of at least three independent experiments. (C) Scatter plots of gene expression profiling between fetal and adult ProEs (F5 and A5). Dashed blue lines indicate the 1.5-fold differential expression cutoff to define the F5-high (increased expression in F5 relative to A5 ProEs) or A5-high (increased expression in A5 relative to F5 ProEs) genes. The numbers of differentially expressed genes are indicated. (D) Numbers of differentially expressed genes conditional on changes in expression levels between A5 and F5 ProEs. (E) The genome-wide distribution of the profiled histone marks and TFs. Total numbers of enriched regions in distal promoters (blue), proximal promoters (red), exons (green), introns (purple), and intergenic regions (light blue) were identified (Experimental Procedures). The graph shows the fraction of enriched regions for each histone mark and TF in fetal and adult ProEs, respectively. See also Figures S1 and S2, and Tables S1 and S2.

Figure 2. Chromatin State Maps and TF Occupancy within the Human β-Globin Gene Cluster

ChIP-seq density plots were generated from raw read data and loaded into the UCSC genome browser as custom tracks. Profiles for histone marks, TFs, and DHS in fetal (green) and adult (red) ProEs are shown at the human β-globin gene cluster. The human β-globin locus is depicted at the bottom containing five β-like globin genes (ε, Gγ, Aγ, δ, and β). Dashed vertical line indicates the location of the upstream DHS (1–5) within the LCR. Solid vertical line indicates the TSS of each globin gene.

Figure 3. Promoter Activities in Fetal and Adult ProEs

(A) Promoters were categorized into active, bivalent, repressed, and null promoters. The fraction of each promoter category is shown for both fetal and adult ProEs. (B) mRNA expression values are shown for each promoter category in fetal and adult ProEs. Boxes show median line and quartiles. Whiskers show the boundary (1.5 times of the inter-quartile range from the first or third quartile) to define outliers (red dots). (C) Venn diagrams show genome-wide overlaps between fetal and adult ProEs for each promoter category. (D) Unsupervised hierarchical clustering of all ChIP-seq datasets within the proximal promoter regions (−2kb to +1kb of TSS) between fetal and adult ProEs. Heatmap depicting the Pearson correlation coefficient of ChIP-seq read densities of indicated TFs and histone marks is shown. (E) ChIP-seq density heatmaps are shown for the profiled histone marks and TFs within each promoter category. See also Figure S3 and Table S3.

Figure 4. Identification of Developmental Stage-Specific Enhancers

(A) Active enhancers were identified as those genomic regions harboring H3K4me1, H3K9ac or H3K27ac, DHS, and absence of H3K27me3. Venn diagram shows the overlap between fetal and adult active enhancers. (B) Representative fetal-only, adult-only, and common enhancers are shown. The putative active enhancers are depicted by shaded lines. (C) ChIP-seq density heatmaps are shown for the profiled histone marks and TFs within active (fetal-only, adult-only, and common) and bivalent (marked by H3K4me1, H3K9ac or H3K27ac, and H3K27me3) enhancers in both fetal (upper panel) and adult (lower panel) ProEs. (D) The distribution of fetal-only, adult-only, and common enhancers around F5-high or A5-high genes, compared to all genes. (E) Predicted enhancers were active in reporter assays in primary erythroid cells. Data are means ± SD of three independent experiments. (F) 3C analysis within the representative fetal-only (COL4A5), adult-only (IRF2), and common (ZFPM1) enhancers. The interaction frequency between the anchoring point (black bar) and distal fragments (shaded bar) is shown. Results are means ± SD. ChIP-seq density plots for H3K4me1, H3K4me3, H3K27ac, and DNase-seq (DHS) density plots are shown. See also Figure S4 and Table S4.

Figure 5. Enhancers Control Developmental Stage-Specific Gene Expression Programs

(A) Target genes of HC-fetal-only, HC-adult-only, and HC-common enhancers were compared with genes differentially expressed between F5 and A5 ProEs. The enrichment scores represent the fold change of the number of overlapped genes between enhancer target genes and differentially expressed genes compared to the number expected at random using all genes as background. P-values are calculated by Fisher’s exact test to quantify the significance of the relative bias towards F5-high or A5-high genes using all differentially expressed genes as background. (B) Fractions of genes associated with HC-fetal-only, HC-adult-only, or HC-common enhancers, conditional on changes in expression levels between F5 and A5 ProEs. (C) Heatmap of DNase-seq intensities within the HC-fetal-only, HC-adult-only, or HC-common enhancers in fetal and adult erythroid progenitors, respectively. (D) TF motifs associated with HC-fetal-only or HC-adult-only enhancers. Top 10 most enriched motifs from mammalian TRANSFAC and JASPAR databases are shown. The ratios represent the fold change of the frequency to observe motif targets in HC-adult-only enhancers compared to HC-fetal-only enhancers. P-values are calculated by hypergeometric distribution to compare the motif presence between HC-adult-only and HC-fetal-only enhancers. (E) Expression of IRF2, IRF6, and MYB proteins in F5 and A5 ProEs. (F,G) IRF2 and IRF6 proteins were depleted by lentiviral shRNAs in adult ProEs, respectively. (H) Boxplots show mRNA expression changes of A5-high or F5-high genes upon depletion of IRF2, IRF6, or MYB in adult ProEs, respectively. Boxplots are constructed as described in Figure 3B. (I) Gene set enrichment analysis (GSEA) of A5-high signature genes using the expression array data of shIRF2, shIRF6, or shMYB relative to controls, respectively. (J) ChIP-seq density plot is shown for IRF2 occupancy within the HC-fetal-only and HC-adult-only enhancers. The average read density in unit of RPKM (number of reads per kilo base-pairs per million mapped reads) is shown for the 4kb region surrounding the enhancer summit. (K) Venn diagrams show genome-wide overlaps between IRF2 peaks, fetal, and adult enhancers. (L) Fractions of IRF2-occupied HC-fetal-only, HC-adult-only, and HC-common enhancers are shown, respectively. (M) Representative IRF2-occupied adult-only (AIG1 and ANKFY1) and common (TCTN3) enhancers. ChIP-seq density plots for H3K4me1, H3K4me3, H3K27ac, GATA1, TAL1, IRF2, and DHS density plots are shown. See also Figures S5, S6, and S7, and Tables S5 and S7.

Figure 6. Combinatorial Regulation of Enhancer Functions

(A) K-means clustering of all TF-associated adult enhancers (see Experimental Procedures). The heatmap on the left shows ChIP-seq read density of all TFs used for clustering in adult ProEs. The heatmap on the right shows a side-by-side comparison of ChIP-seq read densities between fetal and adult ProEs within the same enhancer clusters (AE1 to AE4). (B) Target genes of all adult enhancers and each enhancer cluster were compared with genes differentially expressed between F5 and A5 ProEs. The enrichment scores and P-values were calculated as in Figure 5A. (C) Enrichment of selected JASPAR motifs in each enhancer cluster. The fold enrichment represents the frequency to observe motif targets in the selected enhancer cluster compared to randomly selected regions from human genome. (D) ChIP-seq density plot is shown for IRF2 within each enhancer cluster in adult ProEs. (E) Fraction of IRF2-occupied enhancers within each enhancer cluster is shown.

Figure 7. IRF Proteins Contribute to Adult Enhancer Activity

(A) IRF2 or IRF6 interacting proteins identified by proteomic screen in K562 stable cell lines. The number of peptides obtained from at least two independent experiments is shown. (B) Validation of protein interaction by co-IP experiments in K562 stable cell lines. (C) Validation of protein interaction by co-IP experiments in primary adult ProEs. (D) Venn diagrams show genome-wide overlaps between IRF2, GATA1, and TAL1 ChIP-seq peaks. (E) Venn diagrams show overlaps between IRF2, GATA1, and TAL1 ChIP-seq peaks within adult enhancers. (F) ChIP-seq density plots for H3K4me1, H3K4me3, H3K27ac, GATA1, and TAL1, and DHS density plots are shown for both IRF2 and IRF6 loci in fetal (green) and adult (red) ProEs, respectively. Putative adult-specific enhancers are indicated by shaded vertical lines. (G) ChIP-qPCR analysis of representative GATA1/IRF2 co-bound enhancers and GATA1 only bound enhancers in control and IRF2 lentiviral shRNA (sh1 and sh2) transduced adult ProEs. Results are means ± SD of three independent experiments. (H) Model of the temporal regulation of adult-high genes through gene-distal enhancers. See also Figure S8 and Tables S6 and S7.

Figure 8. Gene Regulatory Networks underlying Human Erythropoiesis

(A) Enhancer-to-Gene network. The edges represent the mapping of enhancers to their target genes. The inner ring of the network represents enhancers, colored according to whether they were identified in fetal (green), adult (red), or both (blue) samples. The outer ring represents differentially expressed genes that are upregulated in fetal (green) or adult (red) ProEs. Edges between enhancers and genes are colored according to their originating enhancers. (B) TF-to-Enhancer network. The outer ring represents the same set of enhancers as shown in the inner ring of panel (A). The inner ring contains the seven profiled TFs. Edges extending from these TFs to an enhancer represent the presence of an identified TF binding site in this enhancer region in fetal (green), adult (red), or both (blue) samples. Known protein-protein interactions (purple edges) between the profiled TFs are also indicated. See also Table S8.