The multi-lineage transcription factor ISL1 controls cardiomyocyte cell fate through interaction with NKX2.5

Abstract

Congenital heart disease often arises from perturbations of transcription factors (TFs) that guide cardiac development. ISLET1 (ISL1) is a TF that influences early cardiac cell fate, as well as differentiation of other cell types including motor neuron progenitors (MNPs) and pancreatic islet cells. While lineage specificity of ISL1 function is likely achieved through combinatorial interactions, its essential cardiac interacting partners are unknown. By assaying ISL1 genomic occupancy in human induced pluripotent stem cell-derived cardiac progenitors (CPs) or MNPs and leveraging the deep learning approach BPNet, we identified motifs of other TFs that predicted ISL1 occupancy in each lineage, with NKX2.5 and GATA motifs being most closely associated to ISL1 in CPs. Experimentally, nearly two-thirds of ISL1-bound loci were co-occupied by NKX2.5 and/or GATA4. Removal of NKX2.5 from CPs led to widespread ISL1 redistribution, and overexpression of NKX2.5 in MNPs led to ISL1 occupancy of CP-specific loci. These results reveal how ISL1 guides lineage choices through a combinatorial code that dictates genomic occupancy and transcription.

Keywords: ISL1; NKX2.5; cardiac development; cardiac progenitor; cell specification; combinatorial code; transcription factor motifs; transcription factors; transcriptional regulation.

Figure 1

ISL1 is required for cardiac cell-fate-specific gene regulation (A) Schematic of CM differentiation from hiPSCs. (B) Western blot of ISL1 expression during hiPSC-CM differentiation (days 0–12). Vinculin is the loading control. (C) Hierarchical UMAP clustering of WT (n = 3 independent experiments) or ISL1^−/− (n = 3 independent experiments) day 8 CPs after scRNA-seq. (D and E) Expression of HEY1 (D) and MYH6 (E) in day 8 CPs, superimposed on the UMAP from (C). (F) Monocle pseudotime analysis of day 8 CPs. Colors represent arbitrary units of pseudotime. (G) FeaturePlot display of the distribution of WT and ISL1^−/− day 8 CPs across clusters depicted in (C). (H) Percentages of each genotype in each cluster. Cell numbers for each genotype are normalized to the total number of cells in each cluster. (I) Heatmap of genes and accompanying biological GO terms significantly dysregulated in scRNA-seq of ISL1^−/− as compared with WT day 8 CPs.

Figure 2

scRNA-seq analyses of ISL1 function in MNPs (A) Schematic of MN differentiation from hiPSCs. (B) Western blot of ISL1 in WT and ISL1^−/− day 18 MNPs. Vinculin served as the loading control. (C) Immunofluorescence of neural factors SMI32 and NKX6.1 in day 18 MNPs. Scale bar, 100 μM. (D) Hierarchical UMAP clustering of WT (n = 3 independent experiments) and ISL1^−/− (n = 3 independent experiments) scRNA-seq data in day 18 iPSCs differentiated toward MNPs. (E and F) Expression levels of the MN-related TF NKX6-1 (E) and the IN-related TF IRX3 (F) in day 18 MNPs, superimposed on the UMAP from (D). (G) Heatmap of genes enriched with accompanying GO terms in each of the MNP or IN clusters depicted in (D). (H) Replicate (Rep) comparison in WT (n = 3; 14,609 cells) or ISL1^−/− (n = 3; 6,372 cells) day 18 MNPs across the clusters depicted in (D). (I) Heatmap clustering of genes significantly dysregulated in scRNA-seq of ISL1^−/−compared with WT day 18 MNPs. (J) Venn diagram of significantly downregulated genes shared between ISL1^−/− day 8 CPs and ISL1^−/− day 18 MNPs. (K) Venn diagram of significantly upregulated genes shared between ISL1^−/− day 8 CPs and ISL1^−/− day 18 MNPs.

Figure 3

Grammar of distinct ISL1 binding patterns in cardiac or MNPs based on nearby motifs using deep learning (A) Venn diagram of ISL1 peaks in CPs or MNPs. Overlap is statistically significant (Regione R) (p = 0.001). (B) Histogram of ISL1-bound peaks based on distance from gene TSSs. (C) Example ChIP-seq tracks showing ISL1-bound loci with CP-specific, MNP-specific, or shared peaks. Corresponding de novo motifs identified by BPNet that predict cell-type-specific ISL1 binding are shown below. Data shown from a single representative replicate. (D) De novo motifs identified by BPNet after simultaneous training of the CP-specific or MNP-specific binding profiles. Displayed for each identified motif is: the total count of motifs and percentages, name based on the most likely TF(s) that binds them, the position weight matrices (PWMs), and average ChIP-seq intensity and average contribution to the ISL1 binding predictions as extracted from BPNet.

Figure 4

Proteomics and BPNet identify NKX2.5 and other TFs that guide ISL1 DNA-binding specificity in CPs (A) Schematic for ISL1 co-IP and subsequent SRM targeted proteomics analyses. (B) Barplot results from targeted affinity purification followed by MS, depicting enrichment of individual peptides corresponding to cardiac factors, positive controls (ISL1, LDB1), and negative controls (LHX3, PHOX2A) (n = 1 independent experiment). (C) Cooperative relationship among all de novo identified motifs binned by distance from ISL1 binding predicted by in silico motif mutation.

Figure 5

NKX2.5 binds a majority of ISL1 peaks and cooperates with ISL1 to regulate gene expression in CPs (A) Venn diagram of NKX2.5 and ISL1 ChIP-seq peaks in WT day 6 CPs. (B) Example ChIP-seq tracks showing loci bound by ISL1, NKX2.5, or both in day 6 CPs. Putative NPPA/NPPB enhancer labeled as related to Figure S5I. Data shown from single representative replicate. (C) Hierarchical UMAP scRNA-seq clustering of WT, ISL1^−/−, and NKX2.5^−/− day 8 CPs (n = 3 independent experiments). (D) Expression of MYL9 in day 8 CPs, superimposed on the UMAP from (C). (E) FeaturePlot display of the distribution of WT and ISL1^−/− day 8 CPs across the clusters depicted in (C). (F) Heatmap of significantly dysregulated genes with biological GO terms in day 8 ISL1^−/− or NKX2.5^−/− CPs compared with WT.

Figure 6

ISL1 depends on NKX2.5 for DNA localization at a subset of cardiac loci (A) Venn diagram of ISL1-bound DNA regions that are lost, maintained or gained in NKX2.5^−/− day 6 CPs. (B) Example of an NKX2.5-dependent ISL1-bound peak. Data shown from a single representative replicate. (C) Venn diagram of ISL1-bound loci from (A) compared with GATA4-bound loci. (D) Quantitative analysis of change in ISL1 ChIP peak intensity in WT CPs compared with ISL1 binding in NKX2.5^−/− day 6 CPs. Regions with statistically significant change using DESeq are in red (p ≤ 0.05); points in gray are not significant. Fold change of ISL1 binding (NKX2.5^−/−/NKX2.5^+/+) is graphed across ISL1 ChIP signal intensity in WT CPs. (E) ChIP-seq track of ISL1 detailing increased intensity of ISL1 binding in the absence of NKX2.5 at the NRP1 locus. Data shown from single representative replicate. (F) Violin plot of NRP1 expression in WT, NKX2.5^−/−, and ISL1^−/− day 8 CPs (n = 3 independent experiments for each). (G) Schematic of NKX2.5 experimental overexpression in MNPs. (H) ChIP-seq tracks of an ISL1 peak gained within the MYH11 locus when NKX2.5 was overexpressed in MNPs. Data shown from a single representative replicate.