Author response: Single-cell expression profiling reveals dynamic flux of cardiac stromal, vascular and immune cells in health and injury (2019)

Article Figures and data Abstract eLife digest Introduction Results Discussion Materials and methods Data availability References Decision letter Author response Article and author information Metrics Abstract Besides cardiomyocytes (CM), the heart contains numerous interstitial cell types which play key roles in heart repair, regeneration and disease, including fibroblast, vascular and immune cells. However, a comprehensive understanding of this interactive cell community is lacking. We performed single-cell RNA-sequencing of the total non-CM fraction and enriched (Pdgfra-GFP+) fibroblast lineage cells from murine hearts at days 3 and 7 post-sham or myocardial infarction (MI) surgery. Clustering of >30,000 single cells identified >30 populations representing nine cell lineages, including a previously undescribed fibroblast lineage trajectory present in both sham and MI hearts leading to a uniquely activated cell state defined in part by a strong anti-WNT transcriptome signature. We also uncovered novel myofibroblast subtypes expressing either pro-fibrotic or anti-fibrotic signatures. Our data highlight non-linear dynamics in myeloid and fibroblast lineages after cardiac injury, and provide an entry point for deeper analysis of cardiac homeostasis, inflammation, fibrosis, repair and regeneration. https://doi.org/10.7554/eLife.43882.001 eLife digest In our bodies, heart attacks lead to cell death and inflammation. This is then followed by a healing phase where the organ repairs itself. There are many types of heart cells, from muscle and pacemaker cells that help to create the beating motion, to so-called fibroblasts that act as a supporting network. Yet, it is still unclear how individual cells participate in the heart's response to injury. All cells possess the same genetic information, but they turn on or off different genes depending on the specific tasks that they need to perform. Spotting which genes are activated in individual cells can therefore provide clues about their exact roles in the body. Until recently, technological limitations meant that this information was difficult to access, because it was only possible to capture the global response of a group of cells in a sample. A new method called single-cell RNA sequencing is now allowing researchers to study the activities of many genes in thousands of individual cells at the same time. Here, Farbehi, Patrick et al. performed single-cell RNA sequencing on over 30,000 individual cells from healthy and injured mouse hearts. Computational approaches were then used to cluster cells into groups according to the activities of their genes. The experiments identified over 30 distinct sub-types of cell, including several that were previously unknown. For example, a group of fibroblasts that express a gene called Wif1 was discovered. Previous genetic studies have shown that Wif1 is essential for the heart's response to injury. Further experiments by Farbehi, Patrick et al. indicated that this new sub-type of cells may control the timing of the different aspects of heart repair after damage. Tens of millions of people around the world suffer from heart attacks and other heart diseases. Knowing how different types of heart cells participate in repair mechanisms may help to find new targets for drugs and other treatments. https://doi.org/10.7554/eLife.43882.002 Introduction Cardiovascular disease including myocardial infarction (MI) remains a leading cause of morbidity and mortality in the Western and developing worlds. After acute MI, millions of cardiomyocytes (CM) are lost by necrosis and apoptosis, and an initially adaptive collagen-rich scar is laid down to preserve chamber geometry and prevent rupture. The mammalian heart is regarded as being poorly regenerative as the long-term sequelae in virtually all etiologies of heart disease involve increased wall stiffness, reduced heart function and progression to heart failure. However, some inbred strains of mice show surprising cardiac reparative abilities (Patterson et al., 2017), and CM renewal and heart regeneration can be stimulated experimentally (D'Uva et al., 2015; Mohamed et al., 2018; Srivastava and DeWitt, 2016; Wang et al., 2018), garnering optimism that heart regeneration can be achieved in humans. Cardiac chamber walls are composed of a complex, interdependent community of interstitial cells, including vascular, fibroblast, immune and neuronal cells, although how they interact in cardiac homeostasis, injury and repair, is relatively unexplored. In regenerative systems, connective tissues play key roles in defining positional information, and organizing tissue architecture and niche environments (Nacu et al., 2013; Chan et al., 2013; Greicius et al., 2018). Cardiac fibroblasts represent ~10% of all cardiac cells (Pinto et al., 2016) and are distributed throughout the cardiac interstitial, perivascular and sub-epicardial spaces, where they are proposed to have sentinel, paracrine, mechanical, extracellular matrix (ECM) and electrical functions (Shinde and Frangogiannis, 2014; Tallquist and Molkentin, 2017). After injury, inflammation is principally executed by poly-functional monocytes (Mo) and macrophages (MΦ), and is necessary to protect against pathogens and autoimmunity, and to coordinate healing. Fibroblasts also participate in inflammation and phagocytosis and are the principal drivers of fibrotic repair (Shinde and Frangogiannis, 2014; Gourdie et al., 2016). In heart repair, timely resolution of inflammation is necessary for limiting fibrosis and enabling tissue replacement, while uncontrolled inflammation leads to increased fibrosis and chamber wall stiffening, poor electro-mechanical coupling, continued loss of CMs and worsening outcomes (Mescher, 2017; Lai et al., 2017; Williams et al., 2018). The general principles of inflammation and fibrosis have been mapped in different organs, and the implementation of specific lineage-tracing tools has provided significant new insights into cardiac leukocyte and fibroblast origins and fate (Tallquist and Molkentin, 2017; Williams et al., 2018; Fu et al., 2018; Ivey et al., 2018; Kanisicak et al., 2016; Moore-Morris et al., 2014; Ensan et al., 2016; Heidt et al., 2014; Molawi et al., 2014; Epelman et al., 2014; Plein et al., 2018). However, controversies persist around nomenclature, defining markers, origins, heterogeneity and plasticity (Tallquist and Molkentin, 2017; Epelman et al., 2015; Swirski and Nahrendorf, 2018). Even the question of whether the transitions from quiescent to activated fibroblast, then to myofibroblast, should be seen as differentiation in the classical sense, degrees of a scalable and reversible continuum governed by the injury environment, or a branched dynamic network, is unresolved (Tallquist and Molkentin, 2017; Ivey and Tallquist, 2016; Travers et al., 2016). One approach to a deeper understanding of cardiac population biology is through single-cell genomics, including single-cell RNA sequencing (scRNA-seq). Single-cell methods have the power to overcome the limitations of bulk cell analyses, where insights into complex cell system dynamics are lost (Tanay and Regev, 2017). The rich data generated by single-cell methods allow new computational frameworks for inferring cell dynamics and causality, unencumbered by strict a priori notions of cell identity, hierarchy, trajectory and markers. Here, we present the first comprehensive analysis of cellular lineage heterogeneity, dynamics and intercellular communication among immune and stromal (non-CM) cells in healthy and injured adult mouse hearts using scRNA-seq. Clustering analysis of >30,000 cells identified over 30 cell populations across the total non-CM fraction and enriched (Pdgfra-GFP+) fibroblast lineage cells. These populations comprised most of the known cell types and their dynamics after injury, as well as novel cell types and their intermediates. We describe a novel population of activated fibroblasts present in both sham and injured hearts expressing a strong anti-Wingless-related integration site (WNT) transcriptome signature, a putative pre-proliferative state, and three novel myofibroblast subtypes expressing pro-fibrotic or anti-fibrotic (including anti-WNT) signatures. We were also able to distinguish the major tissue-resident and infiltrating Mo/MΦ, and numerous minor populations. Overall, our data reveal dynamic, multi-dimensional lineage trajectories in the injured heart. This deep resource will provide novel insights into the inflammatory and fibrotic cascades in the injured mouse heart that may suggest novel molecular or cellular targets for enhancing heart repair and regeneration in man. Results Single-cell RNA-seq of total cardiac interstitial cell population We performed single-cell expression profiling on the total cardiac interstitial cell population (TIP) using the 10x Genomics Chromium platform, from hearts of 8 weeks old male PdgfraGFP/+ mice at days 3 and 7 post-sham or MI surgery. To enrich for cells relevant to cardiac ischemic injury and repair, we isolated TIP cells from dissected ventricles and interventricular septum, excluding cells of the atria, annulus fibrosus and atrioventricular valves (Figure 1—figure supplement 1A). Transcriptional profiles of 13,331 cells were captured after quality control filtering (sham: 5,723; MI-day 3: 3,875; MI-day 7: 3,733). To identify cells with distinct lineage identities and transcriptional states, we performed unbiased clustering on an aggregate of cells using the Seurat R package (Butler et al., 2018), with cell populations visualized in t-SNE dimensionality reduction plots (Materials and methods). For initial analyses, populations expressing markers of multiple lineages (hybrids) were removed; however, select examples are discussed in more detail below. TIP cells were represented by a total of 24 populations and nine distinct cell lineages (Figure 1A–D; Figure 1—figure supplement 1B–F). Major cell types comprised fibroblasts/myofibroblasts (Col1a1+Pdgfra+GFP+), endothelial cells (ECs; Kdr+Pecam1+), mural cells (Cspg4+Pdgfrb+), Mo and MΦ (Cd68+Itgam+), dendritic-like (DC) cells (Cd209a+Itgam+), glial cells (Plp1+Kcna1+), B-cells (Cd79a+Ms4a1+), T-cells (Cd3e+Cd3d+Lef1+) and natural killer cells (NKCs; Klrk1+Ccl5+) (Figure 1A–E; Figure 1—figure supplement 1D; Supplementary file 1). New lineage markers were identified; for example, Vtn, encoding Vitronectin, was specifically expressed in mural cells, whereas Kcna1, encoding the potassium voltage-gated channel subfamily A member 1, was highly specific to glial cells (Figure 1D and Figure 1—figure supplement 1D). Figure 1 with 8 supplements see all Download asset Open asset TIP scRNA-seq. (A) t-SNE plots showing detected lineages and sub-populations in TIP across conditions. (B) t-SNE plot of aggregate TIP cells with identified sub-populations. (C) Dendrogram of sub-populations according to average RNA expression. (D) Expression of select marker genes across TIP cells as visualized on t-SNE plots. (E) Cell population percentages across conditions determined to be significantly modulated according to Differential Proporation Analysis (DPA) (p<0.01). https://doi.org/10.7554/eLife.43882.003 Figure 1—source data 1 Source data for FACS quantifications summarized in Figure 1—figure supplement 6D,E and Figure 1—figure supplement 7B,C. https://doi.org/10.7554/eLife.43882.012 Download elife-43882-fig1-data1-v2.xlsx Within the fibroblast lineage, we identified several sub-populations, each marked by expression of Pdgfra, Pdgfra-GFP, Col1a1 and other canonical fiboblast markers (Tallquist and Molkentin, 2017; Ivey and Tallquist, 2016) (Figure 1A–D; Figure 1—figure supplement 2A). We describe these in more detail after enrichment below. There were three major sub-populations of ECs (EC1, EC2, EC3), comprising vascular and lymphatic lineages (Pecam1+Kdr+Ly6a+) (Figure 1D; Figure 1—figure supplement 2B). The majority EC1 population expressed Ly6a (encoding SCA1) as well as the vascular transcription factor (TF) Sox17, and likely represents microvascular ECs. EC2 expressed canonical arterial endothelial markers such as Bmx, Sema3g and Efnb2, as well as TF genes Sox17 and Hey1 (Figure 1—figure supplement 2B), the latter acting downstream of NOTCH which is required for arterial EC fate. EC3 almost uniquely expressed venous EC marker Nr2f2 (encoding COUPTFII) and Von Willebrand factor (Vwf), and a minority (~3%) expressed Prox1 and Lyve1 (Figure 1—figure supplement 2B), consistent with a lymphatic identity. A small number of ECs were cycling (Cyc; Figure 1A–C; Figure 1—figure supplement 2C). Our EC data are broadly consistent with recently published single-cell data (Zhao et al., 2018). Among lymphocytes, a single B-cell population (BC) expressed Cd79a, Ms4a1 and Ly6d (Figure 1D; Figure 1—figure supplement 2D). T-cell sub-populations TC1-Cd8 (Cd8a+) and TC2-Cd4 (Cd4+Lef1+) likely represent cytotoxic and helper T-cells, respectively. NKCs exclusively expressed Klrk1 and upregulated Lck, Ccl5 and Ctsw (Figure 1—figure supplement 2D). Differential proportion analysis for detecting cell population dynamics We observed major injury-induced cellular responses and flux after MI, including expansion of Mo/MΦ populations at MI-day 3, as well as myofibroblasts and an additional MΦ population at MI-day 7 (Figure 1A; Figure 1—figure supplement 1B,C). To analyze whether changes in the proportion of populations were greater than expected by chance, we developed a novel permutation-based statistical test (differential proportion analysis; DPA) that considered sources of variation which could arise from experimental procedures (such as differing cell numbers and cell-type capture bias) or in silico analysis (cluster assignment accuracy) (Materials and methods; Figure 1—figure supplement 3A–H). DPA identified 12 populations showing significant (p<0.01) flux between conditions (Figure 1E; Supplementary file 2); for example, the fibroblast sub-populations F-SL and F-SH (see below) decreased sharply in proportion at MI-day 3, while M1 and M2 MΦ populations expanded at days 3 and 7 after MI, respectively. Monocyte/macrophage cell identity and dynamics Cardiac tissue-resident MΦ originate from CX3CR1+ progenitors in the yolk sac and Mo from fetal liver and post-natal bone marrow (Ensan et al., 2016), and have roles in immunity, coronary artery and pacemaker development, and heart regeneration (Lavine et al., 2014; Hulsmans et al., 2017). Resident MΦ are long lived and self-renewing (Epelman et al., 2014; Bajpai et al., 2018), although some are supplanted by blood-derived Mo with age or injury (Heidt et al., 2014; Molawi et al., 2014; Dick et al., 2019). MI triggers a biphasic cascade of inflammation and remodeling, with the acute phase involving early influx of neutrophils and CCR2+LY6C2high pro-inflammatory M1 Mo/MΦ, which phagocytose debris and secrete pro-inflammatory factors IL-1β, IL-6 and TNFα to amplify the inflammatory response (Swirski and Nahrendorf, 2018). The repair phase begins around MI-day 3 when non-classical LY6C2-F4/80high M2 MΦ accumulate and secrete anti-inflammatory cytokines such as Il-10 and TGF-β, and stimulate angiogenesis (Epelman et al., 2015; Swirski and Nahrendorf, 2018). In sham hearts, we identified cardiac tissue-resident MΦ with the signature Cx3cr1highAdgre1(F4/80)highH2-Aa(MHC-II)+Itgam(CD11b)lowLy6c2lowCcr2- (MAC-TR; Figure 1A,B; Figure 2A–D) (Ensan et al., 2016; Epelman et al., 2014; Swirski and Nahrendorf, 2018; Lavine et al., 2014; Lavine et al., 2018). Recent work using flow cytometry and scRNA-seq has delineated several subsets of cardiac tissue-resident MΦ, including a pro-regenerative population with the signature TIMD4+LYVE1+MHC-IIlowCCR2-, that self-renew and are not replaced by blood monocytes even after injury (Dick et al., 2019). We could discern this same population at the scRNA-seq level as a subset within MAC-TR, which persisted after injury (Figure 2—figure supplement 1A; Figure 2B). The additional major subset of CCR2- tissue-resident MΦ (Dick et al., 2019) could also be recognised at the scRNA-seq level as the Timd4-Lyve1-H2-Aa(MHC-II)highCcr2- subset of MAC-TR – this population has been shown to have a low monocyte dependence during homeostasis but is almost fully replaced by monocytes after MI (Dick et al., 2019). Figure 2 with 1 supplement see all Download asset Open asset Cardiac Mo/MΦ populations. (A) t-SNE plot showing extracted Mo/MΦ populations. (B) Expression of select immune cell markers as visualized in box-plots and t-SNE plots. Arrows indicate Ccr2low sub-population of M2MΦ. (C) Dot-plot of top upregulated genes for each Mo/MΦ population where color indicates experimental conditions. (D) Heatmap of differentially expressed genes between Mo/MΦ populations with representative significant GO Biological Process terms. https://doi.org/10.7554/eLife.43882.013 Among other minor resident Mo/MΦ populations detected, the most abundant (pale green cells in Figure 1A,B) clustered with the M2 MΦ present at MI-day 7. In fact, all minor Mo/MΦ populations in sham hearts aligned with adult monocyte-derived Mo/MΦ populations which influx after MI (Figure 1A,B), consistent with recent findings (Dick et al., 2019). A prominent B-cell, and minor DC, T- and NK cell populations were also present in sham hearts. These populations may represent a mixture of resident cells and those involved in homeostatic immunosurveillance (Lavine et al., 2018), although we cannot exclude a response to sham operation. At MI-day 3, a major influx population was identified as classical blood-derived M1 Mo, based on the signature Adgre1(F4/80)+Itgam(CD11b)+Fcgr1(CD64)+Ly6c2highCcr2highH2-Aa(MHC-II)low (M1Mo; Figure 2A–C) (Epelman et al., 2015; Swirski and Nahrendorf, 2018). Differentially expressed genes showed Gene Ontology (GO) term over-representation for cell migration, inflammation and T cell activation (Figure 2D). In FACS, Mo are distinguished from MΦ by having lower size and granularity, and lower levels of MΦ markers including Adgre(F4/80), Itgam(CD11b) and H2-Aa(MHC-II) (Bajpai et al., 2018; Hilgendorf et al., 2014). M1Mo identified at MI-day 3 were also low or negative for other MΦ markers Siglec1, Mrc1, Maf, Trem2 and Mertk, the latter involved in phagocytosis (Figure 2—figure supplement 1A) (Bajpai et al., 2018), and C1 complement genes C1qa, b and c (Figure 2C), which are involved (in addition to complement fixation) in recruitment of new inflammatory cells and protection against autoimmunity (Emmens et al., 2017; Thielens et al., 2017). The more abundant population at MI-day 3 was identified as classical Mo-derived M1 MΦ based on the signature Ccr2highAdgre1(F4/80)+Ly6c2+H2-Aa(MHC-II)+ (M1MΦ; Figure 1A,B; Figure 2A–D). This assignment was supported by expression of the additional MΦ markers cited above, including Mertk and C1q, and hierarchical clustering, which showed M1MΦ to be most closely related to M1Mo (Figure 1C), as for human cognates (Bajpai et al., 2018). Differentially expressed genes showed GO term over-representation for leukocyte migration and responses to interleukin-1 (Figure 2D). The most prominent population at MI-day 7 was identified as non-classical M2 MΦ involved in inflammation resolution and repair, with the signature Ccr2highAdgre1(F4/80)+H2-Aa(MHC-II)highLy6c2- (M2MΦ; Figure 1A,B; Figure 2A–D) (Epelman et al., 2015; Swirski and Nahrendorf, 2018). Differentially expressed genes showed GO term over-representation for antigen presentation via MHC class II (Figure 2D). As expected, the non-classical M2 MΦ population increased late during injury repair from <2% of TIP in sham and MI-day 3 hearts, to 16% at MI-day 7. Interestingly, the population dendrogram showed that M2MΦ were most closely related to MAC-TR (Figure 1C), and similarities between resident and subsets of infiltrating Mo/MΦ have been discerned recently using single-cell methods (Dick et al., 2019). Both MAC-TR and M2MΦ expressed Cx3cr1, often used to define tissue-resident MΦ (Figure 2B), and both upregulated pro-regenerative genes Igf1 (Figure 2B) and Pdgfb/c (Figure 2—figure supplement 1A). The majority of M2MΦ were Ccr2high (important for migration); however, a minor sub-populaion was Ccr2low (arrows, Figure 2B) and these expressed the highest levels of Igf1 and lower levels of MHC-II (Figure 2—figure supplement 1B). In this sense they are similar to the CCR2-MHC-IIlow subset of tissue-resident MΦ which appear to be yolk (Dick et al., et al., 2016), and which play a through expression of and in of the fetal coronary vascular et al., 2016). However, whether in the adult heart they represent resident cells or an infiltrating population that has into a MΦ state will lineage with the anti-inflammatory functions of was expressed in only cells in our and may be at the of (Figure 2—figure supplement 1A). Expression of mouse genes encoding cognates of human and previously used to define non-classical and Mo in human blood et al., 2017), not help to the Mo/MΦ new markers recently from analysis et al., 2018). the M2MΦ marker (encoding was more expressed in M2MΦ than in consistent with findings that not M2 cells et al., the M2 MΦ, other myeloid expressed likely the of fibroblasts et al., et al., 2016) analysis to possible changes in major Mo/MΦ populations (Figure 2—figure supplement a continuum of into a trajectory from early infiltrating M1Mo and inflammatory M1MΦ to the late M2MΦ similar to a recent scRNA-seq study (Dick et al., 2019) and consistent with the in which M1 Mo into M2 cells in (Lavine et al., 2018; Hilgendorf et al., 2014). The also the of M2MΦ present at MI-day 7 with tissue-resident MΦ in sham hearts, a in the population dendrogram (Figure see The minor myeloid populations also showed different expression profiles and dynamics (Figure Figure 2A–D). For example, showed of markers including and file with sub-populations expressing markers of neutrophils and (Figure 2—figure supplement 1A). cells showed strong of genes including and (Figure 2C), consistent with GO term analysis responses to and (Figure 2D). These cells appear to arise from MΦ as to monocytes (Dick et al., and likely to the recently inflammatory MΦ that has negative on heart repair after MI through of inflammatory cell and and expression et al., 2017). communication analysis in TIP We communication with expression of and in and (Materials and methods). were from a of human et al., with after to the et al., 2017). on of with than expected by were identified (Figure and MΦ populations M1MΦ and the number of with having the highest ECs by showed the number and of significant fibroblast populations and to exclusively with vascular and and glial cells. In with this analysis of sham and MI-day 3 hearts showed that fibroblasts were observed in or with endothelial cells (Figure supplement Figure 3 with 1 supplement see all Download asset Open asset (A) of significant Arrows and color indicates and indicates the of between populations. (B) of total total across populations. (C) across and for top in (D) plot showing from the cells. to to to and represent populations. involving were by including 1, and that can with expressed in numerous populations (Figure The minority glial cell population expressed canonical neuronal markers, including and (Figure 1—figure supplement and likely cardiac essential for cardiac regeneration in et al., cells also to with the three EC populations and mural cells (Figure consistent with the of in the cardiac et al., In of we detected highly specific to glial cells in of other TIP including and genes and (Figure involved in both and et al., 2017). these suggest the and of between cardiac cell types in homeostasis and injury. populations in TIP We detected minor populations expressing markers of lineages (Figure 1—figure supplement cells may or in that are to the conditions of are also known to a significant number of et al., the of cells ECs are highly and has been to fibroblasts after cardiac injury et al., 2019). cardiac fibroblasts have been observed to to ECs et al., that this has been et al., 2017). The population markers of fibroblasts and and with other fibroblast populations (Figure 1—figure supplement To this we isolated interstitial cells from dissected ventricles of PdgfraGFP/+ mice after sham or MI and whether we could cells using flow cytometry and a that cell (Materials and methods; Figure 1—figure supplement We detected of cells in sham hearts, in MI-day 7 hearts, and in (Figure 1—figure supplement while cells were they not appear to injury. of Mo/MΦ to into cells in different has been in and in and has et al., although a natural plasticity in Mo an endothelial cell fate in not have strong and 2014). The population markers of ECs and M2MΦ and with M2MΦ (Figure 1—figure supplement single cells in 7 hearts, to in MI-day 7 hearts, an in injury (Figure 1—figure supplement Among were signature of M1 and M2 in sham hearts, to in MI-day 7 hearts. is well known that the expression of EC markers on the of cells is to define as although they can be angiogenesis cells and 2014). these data not exclude the of in our scRNA-seq they the of distinct and populations with and different responses to injury. These Single-cell RNA-seq of the cardiac fibroblast lineage A major subset of fibroblasts in the adult murine heart express the cell