Introduction

Fibroblast is one of the most abundant cardiac cell types and plays important roles in normal heart function and pathological heart remodeling at the adult stage1,2. The fibroblast composition in the heart has been analyzed with multiple approaches such as flow cytometry, with a focus on the left ventricle3,4. However, the anatomical location of fibroblasts in the heart along the developmental progression is still unclear. The heart in mice starts to develop into a four-chambered structure at E9.5. The atrial and ventricular chambers are connected by the atrioventricular canal (AVC), a transient structure that develops into the septum and valve cells5,6. Although the main fibroblast population does not develop until E13.57, fibroblast-like cells with the expression of Col1a1 and other fibroblast genes start to develop in the AVC at E9.56.

CM and Vas_EC are the other two major cell types in the heart. Vas_EC develop at about the same time as FB, while CM develop before the appearance of FB at the chambers. FB has been shown to couple with CM through gap junctions and promote their maturation in human three-dimensional micro-tissue consisting of FB, CM, and EC8. In mice, epicardium-derived cells, mostly fibroblasts, have been reported to serve as guidepost cells by secreting the chemokine Slit2, guiding Vas_EC to the correct anatomical locations and promoting their maturation9. However, the lack of direct assessment of FB function in vivo prevents a comprehensive understanding of its role in regulating CM and Vas_EC development.

To study the function of a cell type rather than a gene, toxin genes such as NTR and DTA have been engineered into conditional ablation systems to eliminate specific groups of cells10,11. The DTA system has been used in mice to ablate cardiac progenitor cells and cardiomyocytes at early developmental stages12. Recently, a group of proliferating fibroblasts with the expression of Postn was identified in neonatal hearts. Further elimination of this cell population using the DTA ablation system showed its importance in promoting cardiomyocyte maturation13. Additionally, the DTA system has been used to ablate fibroblasts at the adult stage to assess their function in adverse cardiac remodeling14,15. However, the function of the primary fibroblast population during embryonic and neonatal stages remains unclear.

One of the major functions of FB is to synthesize and secrete extracellular matrix (ECM) proteins, which influence the behavior of other cell lineages2,16. The ECM is heterogeneous and composed of various types, such as collagen, fibronectin, elastin, and laminin, each with different subtypes. For example, collagen has 28 genes in vertebrates17. The ECM plays important roles in heart development and regeneration. For instance, hyaluronan and proteoglycan link protein 1 (Hapln1) is critical in regulating myocardial compaction and heart regeneration in zebrafish18. ECM has also been found to regulate mouse CM proliferation through β1-integrin signaling in a cell co-culture assay19. Moreover, collagen has been shown to be essential for heart regeneration in adult zebrafish and mice20. However, the heterogeneity and function of ECM genes in mouse heart development remain unclear.

In this study, we identified the anatomical patterns of cardiac FB, CM, and Vas_EC during embryonic development. We further examined the ligand-receptor interactions among these cell types and identified collagen as the top signaling pathway. Additionally, we investigated fibroblast function at different embryonic stages using a cell ablation system and analyzed the resulting molecular defects with scRNA-seq. Finally, we performed short-term and long-term fibroblast ablation at the neonatal stage and found that fibroblasts are dispensable for normal heart function development. In summary, this study provides a comprehensive understanding of the role of cardiac fibroblasts in myocardium and coronary vasculature development during embryonic and neonatal stages.

Results

Cellular analysis of cardiac fibroblast anatomical location at different stages

To understand the spatial distribution of cardiac FBs along the developmental progression, we stained the fibroblast marker gene Col1a1 on heart sections at different stages from E11.5 to P3. We observed a strong expression of Col1a1 in epicardial cells at all analyzed stages (Figure 1A). Inside the heart, we found Col1a1-positive cells enriched at the boundary of the atrial and ventricular regions in E11.5 and 12.5 hearts (Figure 1A, S1), indicating their presence in the atrioventricular canal (AVC). By E13.5 and 14.5, the AVC cells had developed into valve structures, and the valve interstitial cells highly expressed Col1a1 (Figure 1A, S1). In the ventricles, a small proportion of Col1a1-positive cells were observed adjacent to the epicardium. In contrast, few Col1a1 signal was found in the atrium at this stage (Figure 1A). At E15.5 and 16.5, both ventricles were filled with Col1a1-positive cells, but the middle of the ventricular septum remained devoid of such cells (Figure 1A, S1). FBs were also found in the atrium, but the signal at the tips of both atrial chambers was sparser (Figure 1A, S1). The distribution of FBs at E17.5 and E18.5 was highly similar to that observed at E16.5 (Figure 1A, S1). Furthermore, we analyzed the Col1a1 signal in P0, P2 and P3 hearts and found dense Col1a1 expression at the valves and in all four chambers (Figure 1A, S1). In summary, this analysis revealed the spatial location of cardiac FBs along the developmental stages. Specifically, we found that the FBs developed first in the AVC and valve, then at the ventricular wall and atrial, and last at the ventricular septum (Figure 1B).

Anatomical location of cardiac fibroblasts and their spatial relationship with cardiomyocytes and endothelial cells at different stages. (A) RNA staining of Col1a1 revealed the spatial pattern of cardiac fibroblasts at different developmental stages. Scale bar=500µm. (B) The development of cardiac fibroblasts can be grouped into four phases. (C) RNA staining analysis of Col1a1, Actn2, and Cd31 revealed the spatial proximity of FB, CM, and EC in E17.5 hearts. (D) Quantification of the number of CMs and ECs that contact with each FB (n=100). Scale bar=500µm and 100µm in the whole heart sections and enlarged sections, respectively.

Next, to analyze the spatial relationship between FB and other cardiac cell types, such as CM and EC, we performed RNA staining for their marker genes Col1a1, Actn2, and Cdh5. Consistent with the previously described pattern, we observed that FB were mainly enriched in valves and began to develop at the edges of ventricular chambers at E13.5. We also found that FB were intermingled with CM populations and spatially adjacent to ECs (Fig. S2A). Next, we stained the three lineage markers on E17.5 and P3 heart sections. We observed consistent results as E13.5 (Fig. 1C, S2B). Given that FB starts emerging at E13.5, we primarily quantified the contacts at E17.5 and P3. The quantification results showed that each FB was in contact with approximately one Vas_EC and four CMs at both stages. (Fig. 1C, D). These results suggest that cardiac FB are physically close to CMs and ECs, likely facilitating signal communications that impact their development.

Molecular analysis of the signaling communications between fibroblast and the other cardiac cell types

To understand the interactions between FB and other cardiac cell types, we analyzed scRNA-seq data from developing CD1 hearts. By quantifying the number and strength of these interactions, we found that each cell type, except blood cells, had a similar number of interactions. However, in terms of interaction strength, FB exhibited the highest values compared to those of other cell types (Fig. 2A). Next, we analyzed the signaling pathways secreted from FB to ventricular cardiomyocytes (Ven_CM) and Vas_EC. Interestingly, we found that the collagen pathway was the predominant pathway in both FB-VenCM and FB-VasEC interactions (Fig. 2B). In addition to collagen, we identified other signaling molecules such as laminin, Ptn, MK, Fn1, and Thbs, which are mostly extracellular matrix (ECM) proteins, as being secreted by FB to function in Ven_CM and Vas_EC (Fig. 2B). Within the collagen pathway, we further analyzed the detailed ligand-receptor interactions and found more interaction pairs in FB-VasEC than in FB-VenCM. Some interactions were shared between the two cell types, but most were distinct (Fig. 2C).

Molecular analysis of fibroblast communications with other cardiac cell types in developing hearts. (A) The number and strength of interactions among different cardiac cell types identified from scRNA-seq data. (B) The top signaling pathways derived from FB that impact Ven_CM and Vas_EC development. (C) The detailed ligand-receptor interactions in the collagen pathway between FB and Ven_CM or Vas_EC. (D) Collagen expression pattern in CD1 mouse hearts at E11.5, E14.5, E16.5, and P3. Scale bar=500µm and 150µm in the whole heart sections and enlarged sections, respectively. (E) The signaling interactions and top pathways between FB and Ven_CM or Vas_EC in primary human hearts.

To better understand the role of ECM genes in FB_VasEC and FB_VenCM interactions, we analyzed their expressions (GO number, 0031012) in the chamber-derived fibroblasts (main_fb) using the CD1 and C57BL/6 scRNA-seq datasets from 18 stages10 (Fig S3A-F). Interestingly, we identified five groups of ECM genes that displayed stage or chamber specific expression in both datasets (Fig S4-6, Supplementary table S2, 3). Specifically, we found one group of genes (G1) that was highly expressed at stages before E14.5. This gene group includes signaling molecules such as Slit2 and Wnt5a, and collagen genes such as Col13a1 and Col26a1 (Fig S2-4). We found the second group of genes (G2) were mainly expressed in the LV and RV from E14.5 onwards with a significant increase at the neonatal stage. This gene group includes collagen genes such as Col5a1, Col6a1, and Col6a2, and metalloproteinase genes such as Adamts13 and Adamts14. Additionally, we identified the third group of genes were highly expressed in the LA and RA at most stages, which includes Col12a1, Adamts1, and Adamts4. Interestingly, Fn1 also exists in this group, suggesting that it has an important function in atrial fibroblast development. Furthermore, we found the fourth group of genes were preferentially expressed in the other three chambers than the RA and were mainly expressed at late neonatal stages after P5 in certain chambers. This group of genes includes Col4a4, Col4a5, Col4a6, and Ntn4. Finally, we identified the fifth group of genes displayed universal expression in all four chambers at most stages. This gene group includes Col1a1, Col1a2, and Postn, which are the ubiquitously expressed extracellular matrix genes (Fig S4-6). In summary, we identified five groups of ECM genes displaying stage- and chamber-specific expression patterns in FB and observed that each gene group includes different collagen pathway members. Additionally, we analyzed the ligand-receptor interactions between FB and Ven_CM or Vas_EC at each stage. We found collagen pathway interactions were also enriched at most analyzed stages (Fig S7, supplementary table S4, 5). These results together suggested the importance of collagen pathway in mediating the FB function in heart development.

Furthermore, we studied collagen deposition dynamics in developing hearts. We used collagen hybridizing peptide (CHP) to stain the total collagen after denaturing it via an antigen retrieval process. Consistent with Col1a1 expression, we observed strong collagen signals in the epicardium at all analyzed stages (Fig. 2D). Interestingly, at E11.5, we observed collagen in both AVC and ventricular chambers. Given that fibroblasts in the chambers have not developed yet at this stage, the observed collagen is likely derived from other cell types, such as endocardial endothelial cells. At E14.5, we observed a strong collagen signal in the ventricular free walls, while the signal in the ventricular septum and atrium was relatively weaker. At E16.5, the collagen signal became brighter and denser throughout the heart, including the ventricular free wall, septum, and atrium. This distribution was slightly broader compared to the fibroblast distribution at the same stage (Fig. 1A). By P3, the collagen signal had further increased and formed stripe patterns (Fig. 2D). In summary, we found that the progression of collagen accumulation largely coincided with the development of FBs.

Lastly, we analyzed the signaling interactions among human cardiac cells using a human embryonic heart scRNA-seq dataset21. Consistent with the mouse analysis results, we identified strong interactions between FB and other cardiac cell types (Fig. 2E). Additionally, we found collagen to be the most enriched pathway in FB-VasEC interactions. However, interestingly, we did not observe the collagen pathway among the top enriched pathways in FB-VenCM interactions. This observation was confirmed in two other human embryonic heart scRNA-seq datasets22,23. These results suggest a partially similar mechanism underlying fibroblast regulation of heart development in both mice and humans.

Functional analysis of fibroblasts in embryonic heart development

To identify the right transgenic mice labeling the developing cardiac fibroblasts, we analyzed the expression of two well-known fibroblast marker genes, Postn and Pdgfra, in the CD1 scRNA-seq data. We observed that Postn was highly expressed in fibroblasts and mural cells. Among fibroblasts, Postn was expressed in all subpopulations (Fig S8Ai). On the other hand, Pdgfra was specifically expressed in fibroblasts, although its expression was relatively lower than Postn, especially in postnatal staged atrial cells (Fig S8Aii). To evaluate their efficiency in labeling fibroblasts, we bred Postn-CreER mice with Rosa-mTmG reporter mice and administered tamoxifen (200µg/g) to pregnant mice at E13.5 (Fig S8Bi). We collected the hearts at E17.5 and observed only a small proportion of GFP-positive cells in the valves and septum (Fig S8Bii), suggesting a low labeling efficiency from the Postn-CreER mice. In contrast, we conducted the same experiments using Pdgfra-CreER mice and observed a significant proportion of GFP-positive cells in both the valves and heart chambers (Fig 3A). These results demonstrate the efficient labeling capability of the Pdgfra-CreER; Rosa-mTmG mice.

Functional analysis of FB in embryonic heart development. (A) Lineage tracing revealed an abundant of Pdgfra-CreER labeled cells were GFP+ in developing mouse hearts. Scale bar = 500µm. (B) The ablated embryos treated with tamoxifen at E13.5 were smaller than control embryos, but the ablated hearts were not obviously different from the control hearts. Scale bar=1mm. (C) Representative CD31 staining images in control and ablated hearts. The zoomed-in images showed a reduction of CD31 positive area in the ablated hearts compared to controls. (D) Quantification of the thickness of compact and trabecular myocardium (19 sections from 7 ablated hearts and 20 sections from 7control hearts), and CD31 positive areas in control and ablated hearts (15 sections from 5 ablated and 15 sections from 5 control hearts). (E) The ablated mice and hearts treated with tamoxifen at E15.5 were smaller than controls. (F) CD31 staining analysis of control and ablated hearts. (G) Quantification of the compact and trabecular myocardium thickness, and CD31 positive areas in control and ablated hearts (12 sections from 4 ablated hearts and 12 sections from 4 control hearts.). Scale bar =500µm in the whole heart sections and 100 µm in the zoomed in images. * represents p<0.05; ** represents p<0.01.

To study the function of fibroblasts in heart development, we utilized a genetically encoded DTA system to ablate them, which involves the encoding of a toxin gene that induces apoptosis and allows for the ablation of target cells. Initially, we crossed floxed Rosa-DTA mice with Pdgfra-CreER mice and administered tamoxifen to pregnant female mice at E10.5. After 3 days, we harvested the embryos and observed that all the embryos with ablated fibroblasts (Pdgfra-CreER+/-; Rosa-DTA+/-) had perished, while the control embryos remained unaffected. This outcome indicates an essential role for fibroblasts in early-stage embryo development.

Subsequently, we administered tamoxifen to the female mice from the breeding pair at E13.5 and collected the embryonic hearts at E16.5 (Fig 3B). Our findings revealed that the ablated embryos exhibited smaller sizes compared to the control embryos, although the size of the ablated hearts did not significantly differ from that of the control hearts (Fig 3B). Additionally, we conducted TUNEL analysis to assess cell death in both control and ablated hearts and observed a significantly higher number of TUNEL-positive dots in the valve and chamber regions of the ablated hearts compared to the control hearts, confirming the successful ablation of fibroblasts in the experimental hearts (Fig S9A, B). Furthermore, we examined the anatomical structure of the hearts but did not identify obvious defects in the ablated groups compared to control groups. However, we observed a reduction in endothelial cell density in the ablated hearts, as indicated by the CD31 staining signal (Fig 3C). Additionally, we measured the thickness of the compact and trabecular myocardium in both the left and right ventricles. Our observations revealed that, compared to the control hearts, the ablated hearts exhibited a thinner left ventricular (LV) compact myocardium and a thicker LV trabecular myocardium. Furthermore, we noted an elevated ratio of LV trabecular to compact myocardium in the ablated hearts compared to the control hearts (Fig 3D). However, no similar differences were observed in the right ventricle (RV). Finally, we examined cell proliferation post-ablation using pHH3 staining. We observed a reduction in pHH3-positive cells in the ventricular region, but not in the atrial region, of the ablated hearts compared to the control hearts (Fig S9B). These results suggest that cell proliferation in ventricle was impaired following fibroblast ablation.

Next, we performed ablation at late embryonic stages by treating the pregnant mice with tamoxifen at E15.5 and harvested the embryos at E18.5 (Fig S10A). We did not observe obvious differences in the size of the embryos and hearts between the ablated group and the control group (Fig S10B). However, the CD31 staining of heart sections revealed large holes in the ablated atrium but not in the ablated ventricles or control hearts (Fig S10C). Additionally, quantification of the compact and trabecular myocardium thickness in the LV and RV revealed a thinner LV compact myocardium (Fig S10D). We also calculated the ratio of trabecular to compact myocardium thickness in the LV and RV and found that the ratio was increased in the RV but not in the LV of the ablated hearts compared to control hearts (Fig S10E). These results indicated that both the LV and RV of the ablated hearts have subtle defects. Moreover, we analyzed cell proliferation by staining pHH3 and found no significant differences between the ablated and control ventricles (Fig S10F).

Given that the animals did not show obvious morphological changes after one dose of tamoxifen treatment, we increased the tamoxifen treatment to three consecutive days, ranging from E15.5 to E17.5 (Fig 3E). We harvested the embryos at E18.5 and found that both the ablated embryos and their hearts were notably smaller in size compared to the controls (Fig 3E). Next, to confirm the ablation efficiency, we analyzed Col1a1 expression with RNA staining and found that the ablated hearts lost almost all of their Col1a1 signal in the four chambers (Fig S10G), indicating a high ablation efficiency with three doses of tamoxifen treatments. By conducting a staining analysis of CD31 in the control and ablated heart sections, we discovered that the ablated hearts exhibited a significantly lower endothelial cell density (Fig 3F, G). Further quantification of the LV compact and trabecular myocardium thickness showed a significant reduction in the ablated hearts compared to control hearts (Fig 3G). However, interestingly, the ratio of trabecular to compact myocardium thickness in the LV did not differ between the control and ablated hearts, indicating a proportional reduction in both compact and trabecular myocardium thickness in the ablated hearts. Finally, we found that the ablated ventricles had less pHH3 signal than the control ventricles (Fig S10F), indicating impaired cell proliferation in the ventricles following fibroblast ablation. These results together indicated that fibroblast plays an important role in embryonic heart development at all the analyzed stages.

Molecular analysis of the defects in the ablated embryonic hearts

Next, we performed scRNA-seq to analyze the defects in the ablated embryonic hearts. We utilized MULTI-seq to pool the control and ablated hearts from E16.5 and E18.5 (Fig. 4A). After filtering out low-quality cells based on our QC standards, we grouped the single cells from two experimental replicates using unsupervised methods and found a high degree of overlap (Fig. 4B, S11A). We then identified cell types and cell cycle phases within the single cells. In total, we identified eight cell types: Atrial_CM, Ven_CM, Vas_EC, Endo_EC, Mural_cell, FB, Epicardial cell, and Immune cell, along with three cell cycle phases (G1, S, G2M) in each cell type (Fig. 4C). Cells from both control and ablated conditions at E16.5 and E18.5 were present in each cell type. We further quantified the percentages of each cell type across conditions and found that the percentage of FB was significantly reduced in the ablated hearts at both stages, indicating effective ablation. Interestingly, we also observed a reduction in mural cells in the ablated hearts at both stages (Fig. 4D). For the other cell types, we did not identify consistent changes across conditions (Fig. 4E). For example, the percentage of Vas_EC decreased in ablated samples at E16.5 but exhibited inconsistent changes at E18.5. Additionally, while the percentage of Ven_CM increased at E18.5, it showed different trends at E16.5 (Fig. 4E). This inconsistency may be attributed to varying cell capture rates in each sample.

Single cell analysis of the embryonic heart defects after fibroblast ablation. (A) Diagram of the MULTI-seq experiments to profile control and ablated hearts at two developmental stages. (B) The scRNA-seq data from two replicates are highly consistent. (C) UMAP plots of scRNA-seq data labeled by cell type, cell cycle phase, and genotypes. (D) Quantification of the FB and mural cell percentages at each conditions. (E) The percentages of other cardiac cell types under each condition. (F) Detailed analysis of the fibroblast population revealed four groups, including a dying fibroblast subpopulation. (G) UMAP plot of Ven_CMs from different conditions. (H) Quantification of the different cell cycle phased Ven_CMs under each condition. (I, J) Heatmap and pathway enrichment of genes that are differentially expressed in control and ablated Ven_CMs at E16.5. (K) UMAP plot of Vas_EC labeled by conditions. (L, M) Heatmap and pathway enrichment of genes that were upregulated in ablated Vas_ECs compared to control Vas_ECs.

Next, we selected FB for further study of its cellular heterogeneity. By analyzing the expression of cardiac lineage genes, cell cycle phases, and sample sources, we identified four FB subpopulations (Fig. 4F). Specifically, we identified a group of valve fibroblasts (valve-fb) expressing valve mesenchymal cell genes such as Hapln1, a group of neural crest-derived fibroblasts (Neural_crest_fb) expressing Sox10, and a group of chamber-derived main fibroblasts (main_fb) expressing Wt1 and Tcf21 (Fig. S11B). All these FB groups contained cells from both control and ablated samples. In contrast, we identified another group of cells exclusively derived from the ablated samples, which expressed cell death-related genes (Fig. S11C). We named this group dying fibroblasts (dying_fb).

In addition to FB, we analyzed Ven_CM and found that cells from control and ablated hearts were highly overlapped on the UMAP plot (Fig. 4G). We quantified the proportions of cells at each cell cycle phase and found that the ablated samples had slightly higher percentages of cells in the G1 phase compared to control samples (Fig. 4H, S11D). We then performed differential gene expression analysis of control and ablated samples at E16.5 and identified a set of differentially expressed genes (Fig. 4I, supplementary table S6, 7). Gene ontology analysis of the upregulated genes in the ablated CMs showed enrichment in various pathways, including heart development, sarcomere organization, heart tube morphogenesis, and cell proliferation. In contrast, the downregulated genes were enriched in pathways related to chaperone-dependent protein folding, mitochondrial membrane permeability, fatty acid import, and mitochondrial ATP transmembrane transport (Fig. 4J). These results indicated that ablated Ven_CMs upregulated genes associated with heart development and cell proliferation while downregulating genes related to cell maturation, including those involved in mitochondrial function and cell metabolism.

Furthermore, we selected Vas_EC for further analysis. We did not observe distinct populations between control and ablated samples on the UMAP plot (Fig. 4K, S11E). However, through further differential expression analysis, we identified a group of genes with upregulated expression in the ablated samples (Fig. 4L, supplementary table S8). Gene ontology analysis of these genes revealed that they were primarily involved in immune response pathways, including type I and II interferon production, interleukin-27-mediated pathways, NLRP3 inflammasome assembly, and cellular responses to exogenous dsRNA (Fig. 4M). These results suggest that Vas_EC exhibits a strong immune response following fibroblast ablation.

The signaling pathway defects in the ablated embryonic hearts

Next, we studied the signaling defects in the ablated hearts. Through the analysis of signaling pathways between main_fb and Vas_EC or Ven_CM, we found that the collagen pathway consistently ranked at the top in both control and ablated samples (Fig. 5A). We then stained for collagen levels using collagen hybridizing peptide and observed that collagen deposition was dramatically reduced in both the septum and ventricular wall of the ablated hearts compared to the control hearts (Fig. 5B).

Identification of signaling defects in ablated hearts. (A) The top signaling pathways between FB and Ven_CM or Vas_EC in control and ablated hearts. (B) Collagen staining revealed a significant reduction in collagen deposition in ablated hearts compared to control hearts. (C, D) Regulatory analysis predicted the ligands that regulate the genes differentially expressed in control and ablated Ven_CMs and Vas_EC. (E) Representative signaling pathways that were downregulated or upregulated in the dying FBs.

Moreover, we identified signaling ligands that potentially regulate the differentially expressed genes in control and ablated Ven_CM. Interestingly, the collagen pathway, which has four different collagen ligands on the list, emerged as a major regulator of these genes (Fig. 5C). In addition to collagens, we also identified other ligands, such as TGF-beta, IGF, and BMP, which were derived from FB and regulate the expression of differentially expressed genes in control and ablated Ven_CMs (Fig. 5C). We also predicted the receptors in Ven_CM for each ligand (Fig. S12A). Furthermore, we predicted the ligands in FB that regulate the differentially expressed genes in Vas_EC and identified several signaling molecules. These included the DTA pathway-related signal Hbegf, the inflammation-associated signal Il33, Wnt ligands (Wnt11, Wnt4, and Wnt5b), and collagen Col4a6 (Fig. 5D, S12B).

Finally, we examined the signals in the dying_fb cells. Compared to normal main_fb and valve_like_fb, the dying_fb downregulated certain signaling pathways, such as the IGF pathway, EphA pathway, Tenascin pathway, and collagen pathway (Fig. 5E, S13A). Interestingly, they also upregulated several signaling pathways, including the EGF pathway, BMP pathway, VEGF pathway, and FGF pathway (Fig. 5E, S13A). Understanding the role of these altered signaling pathways from dying_fb in contributing to heart defects will be an interesting area for future research.

Functional analysis of fibroblasts at neonatal stages

Next, we assessed the fibroblast function at neonatal stage. We treated the mice with tamoxifen at P1 and harvested them at P4 (Fig S14A). We did not observe obvious differences between the ablated and control hearts (Fig S14B). However, CD31 staining of heart sections revealed defects in the atrium, which had a thinner and segmented atrial wall (Fig S14C). We further quantified the compact and trabecular myocardium in both the LV and RV but did not observe differences between control and ablated hearts (Fig S14D). We also did not find differences in the ratio of trabecular to compact myocardium thickness between control and ablated hearts in the LV and RV (Fig S14E). Finally, we analyzed cell proliferation by staining pHH3 and did not observe differences between the two conditions (Fig S14F). Next, we increased the tamoxifen treatments to three times, ranging from P1 to P3, and collected hearts at P4 (Fig 6A). The heart morphology and size did not reveal obvious differences between control and ablated groups (Fig 6B). The CD31 staining revealed a reduction in endothelial cell density in the ablated hearts when compared to the control hearts (Fig 6C). These reductions in both LV and RV were further confirmed through quantifications (Fig 6D). Finally, we quantified the chamber thickness and cell proliferation. Similar to hearts treated with one dose of tamoxifen, we observed no apparent defects and cell proliferation changes in ablated hearts with three days of tamoxifen treatments (Fig 6D, S14G). Furthermore, we stained the P4 hearts for collagen depositions. We observed that the ablated hearts had a sparser collagen signal and lost most of the stripe pattern in the septum compared to the control hearts. However, a significant amount of stripe-patterned collagen was still observed in the ventricular free wall, although the overall signal was sparser (Fig 6E). These collagens were likely derived from fibroblasts before ablation.

Short term and long term functional analysis of fibroblasts at neonatal stage. (A) Diagram of the experiment to ablate fibroblasts with three doses of tamoxifen treatment from P1 to P3. Scale bar=1cm. (B) No obvious size difference was observed between control and ablated hearts. (C) Representative images of control and ablated hearts stained with CD31. (D) Quantification of the compact and trabecular myocardium thickness and CD31 positive areas in ventricular at control and ablated hearts (12 sections from 4 ablated hearts and 12 sections from 4 control hearts.). (E) Collagen accumulation in control and ablated hearts at P4. Scale bar=500µm and 150µm in the whole heart sections and enlarged sections, respectively. (F) Tamoxifen was given from P3 to P5 to ablate the fibroblasts and hearts were collected at P17 for analysis. (G) The percentage of mice died at different days. (H) Representative echo images of the control and ablated hearts. (I) Quantification of the heart rate, body weight, and heart function in control and ablated mice at P17. * represents p<0.05; ** represents p<0.01; *** represents p<0.001.

Lastly, we conducted a long-term analysis by administering tamoxifen to the mice from P3 to P5 and harvesting them at P17 or P18 (Fig. 6F). We observed that approximately two-thirds of the ablated mice had perished before P17 (Fig. 6G), and the surviving mice were significantly smaller than the control mice. Quantification revealed reduced heart rates and body weights, but a similar heart-to-body ratio in the ablated mice compared to controls (Fig. 6I). A detailed examination of the major tissues showed reduced organ sizes in the ablated mice, including smaller hearts, lungs, brains, kidneys, and livers (Fig. 6F, S16A). Additionally, we assessed mouse heart function using echocardiography before sacrificing them. Interestingly, we found that heart function, including ejection fraction (EF) and fractional shortening (FS), was not reduced in the ablated mice compared to controls (Fig. 6H, I). Furthermore, we evaluated collagen levels in the ablated hearts and did not observe a clear reduction or pattern change (Fig. S15). This finding contrasts with our observations at P4 (Fig. 6E) and suggests potential collagen expression compensation at later stages. In contrast, we found that the lungs of the ablated mice exhibited reduced collagen levels and larger empty spaces compared to the controls (Fig. S16B). These findings suggest that fibroblasts play a crucial role in the growth of multiple tissues during the neonatal phase. Additionally, it will be interesting to explore the mechanisms underlying the preservation of heart function when fibroblasts are eliminated in future studies.

Discussion

In this study, we investigated the cellular and molecular interactions between FB and CM or Vas_EC, identifying collagen as the primary signaling molecule. We also analyzed the role of FB in embryonic heart development using a cell ablation system. Additionally, we employed scRNA-seq to identify the molecular defects in CM and Vas_EC within the ablated hearts. Finally, we examined the function of FB in neonatal heart development through both short-term and long-term ablation experiments. Altogether, our study represents the first comprehensive analysis of fibroblast function in heart development. The insights we gained from the study will not only enhance our understanding of the role of this major cardiac cell type in heart development but also shed light on its function in adult cardiac remodeling, as well as other tissue development and regeneration.

Through RNA staining of the three lineage genes, we found that each fibroblast directly interacts with approximately one Vas_EC and four CMs. Given that the quantification was conducted on 2D tissue sections, the number of contacts may be underestimated. Considering that Vas_EC and FB emerged around the same time, it will be important to understand in the future how their connections are established at the single-cell level and what these connections mean for the development of each lineage. Two studies in mice and zebrafish have shown that epicardial cells or epicardium-derived cells regulate Vas_EC development through the expression of Hapln1 or Slit29,18. It would be interesting to test whether collagen or other fibroblast-derived ECM proteins have similar roles. Additionally, in vitro culture experiments have demonstrated that fibroblasts interact with CMs through gap junction proteins or integrins to regulate their proliferation8,19. Future studies should aim to clarify how these interactions are established in vivo and how they evolve during developmental progression.

Collagen, as the major ECM family, consists of 28 members17. In addition to collagen, there are many other types of ECM genes. Our scRNA-seq data revealed the temporal and spatial expression patterns of ECM genes in cardiac FBs. Future research should investigate their differential functions. For example, we found that most ECM genes in group 4 upregulate their expression by postnatal day 7, when mouse hearts begin to lose their regenerative capacity, suggesting a potential role in heart regeneration. This warrants further investigation. Additionally, while we have analyzed the collagen pathway as a whole in this study, our scRNA-seq analysis indicated that collagen genes in different groups exhibit distinct expression patterns. It will be essential to explore their unique roles in heart development in future studies. Furthermore, we found that collagen is the primary signaling molecule in FB_CM and FB_VasEC pairs in mice; however, we only observed it in FB_VasEC pair in humans. This discrepancy may be due to the loss of collagen receptors in human CMs, which merits further exploration.

The DTA-based cell ablation system has been used to ablate cardiac progenitors and cardiomyocytes before. It revealed that a robust regeneration response occurred after cell ablation at early stages12. However, in our study, based on the pHH3 staining, we did not observe an increase in cell proliferation after cell ablation. This could be due to multiple possibilities. First, the Pdgfra-CreER; Rosa-DTA-based ablation may be too severe to initiate the regeneration response. Considering that ablation at E13.5 and ablation with three doses of tamoxifen at E15.5 both led to smaller embryos, the ablation of a large number of cells may compromise the regeneration/compensation mechanisms. However, it is difficult to explain the embryos at E15.5 with one dose of tamoxifen treatment, which did not display obvious defects after ablation and have no increase in cell proliferation. Second, the ablation stages may be relatively late, and fibroblasts have lost their regeneration capability. Third, fibroblasts may be different from other cardiac cells and have a limited regeneration response to their own ablation, although they can respond to CM loss at neonatal and adult stages.

ScRNA-seq results revealed the upregulation of inflammation-related pathways in Vas_EC in the ablated hearts. It will be important to investigate whether these pathways are initiated by the ablation process or the loss of FBs in the future. Given that these pathways were primarily upregulated in Vas_EC but not in CM, this suggests a cell type-specific response. Studying the function of these signals in vascular development will be essential. Moreover, pathway analysis revealed both upregulation and downregulation of signaling pathways associated with dying fibroblasts. Further studies on these pathways are important, as they could offer valuable insights into the cell death process during normal heart development and heart injury and repair.

Fibroblast ablation at the neonatal stage led to two-thirds of the mice experiencing a lethal outcome, while the remaining mice developed a dwarf phenotype. Functional analysis of the mouse hearts showed preserved cardiac function post-ablation, suggesting that the lethal phenotype may be due to defects in other tissues, such as the lungs, which exhibited disrupted structures. Additionally, we conducted functional measurements at P18 under anesthesia in a small cohort of mice and found that the ablated mice demonstrated improved heart function (data not shown). This finding indicates that fibroblasts may play a role in regulating heart function during the neonatal stage, aligning with observations in adult mice. In the adult stage, cardiac function was found to be superior in residential fibroblast-ablated mice compared to control mice following myocardial infarction or angiotensin II/phenylephrine (AngII/PE) treatments14.

Methods

Experimental methods

Mouse strains

The animal experiments have been approved by the University of Pittsburgh Institutional Animal Care and Use Committee (IACUC). CD1 male and female mice were purchased from Charles River Laboratories and bred in our laboratory to generate embryos and neonatal pups at specific stages for RNA staining experiments. The transgenic mice, including Rosa26-mTmG (Strain #:007676)24, Pdgfra-CreERT2 (Strain #:032770)25, Postn-CreER (Strain #:029645)26, and ROSA26-eGFP-DTA (Strain #:032087)10, were ordered from the Jackson Laboratory.

Tamoxifen treatment and mouse dissection

The male and female mice of specific strains were bred together. The default dosage of 200µg of tamoxifen per gram of body weight (200µg/g) was given to pregnant mice through oral gavage, and neonatal mice were given 10µg/g of tamoxifen by direct injection into their stomach to induce Cre activity13. The pregnant mice and neonatal pups were euthanized using CO2 and decapitation-based methods, respectively. The mouse hearts were isolated following the standard procedure described previously and were directly embedded in OCT without fixation for RNA staining or fixed at 4% paraformaldehyde for other experiments, such as immunofluorescence staining.

Proximity ligation in situ hybridization (PLISH)

PLISH was performed following a published protocol27. Specifically, the embryonic or postnatal mouse hearts were embedded in OCT (Sakura, 4583) without fixation. After sectioning at a thickness of 10μm, the tissue sections were treated with post-fix medium (RNase-free PBS with 3.7% formaldehyde and 0.1% DEPC) followed by 0.1mg/ml pepsin treatment (RNase-free H2O with 0.1mg/ml pepsin and 0.1M HCl). After dehydration, the sections were sealed with hybridization chambers (Invitrogen, S24732) and hybridized with H probes (supplementary table S1) in Hybridization Buffer (1 M NaTCA, 5 mM EDTA, 50 mM Tris pH 7.4, 0.2 mg/mL Heparin). Next, after being treated with circularization reaction and rolling cycle amplification, the samples were hybridized with detection probes conjugated with Cy3 or Cy5 fluorophore. Finally, the samples were stained with DAPI (Invitrogen, D1306), mounted with fluoromount-g (SouthernBiotech, OB100-01), and imaged under confocal microscopy (Leica TSC SP8).

RNAscope Multiplex Fluorescent V2 Assay

The RNAscope Multiplex Fluorescent Reagent Kit v2 (Advanced Cell Diagnostics, 323270) was performed according to the manufacturer’s manual. Briefly, the tissue sections were fixed in pre-chilled 4% PFA (Electron Microscopy Sciences,15710) at 4℃ for 1 hour. After progressive dehydration, the sections were sequentially treated with hydrogen peroxide for 10 minutes and protease IV for 15 minutes (for embryonic) or 20 minutes (for postnatal) at room temperature (RT). After that, the samples were hybridized to gene-specific Z probes for 2 hours at 40℃ using the HybEZ II Hybridization system (ACD, 321721). Following further signal amplification, the hybridization signals were detected with TSA Vivid Fluorophores. The samples were then stained with DAPI and mounted with ProLong Gold Antifade Mountant (Invitrogen, P36930), and imaged under confocal microscopy (Leica TSC SP8). The RNAscope probes used in the study are Mm-Col1a1-C2 (319371-C2), Mm-Actn2 (569061), and Mm-Cdh5-C3 (312531-C3).

Quantification of the interacting cells

Based on the RNA scope staining results, the direct contacts between fibroblasts and cardiomyocytes as well as endothelial cells was quantified. Given that the fibroblast starts to emerge at E13.5, we mainly quantified the contacts at E17.5 and P3. Approximately 100 fibroblasts from the left ventricular wall at each stage were selected, and the number of Actn2+ and Cdh5+ cells surrounding each fibroblast was counted respectively.

Immunofluorescence staining

The immunofluorescence staining was performed following a standard procedure. Briefly, mouse hearts were fixed in 4% PFA overnight before embedding in OCT. Afterwards, the samples were sectioned at 10µm and briefly washed in PBS to clean the OCT. The sections were then blocked for 1 hour in blocking buffer (10% goat serum, 1% BSA, 0.1% Tween 20) and incubated with primary antibodies in the primary antibody buffer (1% BSA in PBST) at 4℃ overnight. On the second day, the samples were stained with fluorophore-conjugated secondary antibodies in blocking buffer for 1 hour at room temperature. Finally, the samples were stained with DAPI, mounted with fluoromount-g, and imaged with a confocal microscope. The primary antibodies used in the study include anti-CD31 (BD, # 550274), and anti-pHH3-488 (abcam, #ab197502). The secondary antibodies used are Goat Anti-Rabbit IgG-488 (ThermoFisher, # A-11008) and Goat Anti-Rabbit IgG-647 (ThermoFisher, #A21247).

Collagen staining

Fresh frozen sections from embryonic and postnatal mouse hearts were fixed in 4% PFA in 1X PBS for 15 minutes at room temperature (RT). After that, the sections were soaked in a sodium citrate-based solution (distilled water with 10mM sodium citrate and 0.5% Tween-20, HCl was added to adjust pH to 6.0) for antigen retrieval, which was performed in a steamer for 30 minutes. Subsequently, the sections were incubated with 20 µM biotin-conjugated collagen hybridizing peptide (Advanced Biomatrix, 50-196-0307) in 1X PBS at 4℃ overnight. The next day, streptavidin-cy5 (Invitrogen, SA1011) was used (1:500 in 1X PBS with 1% BSA) for 1 hour at RT. Finally, the samples were stained with DAPI, mounted with FLUOROMOUNT-G, and imaged under a confocal microscope.

TUNEL staining

TUNEL staining was performed following the manufacturer’s protocol (Roche, 11684795910). Briefly, mouse hearts were embedded and sectioned as described in the immunofluorescence staining procedure. After fixation with 4% Paraformaldehyde in PBS for 20 minutes and three washes with PBS, the sections were incubated in permeabilization solution (0.1% Triton X-100, 0.1% sodium citrate) for 2 minutes on ice, followed by two PBS rinses. Subsequently, 60 µl of TUNEL reaction mixture was added to each section, and they were incubated in a humidified atmosphere for 1 hour at 37℃ in the dark. Finally, the sections were washed with PBS three times, stained with DAPI, mounted with fluoromount-g, and imaged with a confocal microscope.

Single cell mRNA-sequencing experiments

ScRNA-seq was performed following the MULTI-seq procedure previously reported28,29. Briefly, tamoxifen-treated E16.5 and E18.5 Pdgfra-CreER; Rosa-DTA mouse embryos were harvested on the same day. Since the ablated mouse embryos at both stages were smaller than the control embryos, we selected one control and one ablated embryo from each stage for heart dissection. The genotype of the selected embryos was further validated using PCR. The hearts, with four chambers, were dissociated into single cells and labeled with MULTI-seq barcodes. The pooled libraries were then loaded onto the 10X Genomics Chromium iX and profiled using the single-cell 3’ V3.1 kit. The generated libraries were sequenced on the Illumina NovaSeq X Plus platform. The experiments were repeated twice.

Echo imaging

Echocardiography was performed using a standard protocol30. Briefly, heart function of awake mice at P17 was assessed using the Vevo 2100 micro imaging platform (FUJIFILM Visual Sonics Inc., Canada). The measurements were conducted at the level of papillary muscle in M mode. Heart rates, left ventricular end-diastolic diameter (LVEDD), and left ventricular end-systolic dimensions (LVESD) were calculated from the 2D short-axis view.

Statistical analysis

To quantify pHH3 signal, compact, and trabecular myocardium thickness, we used two to three sections from each heart and more than three hearts per genotype. The exact number of embryos profiled can be found in each figure legend. Statistical analyses were conducted using ImageJ and Prism 9 software, employing a two-tailed Student’s t-test to compare groups. P-values below 0.05 were considered significant.

Data analysis

Data Availability

The newly generated scRNA-seq datasets have been deposited in the Gene Expression Omnibus (GEO) under the accession number GSE272048. Please use this secure token (sdcxiewmxfyztar) to review the data. The scRNA-seq datasets from CD1 and C57BL/6 mouse hearts were generated in a previous study and can be downloaded from the GEO database using the accession number GSE19334629. The human scRNA-seq data were obtained from GEO or the European Genome-Phenome Archive (EGA) under the accessions GSE106118 and EGAS000010039962123 for CellChat analysis.

ScRNA-seq data analysis

ScRNA-seq data were mapped to the mouse genome mm10 using CellRanger and further de-multiplexed with the R package deMULTIplex28. Cell types in the newly generated scRNA-seq data were annotated using lineage genes previously published by our team29. We performed quality control by removing cells with fewer than 200 genes and those with more than 40% mitochondrial content. Experimental batch integration was conducted using Harmony 1.2.031, followed by unsupervised clustering with Seurat V532. Next, we utilized the Cell Cycle Scoring function in Seurat to annotate the cell cycle phases for each individual cell. Finally, we applied the FindMarkers function in Seurat with default settings to identify genes that were differentially expressed in Ven_CM or Vas_EC under control and ablation conditions.

Extracellular matrix genes expression analysis

We downloaded the extracellular matrix genes from the Jackson Laboratory under the term “extracellular matrix” and gene ontology ID 0031012 (https://www.informatics.jax.org/go/term/GO:0031012; Download date: 2023-07-07). After downloading, we cleaned the gene list by removing duplicates and filtering out genes with null expression, resulting in a total of 440 genes. Next, we assessed the enrichment of each gene in the main population of fibroblasts at each stage and chamber using the R package AUCell33, with heatmaps drawn with the R package ggplot2 (Fig 9, S9, S10).

Ligand-receptor interaction analysis

To identify ligand-receptor interactions across stages and zones for fibroblasts and cardiomyocytes, we analyzed relevant subsets of the CD1 data, including the main population of FB, ventricular CM, and Vas_EC using CellPhoneDB v3.3434. We used a p-value threshold of 0.2 and ran the analysis with 10 threads, while keeping the rest of the parameters at their default settings. The significant mean value of all interactive partners (log2) and enrichment p-values (-log10) obtained from the CellPhoneDB outputs were plotted as dotplots in R.

Cellchat analysis

The Seurat object containing mouse CD1 single-cell mRNA sequencing data was converted into a CellChat object in R. Default settings in CellChat 1.6.135 were used to identify the number and strength of interactions among cell types. FB were designated as sender cells, while Ven_CM and Vas_EC were designated as receiver cells to identify signaling between these related cell types. The same default settings were applied to analyze signaling pathways among human cardiac cell types. For pathway analysis in control and DTA-ablated cardiac cells, their respective scRNA-seq objects were utilized.

NicheNet analysis

Nichenetr 2.1.536 was used for the analysis with default settings. First, we identified differentially expressed genes between control and ablated cells. Due to the varying number of recovered cells in each cell type under different conditions, we selected different staged samples for comparison and assessed their consistency across the remaining samples. For Ven_CMs, we used the E16.5_control_1 and E16.5_DTA_1 samples, while for Vas_EC, we included all four E16.5 samples. Next, we designated the control sample as the background condition and the ablated sample as the experimental condition, using main_fb as the sending cells and Ven_CM or Vas_EC as the receiving cells to predict their regulatory potential.

Acknowledgements

We’d like to thank all the members in the Li laboratory for their insightful discussions of this work. We are thankful for Dr. Wei Feng for his help on the MULTI-seq experiment. This research was supported in part by the University of Pittsburgh Center for Research Computing, RRID:SCR_022735, through the resources provided. Specifically, this work used the HTC cluster, which is supported by NIH award number S10OD028483.

Sources of Funding

This work was supported by R00HL133472 and DP2HL163745 from the NIH and the CMRF grant from the University of Pittsburgh.

Author Contributions

Y.D., Y.H., G.L. designed the experiments; Y.D. analyzed the ablated hearts with different assays; J.X. bred the mice, treated them with tamoxifen, and harvested the hearts from embryonic and neonatal mice. Y.H. performed the RNAScope and collagen staining experiments, Y.H. and G.L. performed the PLISH experiments; Y.D. performed the scRNA-seq experiment. H.T. and G.L. analyzed the scRNA-seq data; M.Z. and J.X. has performed the Echo scanning and data analysis. Y.D., Y.H., H.T., and G.L. prepared the manuscript; All the authors edited the manuscript.

Competing interests

The authors declare no competing interests.