Smad3 regulates smooth muscle cell fate and mediates adverse remodeling and calcification of the atherosclerotic plaque

Atherosclerotic plaques consist mostly of smooth muscle cells (SMC), and genes that influence SMC phenotype can modulate coronary artery disease (CAD) risk. Allelic variation at 15q22.33 has been identified by genome-wide association studies to modify the risk of CAD and is associated with the expression of SMAD3 in SMC. However, the mechanism by which this gene modifies CAD risk remains poorly understood. Here we show that SMC-specific deletion of Smad3 in a murine atherosclerosis model resulted in greater plaque burden, more outward remodelling and increased vascular calcification. Single-cell transcriptomic analyses revealed that loss of Smad3 altered SMC transition cell state toward two fates: a SMC phenotype that governs both vascular remodelling and recruitment of inflammatory cells, as well as a chondromyocyte fate. Together, the findings reveal that Smad3 expression in SMC inhibits the emergence of specific SMC phenotypic transition cells that mediate adverse plaque features, including outward remodelling, monocyte recruitment, and vascular calcification.

D ecades of research and drug development have led to therapies and interventions that have markedly diminished morbidity and mortality from cardiovascular disease 1,2 . However, CAD remains the leading cause of death in the United States and worldwide, with a sharp decline in the rate of improvement in mortality observed over the past decade 3 . Recent clinical studies targeting well-characterized risk factors such as lipids 4,5 and molecular targets related to inflammation 6,7 have had only modest success 8 , suggesting a continued need to identify additional disease modifiers. Over the past decade, genome-wide association studies (GWASs) have identified over 160 loci that contribute to CAD risk 9,10 . Causal variation identified in these loci point primarily to genes and pathways predicted to function in the blood vessel wall to regulate disease risk [11][12][13] .
Specific features of atherosclerotic plaque have been increasingly recognized to offer considerable prognostic value. For instance, it has been noted that cellular composition, such as SMC contribution to the fibrous cap, influences risk of plaque rupture [14][15][16] . With advances in diagnostic imaging, disease features such as outward remodeling and microcalcification have been found to be highly predictive for myocardial infarctions 17,18 . Although counterintuitive, studies using both intravascular ultrasound and longitudinal studies using computed tomography coronary angiography 17,18 have demonstrated that sites of outward remodeling are much more likely to be the culprit site for plaque rupture and myocardial infarction than sites with more luminal narrowing. However, little is known regarding the cellular and molecular mechanisms by which these high-risk features are controlled.
Several post-genomic GWAS follow-up studies in this and other laboratories investigating the genetic disease-related mechanisms of CAD have focused on SMCs and the relationship of SMC state changes to disease risk [19][20][21][22][23] . These studies have indicated that a substantial portion of the CAD-attributable risk is determined by this cell type 13 . Consistent with this notion, recent lineage tracing studies have demonstrated that most cells inside atherosclerotic plaque, including those expressing some inflammatory markers, are oligo-clonal de-differentiated smooth muscle derivatives [24][25][26][27] . Genes that alter SMC behavior are known to influence the composition of atherosclerotic plaque 19,22,23,26,28,29 . CAD-associated genes TCF21 and AHR have been linked to these processes, and various genomic data suggest that additional CAD genes might also regulate SMC phenotype 21,28,[30][31][32][33] . Although these SMC progeny have been previously lumped together as phenotypically modulated SMCs, advances in single-cell RNA profiling have demonstrated the presence of subsets of these cells with distinct transcriptomes and cell fates. For example, medial SMCs have been shown to give rise to fibroblast-like cells termed fibromyocytes (FMCs) 23 as well as cells similar to endochondral bone-forming cells, chondromyocytes (CMCs) 19,31 , among other bioinformatically defined populations 19,34 . However, the functional significance and relative location of these transcriptionally distinct populations remain to be elucidated.
The TGFβ signaling pathway is central to smooth muscle biology during development and disease and is an important modifier of atherosclerosis 35,36 . Canonical TGFβ signaling is thought to be mediated through Smad family proteins, particularly nuclear signaling factors Smad2 and Smad3. Although themselves poor binders to DNA 37 , through their interaction with other transcription factors, the Smad factors are central to key transcriptional programs that regulate cell fate in development and disease 38,39 . Multiple GWASs have identified rs17293632 (ref. 9 ), a single-nucleotide variant that lies within a functional smooth muscle enhancer that regulates SMAD3 expression 20 , as an important modifier of risk for myocardial infarction. However, how smooth muscle SMAD3 expression influences risk of myocardial infarction is unclear. In vitro, SMAD3 appears to modify SMC differentiation and proliferation through its interaction with other transcription factors critical to SMC biology and risk of CAD 30 . However, the exact effects of SMAD3 expression level on SMC plaque biology remain unknown.
Here we demonstrate that SMC-specific deletion of Smad3 influences the fate of de-differentiated SMCs in atherosclerotic plaques in vivo, promoting both a pro-remodeling SMC transition phenotype that expresses remodeling genes such as Mmp3 and inflammatory chemokines such as Cxcl12 as well as an expansion of the SMC-derived CMC population. These cellular changes are associated with increased outward remodeling and plaque calcification that appear to normally be inhibited by Smad3 in conjunction with transcriptional effects of Hox and Sox factors.

Results
Vascular lesion characterization in Smad3 ΔSMC mice. To understand how smooth muscle expression of Smad3 influences atherosclerotic lesions in vivo, we generated a murine model of atherosclerosis with established smooth-muscle-specific Cre (Myh11-Cre) crossed with a conditional knockout allele of Smad3 (refs. [40][41][42], with concurrent lineage tracing provided by conditional tandem dimer tomato (ROSA Tdt ) expression on the ApoE null background (Smad3 ΔSMC ) (Fig. 1a). To limit confounding created by the critical role of Smad3 during development, the Smad3 gene was deleted via a tamoxifen-inducible Cre only after mice had reached maturity (8 weeks old), immediately before initiation of a Western high-fat diet (HFD; Fig. 1a). The Smad3 conditional knockout mice grew to maturity with no substantial change in weight or mortality compared to control (Extended Data Fig. 1a,b). In these studies, evaluation of Smad3 function in atherosclerotic vascular disease was evaluated primarily in the aortic root (Extended Data Fig. 1c), where single-cell RNA sequencing (scRNA-seq) data revealed highly efficient Cre-mediated deletion of the floxed exons encoding the DNA-binding domain of Smad3 (Extended Data Fig. 1d), with corresponding loss of nuclear phospo-Smad3 protein in knockout SMCs verified by lesion immunohistochemistry (Extended Data Fig. 1e,f). Investigation of the phenotype of atherosclerotic disease lesions was conducted according to published recommendations on the design, execution and reporting of animal atherosclerosis studies 43 .
Examination of atherosclerotic lesions in the aortic root after 16 weeks of HFD demonstrated that SMCs from Smad3 ΔSMC mice were able to migrate into the lesion, expand and contribute to the formation of atherosclerotic plaque and the fibrous cap (Fig. 1b). Quantification of atherosclerotic lesions revealed a significant increase in plaque volume in Smad3 ΔSMC compared to control animals (Fig. 1c,d). To further characterize the anatomy of diseased vessels, specifically regarding outward remodeling versus luminal narrowing, we evaluated the area encapsulated by the diseased vessel as well as lumen area. The lumen area in these sections showed no significant change (Fig. 1e), but the area circumscribed by the external elastic lamina was significantly increased (Fig. 1f). These findings are consistent with expansion of atherosclerotic plaque volume in conjunction with outward remodeling. To determine the cellular anatomy associated with the increased plaque size, we quantified the area occupied by fluorescent tdTomato SMC lineage traced cells (Fig. 1g,h) as well as CD68-stained monocytes and macrophages in the lesions (Fig. 1i,j). This analysis revealed a statistically significant increase in area for both SMCderived cells as well as cells of the monocyte/macrophage lineage. Given that the lineage tracing Myh11-Cre transgene is active only in cells emanating from mature SMCs, these findings suggest that the increased lesion growth has both a cell autonomous as well as a non-autonomous cellular component, with the latter reflecting an SMC-mediated effect on monocyte/macrophage lesion recruitment.
Contents of the plaques were further analyzed via Oil Red O, trichrome and Ter119 staining. There was an increase in the absolute area stained with Oil Red O in the Smad3 ΔSMC compared to control animals (Extended Data Fig. 2a,b), but there was not a difference noted when the lipid-stained area was corrected for total plaque area (Extended Data Fig. 2c), suggesting that the plaques were larger but did not contain proportionally more lipid. There was also an increase in absolute acellular area in Smad3 ΔSMC compared to control animals as determined with trichrome staining, and this difference was sustained when corrected for plaque size (Extended Data Fig. 2d,e,f). This increase in acellular area at the base of the plaque was consistent with the higher-risk features seen in outward remodeling.
Vascular lesion scRNA-seq in Smad3 ΔSMC mice. Given the critical role that Smad3 plays in cell fate decisions during development, we hypothesized that alteration in disease-associated SMC phenotype transitions might account for the observed cellular lesion characteristics as well as recruitment of CD68 + cells and outward remodeling. Thus, to better understand how loss of Smad3 expression produced phenotypic changes in lesion-SMC-derived cells, and their interactions with other lesion cell types, we performed single-cell RNA expression profiling of atherosclerotic lesions from Smad3 ΔSMC and control animals. The atherosclerotic tissue was harvested and processed for single-cell encapsulation, RNA capture, reverse transcription and amplification with the 10x Genomics Chromium V3 platform, and cDNA libraries were sequenced as previously described 23,28 . Subsequently, the single-cell expression data were visualized using uniform manifold approximation and projection (UMAP) to create a two-dimensional (2D) projection representing the organization of individual cells to each other in clusters and of the relationship of the organized clusters to each other. After filtering and normalization, a total of 26,219 cells from control mice and 35,518 cells from Smad3 ΔSMC mice were included in the analysis obtained from three and four independent captures of pooled tissue from two mice in each individual capture. Feature plots were employed to visualize the SMC lineage traced cells (Fig. 2a). SMCtraced cells were identified by expression of the tdTomato gene, and the component of this cluster representing mature medial SMCs was identified by expression of SMC lineage markers Myh11 and Cnn1. As we have shown previously, a significant portion of the SMC lineage traced cells did not express these mature SMC markers during disease and, thus, represented SMCs that had undergone phenotypic transition 22,23,28 . We were, thus, able to determine whether there were alterations in the number of mature differentiated versus transition SMCs in Smad3 ΔSMC mice. We investigated whether there was an alteration in the number of differentiated SMCs among the lineage traced (tdTomato + ) cells, defining 'differentiated' cells using classical Myh11 or Cnn1 expression, along with unbiased clustering (SMC versus transition SMC) (Fig. 2a,b). Cells were clustered based on the scRNA-seq data using the Louvain algorithm to a resolution where the SMC-derived cells are split into two distinct clusters, as we and other have characterized 23 . Regardless of the definition of 'mature SMC' , there was a significant decrease in the proportion of mature SMCs in the Smad3 ΔSMC compared to control mice (Fig. 2a,b).
To determine an optimum biologically relevant clustering resolution, we implemented an empirical approach, aimed at determining the minimum clustering resolution required to separate known biologically distinct populations of endothelial cells (ECs) (Fig. 2c,d)that is, blood vessel ECs (VE-cadherin expressing, Prox1 negative) versus lymphatic ECs (Prox1 expressing) 44,45 . The identities of the cellular subpopulations that were created with this approach were then determined based on specifically expressed canonical 'marker genes' in each cellular cluster (Fig. 2c). Previous work in SMC lineage tracing in atherosclerosis has identified three distinct clusters of SMC-derived cells, including mature SMCs, FMCs and procalcific CMCs 19,28,34,46 . Using the cluster resolution identified above, similar subsets of lineage traced SMCs were identified in these data (Fig. 2d), and, in addition, two distinct additional clusters of Myh11-Cre lineage labeled populations were also identified. One cluster was composed of pericytes, and an additional small population of uncharacterized disease-associated SMC-derived cells clustered by itself as a unique transcriptomic phenotype, which we named 'remodeling-SMC' (R-SMC), to be studied in more detail below.
To validate this clustering resolution of the SMC transition cells identified, we conducted a sensitivity analysis by performing clustering of the scRNA-seq data at different resolutions. At the lowest resolution, lesion transition SMCs segregated from mature SMCs, and, with increasing resolution, the transition cells separated into the FMC and CMC disease-related phenotypes. Further increase in resolution of clustering witnessed the appearance of R-SMC and pericyte clusters. An additional increase in resolution resulted in partitioning of the mature SMCs into two separate clusters that we could not easily spatially discern as distinct populations in the lesion, suggesting possible over-clustering. The robustness of this clustering was further confirmed by the ability to produce similar clusters using different clustering algorithms (Louvain and SLM (smart local moving)), with both all lesion cells and SMC-derived cells alone. Taken together, these various approaches demonstrate the appropriateness of the SMC transition clusters that are expected to reflect relevant biological phenotypes.  To determine if loss of Smad3 expression alters the cell fate transitions of disease-associated SMC lineage cells, we ascertained the fraction of de-differentiated cells that contribute to clusters for each of the two known and the R-SMC transition phenotypes. Smad3 ΔSMC transition SMCs showed an increase in relative proportion of CMCs along with the newly defined population at the expense of the intermediate FMC population (Fig. 2e,f). The increase in fraction of CMCs among disease-associated Smad3 ΔSMC cells corresponded with a significant increase in lesion area expressing CMC marker Col2a1 (Fig. 2g,h), suggesting an absolute increase in the CMC population. Functionally, this increase in Col2a1-expressing cells correlated with an increase in lesion calcified area as assayed by von Kossa staining of atherosclerotic lesions (Fig. 2i,j). These cells appeared to result from an increase in proliferation of de-differentiated  SMCs, as evidenced by a higher level and increased proportion of cells expressing proliferation markers such as Mki67, Ccnd1, Ccnb1 and Myc, with limited alteration in cell cycle inhibitor genes Cdkn2a and Cdkn2b (Fig. 2k). To further assess this possibility, we generated a 'chondrocyte proliferation score' that represents a scaled average expression of all genes associated with Gene Ontology (GO) category 'promote chondrocyte proliferation' with each individual cell's transcriptome 47 . Consistent with higher expression of specific cell cycle regulators, Smad3 ΔSMC transition SMCs had a significantly higher chondrocyte proliferation score than control cells (Fig. 2l).
Applying the analysis using a 'mesenchymal proliferation score' based on GO categories of 'promote mesenchymal proliferation' produced equally significant results (Extended Data Fig. 3a), further suggesting that the increase in number of tdTomato + cells in the lesions likely reflects proliferation of SMC derivative cells. The increase in cell number did not result from alterations in apoptosis, because there was no difference in TUNEL staining of control versus Smad3 ΔSMC sections (Extended Data Fig. 3b).

Disease SMCs promote remodeling and leukocyte recruitment.
In addition to increased CMCs in the Smad3 ΔSMC animals, there was also an increase in the proportion of transition SMCs that contributed to the newly identified cell cluster (Figs. 2d and 3a), which we refer to as R-SMCs due to their high expression of genes involved in extracellular matrix (ECM) remodeling. These cells were identified at the junction of CMC, fibroblast and FMC clusters. The vast majority of cells within this cluster were lineage traced at a proportion similar to that seen with CMCs and FMCs (Fig. 3b), confirming that they were SMC derived. In vascular lesions, they constituted 6% of SMC transition cells in Smad3 ΔSMC mice, whereas they represented only 2% of this SMC population in control animals. To investigate the ontogenic relationship of this group of cells to SMCs, FMCs and CMCs, we ported the scRNA-seq data to Slingshot 48 , a lineage inference tool designed to map trajectories involving multiple branching lineages. Applying this algorithm suggested that SMCs give rise to FMCs, which, in turn, serve as a source for both CMCs and R-SMCs (Fig. 3c). Although the R-SMC gene expression pattern is most similar to that identified in CMCs, as shown by their juxtaposition in UMAP space and by hierarchical clustering, they were easily distinguished from CMCs at the transcriptional level (Fig. 3d). Interestingly, the R-SMC population was marked by a particularly high expression of matrix metalloproteinase-3 (Mmp3) (Fig.  3e,f), an enzyme required for outward remodeling of atherosclerotic vessels 49 . Mmp3 was, in fact, among the most upregulated genes in this population compared to other SMC transition groups. Given our observation of increased outward remodeling identified in Smad3 ΔSMC mice, this finding was further investigated. Singlecell transcriptomic data suggested a higher number of Mmp3expressing cells and an overall higher level of Mmp3 expression in Smad3 ΔSMC than control cells ( Supplementary Fig. 4a,c). In addition to increased mRNA levels, an in vitro functional study using a florescent-based Mmp3 activity assay demonstrated that homogenized aortic tissue from Smad3 ΔSMC mice had more Mmp3 activity than control tissue (Extended Data Fig. 4b).
In situ hybridization with RNAscope showed Mmp3-expressing R-SMCs to be a distinct population from CMCs, which were marked by Col2a1 expression (Fig. 3g), indicating that these R-SMC transition cells are juxtaposed but not overlapping in the lesion plaque and further supporting their distinct phenotype. Consistent with the established role for Mmp3 in outward remodeling, its expression was most prominent at the base of the atherosclerotic lesion in cells juxtaposed to the elastic lamina. Strikingly, these Mmp3-expressing cells were commonly associated with areas of disrupted elastic lamina, consistent with their possible invasion through this structure ( Fig. 3g and Extended Data Fig. 4d). Notably, the same rare population of MMP3-expressing cells were also found in human coronary arteries, where they were also associated with regions of disrupted elastic lamina (Fig. 3h). There were more cells expressing Mmp3 in Smad3 ΔSMC compared to controls as per the scRNA-seq data and more prominent staining in Smad3 ΔSMC versus control cells ( Fig. 3f and Extended Data Fig. 4a,c). Mmp3 expression was altered only in the SMC-derived cells, and not significantly different in the non-SMC-derived population, suggesting that altered Mmp3 levels by Smad3 happens in a cell-autonomous manner in SMC lineage transition cells (Extended Data Fig. 4a). By performing TGF-β stimulation of human coronary artery smooth muscle cells (HCASMCs), we were able to show that MMP3 is suppressed by TGF-β (Fig. 3i). SMAD3 knockdown did not significantly increase the expression of MMP3 at baseline but did negate the TGF-β-dependent suppression of expression, suggesting that SMAD3 is required for the TGF-βdependent regulation of MMP3 in atherosclerotic lesions.
To better understand the vascular function of this cell population, we performed pathway analyses of all significantly differentially expressed genes in the R-SMC as compared to other de-differentiated SMC phenotypes (Supplementary Table 1). PANTHER analysis revealed that the top biological processes enriched in this list of genes are related to remodeling of ECM, followed by regulation of chemotaxis and inflammation (Fig. 3j). Analysis of genes contributing to the enrichment of these pathways reveal that this population of cells express higher levels of multiple chemoattractants whose receptors are primarily restricted to macrophages (Fig. 3k). For example, this population is the main SMC source of Cxcl12 in the lesion, with the Cxcr4 receptor for this ligand being expressed exclusively by lesion macrophages. In addition, R-SMCs differentially expressed cytokine genes, including Saa3 and Ccl2, and all have well-established roles in monocyte recruitment, with their respective receptor genes Tlr1 and Ccr2 expressed mostly by inflammatory cells in the lesion (Fig. 3k). There is high overlap between the Mmp3-high cells and cytokine-expressing cells. In fact, if a marker-based definition of a pro-remodeling and recruitment population was defined as cells with elevation of Mmp3 and Cxcl12, one can recapitulate the unbiased clustering results (Extended Data Fig. 4e). Consistent with the hypothesis that SMAD3 suppresses the SMC chemotaxis signal to inflammatory cells, in vitro transwell migration studies demonstrated augmentation of HCASMC ability to recruit human monocytes with SMAD3 knockdown (Extended Data Fig. 4f). Taken together, these findings suggest that, beyond its role in outward remodeling, the R-SMC population of cells may orchestrate monocyte recruitment and underlie the increase in CD68 + cells in lesions.

Smad3 binds other transcription factors to regulate disease.
Given that Smad3 is likely affecting all SMC lineage phenotypes beyond just the CMC and R-SMC populations described above, we characterized the alteration in gene expression of all de-differentiated SMC-derived cells in the atherosclerotic lesions. To rigorously identify differentially regulated genes, we employed a Wilcoxon rank-sum-based analysis comparing Smad3 ΔSMC versus control mouse data for FMC plus CMC and R-SMC cellular clusters, and 83 genes were identified (Supplementary Table 2). Pathway analysis performed using DAVID demonstrated significant enrichment in TGFβ signaling, which is consistent with Smad3's role in the TGFβ pathway (Fig. 4a). In addition, biological processes broadly related to ECM remodeling, SMC differentiation and chemotaxis were also highly enriched. Using GREAT 50 , we compared the list of 83 identified Smad3 ΔSMC marker genes to the top 1,000 genes expressed by de-differentiated SMCs in diseased murine vessels. The analysis revealed enrichment for genes linked to human atherosclerosis and arterial dilatation (Fig. 4b). In addition to Mmp3, and consistent with the observed phenotype of increased outward remodeling in Smad3 ΔSMC mice, Lox and Mfap5 were also among the top significantly downregulated genes in de-differentiated Smad3 ΔSMC SMC   Table 2). Pathogenic loss of function in both of these genes has been linked to human vascular syndromes, including aortic aneurysms [51][52][53] . To investigate whether these genes are directly regulated by TGFβ, we stimulated HCASMCs with TGFβ in the presence and absence of SMAD3. A subset of genes, such as LOX, appeared to be directly regulated by TGFβ in HCASMCs, and this effect was abrogated by SMAD3 silencing (Fig. 4c), similarly to the regulation of MMP3. However, a subset of genes, such as MFAP5, did not appear to be TGFβ responsive but were sensitive to loss of SMAD3 (Fig. 4c). These findings suggested involvement of additional TGFβ-independent co-regulatory factors.
To better understand the context in which Smad3 regulates this complex transcriptional program and to identify possible coregulatory factors, we analyzed the 5′ regulatory elements of the 83 Smad3 differentially regulated genes to look for enriched transcription factor motifs. Motif analysis conducted with HOMER revealed Regulation of transforming growth factor beta receptor signaling pathway Regulation of cellular response to transforming growth factor beta stimulus   3). g, proximity ligation assay demonstrates close proximity of HOXB2-SMAD3 and SOX9-SMAD3 in the nucleus of HCASMCs. h, relative luciferase activity of MFAP5-Luc in HEK-293 cells transfected with HOXB2 and SOX9 expression vectors in the absence (left) or presence (right) of SMAD3 sirNA (SMAD3-KD) or control sirNA (control), two-tailed t-test. All error bars represent 95% confidence interval of mean. HEK-293 authenticity was validated upon receipt from the ATCC and assessed periodically during the course of this study by cell morphometry and pCr evaluation of lineage markers. Abd, abdominal; NS, not significant; TF, transcription factor. enrichment for Hox/homeobox motifs as well as Sox9/10-related motifs (Fig. 4d). Hox and Sox genes were not dysregulated in the Smad3 ΔSMC SMCs (Extended Data Fig. 5i), suggesting that these factors are likely directly interacting with Smad3 to regulate a joint transcriptional program. Consistent with this hypothesis, Sox9 is known to be a critical transcription factor regulating calcification and has been shown in chondrosarcoma cells to selectively interact with Smad3, but not with Smad2, to modulate an endochondral ossification program 54 . Because Hoxb2 and Sox9 were found to be expressed in transition SMCs (Extended Data Fig. 5g,h), we further investigated the possible interaction of these factors with SMAD3. Knocking down SOX9 and HOXB2, the most highly expressed HOX gene in human coronary SMCs, recapitulated the changes in expression of key vascular remodeling genes MMP3 and MFAP5 observed in vivo in mouse (Fig. 4e). Thus, we hypothesize that HOX factors, whose motifs are also enriched, can directly interact with SMAD3. To test this hypothesis, we performed nuclear co-immunoprecipitation experiments to test their interaction. In cells over-expressing His-HOXB2 or Flag-SOX9, SMAD3 protein was co-immunoprecipitated with anti-His or anti-Flag antibody, respectively ( Fig. 4f and Extended Data Fig. 6). Furthermore, proximity ligation assays demonstrated that endogenous SMAD3 is localized in proximity (<40 nm) to HOXB2 and SOX9 to suggest their involvement in multi-protein complexes in the nucleus of HCASMCs (Fig. 4g). Taken together, these data suggest that HOX family proteins such as HOXB2, as well as SOX family member SOX9, physically interact in the nucleus to regulate transcription of targeted genes.
To test the functional significance and epistatic relationship of these findings, we cloned an evolutionarily conserved regulatory region near the 5′ end of the human MFAP5 gene, which contains conserved putative HOX, SOX and SMAD binding elements, into a luciferase vector and tested its ability to respond to SOX9 and HOXB2 binding. SOX9 and HOXB2 efficiently activated this enhancer element but not the control luciferase vector, further establishing a role in this transcriptional program ( Fig. 4h and Extended Data Fig. 5j). Knocking down SMAD3 diminished the SOX9-dependent and HOXB2-dependent activation of the luciferase construct, suggesting that interaction with SMAD3 is required for regulation of MFAP5 (Fig. 4h).

Discussion
High-resolution scRNA-seq data have allowed us to describe a small but distinct subset of SMC transition cells that expresses proremodeling enzymes. Mmp3, a metalloprotease required for outward remodeling, is expressed most prominently in the R-SMCs that are located in the basal plaque where they appear to migrate into the media toward the adventitia, potentially related to the previously described SMC lineage traced population in the adventitia 55 . Strikingly, in human coronary arteries, this MMP3-expressing proremodeling population resides in the same regions of the artery with observed disruption of the nearby elastic lamina, consistent with a role promoting outward remodeling in human coronary disease. These findings suggest that expansion of this population may be a contributing factor to outward remodeling and other adverse features. The importance of this MMP3-expressing SMC transition population in modifying plaque rupture risk may explain the seemingly paradoxical observation that a single-nucleotide polymorphism (SNP) associated with lower expression of MMP3 is associated with greater luminal coronary artery stenosis on cardiac catheterization, but the higher-expressing variant is associated with more myocardial infarctions 56 . These genetic observations provide further evidence that modulation of R-SMCs alters plaque features and CAD risk in human. Although Mmp3 was previously thought to be expressed by macrophages in the lesions 56 , we found no evidence that Mmp3 was expressed by cells in the macrophage cluster. Previous work has shown that IL-1 drives outward remodeling of plaques, and this function was completely reversed by deletion of Mmp3 (ref. 49 ). This suggests that SMCs are the primary mechanism for the high-risk plaque feature of outward remodeling, as seen in the Smad3 ΔSMC mice, and may account, in part, for the beneficial effect of IL-1 blockade on the risk of plaque rupture in humans 6 . Mmp3 also appears differentially regulated among restricted populations of SMC progeny in carotid artery plaque 19 , suggesting that R-SMCs exist in other atherosclerotic beds.
Interestingly, this population of R-SMCs also express several chemokines, including Cxcl12, Ccl2 and Saa3, whose main receptors are expressed solely on the monocyte/macrophage lineage. This suggests that the R-SMC population likely plays a role in regulating the inflammatory response to the lesion, contributing to the increase in observed monocyte/macrophage population detected in the lesions. The combination of remodeling and inflammatory cell recruitment, both factors that determine plaque stability, highlights the critical role that this specific subpopulation of cells may play in modulating human disease risk. This contributes to existing literature 57 suggesting that plaque SMCs play a central role in regulating inflammatory cell recruitment and retention in atherosclerotic plaque and identifies a specific subpopulation of transition SMCs critical for high-risk plaque features.
There was also an increase in CMCs in Smad3 ΔSMC mice, suggesting that Smad3 actively inhibits differentiation to this phenotype or inhibits their proliferation. The transcriptomic, topological and lineage inference data presented here indicate that they are a distinct population from Mmp3-expressing R-SMCs. The CMCs exhibit a chondrogenic transcriptomic program 28 , expressing Col2a1, Acan and Sox9, with similarities to chondrogenic progenitors in endochondral bone formation and repair. Smad3 has been shown to regulate Sox factor transcriptional activity in a TGFβ-independent manner through physical interactions 54,58 . The observed expansion of the CMCs also provides an interesting parallel to established findings of accelerated bone and wound healing in Smad3 knockout mice 59,60 . The concomitant increase in Col2a1-expressing cells in the plaque and increased vascular calcification suggests that this cell type is at least partially responsible for coronary calcifications seen in human CAD. It remains to be determined whether increased calcification is harmful or protective in terms of plaque rupture risk, because conflicting observational data exist in humans. Although increased coronary calcification is correlated with increased risk of myocardial infarction 61 , local calcification appears to be protective against plaque rupture [62][63][64] , and interventions that lower risk of plaque rupture increase calcification 65,66 . Given the multiple populations of SMC-derived transition cells observed in our studies, their relative ratios could possibly determine the quality of calcification as well, which is also considered to confer differential risk of plaque stability 67,68 .
Identification of the R-SMC transition phenotype adds to the complexity of SMC cell state changes that are associated with vascular disease development. We have previously identified fibroblastlike transition phenotype cells that we have termed FMCs 23 , as well as the noted CMCs. Through our studies of the mechanisms by which CAD GWAS genes modulate disease risk, we have shown that transition of SMCs to FMCs is regulated by Tcf21 (ref. 23 ) and Zeb2 (ref. 29 ), and the transition to CMC is regulated by Ahr 28 . We have now also implicated Smad3 in this transition. Pseudotimeassisted trajectory analyses have suggested that FMC may represent an intermediate transition phenotype that gives rise to the terminal CMC phenotype, and this possibility has been validated by other groups 19 . Trajectory analyses presented here suggest that R-SMCs may represent a distinct terminal phenotype. Published studies by other groups investigating the cell state changes by SMCs in the disease setting in murine models have also identified FMC and CMC phenotype cells 19,34,46 . At least one publication has identified fibromyocytes in human carotid tissue samples 69 . Another murine study has clearly identified groups of cells expressing Mmp3 and Cxcl12 as independent transition clusters, but it did not study the biology of these cells 19 . Finally, it is worth noting that single-cell transcriptomic data have been interpreted to show trans-differentiation of SMCs to the macrophage lineage 34 . These findings are not observed in our study or other recently published work 28,29,70 .
Recent work by Chen et al 46 . has employed single-cell studies to investigate the role of Tgfβ signaling in vascular disease, employing a combined Marfan II/Loeys-Dietz and hypercholesterolemia mouse model. A key finding by these investigators was evidence for an SMC-derived mesenchymal stem cell that gives rise to many cell types, including adipocytes, osteoblasts/chondrocytes and macrophages, in the context of Tgf br2 knockout and HFD. In their study, increased plaque inflammation was due, in large part, to increased SMC-derived macrophage number in the vessel wall through this process. In studies reported here, we did not find evidence for an SMC-derived stem cell that mediates this effect and no evidence that SMCs can transition into adipocytes or macrophages in control or Smad3 ΔSMC mice. By contrast, we found that control SMCs give rise to fibromyocyte phenotype cells that subsequently give rise to CMCs and, with Smad3 knockout, the unique cluster of R-SMC-derived cells. The observed differences between the Tgf br2 knockout and Smad3 knockout suggest a fundamentally different mechanism by which SMCs transition in response to signaling through these two different molecules.
Beyond the changes in proportions of the different SMC derivatives, loss of Smad3 also resulted in alterations in SMC transition phenotype transcriptomes as a whole. Gene knockout downregulated several important ECM genes, including Lox, Mfap5 and Eln, whose loss-of-function mutations have been associated with Mendelian aortopathies. These findings suggest that global transcriptomic changes associated with Smad3 loss weaken the vascular wall and, thus, further promote outward vascular remodeling. These findings may also have implications in non-atherosclerotic vasculopathies, such as Marfan and Loeys-Dietz syndromes. Aortopathies such as Marfan's syndrome have been shown to produce aberrant SMCderived populations that contribute to pathogenesis 71 . In fact, our previous scRNA-seq studies of a murine Marfan model also demonstrated an increase in Mmp3 expressing SMC progeny and lower Mfap5 expression, suggesting that our findings here may extend beyond atherosclerotic disease.
Human genetics data have previously suggested the lead CADassociated SNP rs56062135 at 15q22 is in linkage disequilibrium (LD) with SNP rs17293632 that appears to promote AP-1 binding and increase SMAD3 expression in vitro and in vivo 20 , suggesting that higher SMAD3 expression may be associated with risk of myocardial infarction 9,20,30 . This is contradictory to our finding that complete loss of Smad3 increases plaque size in our murine atherosclerosis model, and there are several possible explanations for this disparity. It is possible that alternate SNPs in LD with rs56062135 have the opposite effect on SMAD3 expression in the context of certain types of cellular stimulation-that is, serve as response quantitative trait loci (QTLs). In this case, the response QTLs may have a greater effect on SMAD3 expression, and the integrative effect of the entire haploblock on SMAD3 expression would be opposite and greater than the effect of rs17294632. Alternatively, it is possible that cell fate changes identified in Smad3 ΔSMC mice overall stabilize the human lesion and, thereby, protect it from plaque rupture, despite there being larger lesion size and plaque burden. This paradoxical effect has previously been observed. For example, the IL-1βblocking antibody canakinumab decreased the risk of myocardial infarction in human trials, but Il1r blockade/knockout increased the lumen obstruction and plaque size in mouse models 49,72,73 . Notably, Il1r1 loss substantially changed plaque composition, suggesting that SMC cell fate in plaques may be a stronger determinant for plaque rupture than plaque size alone. Indeed, it has been observed that the largest plaques seen on coronary angiogram are usually not the ones that rupture and cause myocardial infarction 74 . Recent human epidemiological and clinical data also suggest that the quality of calcification is critical and that some types of more calcified plaques are less likely to rupture 68 . It is possible that the increased calcification seen in Smad3 ΔSMC mice is correlated with lesions in humans that are protected against myocardial infarction, which then contributes to the protective genetic signal.

Methods
Mouse strains. SMC-specific lineage tracing and Smad3 knockout were generated by a well-characterized BAC transgene that expresses a tamoxifen-inducible Cre recombinase driven by the SMC-specific Myh11 promoter (Tg Myh11-CreERT2 ; 019079; Jackson Laboratory). These mice were bred with a floxed-stop-flox tdTomato fluorescent reporter line (B6.Cg-Gt(ROSA)26Sor tm14(CAGtdTomato)Hze /J; 007914; Jackson Laboratory) to allow SMC-specific lineage tracing. Smad3 conditional knockout were obtained from the Matzuk laboratory from the University of Texas Southwestern Medical Center 40,41 with LoxP sites flanking exons 2 and 3, which contain Smad3 DNA-binding domain and create a non-functioning frameshift mutation after deletion 42 . All mice were back-crossed onto the C56BL/6 ApoE −/− background. As the Cre-expressing BAC was integrated into the Y chromosome, all lineage tracing mice in the study were male. The animal study protocol was approved by the Administrative Panel on Laboratory Animal Care at Stanford University.

Induction of lineage marker and Smad3 knockout by Cre recombinase.
For all experiments, the tamoxifen gavage schedule was as follows: two doses of tamoxifen, at 0.2 mg g −1 of body weight, were administered by oral gavage at 7-8 weeks of age, with each dose separated by 72-96 hours. HFD was started (101511, Dyets; 21% anhydrous milk fat, 19% casein and 0.15% cholesterol) after the second gavage.
Mouse aortic root/ascending aorta cell dissociation. Immediately after sacrifice, mice were perfused with PBS. The aortic root and ascending aorta were excised, up to the level of the brachiocephalic artery. Tissue was washed three times in PBS, placed into an enzymatic dissociation cocktail (2 U ml −1 , Liberase (5401127001, Sigma-Aldrich) and 2 U ml −1 elastase (LS002279, Worthington) in Hank's Balanced Salt Solution (HBSS)) and minced. After incubation at 37 °C for 1 hour, the cell suspension was strained and then pelleted by centrifugation at 500g for 5 minutes. The enzyme solution was then discarded, and cells were resuspended in fresh HBSS. To increase biological replication, multiple mice were used to obtain singlecell suspensions at each time point. For each scRNA capture, two mice were used. Four separate pairs of isolation were performed for control and Smad3 ΔSMC , but one control 10x capture unexpectedly failed, resulting in a final of three captures of control and four captures from conditional knockout that were included in the analysis. Cells were sorted by fluorescence-activated cell sorting (FACS) based off tdTomato expression. tdT + cells (considered to be of SMC lineage) and tdT − cells were then captured on separate but parallel runs of the same scRNA-seq workflow (gating strategy and threshold identical to those published in previous work by Wirka et al. 24 ), and datasets were later combined for all subsequent analyses.
Single-cell capture and library preparation and sequencing. All single-cell capture and library preparation were performed at the Stanford Functional Genomics Facility and Stanford Genomic Sequencing Service Center. Cells were loaded into a 10x Genomics microfluidics chip and encapsulated with barcoded oligo-dT-containing gel beads using the 10x Genomics Chromium controller, according to the manufacturer's instructions. Single-cell libraries were then constructed according to the manufacturer's instructions (Illumina). Libraries from individual samples were multiplexed into one lane before sequencing on an Illumina platform with targeted depth of 50,000 reads per cell.
Human coronary artery cell acquisition. Human coronary arteries used in this study were dissected from explanted hearts of transplant recipients and were obtained from the Human Biorepository Tissue Research Bank under the Department of Cardiothoracic Surgery from consenting patients with approval from the Stanford University institutional review board, as previously described.
Generation and studies of aortic root sections. Immediately after sacrifice, mice were perfused with 0.4% paraformaldehyde (PFA). The mouse aortic root and proximal ascending aorta, along with the base of the heart, were excised and immersed in 4% PFA at 4 °C for 24 hours. After passing through a sucrose gradient, tissue was frozen in optimal cutting temperature compound to make blocks. Blocks were cut into 7-µm-thick sections for further analysis.
Immunohistochemistry was performed according to standard protocol. Primary antibodies included: anti-SM22alpha rabbit polyclonal primary antibody (1:300 dilution, ab14106, Abcam), an Mmp3 rabbit monoclonal antibody (1:200 dilution, Abcam, 52915) or a CD68 rabbit polyclonal antibody (1:400 dilution, ab125212, Abcam). Secondary antibody was rabbit-on-rodent HRP polymer (RMR622, Biocare Medical). The processed sections were visualized using a Leica DM5500 microscope objective magnification, and images were obtained using Leica Application Suite X software. Von Kossa staining was performed using the Abcam 150687 kit with manufacturer's recommended protocol, with 90-minute development time. All biological replicates for each staining were performed simultaneously on position-matched aortic root sections to limit intraexperimental variance. Folded sections that were uninterpretable after processing were removed.
Sections obtained at equal distance measured from the superior margin of the aortic sinus were used for comparison of lesion features. Lesion size was defined as the area encompassing the luminal edge of the lesion to the border of Tagln + intimal-medial junction-that is, area circumscribed by intimal border of Tagln staining minus the lumen area. Vessel media area was defined as that circumscribed by the outer edge of Tagln staining-that is, the external elastic lamina minus the area circumscribed by inner edge of the Tagln staining-that is, the internal elastic lamina. Areas of interest were quantified in a blinded fashion using ImageJ (National Institutes of Health) software and compared using a twosided t-test.
RNAscope assays. Slides were processed according to the manufacturer's instructions, and all reagents were obtained from ACD Bio. In short, slides were washed once in PBS and then immersed in 1× Target Retrieval reagent at 100 °C for 5 minutes. Slides were washed twice in deionized water, immersed in 100% ethanol and air dried, and sections were encircled with a liquid-blocking pen. Sections were incubated with Protease III reagent for 30 minutes at 40 °C and then washed twice with deionized water. Sections were incubated with commercially available probes against mouse Mmp3, Col2a1 and human MMP3 or a negative control probe for 2 hours at 40 °C. Colorimetric assays were performed per the manufacturer's instructions.
Analysis of scRNA-seq data. FASTQ files from each experimental time point and mouse genotype were aligned to the reference genome (mm10) individually using Cell Ranger software (10x Genomics). Individual datasets were aggregated using the Cell Ranger aggr command without subsampling normalization. The aggregated dataset was then analyzed using the R package Seurat 75 . The dataset was trimmed of cells expressing fewer than 750 genes and genes expressed in fewer than 50 cells. The number of genes, number of unique molecular identifiers and percentage of mitochondrial genes were examined to identify outliers. As an unusually high number of genes can result from a 'doublet' event, in which two different cell types are captured together with the same barcoded bead, cells with more than 6,000 genes were discarded. Cells containing more than 7.5% mitochondrial genes were presumed to be of poor quality and were also discarded. The gene expression values then underwent library size normalization and normalized using the established Single Cell Transform function in Seurat. Principal component analysis was used for dimensionality reduction, followed by clustering in principal component analysis space using a graph-based clustering approach via the Louvain algorithm. UMAP was then used for 2D visualization of the resulting clusters. Lineage inference was performed using Slingshot with available Slingshot software in R, using the converted Seurat object into single-cell experiment objects. Analysis, visualization and quantification of gene expression and generation of gene module scores were performed using Seurat's built-in function such as FeaturePlot, VlnPlot, AddModuleScore and FindMarker. Lists of genes associated with each GO category were obtained from http://geneontology.org/. Panther/DAVID/GO/GREAT analysis was performed using web-based platforms at http://geneontology.org/, http://great.stanford.edu/public/html/ and https://david.ncifcrf.gov/. The top 1,000 genes expressed in modulated SMCs were defined by the highest expressing 1,000 transcripts (based on average expression) from scRNA data in all de-differentiated lineage traced cells. Promoter/5′-regulatory regions of genes were extracted using the UCSC table browser based off 1 kb upstream of the transcription start site of transcripts. Motif analysis was performed using freely available HOMER software 76 with the findMotifGenome function. The regulatory region of the top 1,000 genes was used as the background as bases for motif enrichment.
HCASMC culture/experiments. Human coronary artery SMCs were obtained from Lonza and cultured in smooth muscle growth medium (Lonza, CC-3182) supplemented with human epidermal growth factor, insulin, human basic fibroblast growth factor and 5% FBS, according to the manufacturer's instructions. All HCASMC lines were used at passages 4-8. Small interfering RNA (siRNA) knockdowns were performed using Lipofectamine RNAiMax (Life Technologies) using the manufacturer's recommended protocol at 50 pg of siRNA per 100,000 cells. Cells were allowed to recover in SMC growth medium (with or without additional growth factor) for 36 hours before RNA harvest. Recombinant TGFβ (PeproTech 100-21) concentration used in stimulation was 10 ng ml −1 . Phenotype of the HCASMCs was surveyed by checking expression of lineage-specific markers Myh11, Tagln and Cnn1.

Co-immunoprecipitation experiment.
Myc-Flag-tagged SOX9 (PS100016, OriGene) and 6×-HIS-tagged-HOXB2 (Addgene, 8522) were obtained from commercial vendors and cloned into pCMV6 vector and transfected into HEK cells. The cells were allowed to recover for 36 hours after media change, and nuclear complex co-immunoprecipitation was performed using the commercially available Nuclear-Complex Co-IP Kit from Active Motif (54001) with manufacturer's recommended protocol using mouse anti-Flag (Sigma-Aldrich, F3165) or mouse anti-His (Abcam 18184) antibody for immunoprecipitation, followed by blotting using rabbit anti-Smad3 antibody (Cell Signaling Technology, 9523S), followed by anti-rabbit HRP (Cell Signaling Technology, 7074S) and detected by Luminata Forte Western HRP substrate (Millipore).

Luciferase experiments.
For luciferase experiments, an evolutionary conserved region of human MFAP5 regulatory region (chr12:8815212-8815569) was cloned from human genomic DNA and placed into pLuc-MCS vector, whereas inert/ scramble similar length spacer was cloned into baseline pLuc-MCS as control. pCMV6-empty and cloned pCMV6-Flag-SOX9 or pCMV6-His-HOXB2 were transfected into cells via Lipofectamine 2000 along with respective luciferase and Renilla vector. Media were changed after 6 hours, and dual-luciferase activity (Promega) was recorded after 24 hours using a SpectraMax L luminometer (Molecular Devices). Relative luciferase activity (firefly/Renilla luciferase ratio) is expressed as the fold change over control conditions. Mmp3 activity assay. Mmp3 activity was measured with the Abcam MMP-3 Activity Assay Kit (ab118972) following their tissue-based activity measurement protocol. Next, 0.5 cm of dissected thoracic aorta from identical locations of control and Smad3 ΔSMC mice was placed in the tissue homogenizer for 10 seconds on ice in chilled assay buffer. After centrifugation, the supernatants were then directedly assayed. Mmp3 activity was measured at 10 minutes and 30 minutes after initiation of the reaction, using a SpectraMax luminometer at Ex/Em = 325/393 nm, with exposure of 600 ms. Two biological replicates with three separate segments of thoracic aorta were used for the assay.
Transwell assay. HCASMCs were grown in 24-well plates at low density at 10,000 cells per well. SMAD3 knockdown was performed 24 hours before initiation of the migration experiment. HCASMCs were washed with PBS and cultured in serumfree HCASMC media. Next, 8-µm cell culture inserts (Corning, 353097) were placed into the wells. THP-1 cells (American Type Culture Collection (ATCC)) were grown in standard culture conditions and then spun down, resuspended in serum-free HCASMC media and placed in the top chamber for 3 hours at 37 °C. After 3 hours, THP-1 cells in the bottom chamber in suspension were quantified. TUNEL assay. TUNEL assay was performed on cryopreserved aortic root sections using commercially available chromogenic TUNEL assay kits. Quantification was performed on ×40 magnification at one random point on each cusp, and TUNEL + nuclei were counted manually in a blinded manner. Statistical methods. Differentially expressed genes in the scRNA-seq data were identified using a Wilcoxon rank-sum test. Distribution of cells within defined populations was tested by the chi-square test. Significance determination of histological measurements, luciferase studies, qPCR results and composite gene scores were done by two-tailed t-test, which was performed after verifying normal distribution of the data using the D' Agostino-Pearson test, performed using Prism software. Multiple comparisons were corrected by Bonferroni correction when necessary.
Reporting Summary. Further information on research design is available in the Nature Research Reporting Summary linked to this article.

Data availability
Primary and processed data, along with all relevant metadata, have been deposited to the National Center of Biotechnology Information Gene Expression Omnibus under accession project number PRJNA794806.