Cell-Type-Specific Gene Regulatory Networks Underlying Murine Neonatal Heart Regeneration at Single-Cell Resolution

Abstract

The adult mammalian heart has limited capacity for regeneration following injury, whereas the neonatal heart can readily regenerate within a short period after birth. Neonatal heart regeneration is orchestrated by multiple cell types intrinsic to the heart, as well as immune cells that infiltrate the heart after injury. To elucidate the transcriptional responses of the different cellular components of the mouse heart following injury, we perform single-cell RNA sequencing on neonatal hearts at various time points following myocardial infarction and couple the results with bulk tissue RNA-sequencing data collected at the same time points. Concomitant single-cell ATAC sequencing exposes underlying dynamics of open chromatin landscapes and regenerative gene regulatory networks of diverse cardiac cell types and reveals extracellular mediators of cardiomyocyte proliferation, angiogenesis, and fibroblast activation. Together, our data provide a transcriptional basis for neonatal heart regeneration at single-cell resolution and suggest strategies for enhancing cardiac function after injury.

Keywords: angiogenesis; cardiomyocyte proliferation; epicardial cells; fibrosis; injury response; myocardial infarction; open chromatin landscape; secreted factors; single cell ATAC-seq; single cell RNA-seq.

Figure 1.. Single-Cell RNA Sequencing (scRNA-Seq) Reveals Cellular Heterogeneity in Neonatal Hearts

(A) Schematic of experimental design for the single-cell analyses. MI or sham surgeries were performed on P1 and P8 hearts. Ventricles below the ligation plane (indicated by the dashed lines) were collected at 1 or 3 days post-surgery for scRNA-seq, and 3-day post-surgery samples were collected for scATAC-seq. n = 8–12 animals were used for tissue collection. n = 1 sequencing library was generated for each time point and condition. (B and C) UMAP plots showing single-cell transcriptomes analyzed in the study, color coded for cell clusters (B) or sample origin (C). For each time point, one scRNA-seq library was generated using pooled tissues dissected from 8–12 individual animals to control the differences among individual animals and dissection variations. (D) Stacked violin plots showing expression of marker genes for each cluster. Cell clusters are color coded according to the UMAP plot in (B). (E) Percentage of each cell type over all nonmyocytes in each scRNA-seq sample. Endothelial cells (ECs) contain Art.EC, VEC1, VEC2, VEC3, Pro.EC, and Endo clusters. Fibroblasts (FBs) contain FB1, FB2, FB3, FB4, and Pro.FB clusters. Immune cells contain macrophage, DC-like, monocyte, and Gra clusters. Art.EC, arterial endothelial cell; CM, cardiomyocyte; dpi, days postinjury; dps, days post sham; Endo, endocardial cell; Epi, epicardial cell; Gra, granulocyte; Pro.EC, proliferating endothelial cells; Pro.FB, proliferating FB; SMC, smooth muscle cell; VEC, vascular endothelial cell. See also Figure S1 and Tables S1 and S2.

Figure 2.. scATAC-Seq Reveals Open Chromatin Landscapes of Single Cells in Neonatal Hearts

(A and B) UMAP plots showing single-cell open chromatin profiles analyzed in the study, color coded for cell clusters (A) or sample origins (B). For each time point, one scATAC-seq library was generated using pooled tissues dissected from 8–12 individual animals to control the differences among individual animals and dissection variations. (C) Heatmap showing activity of top enriched open chromatin peaks and their associated genomic coordinates for each cell cluster. (D) scATAC-seq tracks showing open chromatin peaks associated with cell-type-specific genes, including Myh6 and Myh7 (CM markers), Pecam1 (EC marker), Postn (FB marker), Myh11 (SMC marker), Pdgfrb (pericyte marker), Msln (epicardial cell marker), Fcgr1 (macrophage marker), Cd3g and Cd3d (T cell markers), and Ms4a1 (B cell marker), across different cell clusters. (E) Transcription factor motif enrichment (upper row), gene accessibility (middle row), and gene expression profiles (bottom row) for lineage-specific transcription factors Nkx2–5 (CM specific, first column), Ets1 (EC enriched, second column), Tcf21 (FB specific, third column), and Ebf1 (SMC/Pericyte enriched, fourth column). Art.EC, arterial endothelial cell; CM, cardiomyocyte; dpi, days postinjury; dps, days post sham; Endo, endocardial cell; Epi, epicardial cell; FB, fibroblast; SMC, smooth muscle cell; VEC, vascular endothelial cell. See also Figures S2 and S3 and Tables S1 and S2.

Figure 3.. Decoding Cell-Type-Specific Gene Regulatory Networks in Regenerative and Non-regenerative Injury Responses

(A) Bar graph showing the number of MI-induced open chromatin regions in CM, EC, and FB at 3 days after P1 or P8 MI. (B) Bar graph showing the number of MI-repressed open chromatin regions in CM, EC, and FB at 3 days after P1 or P8 MI. (C) Differential accessibility analysis of scATAC-seq data revealed MI-induced open chromatin peaks from CMs (left column), ECs (middle column), and FBs (right column). Upper row: Venn diagrams showing the number of open chromatin regions that are induced upon MI and specific to P1 hearts, or P8 hearts, or shared, in CM, EC, or FB cell clusters. Middle row: top enriched motifs identified from open chromatin regions that are induced upon MI at P1 (left) or P8 (right). Bottom row: examples of injury-induced open chromatin peaks that are preferentially accessible in P1 (left) or P8 hearts (right) for CM, EC, or FB cell clusters. See also Figure S3.

Figure 4.. Distinctive Epicardial Response during Neonatal Heart Regeneration

(A) Schematic illustration for the bioinformatics pipeline used to deconvolute bulk RNA-seq data using scRNA-seq data as a reference. (B) Heatmap showing expression of cell-type-specific, injury-induced genes across different cell types in P1+1 dpi scRNA-seq data. (C) Venn diagram showing numbers of epicardial injury-induced genes that are specific to or shared among P1+1 dpi, P1+3 dpi, P8+1 dpi, and P8+3 dpi time points. (D) Table summarizing TF motifs enriched in open chromatin regions that are associated with epicardial-specific, MI-induced genes at P1+1 dpi. (E) Injury-induced ligand-receptor interaction networks at P1+1 dpi (2,814 total interaction pairs), P1+3 dpi (474 total interaction pairs), P8+1 dpi (2,103 total interaction pairs), and P8+3 dpi (3,800 total interaction pairs). Size of the nodes denotes the total number of predicted ligand-receptor interactions initiated from the indicated cell types. The color and width of the lines connecting two nodes denote the number of predicted ligand-receptor interactions between connected cell types. The network illustrates the potentiality of cellular crosstalk events mediated by injury-induced ligand-receptor pairs but does not account for the anatomic position of different cell types. See also Figure S4 and Tables S4 and S5.

Figure 5.. Endothelial Cell Heterogeneity in Neonatal Hearts

(A) UMAP representation of different EC sub-populations analyzed. (B) Heatmap showing expression of top enriched genes for each EC sub-population. (C–F) UMAP plots showing expression of Art.EC marker gene Fbln5 (C), endocardial marker gene Npr3 (D), cell-cycle-related gene Mki67 (E), and capillary EC marker gene Gpihbp1 (F) in ECs. (G) UMAP plot showing expression of Rspo1 among all cells analyzed, which is restricted to the epicardial cell population. (H) Heatmap showing fold induction of Rspo1 expression at various time points after P1 or P8 MI detected by bulk RNA-seq. (I) scATAC-seq tracks showing open chromatin landscape of Rspo1 genomic region across various cell types. Positions of predicted KLF14 (blue), TFAP2A (orange), EGR1 (green), and E2F4 (red) binding sites from previous analysis (Figure 4D) are also shown. (J) Representative images showing in vitro human umbilical cord VEC (HUVEC) tube formation after 8 h of 200 ng/mL bovine serum albumin (BSA) (negative control), 10 or 100 ng/mL recombinant RSPO1, or 10 or 100 ng/mL recombinant VEGF treatment, cultured in Matrigel, with quantifications showing the number of branching points per field under each treatment (n = 15 per each group; ****p < 0.0001). Scale bar, 500 μm. See also Figure S5.

Figure 6.. Heterogeneous Fibroblast Populations in Injured and Uninjured Neonatal Hearts

(A) UMAP representation of different FB sub-populations analyzed. (B) Heatmap showing expression of top enriched genes for each FB sub-population. (C) Percentage of each FB sub-population over total nonmyocytes within each sample. (D–H) UMAP plots showing expression of Pro.FB-enriched gene Hmgb2 (D), FBI-enriched gene Cxcl1 (E), FB2-enriched gene Dlk1 (F), FB3-enriched gene Nov (G), and FB4-enriched gene Fbln5 (H). (I) Heatmap showing relative fold induction (Z score) of Serpinb2, Wnt5a, and Ltbp3 expression at various time points after P1 or P8 MI detected by bulk RNA-seq. (J) EdU incorporation (magenta) and vimentin immunofluorescent staining (green) of NRCF cells treated with 200 ng/mL BSA (negative control), 20 ng/mL recombinant SERPINB2, 100 ng/mL recombinant LTBP3, or 100 ng/mL recombinant WNT5A (positive control), with quantification showing the proportion of EdU-positive cells among vimentin-positive cells (fibroblasts) (n = 4 per each group; ****p < 0.0001). Scale bar, 100 μm. See also Figure S6.

Figure 7.. Heterogeneity and Distinct Injury Response of Immune Cells

(A) UMAP representation of different myeloid cell populations analyzed. (B) Heatmap showing expression of top enriched genes for each myeloid cell sub-population. (C–G) UMAP plots showing expression of Cd86, an M1-macrophage marker gene (C); Arg1, an M2-macrophage marker gene (D); Plac8, a monocyte-enriched gene (E); Ly6c2, an M1 monocyte-enriched gene (F); and Naaa, a DC-like cell enriched gene (G). (H) Percentage of each immune cell population over total nonmyocytes within each sample. (I) UMAP plot showing expression of Clcf1 among different cell types present in the P1+1 dpi and dps scRNA-seq data. (J) Representative images showing EdU incorporation (magenta) and cardiac troponin T (cTnT) immunofluorescent staining (green) of NRVM cells treated with 200 ng/mL BSA (negative control), 20 ng/mL recombinant insulin-like growth factor 2 (IGF2) (positive control), or different concentrations of recombinant CLCF1 (5, 10, 25, and 50 ng/mL). Scale bar, 100 μm. (K) Quantification showing the proportion of EdU-positive cells among cTnT-positive cells (CM) after treatment of BSA, CLCF1 recombinant protein, or recombinant IGF2 (n = 4 per each group; **p < 0.01, ****p < 0.0001). See also Figure S7.