FOXA2 Is Required for Enhancer Priming during Pancreatic Differentiation

Abstract

Transcriptional regulatory mechanisms of lineage priming in embryonic development are largely uncharacterized because of the difficulty of isolating transient progenitor populations. Directed differentiation of human pluripotent stem cells (hPSCs) combined with gene editing provides a powerful system to define precise temporal gene requirements for progressive chromatin changes during cell fate transitions. Here, we map the dynamic chromatin landscape associated with sequential stages of pancreatic differentiation from hPSCs. Our analysis of chromatin accessibility dynamics led us to uncover a requirement for FOXA2, known as a pioneer factor, in human pancreas specification not previously shown from mouse knockout studies. FOXA2 knockout hPSCs formed reduced numbers of pancreatic progenitors accompanied by impaired recruitment of GATA6 to pancreatic enhancers. Furthermore, FOXA2 is required for proper chromatin remodeling and H3K4me1 deposition during enhancer priming. This work highlights the power of combining hPSC differentiation, genome editing, and computational genomics for discovering transcriptional mechanisms during development.

Keywords: ATAC-seq and chromatin accessibility; FOXA1; FOXA2; GATA6; PDX1; enhancer priming and activation; hESCs; hPSCs; human embryonic stem cells; human pluripotent stem cells; nucleosome remodeling; pancreatic progenitors.

Figure 1.. Different Stages during Pancreatic Differentiation Are Associated with Distinct TF Motifs

(A) Pancreatic differentiation efficiency was verified by examining lineage markers using flow cytometry. ES, embryonic stem cells; DE, definitive endoderm; FG, posterior foregut; PP1, primary pancreatic progenitor; d, day(s). See STAR Methods for a detailed differentiation protocol. (B) Principal-component analysis (PCA) of ATAC sequencing (ATAC-seq) for the top 3,000 variable peaks. (C) TF motif enrichment during stage transitions. The numbers of opened and closed peaks during stage transitions are indicated. Opened or closed peaks between successive stages were compared with the total atlas to examine the TF motif enrichments using the one-sided Kolmogorov-Smirnov (KS) test. The KS test effect size is shown on the y axis, and the proportion of peaks associated with the TF motif is plotted on the x axis. The size of each circle represents the odds ratio (OR), which was defined as the frequency of the TF in an opened or closed group divided by its frequency in the entire atlas. TF motifs with a KS test effect size ≥ 0.05 (indicated by the dashed lines) and OR ≥ 1.2 are shown; when there were more than 10 such motifs, only the top 10 are shown. TF motifs derived from species other than Homo sapiens are marked with brackets. GATA motifs (GATA2, GATA3, GATA4, and GATA6) are marked in red, the HNF1B motif is in blue, and FOXA motifs (FOXA1, FOXA2, and FOXA3) are in green. See also Table S1. (D) Proportion of ATAC-seq peaks with the GATA6, FOXA2, and HNF1B motif are shown. (E) HOMER motif analysis was performed during the stage transition. The top 5 known motifs are shown after removing redundant motifs.

Figure 2.. Pancreatic Progenitor Stage-Specific Accessible Chromatin Regions Are Enriched with FOXA Motifs

(A) ATAC-seq stage-specific groups; the chromatin accessibility was compared among 4 stages (ES, DE, FG, and PP1), and stage-specific groups were defined based on their patterns of chromatin status across stages. The number of peaks for each group are annotated. ChIP-seq profiles are shown in the same order as ATAC-seq peaks. See also Table S2. (B) Composition of accessible peaks at the PP1 stage. The number of peaks for each group is annotated. (C) TF motif enrichment in the PP1-specific group. The hypergeometric test was used to compare the enrichment of proportions of TF motifs for the PP1-specific group (foreground ratio) versus those for the total atlas (background ratio). The number of peaks containing TF motifs in the PP1-specific group and the total atlas, as well as the size of the PP1-specific group and the total atlas, are supplied. The 15 most highly enriched TFs are shown with the foreground ratio to the right; the horizontal axis shows the binomial Z score, representing the number of SDs between the observed count of PP1-specific peaks containing a TF motif and the expected count based on the background ratio. The Bonferroni-corrected hypergeometric p value of the bottom-most TF is shown with the red vertical line. TF motifs derived from species other than Homo sapiens are marked with brackets. GATA motifs (GATA2, GATA3, GATA4, and GATA6) are marked in red, the HNF1B motif is in blue, and FOXA motifs (FOXA1, FOXA2, and FOXA3) are in green. See Figure S1B and Table S3 regarding other stage-specific groups. (D) TF motif enrichment results using the HOMER algorithm among DE-containing groups; DE specific, DE-FG, and DE-PP1. The top 5 known motifs are shown after removing redundant motifs. (E) Summary of motif enrichment for each stage-specific group. See also Figure S1.

Figure 3.. Requirements for FOXA2 in Pancreatic Progenitor Cell Specification

(A) HUES8 FOXA2 KO lines were generated at the ES stage and analyzed by RNA-seq, ATAC-seq, and ChIP-seq throughout pancreatic differentiation. d, day(s). See also Table S4. (B) Loss of protein from HUES8 FOXA2 KO cells was verified by western blotting at the FG stage. (C) Expression of stage-specific lineage markers in HUES8 FOXA2 KO was examined by flow cytometry. SOX17 expression at the DE stage, HNF1B expression at the FG stage, and PDX1 and NKX6–1 expression at the PP1 and PP2 stages are shown. The differentiation experiments were repeated four times with three independent clonal lines per genotype. Results from clonal lines of the same genotype were combined (n = 12), and results from each experiment are shown, together with mean ± SD. Statistical analysis was performed by unpaired two-tailed Student’s t test. ^****p < 0.0001; ns, not significant (p ≥ 0.05). Results from individual lines are plotted in Figure S2D. (D) Representative plot of flow cytometry are shown for PDX1⁺ NKX6–1⁺ expressions in WT and FOXA2 KO at the PP1 and PP2 stages. (E) PDX1⁺ NKX6–1⁺ expressions at the PP1 and PP2 stages were examined by immunostaining. Scale bar indicates 100 μm. (F) Principal component analysis (PCA) of the top 890 variable genes between WT and FOXA2 KO when controlling for the stage effect by two-factor modeling. n = 3, independent clonal lines. See also Table S5. (G) Hierarchical clustering of the top 890 variable genes from (F). See also Figure S2.

Figure 4.. Pancreatic TF Bindings Are Enriched in Less Accessible Regions in FOXA2 KO

(A) Cumulative plot of opening (19,363) and closing (21,925) peaks from the ES-to-PP1 transition. (B) Tornado plot of less accessible peaks in FOXA2 KO shown for the DE, FG, and PP1 stages. FOXA2 and GATA6 ChIP-seq in less accessible peaks is shown in the same order as ATAC-seq. See also Table S6. (C) Venn diagram showing the overlap between FOXA2 and GATA6 binding sites at the PP1 stage in WT cells. (D) Average of normalized counts visualized as a metapeak. Metapeaks show FOXA2 binding in WT cells and GATA6 binding in WT and FOXA2 KO cells at less accessible peaks. The maximum value of each y axis is annotated in TPM. See also Figures S3 and S4.

Figure 5.. FOXA2 Is Required to Acquire Proper Chromatin Accessibility during Pancreatic Differentiation

(A) FOXA2 KO effects on ATAC-seq stage-specific groups are shown. The average of normalized counts is visualized as a metapeak, and the maximum value of y axis is indicated in the plot as TPM. The WT signal is marked in black and the FOXA2 KO signal is in red in the tornado plots. (B) ATAC-seq insert size of WT and FOXA2 KO for DE-PP1-, FG-PP1-, and PP1-specific groups is visualized. The density of insert size is shown by aggregating based on the distance from the ATAC-seq peak summit. (C) ATAC-seq read pairs are visualized as tornado plots, divided by the insert size: nucleosome-free regions < 150 bp and nucleosome-associated regions ≥ 150 bp. See also Figure S5.

Figure 6.. H3K4me1 Deposition during Pancreatic Differentiation

(A) Tornado plot and metapeak of the ATAC-seq stage-specific group show the profile of H3K4me1 and H3K27ac histone modifications in WT and FOXA2 KO cells. The maximum value (TPM) of the y axis is indicated in the plot. The WT signal is marked in black and the FOXA2 KO signal is in red in the tornado plots. (B) Primed PP1-specific ATAC-seq group was defined by excluding peaks that contained the H3K27ac signal at the DE or the FG stage (n = 1,136). See also Table S7. Metapeak analysis shows FOXA2, H3K4me1, ATAC-seq, and the H3K27ac signal of the primed PP1-specific group. WT is marked in black, and FOXA2 KO is in red. The maximum value (TPM) of the y axis is annotated in the plot. (C) Beeswarm plot to show the expression of 686 genes associated with the primed ATAC-seq group. Z scores of normalized counts are plotted, with each cell type averaged over replicates. Five WT samples, PP1 samples, and PP2 samples are grouped by themselves for comparison. See Figure S6C for the beeswarm plot of all cell types. The graph shows p values for significant increases in gene expression levels when comparing PP1 versus DE and PP2 versus DE using one-sided unpaired Wilcoxon rank-sum tests. There are also significant increases when comparing PP1 versus FG (p = 4.639 × 10⁻¹⁸), PP2 versus FG (p = 1.603 × 10⁻²⁷), and PP2 versus PP1 (p = 6.329 × 10⁻⁶). See also Figure S6.

Figure 7.. FOXA2 Is Required for the Establishment of Primed Enhancer Mark H3K4me1 Deposition

(A) Comparison of normalized levels of H3K4me1 at the primed PP1-specific groups. FOXA2 binding regions were defined by the location of reproducible FOXA2 ChIP-seq peaks (irreproducible discovery rate [IDR] ≤ 0.01). The values of histone modification in WT and FOXA2 KO cells were compared by a one-sided unpaired Wilcoxon signed rank test with Bonferroni correction. The center line in the boxplots shows the median; the box limits are the upper and the lower quartiles, respectively; and whiskers are defined as 1.5 times the interquartile range above and below the box limits. All outliers outside the whiskers are shown as points. (B) Fold change of the H3K4me1 signal in FOXA2 KO compared with WT for the primed PP1-specific group. Peaks were grouped by signal strength quintile of the reproducible FOXA2 ChIP-seq peaks at each stage. At all 3 stages, stronger FOXA2 binding showed a more negative shift in the H3K4me1 signal. The H3K4me1 signals at the top 20% and bottom 20% of FOXA2 peaks were compared at each stage using a one-sided unpaired Wilcoxon signed rank test (p values as annotated). (C) Model of the role of FOXA2 in enhancer priming during pancreatic differentiation. During enhancer priming, FOXA2 recruitment facilitates H3K4me1 deposition. This is followed by recruitment of additional lineage TFs, such as GATA factors (GATA6 and GATA4), and acquisition of the H3K27ac activation mark.