Maturation of heart valve cell populations during postnatal remodeling

Mature heart valves in vertebrates consist of valve interstitial cells (VICs) and hematopoietic-derived cells (Hulin et al., 2018) embedded in a stratified extracellular matrix (ECM) surrounded by valve endothelial cells (VECs) (Hinton et al., 2006). At birth, heart valves are not yet fully mature, and postnatal ECM remodeling in humans and mice generates distinct ECM layers (Aikawa et al., 2006; Hinton et al., 2006). Heart valve remodeling includes elongation of valve leaflets with increased collagen and elastin production, leading to formation of collagen-, proteoglycan- and elastin-rich layers by postnatal day (P)30 in mice and during childhood in humans (Aikawa et al., 2006; Hinton et al., 2006). Although embryonic heart valve development is well studied, regulatory mechanisms of postnatal valve ECM remodeling are relatively unknown. Moreover, postnatal VICs transition from minimally proliferative to increasingly quiescent at the transcriptional level in healthy adults (Hinton et al., 2006). In addition to VICs, other cell lineages are present in developing and mature valves, including multiple hematopoietic cell lineages that change after birth (Hulin et al., 2018). Nevertheless, a comprehensive analysis of cell composition and heterogeneity in postnatal heart valve remodeling has not previously been reported.

Endothelial, interstitial, and immune cell lineages have been identified in heart valves, but diversity and subpopulations within these lineages are not well characterized, and VIC differentiation profiles in remodeling heart valves after birth have not been defined. In mature pig aortic valves, VECs on aortic and ventricular surfaces have distinct gene expression profiles, but it is not known when these differences arise during development or how heterogeneous VECs are in actively remodeling valves (Simmons et al., 2005). Diverse populations of hematopoietic-derived cell types, including macrophages and dendritic cells, are present in developing valves, but the full complement of infiltrating cells has not been previously determined (Anstine et al., 2017; Hulin et al., 2018; Visconti et al., 2006). Pigmented cells with characteristics of melanocytes also are evident in mature and diseased mouse valves, but definitive gene expression profiles of these cells have not previously been reported (Balani et al., 2009; Hulin et al., 2017; Mjaatvedt et al., 2005). Finally, it has long been suspected that VICs, the primary ECM-producing cells within heart valves, are heterogenous, with distinct subpopulations producing specific types of ECM in the developing valve strata (Liu et al., 2007), but data on VIC subpopulations and diversity in gene expression profiles in vivo are still missing. Specifically, it is not known whether distinct VIC cell-types are responsible for production of collagen-, proteoglycan- and elastin-rich ECM layers during postnatal heart valve remodeling. Moreover, it is likely that additional cell types are present in remodeling heart valves, necessitating a full unbiased characterization of cell types based on gene expression at the single cell level during the postnatal period.

Using droplet sequencing (DropSeq), we performed single cell RNA sequencing (scRNA-Seq) on heart valves at two distinct developmental stages, comparing premature valve primordia at P7 with fully stratified leaflets at P30, to define an atlas of heart valve cell diversity, and to identify key cell subsets involved in postnatal collagen production and ECM remodeling. Our study confirms and uncovers distinct subpopulations of endothelial, immune and melanocyte cells present both at P7 and P30. For the first time, postnatal differentiation of VICs is shown, with identification of collagen- and proteoglycan-expressing VICs at P7 that are transcriptionally distinct from P30 VICs, which include multiple different subpopulations.

To study heart valve single cell transcriptomes and subpopulations during postnatal valve maturation, two developmental stages were chosen: P7 and P30. At P7, the formation of a collagen layer starts to be detected in valve leaflets (Fig. 1A, black arrowheads), in contrast to its low expression at P1 (Amofa et al., 2017), but the primitive leaflets remain thickened (Fig. 1A). At P30, valve leaflets are elongated with regionalized distribution of fibrillar collagen and proteoglycan (Fig. 1A) (Amofa et al., 2017). Single cell isolations from aortic and mitral valve leaflets were obtained from valves dissected from P7 and P30 mouse pups, which were pooled and dissociated to obtain single cell suspensions. DropSeq was then performed to generate gene expression profiles from individual cells. After initial sequence analysis and exclusion of outliers and low expressing cells, 18,702 transcripts from 594 P7 cells and 2246 P30 cells were analyzed (Fig. 1B). The scRNA-Seq data were subjected to the Iterative Clustering and Guide-gene Selection (ICGS) algorithm from the open-source software AltAnalyze (Fig. 1A) (Magella et al., 2018; Olsson et al., 2016; Salomonis et al., 2009), which allows identification of cell populations based on highly intra-correlated genes in each cell cluster defined by de novo guide genes. ICGS of the entire scRNA-Seq dataset resolved nine clustered cell populations (Fig. 1B). P7 cells are present within five cell clusters and P30 cells are present within seven cell clusters (Fig. 1B). Although cell clusters 1 to 3 are clearly transcriptionally distinct, cell clusters 4 to 9 display some similar guide-genes but with diverse expression levels (Fig. 1B). This initial clustering of cells clearly separates immune, melanocyte, endothelial and interstitial subpopulations of cells, based on guide gene identity in remodeling valves.

The transcriptional separation of different cell clusters is further validated by a t-Distributed Stochastic Neighbor Embedding (t-SNE) dimensional reduction that indicates similarities between individual cells (Fig. 1C). This t-SNE analysis indicates four major distinct cell populations (Fig. 1C). By superimposing expression of known heart valve cell markers in the t-SNE plot, clusters can be distinguished as immune (Cx3cr1, cluster 1), melanocyte (Dct, cluster 2), endothelial (Pecam1, cluster 3) and interstitial cells (Col1a1, clusters 4 to 9) (Fig. 1C,D). We used the MarkerFinder function in AltAnalyze (Fig. 2A) to identify the top enriched genes for each heart valve cell cluster (Fig. 2B). Enriched genes are analyzed with GO-Elite to identify enriched biological processes (Fig. 2A). As expected, the immune cell cluster exhibits functional processes in defense and immune response, whereas endothelial cells are involved in angiogenesis and locomotion (Fig. 2A). Melanocytes are a small but distinct population characterized by melanin synthesis and pigmentation through expression of Dct, Mitf, Tyr and Pax3 (Fig. 2A,B). During postnatal development, melanocytes undergo differentiation illustrated by decreased Dct expression (Fig. S1A-D), while melanin production is increased from P7 to P30 in AoV and MV leaflets (Fig. S1E-H). Interstitial cells of clusters 4 to 9 express the ECM genes Col1a1, Vim, Col3a1, Vcan, Bgn and Lum, which is indicative of enriched processes of collagen fibril organization and cell adhesion (Fig. 2A). Although closely related when compared with endothelial cells, immune cells and melanocytes (Figs 1C and 2A), interstitial cells display increased heterogeneity relative to other cell types. Interestingly, heart valve cells cluster in the four major cell populations independently of their diverse developmental origins (Fig. 2C). The relative proportions of melanocytes, interstitial, endothelial and immune cells are similar at P7 and P30; however, the percentage of interstitial cells is increased at P30, which is reflected in decreased proportions of the other cell types (Fig. 2D).

As our single cell preparation is a mixture of female and male heart valve cells, we assessed possible sexual dimorphism in gene expression. Cells were segregated as female if expressing Xist or as male if expressing Ddx3y, Eif2s3y, Erdr1, Kdm5d or Uty. Unspecified cells were classified based on expression of neither or both classes of genes (Skelly et al., 2018). Female and male cells are equally distributed in cell clusters (Fig. S2A,B), and similar clustering patterns are observed at P7 (Fig. S2C) and P30 (Fig. S2D) between male and female cells. When comparing overall gene expression between female and male for each developmental stage, the vast majority of differentially expressed genes in males versus females are below a twofold change, with a median relative fold change close to 1.1 (Fig. S2E,F). Thus, heart valve cell clustering is not significantly altered by sex. Together, DropSeq and unbiased ICGS identify major heart valve cell populations as immune, melanocyte, endothelial and interstitial cells with relatively more heterogeneity and developmental maturation of interstitial cells.

To explore the cell clusters in more detail, immune, endothelial and interstitial cell clusters were further analyzed. ICGS of immune cells, illustrated by Ptprc expression (Fig. 3A), highlights five distinct immune cell subpopulations (Fig. S3A), which display different gene signatures (Fig. 3B,C) and enriched biological processes (Fig. S3A). The five populations have gene expression characteristics of T cells (Cd24a, Trbc1, Cd3d, Cd8a), dendritic cells (DCs) (Klrb1, Xcr1), mast cells (Kit, Cpa3, Tpsab1) and macrophages, including two distinct subsets: Mφ1 (Mrc1) and Mφ2 (Ccl2, Ccl4) (Fig. 3C,D). Mφ2 exhibit higher expression of the chemokines Tnf and Cd74 (Fig. 3D), whereas Mφ1 express genes involved in macrophage differentiation (e.g. Dab2, Maf, Csf1r) and CD206 (Mrc1) (Fig. 3D). During postnatal heart valve development, mast cells and T-cell populations decrease from P7 to P30, with a switch from naïve towards effector memory T cells (Fig. 3E,F; Fig. S3C, Table S1). While the relative number of dendritic cells remains unchanged, the Mφ1 population is decreased, whereas Mφ2 represent almost 75% of immune cells, in P30 heart valve leaflets (Fig. 3E,F). During development, immune cell function matures toward antigen presentation (Fig. S3D, Table S2), immune and defense response (Fig. S3E,F, Tables S3 and S4). Therefore, scRNA-Seq of developing heart valves indicates that T cells, mast cells, dendritic cells and macrophages account for immune cells found in heart valves.

ICGS analysis was performed on the endothelial cell cluster, exemplified by Cdh5 expression (Fig. 4A). As a result, endothelial cells can be separated into three endothelial subclusters (Fig. S4A) shown in higher magnification of t-SNE plots (Fig. 4B). The lymph-VECs are characterized by a gene signature including Prox1, which is involved in lymph vessel development and is specific to lymphatic VECs (Fig. 4C, Fig. S4A). Cells expressing Prox1, as detected by immunostaining, are located on the fibrosa side, away from blood flow, in both AoV and MV leaflets, with 22.2±13.7% of endothelial cells expressing Prox1 (Fig. 4D,G; Fig. S4B,E). The second cluster was characterized as VECs due to high expression of typical endothelial cells markers, e.g. Emcn, Edn1, Vwf (Fig. 4B,C; Fig. S4A). Endomucin-expressing cells are located close to Prox1-expressing lymph-VECs in the AoV, as well as in the proximal part of the atrialis side of the MV at P30, with 59.8±13.1% of endothelial cells expressing Emcn (Fig. 4E,H, Fig. S4C,F). These VECs also express several growth factors (Bmp4 and Tgfb2), as well as ECM genes (Fig. 4J). The third endothelial cell subcluster displays a gene signature (Fig. 4C) related to cell migration and actin filament polymerization (Fig. S4A), and is characterized by expression of Hapln1 and genes related to Wnt-β catenin signaling, such as Wnt9b and Dkk2 (Fig. 4B). Interestingly, the cells identified by Hapln1 immunostaining are in the area where leaflet coaptation occurs, with 29.1±5.6% of endothelial cells expressing Hapln1 (Fig. 4F,I, Fig. S4D,G), and are therefore named coapt-VECs. Similarly, the three endothelial cell subsets can be found and exhibit similar location in adult murine (Fig. S5A-F), porcine (Fig. S5G-I) and human (Fig. S5J-L) AoVs. Diversification of endothelial cells in distinct regions of the valves occurs by P7, as the percentage of cells in each subset is similar between P7 and P30 (Fig. 4K,L). During development, some differentially expressed genes (DEGs) are observed for each endothelial subcluster, mostly involved in response to calcium and cAMP (Fig. S4H-J). Specifically, P30 lymph-VECs display genes characteristic of increased interaction with immune cells (Fig. S4H, Table S5), P30 VECs have increased levels of genes required for ribogenesis (Fig. S4I, Table S6) and P30 coapt-VECs have increased expression of genes involved in response to extracellular stimulus (Fig. S4J, Table S7). Altogether, three endothelial cell populations that display specific spatial location are identified in P7 and P30 heart valve leaflets.

More than 75% of heart valve cells express Col1a1 (Fig. 2D, Fig. 5A) at P7 (76.7%) and P30 (84%), and thus are classified as interstitial cells. Unlike endothelial and immune cells, interstitial cells from P7 and P30 valves do not cluster together (Fig. 5B), indicative of ongoing interstitial cell differentiation during the first month after birth. P7 VICs display higher synthetic activity, with enriched biological processes involved in collagen fibril organization, peptide cross-linking and bone formation (Fig. 5C, Table S8). Gene expression profiles of interstitial cell populations at P7 were analyzed to identify cell populations responsible for collagen production and ECM remodeling active at this stage. Based on gene expression, P7 VICs are clustered into two subsets (Fig. S6A). The Col-VIC subpopulation highly expresses Col1a1 (Fig. S6A; Fig. 5D,E) and represents about 50% of P7 VICs (Fig. 5F). DEGs between the two clusters (Table S9) reaffirm that Col-VICs not only express fibrillar collagen genes, but also genes involved in fibril organization and ECM maturation, e.g. Prelp, Fmod, Ctgf, Dcn (Fig. 5G). The other P7 interstitial subset is, in contrast to Col-VIC, named GAG-VIC, as illustrated by glycosaminoglycan (GAG)-related Lum, Fbln2 and Vcan expression (Fig. 5D,G). Interestingly, Postn and Sox9 expression, which is associated with the development of collagen-rich and proteoglycan-rich tissues, does not differ between the two clusters (Fig. S6B). In heart valve leaflets, Col-VICs are in the fibrosa, as illustrated by Col1a1 immunoreactivity (55.6% of leaflet area), whereas GAG-VICs expressing versican are mostly found in the tip of leaflet (44.4% of leaflet area) and are absent in the hinge (Fig. 5H-K). A specific elastin-producing VIC cluster was not identified, although expression of Eln (1.49-fold change) is higher in P7 VICs than in P30 VICs (Fig. S6B,C). Still, Eln is expressed in both P7 Col- and GAG-VICs (Fig. S6B). These analyses define two major heterogeneous subpopulations of VICs that are responsible for compartmentalized ECM stratification during active valve remodeling at P7.

The major events of heart valve stratification and morphogenesis are complete at P30, when interstitial cells are subdivided into four subsets with gene expression distinct from P7 VICs (Fig. S6C). Matrifibrocytes, a recently described fibroblast subpopulation present after cardiac injury (Fu et al., 2018), express numerous interstitial cell markers with characteristic higher expression of Comp, Fmod and Chad (Fig. 6A,B; Fig. S6C), which are related to collagen fibril organization, wound healing and endochondral bone formation (Fig. S6D). Nearly 50% of VICs are characterized as fibrosa-VICs, which have typical fibroblast characteristics and express established interstitial cells markers such as Postn and Vim, and additional genes involved in ECM and collagen organization (Fig. S6D). However, these cells have lower levels of osteogenic gene expression than matrifibrocytes. An additional previously unreported cell type (Tcf21-VIC) expresses multiple complement genes involved in defense response, C4b, C4a and Cp, as well as Tcf21 (Fig. 6A,B; Fig. S6C,D). Finally, antigen-presenting VICs express genes involved in antigen processing and presentation via MHCII, such as Cd74 (Fig. 6A,B; Fig. S6D). These two VIC clusters also express Sox9 and Postn, but not Scx (Fig. S6C), which is a marker of early postnatal interstitial cells (Levay et al., 2008) expressed in almost all P7 VICs (Fig. S6B). Together, these analyses demonstrate the presence of at least four subpopulations of VICs at P30.

Immunohistochemistry and cell lineage analysis were used to localize and determine MV versus AoV specificity of VIC subpopulations in P30 valves. Localization of antigen-presenting VICs is confirmed by expression of MHCII (Fig. 6D′,E′) in CD45-negative (Fig. 6D″,E″) cells, although at lower percentages, in the AoV (2.6±1.1%) and MV (1.4±0.5%) leaflets (Fig. 6D,E). Tcf21-VICs, detected using a tamoxifen-inducible Tcf21Cre-dependent GFP reporter (Acharya et al., 2011a), are found mostly in the non-collagenous part of MV (33.1±5.5%) (Fig. 6F), while only few complement-VICs are detected in the AoV (3.7±1.1%) (Fig. 6E). Conversely, matrifibrocytes, as illustrated by co-expression of Fmod and Chad, are found predominantly in the hinge region only in AoV (Fig. 6G,H, white arrowheads) but not in MV (Fig. 6I,J). Finally, fibrosa-VICs, illustrated by Fmod staining, but lacking Chad expression, can be found in the fibrosa in both the AoV and MV (Fig. 6G,I, red arrowheads). Interestingly, Postn and Sox9 are not definitive gene markers for specific interstitial cell subpopulations, but are found in multiple subpopulations (Fig. S6C,D). Thus, multiple interstitial subpopulations are present in P30 valves with distinct functional processes and predominance in the AoV versus MV. Altogether, this study depicts heart valve cell populations found in aortic and mitral valve leaflets, including differentiation of interstitial cell populations during postnatal heart valve development (Fig. 7).

Although heart valve cells are responsible for production and remodeling of ECM during development, as well as for homeostasis during adult life, few studies have focused on the cellular heterogeneity of valve cells and their putative functions. We took advantage of DropSeq, a technology for the determination of transcriptomes of individual valve cells, to investigate the gene expression profiles of distinct cell populations in postnatal remodeling heart valve leaflets. Interestingly, subpopulations of endothelial and immune cells are relatively constant during postnatal development, but interstitial cell subpopulations undergo significant changes in gene expression and cellular functions in primordial and mature valves. Moreover, the diversity in interstitial cell subpopulations demonstrates an unappreciated heterogeneity of cell types in valve ECM production and homeostasis after birth. Novel matrifibrocyte and complement-expressing populations in mature valves are suggestive of specific cell contributions to valve homeostasis and calcific disease mechanisms. Together, this study is the first unbiased transcriptomic analysis of heart valve cells to highlight heart valve cell heterogeneity and transitions in interstitial cells in heart valve remodeling after birth (Fig. 7).

Here, we confirm that CD45-expressing cells of hematopoietic origin are present in developing valves (Hulin et al., 2018) and demonstrate the presence of macrophage subpopulations, dendritic cells, T cells and mast cells in remodeling valves. While infiltrating immune cell lineages have been noted in mature heart valves for many years, their functions in normal valve homeostasis and disease have not yet been determined (Hulin et al., 2018, 2017). The presence of T cells and mast cells in remodeling valves has not been reported previously, and they may have important, yet uncharacterized, roles in valve maturation. Similarly, pigmented cells have been noted previously in developing and diseased heart valves, but their specific origins and functions are not known (Mjaatvedt et al., 2005). Here, we show that these cells express definitive melanocyte markers, including Dct, confirming their cell identity as melanocytes. As these cells are often misidentified as calcification by von Kossa staining, investigators should be aware of the presence of melanocytes in mouse models of valve disease. Together, our detailed single cell analysis demonstrates multiple cell types and lineages present in remodeling valves beyond multiple VEC and VIC subpopulations.

Although previous studies have shown differences in valve endothelial cells on the fibrosa or elastin-rich sides of valve leaflets (Simmons et al., 2005), our study indicates that endothelial heterogeneity extends beyond side specificity. The VEC subset, which is characterized by classic endothelial markers, seems to correspond to endothelial cells described on the aortic side of porcine aortic valves described with higher Vwf and Selp expression, but we show that these cells can be found on both sides of the valve leaflets. Lymph-VECs are exclusively located away from pulsatile blood flow in a similar proportion to that reported in adult porcine aortic valves (Blancas et al., 2016). The presence of VICs with a lymphatic gene expression signature in developing valves has not been previously reported, and the specific functions of these cells in valve remodeling and homeostasis is unknown. The coaptation-VECs, which display enriched gene expression involved in cytoskeletal organization, are present in regions of high mechanical stress, where valve leaflets are exposed to high shear when open and compressive forces when closed (Thubrikar et al., 1980). Likewise, a recent study demonstrated restricted Wnt9b expression in endothelial layer of developing zebrafish heart valves in response to fluid forces (Goddard et al., 2017). These findings highlight potential contributions of biomechanical forces in valve biology, which could be important in establishing and maintaining VEC diversity in different parts of the valves.

By studying developing heart valves at P7 compared with P30, we have identified multiple subpopulations of ECM-expressing VICs, and these populations change from P7 to P30. During development, VICs arise from distinct embryonic origins. The majority of VICs are derived from the endocardial cushions, which arise by endothelial-to-mesenchymal transition at the AV canal and outflow tract (Gomez-Stallons et al., 2016). Additional valve interstitial cells of the semilunar valves are of neural crest origin and contribute to VICs (Gomez-Stallons et al., 2016) and melanocytes. Differences in neural crest versus endocardial cushion-derived VICs have not been detected, and it is not known whether the specific subpopulations of VICs detected in the current report have distinct embryonic origins. However, one population of VICs in the MV expresses Tcf21 and is likely derived from Tcf21-expressing epicardial cells that migrate into the parietal leaflets of the AV valves during fetal development (Lockhart et al., 2014). Interestingly, these cells also express multiple complement and enriched defense biological process genes at P30, with unknown implications for valve homeostasis or disease. Similarly, transcriptome data from global murine valve tissue showed enrichment of the defense response in mature heart valves (Nordquist et al., 2018), and complement factors have been detected in humans with calcific aortic valve disease (CAVD) (Dikhoff et al., 2015; Schlotter et al., 2018). Still, the roles of Tcf21-VICs and valvular defense response are completely unknown.

We further sought to identify cell populations involved in collagen production and stratification of the ECM. Our analysis distinguished a specific VIC subset producing fibrillar collagens, in addition to genes related to collagen fiber assembly necessary to form the fibrosa layer at P7. As valve maturation and stratification proceeds, interstitial cells differentiate, collagen production decreases and VICs are become quiescent. At P30, fibrosa-VICs might differentiate from collagen-VICs present at P7. These cell populations display similar gene expression, although with lower levels and they are found in the fibrosa layer, consistent with Scx cell lineage-tracing data (Levay et al., 2008). These collagen-producing cells are likely to play crucial roles in valve matrix organization and morphogenesis. Future investigations are needed to determine the onset and molecular mechanisms recruiting VIC progenitors to a collagen-producing subpopulation and inducing differentiation into mature fibrosa-VICs. Postn and Sox9 are expressed in multiple VIC subpopulations, consistent with their expression in embryonic valve progenitors (Lincoln et al., 2007; Snider et al., 2008), and thus may be less restricted than VIC subcluster genes induced during valve remodeling after birth. Still, Postn and Sox9 might display different role in VIC subpopulations upon post-transcriptional modifications and binding partners.

A subpopulation of VICs with properties of matrifibrocytes is localized in the hinge regions of AoVs, but is not detected in MVs. Interestingly, the hinge region of mouse AoVs is prone to calcification, as observed in the Klotho mouse model of premature aging (Cheek et al., 2012). The valve matrifibrocytes display a gene signature similar to fibrosa-VICs, as well as chondrocytes and osteoblasts, notably involved in maintaining collagenous environment. Matrifibrocytes represent a specialized subpopulation of cardiac fibroblasts that remain within the mature scar of an injured heart (Fu et al., 2018). Thus, matrifibrocytes in heart valves may have a role in pathologic calcification, supported by increased expression of the matrifibrocyte marker, Chad, in human bicuspid AoVs, which are susceptible to calcification (Padang et al., 2015). Together, our findings shown that VICs exhibit a complexity and heterogeneity that continues to mature during postnatal heart valve remodeling. Our study focusses on healthy developing aortic and mitral leaflets from P7 to P30, indicating that postnatal valve leaflets are not yet quiescent. Interstitial cell subpopulations and specific biological processes, including wounding and defense response, likely contribute to valve postnatal remodeling and homeostasis. Further exploration will define their specific functions in valve development, physiology and potentially pathophysiology, as activation of developmental pathways has been found to be reactivated during heart valve diseases.

Some potential limitations of our single cell transcriptomic analysis exist. One possible limitation is that the cell isolation protocol underrepresents or over-represents specific cell populations. We acknowledge that our approach might miss putative small populations due to developmental stage, low gene expression or low cell number. Known heart valve cell populations were detected in this study, and even small populations of T cells and melanocytes in P30 heart valve leaflets were present in our cell isolations. Moreover, comparative relative percentages of CD31⁺ endothelial and CD45⁺ hematopoietic cells were detected by flow cytometry and confirmed by immunohistochemistry using similar isolation protocols (Hulin et al. 2018). In both studies, ∼80-85% of cells were CD31⁻CD45⁻, representing interstitial cells, with 5-7% endothelial cells and ∼10% immune cells. These relative proportions of cell types were maintained during postnatal heart valve remodeling. Moreover, quantification of immunostaining is consistent with the percentages obtained with our single cell RNA sequencing. It is interesting to note our study did not identify a specific VIC subcluster responsible for producing the Eln-rich layer. Eln expression is higher in P7 VICs, when compared with P30 VICs, but a specific Eln-expressing VIC cluster was not identified. To identify such cluster, additional studies need to be performed in larger animals that have a more predominant elastin layer than murine heart valves (Hinton et al., 2006). In addition, timing of elastin deposition during heart valve development is not clear and might precede P7. Pooling aortic and mitral valves in this study allowed us to obtain sufficient cell numbers but hindered detection of the putative differences between leaflets. Immunostaining supported the hypothesis that VIC subpopulations are found differentially in valve leaflets. Although Tcf21-VICs are located mainly in MV, matrifibrocytes are detected exclusively in aortic valve leaflets. Future studies will decipher whether VIC clusters and their associated enriched biological processes might predispose leaflets to calcification versus myxomatous degeneration.

In summary, this study highlights the diversity of cells present in remodeling heart valves and confirms that postnatal leaflets are still undergoing differentiation supported by interstitial cell plasticity and collagen layer formation after birth. The specific genetic markers and functionally distinct cell subpopulations identified here will be useful tools for heart valve tissue engineering optimization and open novel avenues to investigate molecular and cellular mechanisms of heart valve homeostasis and disease.

All mouse experiments conform to NIH guidelines (Guide for the Care and Use of Laboratory Animals) and were performed with protocols approved by the Institutional Animal Care and Use Committee at Cincinnati Children's Hospital Research Foundation. C57BL/6J mice were used and were obtained from Jackson Laboratory. Tcf21^{MerCreMer(MCM)} mice (Acharya et al., 2011a) were crossbred with a R26^EGFP reporter line to obtain Tcf21^MCM/+;R26^EGFP mice. Perinatal induction of Cre recombinase was achieved in Tcf21^MCM/+;R26^EGFP mice by intragastric injections of 50 µl of tamoxifen (Sigma-Aldrich) in corn oil at 1 mg/ml at P0, P1 and P2. Mice were sacrificed under isoflurane inhalation followed by cervical dislocation at the specific time point for each experiment. Male and female mice were used together for all analyses.

Murine aortic and mitral valve leaflets were isolated and pooled from 13 (P7) and 10 (P30) animals, and collected in PBS on ice. Tissues were then digested using digestion buffer a previously described (Hulin et al., 2018). Leaflets were incubated in digestion buffer at 37°C for 3×20 min with gentle rocking. Cell suspensions were filtered using a 35 µm filter and washed with 0.1% BSA in 1× PBS. The concentration of cell suspension was adjusted to 10 cells/µl using a hemocytometer. Droplet sequencing was performed as previously described (Macosko et al., 2015). Dissection and droplet sequencing for each developmental stage were performed on separate days.

The data analysis was performed similarly to a previous report (Magella et al., 2018). Paired-end sequences were tagged for the unique molecular identifier (UMI) and cell/bead barcode. Reads were aligned to the mm10 mouse genome using Bowtie2 v2.27. Genes that successfully aligned to an exon were tagged with the corresponding gene name. Reads were sorted by the cell barcode and UMIs for each gene were counted. A digital expression matrix of raw reads for each gene was created from the selected cells. A QC step consisted of filtering cells and retaining cells that expressed a minimum of 500 genes. These raw reads were normalized by dividing each gene count by the total number of reads from that cell and multiplying by 10,000, and the matrices from P7 and P30 were merged. The resulting data have been deposited in GEO under accession number GSE117011. To identify and annotate populations of cells from the P30 valve, we performed an unsupervised analysis of the gene-level UMI counts with the algorithm Iterative Clustering and Guide-gene Selection (ICGS) the software AltAnalyze version 2.1.1., which utilizes the pairwise correlation of dynamically expressed genes and iterative clustering with pattern-specific guide genes to delineate coherent gene-expression patterns, was used to predict de novo cell populations, as previously described (Olsson et al., 2016). ICGS was performed using the default options for Drop-Seq data (rho>0.3, restrict to protein-genes, at least four cells with a fourfold difference in expression, exclude cell-cycle effects). Gene expression clusters were generated in AltAnalyze using the hierarchical-ordered partitioning and collapsing hybrid (HOPACH) algorithm. Following the initial identification of nine cell clusters, ICGS was re-run with the same options on each of the major heart valve cell population, endothelial, immune and interstitial cells, to identify additional subpopulations. From all ICGS analysis, cell subpopulation restricted genes were identified using MarkerFinder, which generates 60 specific genes for each cell cluster determined by Pearson correlation coefficient. The top 20 correlated genes determined with MarkerFinder with at least a 0.3 Pearson correlation coefficient were selected as the gene signature for each cell cluster. In order to name our immune cell clusters, immune cell gene signatures were compared with the mouse immune cell database Immunological Genome Project (ImmGen). t-SNE plots were generated through AltAnalyze. Gene Ontology terms associated with differentially expressed genes were evaluated using GO-Elite in AltAnalyze. Similarly, genes that were differentially expressed between developmental stages for each cell clusters were determined and analyzed using AltAnalyze by indicating corresponding development origin for each cell prior to analysis (empirical Bayes moderated t-test P<0.05 and fold>2). Cells were classified as of female or male origin using expression of Xist (female) or five Y chromosome genes, Ddx3y, Eif2s3y, Erdr1, Kdm5d or Ut (male) in the dataset as previously described (Skelly et al., 2018). Female cells express Xist and no Y chromosome genes, male cells express Y chromosome genes, but not Xist, and unspecified cells express neither or both classes of genes. Evaluation of genes showing sexually dimorphic gene expressing was performed using AltAnalyze by indicating corresponding sex origin, and differentially expressed genes were determined using the empirical Bayes moderated t-test. UMAP plots were generated through AltAnalyze. In heat maps, yellow indicates high relative gene expression, and blue or black indicate low or no gene expression, respectively, of the associated genes (rows).

Murine hearts were harvested at indicated time points. Porcine aortic valve leaflets were isolated from pigs aged from 6 months to 3 years. Human aortic valves from donors with no history of cardiovascular medical conditions were obtained postmortem from the National Disease Research Interchange (NDRI). Human studies were approved by the Institutional Review Boards at Cincinnati Children's Hospital Medical Center. Hearts and aortic valve leaflets were processed for paraffin wax embedding and sectioned at 5 or 7 µm, as described previously (Hulin et al., 2018). Tissues were fixed in 4% paraformaldehyde (PFA) for 2 h at room temperature or overnight at 4°C, dehydrated through a graded ethanol series, cleared in xylenes and embedded in paraffin wax. For immunofluorescence studies, slides were incubated at 60°C for 1 h and placed in distilled H₂O. All slides were treated either with a citric acid antigen retrieval reagent (Vector Laboratories) in a pressure cooker for 5 min or with Proteinase K (Thermo Scientific) for 15 min. For Vcan and Col1a1, the slides were treated for 45 min with 0.2 U/ml chondroitinase ABC (Sigma-Aldrich) at 37°C. For Prox1 and Hapln1 staining, slides were pre-treated with 0.025% Triton X-100 in 1× TBS for 10 min. The primary antibodies used were against: Prox1 (Abcam, ab199359; 1:500), Hapln1 (DSHB, 9/30/8-A-4; 1:75), Endomucin (eBioscience, 14-5851-82; 1:250), Vcan (Abcam, ab1033; 1:400), Col1a1 (Abcam, ab34710; 1:25), Dct (Abcam, ab74073; 1:100), pHH3 (Abcam, ab14955; 1:100), Chad (Sigma-Aldrich, HPA018241; 1:100) and Fmod (Kerafast, ENH085; 1:50). DAB staining was performed according to the Ultra-Sensitive mouse or rabbit ABC Peroxidase Staining Kit user guides (32052, 32054, Thermo Fisher Scientific). For fluorescent detection of antibody staining, Alexa Fluor-488-, AlexaFluor-568- and Alexa Fluor-647-conjugated secondary antibodies (Abcam, Life Technologies) were used at 1:200 dilution. Nuclei were counterstained with 4′, 6-diamidino-2-phenylindole, DAPI (Life Technologies, 1:10,000). Movat's pentachrome staining was performed according to the manufacturer's protocols (American MasterTech). Images were captured using a Nikon A1-R confocal system with NIS-Elements D 3.2 software. The percentage of cells positive for the specified staining was determined as the number of positively stained cells/total number of nuclei per valve leaflet section determined by DAPI staining. The percentage of areas stained were calculated as the positively stained area/total area of the valve leaflet section using NIS Elements Basic Research software (Nikon). A minimum of four mice per timepoint were analyzed.

We thank Andrew Kim and Christina Alfieri for technical assistance. Sequencing was performed by the CCHMC DNA Sequencing Core.