Eight-week-old C57BL/6 male mice were procured from Beijing Vital River Laboratory Animal Technology Co., Ltd. (Beijing, China). The Pdgfra${}^{\mathrm{\text{DreER}}}$ mouse line was described previously [14]. Pdgfra${}^{\mathrm{\text{DreER}}}$-tdTomato male mice were generated by crossing the Pdgfra${}^{\mathrm{\text{DreER}}}$ and R26R-rox-tdTomato reporter lines. All mice were housed in specific pathogen-free conditions in an animal room with a 12/12-hour day/night cycle and with ample access to water and food. Animal experiments in this study were conducted in accordance with the National Institutes of Health Guide for the Care and Use of Laboratory Animals. Approval was obtained from the Ethics Committee on Animal Experiments, Shanghai Jiao Tong University School of Medicine (SYXK-2003-0026).

To induce vascular injury, mice underwent carotid artery ligation as previously described [15]. In brief, mice were anesthetized using isoflurane, and the left carotid artery was completely ligated with 3-0 absorbable suture just proximal to the carotid bifurcation.

Tamoxifen was administered to Pdgfra${}^{\text{DreER}}$-tdTomato mice at 5 weeks of age for 5 consecutive days, followed by a 2-week washout period and subsequent carotid ligation. The carotid arteries were then harvested, rinsed with PBS, fixed in 4% paraformaldehyde overnight, embedded in paraffin, and then sliced into 5-µm thick sections. Immunofluorescence staining was conducted as described previously [15]. The primary antibodies used were Pecam1 antibody (Santa Cruze Biotechnology, Dallas, TX, USA, sc-1506-R, 1:100 dilution), Red Fluorescent Proteins (RFP) antibody (Proteintech, Wuhan, China, 5F8, 1:200 dilution), and Pdgfra antibody (R&D, Minneapolis, MN, USA, AF1062, 10 µg/mL). The Alexa fluor-conjugated secondary antibodies employed were donkey anti-rabbit IgG (Invitrogen, Carlsbad, CA, USA, A21202 for Alexa Fluor 488), goat anti-rat IgG (Invitrogen, A11006 for Alexa Fluor 594), and donkey anti-goat IgG (Invitrogen, A21447 for Alexa Fluor 647). Images were obtained using a Zeiss fluorescence microscope (Axio Imager M2, Zeiss, Göettingen, Germany).

Carotid arteries from six C57BL/6 mice were harvested at 2-weeks post-ligation (Supplementary Fig. 1A). The vessels were enzymatically digested to obtain a single-cell suspension. scRNA-seq libraries were then constructed using the 10x Genomics platform. Raw, single-cell RNA-sequencing data was processed using Cell Ranger software (version 6.0, https://github.com/10XGenomics/cellranger). Briefly, the “Cellranger mkfastq” command was employed to de-multiplex raw data and generate FASTQ files. These files were further processed by “Cellranger count” to align reads to the mouse reference, count the number of barcodes and unique molecular identifiers (UMIs), and generate feature-barcode matrices. The raw gene-expression matrix was converted into a Seurat object using the Seurat package (v4.0.5) of R (v4.1.1, University of Auckland, Auckland, New Zealand) [16]. Harmony was utilized to mitigate batch effects across different groups (Supplementary Fig. 1B) [17]. For quality control, only gene features expressed in at least three cells, and cells with at least 200 detected genes, were retained. Subsequently, cells displaying 200 gene features or having $>$15% mitochondrial counts were filtered out. Doublet cells were identified and excluded from the remaining cells using the R package DoubletFinder (v2.0.3) [18]. A total of 17,605 cells remained after filtering. Data normalization was achieved through “LogNormalize”. The top 2000 highly variable genes were selected and scaled using the “ScaleData” function. Principle Component Analysis (PCA) was performed on these selected, highly variable genes. The first 15 principal components with a resolution of 0.4 were employed for cell clustering and visualization via t-distributed stochastic neighbor embedding (t-SNE), or uniform manifold approximation and projection (UMAP). Marker genes upregulated in each cluster were identified through “FindAllMarkers” (min. pct = 0.25, logfc. threshold = 0.25). For sub-clustering of the major cell populations, the aforementioned procedure was repeated. To ascertain similarities between different clusters, correlation matrix heatmaps for all clusters were generated by Scpub [19].

Genes analyzed by Metascape (https://metascape.org/) were obtained by the “FindAllMarkers” method and filtered using the criteria: log${}_2$FC $>$0.25 and p 0.05 (Supplementary Material) [20]. To investigate whether vSMCs-3 dedifferentiate into fibro-like vSMCs, or ECs-4 to ECs-2, related genes from different modules calculated by monocle were subjected to Gene Ontology (GO) function-annotation analysis [21].

The pySCENIC protocol (version 0.11.0, https://github.com/aertslab/pySCENIC) was used as part of a Python package set to default configurations in order to deduce the activity of transcription factors (TFs) and to identify their associated target genes within all sampled cells [22]. This process was segmented into three primary phases. Initially, we refined the single-cell gene expression dataset by removing genes that appeared in 10 cells. The remainder were kept for the construction of co-expression modules with the GRNBoost2 algorithm, which operates on a regression per-target basis. Subsequently, these modules were refined by eliminating indirect targets via the implementation of cis-regulatory motif analysis (cisTarget). Finally, the operational status of these regulons was determined by using the AUCell method (version 3.18, https://bioconductor.org/packages/release/bioc/html/AUCell.html) to calculate an enrichment score for the genes targeted by each regulon.

The “AddModuleScore” function built in the Seurat package was employed to quantify the activity of a specific gene set across different clusters.

Monocle3 (version 1.0.0, https://cole-trapnell-lab.github.io/monocle3/) was used to deduce potential differentiation pathways among various EC subtypes, without any prior assumptions about the sequence or orientation of differentiation [23]. For the analysis of RNA velocity within the subsets of fibroblasts and vSMCs-3, we applied the scvelo tool (version 0.2.4, https://github.com/theislab/scvelo) as recommended by its creators [24]. This allowed investigation of the developmental trajectories from vSMCs-3 to fibro-like vSMCs using the trajectory analysis capabilities of the Monocle algorithm [25].

NicheNet (version 2.0.4, https://github.com/saeyslab/nichenetr) was used for the identification of potential ligands secreted by fibro-like vSMCs. This is an advanced approach designed to predict ligands that may cause changes in the target cell transcriptome [26]. Our analysis also incorporated CellChat (version 1.5.0, https://github.com/sqjin/CellChat) by using its comprehensive database, CellChatDB.mouse. This includes data on 1211 secreted signaling pairs, 432 extracellular matrix (ECM) receptor pairs, 378 cell-cell contact pairs, 229 families of signaling pathways (e.g., WNT, transforming growth factor-$\beta$ (TGF$\beta$), bone morphogenetic protein (BMP), epidermal growth factor (EGF), TGF$\alpha$, fibroblast growth factor (FGF), platelet-derived growth factor (PDGF), vascular endothelial growth factors (VEGF), insulin-like growth factor (IGF)), and chemokine and cytokine signaling pathways (chemoattractant cytokine ligand (CCL), C-X-C motif chemokine ligand (CXCL), interleukins (IL), interferon (IFN)). CellChat assigns a probability of communication between ligand-receptor pairs based on the principle of mass action. It considers the average expression level of ligands in one cell population, and of receptors in another. The significance of these communication probabilities was determined using permutation tests, with p-values 0.05 deemed significant. This analytical framework facilitates the assessment of interaction strengths between various cell clusters identified from the scRNA-seq data, with the results represented graphically as either circle or dot plots.

To create a detailed cellular atlas of NIH at the single-cell level, we performed single-cell sequencing for both self-control and ligation-injured carotid arteries from mice. After quality control and the removal of batch effects between the two groups, 17,605 cells were further analyzed. Unsupervised clustering algorithm yielded 14 distinct clusters, including seven major cell types, and dispersion into two groups (Fig. 1A, Supplementary Fig. 1B). The 14 cell clusters were comprised of 7 cell types annotated with canonical markers: fibroblasts (3 clusters; Pdgfra, Pi16), vSMCs (4 clusters; Myh11, Tagln, Myl9, Acta2), ECs (Pecam1, Cdh5), lymphatic ECs (Lyve1, Pecam1), T cells (2 clusters; Ptprc, Cd3d, Cd3e), monocytes/macrophages (2 clusters; Lyz2, S100a8), and Schwann cells (Mbp, Sox10) (Fig. 1B). The signature genes of all clusters are shown in Fig. 1C.

Fig. 1.

The lineage and characteristics of all cells in healthy (Sham) and ligated (Ligation) carotid artery. (A) Annotated UMAP of mouse carotid artery single-cell RNA sequencing (scRNA-Seq). (B) Relative expression of marker genes in all cells from all samples projected onto a Dotplot. (C) Heatmap showing the mean expression of highest specificity gene markers of each clusters. Color represents average normalized expression level within the denoted cluster. (D) The correlation between different clusters. vSMCs, vascular Smooth Muscle Cells; ECs, Endothelial Cells; LECs, Lymphatic Endothelial Cells; SWs, Schwann Cells; UMAP, Uniform Manifold Approximation and Projection.

Tissue inflammation occurs in response to vascular injury and involves the infiltration and activation of immune cells [27]. Two subgroups of macrophages and two subgroups of T cells are present in the carotid artery, and these cells are mainly found in NIH vessels. The elevated transcription levels of Foxp3 and Satb1 in TCs1 indicate this subpopulation possessed Treg (regulatory T cell) properties (Fig. 1C, Supplementary Fig. 1C) [28]. TCs2 exhibit proliferative features due to high transcription levels of the MKI67 and TOP2A genes [29]. The increased number of lymphatic endothelial cells (LECs) suggests that lymphangiogenesis occurs in NIH induced by vascular ligation. SWs are found in both healthy and diseased carotid arteries, which may be attributed to nerve fibers on the vascular wall. Additionally, subpopulations with different annotations but belonging to the same cell type showed high correlation (Fig. 1D).

These results demonstrate the heterogeneity of fibroblasts, ECs and vSMCs, which are the predominant cell populations within the artery. Our study therefore focused on these three cell types to gain a better understanding of their role in NIH.

To further characterize the heterogeneity of ECs in NIH, we performed a re-clustering of the EC cohort (C8). This identified four distinct subpopulations (ECs1, ECs2, ECs3, and ECs4) based on their transcriptomic profiles. Interestingly, the ECs4 subpopulation was found primarily in the sham group, whereas the frequency of cells from the other subpopulations was notably higher in the ligation group (Fig. 2A). A novel finding was the exclusive presence of ECs3 in the ligation group, characterized by the expression of genes typically associated with vascular fibroblasts, including Dcn and Col3a1 (Supplementary Fig. 2A). These cells express both Pdgfra and Pecam1 and are termed fibro-like ECs (Fig. 2B). Comparative analysis revealed upregulation of mesenchymal markers (Prss23, Cd34, and Ly6a) in ECs1 and ECs2 relative to ECs4. This suggests that ECs1 and ECs2 have a mesenchymal-like phenotype, in contrast to the more endothelial-like phenotype of ECs4 (Fig. 2C).

Fig. 2.

Description of phenotype and function of ECs in NIH. (A) Re-clustering of ECs identified in carotid tissue before and after acute injury. (B) Gene expression of Pdgfra and Pecam1 for each clusters. (C) Violin plots on the expression of selected mesenchymal markers of ECs1, ECs2 and ECs4. (D) UMAP visualization showing the pseudo-time trajectory of the ECs1, ECs2 and ECs4. (E) Heatmap hierarchical clustering showing differentially expressed genes along with the pseudotime curve (left). Gene Ontology biological processes of each DEGs clusters (right). Color key from blue to red indicates relative expression levels from low to high. (F) Top-ranked enrichment analysis of feature genes of ECs3 on Metascape platform. Bars are discretely colored to encode p-values of increasing statistical significance. (G) The immunofluorescent staining of Pecam1 (green), tdTomato (red) and Pdgfra (purple, autofluorescence) in carotid artery from Pdgfra${}^{\text{DreER}}$-tdTomato mice. Nuclei were counterstained with 4${}^{\prime }$,6-diamidino-2-phenylindole (DAPI) (blue). Scale bar = 50 µm (n = 3). A, adventitia; M, media; N, neointimal hyperplasia; NIH, Neointimal Hyperplasia; DEGs, Differentially Expressed Genes; IGF, insulin like growth factor; ECM, extracellular matrix; MAPK, mitogen-activated protein kinase.

Further exploration using the Monocle3 R package revealed a positional dichotomy between ECs2 and ECs4. ECs2 was centrally located, indicating a potential differentiation pathway from ECs4 to ECs1 and ECs2, but not to ECs3 (Fig. 2D). Analysis of an RNA-seq dataset (GSE220512) revealed upregulation of feature genes for ECs1 and ECs2 in injured blood vessels at days 7 and 14 after vascular injury (Supplementary Fig. 2B). Dynamic expression of Teme252, Flt1, Nrp2 and Lrg1 showed an increasing trend along the pseudo-time trajectory (Supplementary Fig. 2C). Our analysis identified 5705 genes displaying four distinct expression patterns across the trajectory. Functional enrichment analysis found that genes associated with cell-substrate adhesion were highly expressed in ECs4, but were downregulated upon transition to ECs1 and ECs2. Conversely, genes implicated in cell proliferation and activation of the inflammasome complex were significantly upregulated in the later stages of pseudo-time (Fig. 2E).

Pathway analysis of ECs3 signature genes using Metascape provided insights into their biological functions, including ECM organization, elastic fiber formation, and cell motility. This suggests that ECs3 have ECM remodeling and cell migration activities in NIH (Fig. 2F). Previous studies have reported trans-differentiation of human fibroblasts into ECs [30]. To investigate whether ECs3 represents a novel subpopulation derived from Pdgfra${}^{\mathrm{+}}$vascular fibroblasts, we utilized Pdgfra${}^{\mathrm{\text{DreER}}}$-tdTomato mice to trace fibroblast lineage cells within the carotid artery (Supplementary Fig. 3A) [31]. Tamoxifen-induced Dre-rox recombination facilitated effective tdTomato labeling of vascular fibroblasts (Supplementary Fig. 3B). This was followed by vessel ligation to induce NIH in these mice (Supplementary Fig. 3C). Of note, intense tdTomato fluorescence within the neointima signified fibroblast migration from the vascular adventitia to the intima (Fig. 2G), with a subset of tdTomato${}^{+}$ cells co-expressing Pecam1 and Pdgfra (Fig. 2G). Negative control observations are documented in Supplementary Fig. 3D.

Together, these findings revealed two subsets of ECs in carotid artery with NIH. These showed mesenchymal characteristics, markers for proliferation and the inflammasome complex, and existed alongside a subpopulation of fibro-like ECs.

Re-clustering of the fibroblasts and vSMCs shown in Fig. 1A revealed four subpopulations of fibroblasts (Fibs-1, Fibs-2, Fibs-3 and Fibs-4), three subpopulations of vSMCs (vSMCs-1, vSMCs-2 and vSMCs-3) and a special subpopulation named fibro-like vSMCs [32]. These were annotated according to their expression of Pdgfra and Myh11 (Fig. 3A,B), and by the top 10 defining genes for each cluster (Supplementary Fig. 4A). vSMCs-2 were mainly distributed in injured artery and expressed osteoblast differentiation-associated genes, including Mef2c, Tnc, Junb, and Tpm4, indicating chondrogenic features (Supplementary Fig. 4B).

Fig. 3.

Functional analysis of different subtypes of fibroblasts. (A) A t-SNE plot of all fibroblasts and vSMCs colored according to cluster (left) and group (right). (B) Relative expression of Myh11 and Pdgfra in fibroblasts and vSMCs projected onto a t-SNE plot. (C) Chemokines, Collagens and Proteinase score calculated by the AddModuleScore function. (D) Heatmap showing chemokine genes, collagen genes, and proteinase genes in eight clusters. (E) The cell‒cell junction scores between any two clusters (left) and cell matrix scores of each clusters (right).

We then analyzed the expression level of all chemokines, collagens and proteinase in all clusters. Of note, Fibs-1 and fibro-like vSMCs exhibited elevated levels of inflammatory chemokines (Cxcl14, Ccl11, Cxcl9, Ccl7 for Fibs-1 and Ccl4, Ccl22, Ccl6, Ccl25), while the latter cells also showed moderate collagen expression and elevated proteinase activity (Fig. 3C,D). Moreover, vSMCs-1 exhibited the highest scores in cell-cell and cell-ECM junctions, whereas fibro-like vSMCs had the lowest score for cell-ECM junctions, hinting at their migratory capability (Fig. 3E).

Fibs-1 emerge as pro-inflammatory, Fibs-3 as collagen-producing, and Fibs-4 as indicative of fibroblasts in homeostasis. This raises the question of the role of Fibs-2 in NIH. By conducting monocle analysis on Fibs-1, Fibs-2, Fibs-3 and Fibs-4, it was discovered that Fibs-3 likely originates from both Fibs-2 and Fibs-4, with Fibs-1 being derived from Fibs-2 (Supplementary Fig. 5A). We also examined the transcriptional level of stemness genes in different fibroblast subpopulations. Significant expression of the mesenchymal marker Gli1 was observed in Fibs-4, indicating this subpopulation was MSC-like fibroblasts (Supplementary Fig. 5B). The highest level of Sca1 transcription was seen in Fibs-2 (Supplementary Fig. 5B). The CytoTRACE algorithm revealed that Fibs-2 and Fibs-4 cells exhibited the highest level of stemness, followed by Fibs-1, and with Fibs-3 showing the least stemness (Supplementary Fig. 5C). The integration of monocle2 and CytoTRACE helps to differentiate the various fibroblast subpopulations (Supplementary Fig. 5D).

The above findings highlight the cellular heterogeneity present in NIH. The results identify fibro-like vSMCs as key players due to their robust chemokine secretion, moderate ECM production, relatively elevated proteinase levels, and diminished adhesion, thereby making them potential targets for NIH therapy. In addition, these results suggest that Fibs-2 are a minor, yet critical cell population within healthy blood vessels, and are capable of proliferation and differentiation into pro-inflammatory fibroblasts following vascular injury.

Previous scRNA-seq analyses conducted on atherosclerotic lesions from murine and human arteries have shown that vSMCs tend to shift toward genuine fibroblasts, while remaining distinct from these [13]. We hypothesized that this transdifferentiated phenotype also exists in NIH. The transcriptional kinetics from vSMCs-3 to fibro-like vSMCs was revealed by RNA velocity, which is a high-dimensional vector that predicts the future state of individual cells over several hours (Fig. 4A) [33]. The R package Monocle was used to construct the comprehensive lineage differentiation trajectory for vSMCs-3 to fibro-like vSMCs, highlighting the overlap between the two clusters on this trajectory (Fig. 4B). Subsequent analysis of gene expression patterns during the cell state transitions along this trajectory identified 6705 dynamically expressed genes, categorizing them into four distinct patterns. Functional enrichment analysis was carried out on three of these patterns (Fig. 4C). Canonical vSMC genes, such as Myh11 and Cnn1, were down-regulated, whereas fibroblast-specific genes such as Pdgfra and Dcn were up-regulated (Fig. 4D). Furthermore, SCENIC analysis found that Dlx5, Ets2, Rxrb, Nfya and Batf were the top five most active regulons in individual cells (Fig. 4E,F).

Fig. 4.

The transdifferentiation of vSMCs-3 into fibro-like vSMCs. (A) RNA velocity overlaid on tSNE of the extracted vSMC and fibroblast subset. (B) The developmental trajectory of selected clusters, coloured-coded by the groups (upper) and associated cell subpopulations (lower). (C) Heatmap hierarchical clustering showing differentially expressed genes along with the pseudotime curve (left). Gene Ontology biological processes of each DEGs clusters (right). Color key from green to purple indicates relative expression levels from low to high. (D) Vlnplot showing the expression levels of Myh11, Cnn1, Pdgfra and Dcn in fibro-like vSMCs and vSMCs-3. (E) The regulatory activity of all clusters. (F) Scatter plot showing the specificity scores of regulons of fibro-like vSMCs. The top 6 regulons are highlighted.

Taken together, the above results suggest that the key pathogenic cell subpopulation, fibro-like vSMCs, may be derived from vSMCs-3.

CellChat was employed to build a cell–cell-communication network based on known ligand-receptor pairs and their cofactors [34]. Fibro-like vSMCs had distinct effects on different cell types, and interact mainly with ECs (Fig. 5A,B). NicheNet analysis [26] revealed that Vegfa was secreted predominantly by fibro-like vSMCs (Fig. 5C) and may drive biological functions in other cell types. In NIH, Vegfa influences ECs via Vegfr2 (Fig. 5D). Additionally, Vegfa expression was significantly higher in wire-induced or ligation-induced vessels, as shown by the RNA-seq datasets GSE70410 and GSE220512 (Fig. 5E).

Fig. 5.

Cell–cell communications of fibro-like vSMCs with other cells. (A) Arrow plot showing the interaction strength between fibro-like vSMCs with other cell populations, based on the interaction scores. (B) Heatmap showing the number of interactions between different structual cells in vessel. Blue block indicate that the displayed communication is decreased in ligated carotid artery (CA), whereas red cblock indicate that the displayed communication is increased in ligated CA compared with healthy carotid artery. The darker the color, the greater the communication probability between the two cell type. (C) Genes are ligand that predicted by Nichenet, secreted by fibro-like vSMCs, ranked by prioritized liand activity. The darker the color, the higher the ranking. Dot plot showing average expression of ligand in all clusters. (D) Dot plot showing vascular endothelial growth factor (VEGF) related chemokine-receptor interaction pairs between fibro-like vSMCs and three ECs subsets or LECs, based on the communication probability. (E) Barplots showing expression level of Vegfa in different group of different database. **, p $\le$ 0.05; *, p $\le$ 0.05; ns = not significant. Student’s t test (left) and paired Student’s t test (right).