https://orcid.org/0000-0002-2766-8882
Abstract
The role of vascular smooth muscle cells (SMCs) in atherosclerosis has evolved to indicate causal genetic links with the disease. Single cell RNA sequencing (scRNAseq) studies have identified multiple cell populations of mesenchymal origin within atherosclerotic lesions, including various SMC sub-phenotypes, but it is unknown how they relate to patient clinical parameters and genetics. Here, mesenchymal cell populations in atherosclerotic plaques were correlated with major coronary artery disease (CAD) genetic variants and functional analyses performed to identify SMC markers involved in the disease.
Bioinformatic deconvolution was done on bulk microarrays from carotid plaques in the Biobank of Karolinska Endarterectomies (BiKE, n = 125) using public plaque scRNAseq data and associated with patient clinical data and follow-up information. BiKE patients were clustered based on the deconvoluted cell fractions. Quantitative trait loci (QTLs) analyses were performed to predict the effect of CAD associated genetic variants on mesenchymal cell fractions (cfQTLs) and gene expression (eQTLs) in plaques.
Lesions from symptomatic patients had higher fractions of Type 1 macrophages and pericytes, but lower fractions of classical and modulated SMCs compared with asymptomatic ones, particularly females. Presence of diabetes or statin treatment did not affect the cell fraction distribution. Clustering based on plaque cell fractions, revealed three patient groups, with relative differences in their stability profiles and associations to stroke, even during long-term follow-up. Several single nucleotide polymorphisms associated with plaque mesenchymal cell fractions, upstream of the circadian rhythm gene ARNTL were identified. In vitro silencing of ARNTL in human carotid SMCs increased the expression of contractile markers and attenuated cell proliferation.
This study shows the potential of combining scRNAseq data with vertically integrated clinical, genetic, and transcriptomic data from a large biobank of human plaques, for refinement of patient vulnerability and risk prediction stratification. The study revealed novel CAD-associated variants that may be functionally linked to SMCs in atherosclerotic plaques. Specifically, variants in the ARNTL gene may influence SMC ratios and function, and its role as a regulator of SMC proliferation should be further investigated.
An overview of the scientific questions, workflow design, and main findings of the study. CAD, coronary artery disease; GWAS, genome-wide association study; scRNAseq, single cell RNA sequencing; SMC, smooth muscle cell; SNP, single nucleotide polymorphism.
SMC-derived cells are the dominant cell type in human atherosclerotic plaques and identifying the genes associated with various SMC sub-phenotypes offers a potential for therapeutically targeting SMC phenotypic switching in atherosclerosis. This study integrates scRNAseq data with per patient-matched clinical, genetic, and transcriptomic data from a large biobank to stratify patient risk prediction and to identify genes associated with SMC sub-phenotypes in human carotid plaques. Notably, several variants in ARNTL gene associated with circadian rhythm were discovered, which may regulate SMC proliferation in atherosclerosis.
This editorial refers to ‘Dissecting the effect of genetic variants on atherosclerosis: integrating bulk and single-cell transcriptomics’, by E. Duregotti and M. Mayr, https://doi.org/10.1093/eurheartj/ehae489.
Introduction
Smooth muscle cells (SMCs) are an integral part of the arterial wall that controls vascular tone and blood pressure.1 In adult vessels, SMCs are not terminally differentiated and in response to exogenous stimuli, can undergo phenotypic switching,2 attaining a proliferative and migratory state. In atherosclerotic plaques, transdifferentiation of multiple cell types of mesenchymal origin, leads to a variety of sub-phenotypes with both benefit3 and/or detriment4 to the progression of the disease. For example, adventitial fibroblasts contribute to vascular remodelling by differentiating into myofibroblasts, which have high proliferative capacity, coupled with expression of markers such as α-smooth muscle actin (α-SMA) similar to SMCs and secretion of inflammatory cytokines.5,6 Pericytes, typically found in the vasa vasorum lining neovessels, are also highly plastic and difficult to distinguish from SMCs in plaques solely by the expression of classical markers. Thus, SMCs, fibroblasts and pericytes collectively represent mesenchymal cells in the vascular wall, all reported to contribute to the healing response after arterial injury and in atherosclerotic plaque formation.7,8 However, whether patient genotypes and phenotypes are associated with mesenchymal cell ratios in late-stage atherosclerotic plaques, have not been demonstrated.
In patients with coronary artery disease (CAD), genome-wide association studies (GWAS) have identified nearly 300 genetic variants with causal links to the disease.9–11 Of these loci, about 30% have been mapped near the genes regulating lipoprotein metabolism and blood pressure, two of the most common risk factors for cardiovascular disease (CVD).12,13 Recent advances in genetics and genomics have also provided links between these loci and different cell types in the vasculature,14 particularly genes with vital functions in SMCs.3,15–17 Nevertheless, full knowledge of how CVD genetics influences the contribution of different SMC sub-types in atherosclerosis is still missing, as well as which specific genes may mediate this effect, ultimately influencing overall plaque vulnerability.
Single cell RNA sequencing (scRNAseq) has improved the characterisation of different cell populations in atherosclerosis, such as the identification of various sub-populations of T-cells and macrophages,18–20 as well as understanding the role of genetic variants in the disease.3,21,22 However, due to large costs,23 current scRNAseq studies of plaques have been restricted to low-sample numbers, failing to capture the full heterogeneity of human disease. Alternatively, the proportion of different cell types can be estimated by deconvolution of bulk gene expression data from many samples,24 using dedicated software such as CIBERSORTx25,26 that has been experimentally validated also for atherosclerotic plaques.26,27 Such strategies can take advantage of publicly available scRNAseq data, scaling up single cell analyses from a few patients to large human disease cohorts. Importantly, these approaches can result in identifying novel atherosclerosis targets relevant for the contribution of mesenchymal cells in plaque instability, with increased translational potential.
In this study, we deconvoluted microarray data from 125 carotid plaques in the Biobank of Karolinska Endarterectomies (BiKE)28 using publicly available scRNAseq data from human coronary lesions.3 These data were then investigated for associations with patient phenotypes, i.e. symptomatology, medication and other clinical parameters. Moreover, we analysed associations between a large GWAS catalogue of published CAD genetic variants and mesenchymal cell fractions in plaques, which yielded novel cell type—specific signatures in symptomatic atherosclerotic patients, studied further using functional in vitro assays.
Methods
Human carotid atherosclerosis cohort
Patients undergoing carotid endarterectomy at the Department of Vascular Surgery, Karolinska University Hospital, Stockholm, Sweden were consecutively included in BiKE and clinical data recorded upon admission. Individual data for death, major adverse cardio- and cerebro-vascular events (MACCEs), comorbidities and blood parameters were extracted from the national registries, patient charts, and BiKE database. Carotid endarterectomies and blood were collected during surgery and retained within BiKE. Plaques were divided transversally at the most stenotic part and the proximal half of plaques was used for transcriptomic analysis. The distal half was fixed in 4% zinc-formaldehyde and processed for histology. Microarray analysis was performed on RNA extracted from plaques (n = 127; symptomatic 87+ asymptomatic 40), normal arteries (n = 10) and peripheral blood mononuclear cells. Normal, non-atherosclerotic arteries were nine disease-free iliac arteries and one aorta from organ donors without a history of CVD. The vessels were dissected and the intima and media used for RNA isolation and arrays. One normal carotid artery was fixed in 4% zinc-formaldehyde and processed for histology. Symptomatic patients were categorized as minor stroke, transient ischemic attack, or amaurosis fugax and patients without qualifying symptoms within 6 months prior to surgery were considered asymptomatic. The use of human samples from BiKE is approved by the Ethical Committee of North Stockholm and follows the guidelines of the Declaration of Helsinki. All human samples were collected with informed consent from patients or organ donors’ guardians.
Additional methods can be found in the Supplementary Material.
Results
Cell fractions obtained from plaque deconvolution reveal differences in lesion composition associated with patient phenotype and survival
Deconvolution of the BiKE microarray data from plaques (n = 127) and normal arteries (n = 10) was initially performed using the signature matrix from scRNAseq data of human coronary plaques3 (n = 5, workflow shown in Figure 1A, raw deconvolution data in Supplementary data online, Tables S1–S3). Comparison of cell fractions from normal arteries and BiKE plaques confirmed that plaques had significantly higher proportion of immune cells, i.e. Type 1 and 2 macrophages, T-cells and natural killer (NK) cells, while fractions of both classical and non-typical SMCs were decreased (see Supplementary data online, Figure S1A). Cell fractions were then stratified based on clinical parameters (sex, medications, and symptoms). Smoking, diabetes, or statin therapy did not appear to influence the plaque cellular composition significantly (Figure 1B, Supplementary data online, Figure S1B). Major differences were observed in the cell proportions when stratifying according to patient sex. Plaques from women had more Type 1 macrophages and pericytes, while relatively less classical SMCs even when age was considered (Figure 1B, Supplementary data online, Figure S1C). Plaques from symptomatic patients had a higher fraction of Type 1 macrophages and pericytes (Figure 1B). The most striking difference was that plaques from symptomatic patients had significantly less classical (SMC1) and modulated SMCs than those from asymptomatic patients (Figure 1B). This indicated that alterations in mesenchymal sub-phenotype ratios might have influenced plaque vulnerability and symptomatology in atherosclerotic patients.
Deconvolution of Biobank of Karolinska Endarterectomies (BiKE) microarrays reveals differences in plaque cell composition with respect to patient symptoms and sex. (A) An illustration of the deconvolution workflow: CIBERSORT was used to deconvolve cell fractions from BiKE plaque microarrays (n = 127) using the publicly available gene expression signature matrix from the scRNA sequencing data of coronary arteries (n = 5). (B) The cell fractions from n = 125 plaques (due to the lack of clinical data for two patients) were stratified by clinical parameters such as symptoms, sex, and statin treatment. Statistical differences were calculated using Student’s t-test and did not survive a correction for multiple comparisons using FDR = 0.05. Data represented as mean ± SD Correlations among cell fractions were performed in plaques from symptomatic (C) and asymptomatic (D) patients. Asymptomatic patients n = 40, symptomatic n = 87. Correlations were calculated using Pearson correlation, coefficient indicated in the legend below the plots, *P < .05; **P < .01; ****P < .0001. (E) Kaplan—Meier plots showing the risk of future MACCEs in patients stratified by higher (75%) and lower (25%) quartile of deconvoluted Type 1 macrophage cell fractions. Statistical differences between survival curves were computed using log-rank test
To cross validate the distribution of cell fractions, we also deconvolved two additional public bulk transcriptomic datasets, one from carotid lesions29 and one from coronary lesions,30 using the scRNAseq data from coronary plaques. As shown in Supplementary data online, Figure S2A, the proportions of cell fractions deconvoluted from the carotid bulk dataset were similar to the ones from BiKE (Figure 1B, Supplementary data online, Figure S1A). Interestingly, there was a marked difference in the proportion of cell fractions deconvoluted from coronary bulk dataset (see Supplementary data online, Figure S2B), highlighting previously reported differences between the carotid and coronary vascular beds.31 Since these biological differences could contribute to the skewing of cell fractions in our study, we later aimed to validate our findings using available scRNAseq data from carotid plaques.31 Thus, BiKE carotid plaque microarrays were deconvoluted again using carotid scRNAseq data, which confirmed that the cell fraction distribution in BiKE plaques was indeed dominated by mesenchymal and macrophage/monocyte fractions (see Supplementary data online, Figure S2C). Interestingly, among the mesenchymal cell types there was a higher proportion of modulated SMC signature (fibromyocyte and fibrochondrocyte) in plaques. Overall, our results from both coronary-to-carotid and carotid-to-carotid plaque cell fraction deconvolution were in good agreement.
Plaque cell fractions deconvoluted by coronary scRNAseq data were then correlated with each other, and correlations stratified based on patient symptoms. Mesenchymal cell fractions were negatively correlated with macrophage, T-cell, and NK cell fractions. This tendency was expected and there was no difference between symptomatic and asymptomatic patients (Figure 1C and D). Interestingly, endothelial cell (EC) fraction 1 showed positive correlations with immune cells such as Type 1 macrophages, B-, T-, and NK-cells and negatively correlated with all SMC cell fractions, specifically in symptomatic patients (Figure 1C and D). These results could indicate the presence of typical plaque vulnerability features in symptomatic patients, such as increase in neovessel formation and immune cell infiltration.
To find out if the different cell clusters might contribute to the onset of symptomatic events, we performed a survival analysis for ischemic stroke, myocardial infarction (MI), MACCE, or death between patients with higher quartile (75%) and lower quartile (25%) of cell clusters in their plaques. Our analysis showed that patients with higher fraction of Type 1 macrophages had a significantly higher risk of MACCEs during long-term follow-up after carotid endarterectomy (Figure 1E and Supplementary data online, Figure S3). Age and sex did not affect the incidence of MACCEs as shown by the Cox proportional hazard analysis (see Supplementary data online, Table S4A). Of note, when the Macrophage 1 cluster was considered as a continuous variable (as opposed to categorical shown in Figure 1E), the difference in the incidence of MACCE lost significance (see Supplementary data online, Table S4B). Taken together, these results show that differences in plaque composition, in particular the ratio between SMCs and Type 1 macrophages, may be responsible not only for symptom onset prior to surgery, but also influence the association to future poor outcomes in patients, corroborating recent findings that changes in plaque SMCs are an important factor associated with poor outcomes after surgery in men.32
Clustering based on cell fractions identifies vulnerable patients
To identify molecular plaque phenotypes and their association with patient clinical data, an unbiased clustering of patients using plaque cell fractions was performed, identifying three distinct clusters with minimal overlap (Figure 2A). Cluster 2 was slightly smaller in number of patients (n = 23) compared with Cluster 1 (n = 53) and Cluster 3 (n = 49). As shown in Figure 2B, the clusters varied significantly in their cell fraction content, providing insights into the differences of the underlying plaque compositions. Cluster 1 was characterized by higher fraction of Type 1 and Type 2 macrophages, Cluster 2 was dominated by SMC1, SMC2, and modulated SMCs, while Cluster 3 contained significantly higher proportions of all immune cells, i.e. macrophages, T cells, B cells, and NK cells combined with an increased proportion of pericytes and typical ECs (Figure 2B). These results implied that plaques from patients from Cluster 2 may have a relatively stable phenotype, while plaques from Clusters 1 and 3 represent an inflammatory phenotype, with Cluster 3 further showing features of rupture prone lesions such as angiogenic vessels.
Clustering of the patients based on plaque cell fractions reveals distinct phenotypes associated with clinical symptoms. (A) Biobank of Karolinska Endarterectomies patients were clustered based on plaque cell fractions using a hierarchical K-means clustering algorithm. Middle Cluster in green represents Cluster 1 (n = 53 patients), blue represents left Cluster 2 (n = 23) and red represents right Cluster 3 (n = 49). (B) Plaque cell fractions were stratified according to the 3 patient clusters. Statistical differences were calculated using one-way ANOVA and corrected for multiple comparisons using FDR = 0.05. Corrected P-values are displayed in each plot. Tukey’s post hoc test was performed to identify statistical differences between the clusters, *P < .05. (C) Kaplan—Meier plot showing the long-term risk of future ischemic stroke after carotid endarterectomy in patients stratified by three clusters (number of stroke events—Cluster 1: 50, Cluster 2: 22, Cluster 3: 47). Statistical differences between the survival curves were computed using log-rank test
A detailed analysis of differences in blood parameters and comorbidities among the three clusters (see Supplementary data online, Table S5), showed that there was a predominance of chronic kidney disease among patients in Cluster 2 (45%). Importantly, we found that the percentage of patients with an incidence of ischemic stroke prior to surgery was significantly higher in Cluster 3 compared with other two clusters (see Supplementary data online, Table S5). Moreover, long-term survival analysis after carotid endarterectomy also showed a tendency towards higher risk for ischemic stroke in Cluster 3 (Figure 2C), confirming that patients from Cluster 3 harboured vulnerable plaques. Of note, survival analyses for other cardiovascular endpoints such as MACCE, MI, and death did not show significance (see Supplementary data online, Figure S4).
To minimize the vascular bed bias in risk prediction analyses, we then performed clustering of BiKE patients based on the cell fractions deconvoluted from carotid scRNAseq data. The patients again grouped into three clusters (see Supplementary data online, Figure S5A), among which the proportion of patients with an incidence of MI was significantly higher in Cluster 1 compared with Clusters 2 and 3 (see Supplementary data online, Figure S5B). Patients in Cluster 1 had a higher presence of inflammatory cells such as monocytes, dendritic cells, foamy macrophages, mast cells and activated NK cells, angiogenic ECs and pericytes, and relatively lower presence of mesenchymal sub-clusters including fibromyocytes (see Supplementary data online, Figure S5C), shown to confer stability in plaques.3
Taken together, these analyses showed that plaque cell fractions strongly associate with patient clinical data and can be used to refine the patient vulnerability profile, risk stratification, and even long-term prognosis.
Mesenchymal cell—specific CAD risk loci link with gene expression in plaques and patient symptoms
To investigate how CAD genetics may link specifically with mesenchymal cell fractions in plaques, we developed a novel strategy where a catalogue of CAD genetic variants was assembled from literature and GWAS studies9–12 and filtered based on their association with mesenchymal cell fractions (cfQTLs) in BiKE (Figure 3, Supplementary data online, Tables S5 and S6). Concurrently, we also associated CAD variants directly to patient symptoms, however as shown in Supplementary data online, Table S8, only three variants were significantly associated with symptoms in our cohort. As an alternative strategy, we sought to identify genes in mesenchymal cells that may be modified by CVD genetics by performing eQTL analysis on the cfQTL variants and their proxies (3152 in total), with and without symptom as a covariate. In this way, we obtained variants associated with gene expression in mesenchymal cell fractions and indirectly influenced by the presence of symptoms (see Supplementary data online, Table S9).
Strategy for obtaining mesenchymal cell-specific CAD risk loci associated with target gene plaque expression and patient symptoms. A schematic representation of the strategy for filtering CAD-associated GWAS loci using their association with mesenchymal cell fractions (cfQTL), gene expression in Biobank of Karolinska Endarterectomies plaques (BiKE eQTL), and with respect to patient symptom as a covariate. The resultant mesenchymal cell-specific, symptom-specific BiKE eQTLs were used for further functional studies
Interestingly, many of the top SNPs from this analysis were linked with genes previously shown to be relevant for atherosclerosis. The top identified SNP rs9595908, an intronic variant in the gene PDS5 cohesin adhesion factor 5 (PDS5B), was a BiKE eQTL for the neighbouring gene NEDD4 Binding Protein 2 Like 2 (N4BP2L2). CAD-GWAS SNP rs9513116 was a BiKE eQTL for FLT1, the role of which in SMC function in the context of atherosclerosis has been reported before.33 Another interesting SNP was rs8141797, a missense variant in the Sushi domain containing 2 (SUSD2) and BiKE eQTL for matrix metalloproteinase-11 (MMP11). The presence of MMP11 in atheroma regulated via CD40-CD40L is well-known.34 One more missense variant rs8141797 in SUSD2 was a BiKE eQTL for Macrophage migration inhibitory factor. The importance of this gene in atherosclerosis and its abundance in various plaque cell types has already been investigated.35 Although these genes have been studied in atherosclerosis before, mostly in macrophages, our analysis shows that they may have yet unknown functions for mesenchymal cell biology too.
Most notably, multiple SNPs on top of this list were found to be BiKE eQTLs for Aryl Hydrocarbon Receptor Nuclear Translocator Like (ARNTL), suggesting that this gene could play a previously unexplored role in the context of atherosclerosis and mesenchymal cells in plaques. ARNTL is part of the core complex of Clock genes regulating circadian rhythm in mammals,36 a mechanism shown to be of key importance in patients with atherosclerosis.37
Risk alleles of top mesenchymal—specific variants in the ARNTL locus associate with decreased expression of ARNTL in carotid plaques
By inspecting the UCSC genome browser, all six eQTLs for ARNTL were located at the intergenic and regulatory regions of this gene38 (Figure 4A). We found that the minor alleles of all these variants were associated with significantly lower expression levels of ARNTL in carotid plaques (Figure 4B), but also with increased presence of SMC2 and modulated SMC cell fractions (see Supplementary data online, Figure S6). Further evidence from publicly available phenome-wide analysis studies from the Open Target Genetics webtool,39 showed that the most significant variant rs4757138 was also significantly associated with CAD and other CVD traits such as blood pressure and hypertension (see Supplementary data online, Table S10). Thus, our analysis identified several variants in ARNTL gene as mesenchymal cell-specific eQTLs in atherosclerotic plaques.
Several top mesenchymal cell—specific variants are eQTLs for ARNTL gene expression in plaques. (A) Snapshot view of the top ARNTL variants on the UCSC genome browser. The view shows that the variants are located at the intergenic (rs900145, rs2403661, rs11022742, and rs11605518) and regulatory regions (rs998089, rs4757138) of ARNTL (blue lines indicate each of the six variants). (B) Gene expression QTL in plaques was calculated for the top mesenchymal cell—specific variants in ARNTL—rs4757138, rs11022742, rs998089, rs900145, rs2403661, and rs11605518. eQTL—expression quantitative trait locus. Statistical differences were calculated using one-way ANOVA and Tukey’s post hoc test was performed to identify statistical differences between individual alleles, *P < .05. (C) Association of ARNTL variants rs2403661, rs4757138, and rs900145 with sub-clinical atherosclerosis phenotypes, i.e. the presence of plaques and growing carotid intima-media thickness, was tested in the IMPROVE study cohort
To verify the importance of ARNTL eQTL variants from BiKE, we investigated them further using two independent large vascular cohorts, namely the IMROVE40 and ASAP studies.41 In the IMPROVE cohort, consisting of individuals with at least 3 CV risk factors, the 3 ARNTL variants (rs4757138, rs2403661, and rs900145) associated with a significantly higher odds ratio for the presence of carotid atherosclerotic plaques, defined as carotid intima media thickness ≥1.5 mm (Figure 4C). In the ASAP cohort, which enrols patients with aortic valve and ascending aortic disease undergoing elective open-heart surgery, eQTL analysis was performed between the variants and ARNTL gene expression separately in the aortic intima-media and adventitia of patients with bicuspid aortic valve. Here, patients with the risk alleles of the same 3 variants showed a significant increase in ARNTL expression levels compared with that of the heterozygous allele in adventitia (see Supplementary data online, Figure S7A). In media, there was no significant association between ARNTL levels and the variants, although, there was a positive trend (see Supplementary data online, Figure S7A). Together, the addition of results from IMPROVE and ASAP cohorts validate and extend BiKE data and suggest a complex relationship between ARNTL genetics and vascular pathologies in general.
ARNTL is expressed by SMA+ mesenchymal cells in carotid plaques
Next, we sought to investigate the overall expression, localisation, and role of ARNTL in plaques in more detail. As shown in Figure 5A, the gene expression of ARNTL was significantly reduced in plaques compared with normal arteries in bulk microarray data. Despite this overall down-regulation, it is noteworthy that residual ARNTL expression based on scRNAseq data was relatively broad among various cell types in plaques, including immune cells, SMCs, fibroblasts, and pericytes (see Supplementary data online, Figure S7B). Indeed, BMAL1 protein encoded by the ARNTL gene could be localized in SMA + SMCs in the normal aortic media as well as the carotid plaque fibrous cap and in the sub-intima of plaques (Figure 5B). Nevertheless, overall ARNTL gene expression in BiKE plaques was negatively correlated with typical contractile SMC markers such as SMA alpha 2 (ACTA2), myosin heavy chain 11 (MYH11), smoothelin (SMTN), calponin 1 (CNN1), leiomodin 1 (LMOD1) (Figure 5C), as well as SMC cell fractions (SMC1, SMC2, and modulated SMC) (Figure 5D).
ARNTL is down-regulated in carotid plaques and negatively correlated with typical SMC markers. (A) Boxplot of ARNTL gene expression in Biobank of Karolinska Endarterectomies (BiKE) normal arteries and carotid plaques. Statistical difference was calculated using student t-test, **P < .01. (B) Representative images of immunohistochemical staining of BMAL1 (ARNTL) and α-SMA in normal arteries, as well as remnants of media and fibrous cap regions of human carotid plaques. Scale bar for original image = 500 µm; scale bar for inset image = 100 µm. (C) ARNTL expression was negatively correlated with SMC contractile markers (ACTA2, MYH11, SMTN, CNN1, and LMOD1) in BiKE plaques. (D) ARNTL expression in plaques was negatively correlated with mesenchymal cell fractions, i.e. SMC1, SMC2, modulated SMC, but positively with pericytes and marginally also with fibroblasts. Correlations were computed using Spearman rho method and P-values after correction for multiple testing (FDR = 0.05) are displayed within the correlation plots. (E) Box plot of ARNTL gene expression in each of the three patient clusters. Statistical differences were calculated using one-way ANOVA and Tukey’s post hoc test was performed to identify statistical differences between clusters, ***P < .001. (F) Global correlations of ARNTL gene expression with all genes in BiKE plaques were performed, followed by GSEA. Highly enriched Hallmark pathways from the GSEA analysis are shown in a dot plot. ACTA2, α-SMA; CNN1, calponin1; MYH11, myosin heavy chain-11; LMOD1, leiomodin1; SMTN, smoothelin
ARNTL gene expression also segregated among the three patient clusters defined by deconvoluted cell fractions, where Cluster 2 had significantly reduced ARNTL levels compared with Clusters 1 and 3, consistent with its negative correlations to SMC cell fractions and the distribution of SMC markers in these clusters (Figure 5E, Supplementary data online, Figure S7C).
Furthermore, we probed for the expression of Arntl in rat carotid arteries subjected to balloon catheter injury, since the dominance and functional relevance of SMCs has been well-studied in this response-to-injury model.42,43 When comparing the gene expression levels of Arntl with various SMC markers throughout the longitudinal healing process (see Supplementary data online, Figure S8A), it was evident that Arntl expression increased from 0 to 2 h and from 2 to 5 days post-injury, while the expression of classical cytoskeletal SMC markers (except transcription factor Myocd) decreased through the same timepoints, in agreement with the negative correlations obtained also in human plaques. At other time points the expression of Arntl followed a similar pattern to that of classical SMC markers. In addition, immunostaining for Bmal1 in the injured vessels confirmed its localisation in SMA+ cells on protein level (see Supplementary data online, Figure S8B).
ARNTL is associated with regulation of cell growth processes in plaques
To investigate the possible influence of ARNTL on the overall plaque phenotype, we performed global gene correlations with ARNTL in BiKE microarrays (see Supplementary data online, Table S11) coupled to gene set enrichment analysis (GSEA) of the highly correlated genes. The enriched pathways were related to adipogenesis, heme metabolism, glycolysis, mTORC1 signalling, and MYC targets (Figure 5F and Supplementary data online, Table S12). Most of these pathways are linked with cell growth,44–47 suggesting a mechanistic role for ARNTL in atherosclerosis.
To gain a broader understanding of the genes that may be mechanistically associated with ARNTL, a protein–protein interaction network was created using the webtool Stringdb (see Supplementary data online, Figure S9). ARNTL was found to be connected with factors responsible for muscle cell differentiation (MYOD1) and proliferation (CTNNB1, MDM2, and MYB). Apart from coordinating the transcription of Clock genes (CRY1, CRY2, and PER1), ARNTL was also directly connected to the various epigenetic modulators (EP300, HDAC1, and SIRT1), while indirectly to genes controlling hypoxia (HIF1A), inflammation (STAT1, NFKB1, and JUN) and steroid hormone signalling (AR and ESR1).
ARNTL silencing inhibits SMC proliferation
The next aim was to characterize the functional relevance of ARNTL in human primary carotid SMCs in vitro. In response to PDGF-BB stimulation that has been associated with ARNTL activation in previous studies, ARNTL gene expression increased significantly (Figure 6A). ARNTL silencing, coupled to significantly decreased ARNTL gene expression and a modest decrease in BMAL1 protein (see Supplementary data online, Figures S10A, B and S11), led to a decrease in the proliferative capacity of these cells (Figure 6B), with a concordant increase in gene expression of SMC contractile markers (ACTA2, CNN1, LMOD1, and SMTN) compared with controls. The protein levels of α-SMA did not significantly increase with ARNTL silencing, (Figure 6C and D, Supplementary data online, Figure S11).
Silencing ARNTL in human carotid SMCs inhibits proliferation and induces cellular senescence. Human carotid SMCs were treated with siRNA against ARNTL followed by functional assays. SMCs were also treated with or without PDGF-BB after ARNTL silencing and ARNTL gene expression levels were measured (A). (B) Proliferation was measured as a function of cell confluence using Incucyte Zoom system (Essen BioScience) (n = 4). (C) Gene expression levels of contractile SMC markers were measured using qRT-PCR (n = 4) and (D) protein levels of contractile protein α-SMA (42 kDa). Protein levels were quantified using immunoblots from 3 experiments and representative blots are shown. ACTA2, α-SMA; CNN1, calponin1; MYH11, myosin heavy chain-11; LMOD1, leiomodin1; SMTN, smoothelin; α-SMA. Data points are represented as mean ± SEM. Statistical differences were calculated using student t-test. *P < .05, **P < .01, ***P < .001
The reciprocal relationship between Clock genes and cell cycle has been well studied,48 while a recent publication also showed that ARNTL non-transcriptionally regulates cellular senescence in mesenchymal progenitor cells.49 We detected an increase in expression levels of genes involved in mTORC1 and glycolysis pathways (see Supplementary data online, Figure S10C and D), concordant with publications showing that higher glycolytic activity is typical for senescent cells, but also with pathways enriched in association with ARNTL in BiKE plaques. Given such evidence, we investigated if the decrease in cell proliferation was due to ARNTL down-regulation-induced cellular senescence. Indeed, there was a significant increase in mRNA expression for key senescence markers such as CDKN2A, CDKN1A, RB1, and TP53 (see Supplementary data online, Figure S10E), although on protein level the increase in p16-INK4A was not significant (see Supplementary data online, Figures S10F and S11). Additionally, we also detected higher β-galactosidase activity in ARNTL silenced cells compared with scrambled controls (see Supplementary data online, Figures S10G and H).
Discussion
In this study, bioinformatic deconvolution approaches were applied to obtain relative cell fractions in human carotid plaques, which showed that classical and modulated SMC fractions were significantly decreased in plaques from symptomatic patients, while higher Type 1 macrophages associated with a poorer survival after surgery. Clustering based on plaque cell fractions, revealed three distinct patient groups, with relative differences in their plaque stability profiles and associations to stroke, even during long-term follow-up. A novel strategy was also devised to obtain mesenchymal cell-specific genetic variants in symptomatic patients, which identified several variants in the regulatory region of the gene ARNTL as top hits. Upon further investigation of the role of ARNTL, we show that: (ⅰ) the risk allele of these mesenchymal cell-specific variants associates with decreased ARNTL gene expression in carotid plaques, (ⅱ) overall ARNTL expression is reduced in atherosclerotic plaques compared with normal arteries and negatively correlated with contractile SMC markers, and (ⅲ) silencing of ARNTL in primary human carotid SMCs leads to a decrease in cell proliferation (Figure 7).
This study revealed novel mesenchymal cell-specific, CAD-related plaque variants (rs4757138 shown here as an example) in the ARNTL gene locus. Presence of the rs4757138 minor allele (A/A) (right panel) associated with a significant decrease in ARNTL expression levels relative to that of the major G/G allele (left panel) in carotid plaques, presumably caused by mesenchymal cells. Indeed, silencing of ARNTL gene expression in primary human SMCs in vitro led to a decrease in proliferation
Deconvolution of different plaque datasets yielded similar predominance of SMC and macrophage cell fractions, validating the reproducibility of our bioinformatic approach (Structured Graphical Abstract).
Considering the importance of SMCs in conferring plaque stability,2 it is not surprising that plaques from symptomatic patients had reduced SMC fractions compared with those from asymptomatic patients. It may be suggested that higher fractions of Type 1 macrophages and pericytes coupled to a positive correlation found between immune and EC fractions exclusively in plaques from symptomatic patients, could indicate neovessel formation that also facilitates the homing of inflammatory macrophages into the plaque,50 altogether promoting plaque instability. Similarly, plaques from women had higher fraction of Type 1 macrophages and pericytes and lower fraction of classical SMCs compared with those of men, even when age stratification was considered, suggesting again a more vulnerable profile.
While it is well-known that men develop atherosclerosis earlier in life and have higher plaque burden, age stratifications have shown that this trend reverses with years and the disease in females worsens with old age.51 Since BiKE is a biobank of end-stage atherosclerosis, where the average age of patients is >65,28 the differences in plaque cell composition between men and women are not surprising. Given the impact of smoking and diabetes on the pathogenesis of atherosclerosis52,53 and the broad anti-inflammatory effects of LDL-lowering therapies such as statins,54 one might expect these factors to reflect on the plaque composition.52,53 Intriguingly, neither the presence of diabetes or smoking status, nor statin treatment, showed any association with the cellular composition of the plaques, although other studies have shown that they profoundly impact the gene expression levels and pathways ongoing in these lesions.28
Clustering of plaques according to their cellular composition, identified that Cluster 3 had particularly vulnerable features with higher immune and EC content and a corresponding increased incidence of ischemic stroke both at surgery and during follow-up compared with Clusters 1 and 2. Cluster 2 exhibited features of a stable plaque, including the higher presence of SMC fractions and lower immune cell fractions, however with a strong association to CKD. More research would be needed to explore this association, but the mechanisms for this could be attributed to kidney fibrosis with proliferation of glomerular mesangial cells and production of extracellular matrix. Although previous studies have performed clustering based on the bulk transcriptome of plaques,55 we have shown here that clustering at single cell resolution and coupling to extensive clinical and follow-up data in a large cohort, can provide further refinement and increase the precision of long-term risk prediction in atherosclerotic patients. To validate this, we have performed deconvolution and clustering of BiKE patients using scRNAseq data from both carotid and coronary plaques with similar results. Whether this bioinformatic approach has a practical translational value remains to be tested, but it is tempting to speculate that some molecular species, which characterize these clusters could be readily detected at systemic plasma level and therefore employed as biomarkers.
Recent studies have also applied single cell data from plaques combined with pseudotime tracking of cellular sub-phenotypes, genetic,22,56 and transcriptomic data55 to decipher the contribution of various cell types to plaque biology. Our integrated approach, encompassing many patients and combining deconvoluted plaque transcriptomic data with genetic and clinical data to prioritize CAD-GWAS loci, provided novel insights into the variants linked specifically with different mesenchymal cell types. Our analysis identified three variants that significantly associate with symptomatology of the patients. While one of these variants was located in a long non-coding RNA region (rs7642179), the other two were located in coding regions (rs7625059, rs12298484). The intronic variant rs7625059 in the gene Ventricular Zone Expressed PH Domain Containing 1 (VEPH1) had the highest correlation with patient symptoms. VEPH1 has recently been shown to modulate phenotypic switching of aortic SMCs in the context of abdominal aortic aneurysms,57 but its role in atherosclerosis has not been explored yet. Another intronic variant rs12298484 is located in the gene Dynein Axonemal Heavy Chain 10, the role of which is unknown in atherosclerosis.
BiKE analysis also identified many variants in genes that may be relevant for SMC function, most significantly the circadian rhythm gene ARNTL. We showed that minor alleles of several mesenchymal-specific variants in ARNTL associate with decreased gene ARNTL expression in bulk plaque transcriptomic data. Investigation of these variants in other independent cohorts showed their association with the presence of atherosclerotic plaques in high CVD risk individuals (IMPROVE cohort), as well as the expression of ARNTL in aortic adventitia from bicuspid valve patients (ASAP cohort). Notably, the trends of ARNTL expression levels samples from in BiKE and ASAP patients in association with the six variants were opposite, likely indicating differences in the biology of these diseases (atherosclerosis vs. aneurysm), location of the diseases (intima vs. media vs. adventitia), as well as pointing to the role of SNPs in regulating the expression levels of ARNTL in different vascular beds (carotid vs. aortic). Nevertheless, these results collectively validate and extend BiKE data and suggest a complex relationship between ARNTL genetics and vascular pathologies in general.
Also, overall expression levels of ARNTL were significantly decreased in plaques compared with normal arteries. However, it is difficult to conclude that this effect is due to ARNTL expression from mesenchymal cells only, although it is likely dominated by it. The reason is that scRNAseq data clearly shows broad ARNTL expression within plaques, with higher levels in immune (T cells and NK cells) and mesenchymal compartments, confirmed by the prevalent protein localisation in SMA+ regions of the plaque fibrous cap and remnants of the media. Of note, while a previous report showed ARNTL to be reduced in symptomatic compared with asymptomatic plaques,58 we could not confirm this in our cohort, possibly due to different criteria for categorising plaques, which in BiKE are driven by the patient classification based on clinical presentation and not plaque features.
ARNTL is important in regulation of the core circadian rhythm genes in mammals, deeply interconnected with atheroprogression, since it has been shown that disruption of circadian rhythm in humans and animal models results in atherosclerosis.59–61 In particular, a myeloid-specific Arntl deletion in an Apoe−/− mouse model increased the size of atherosclerotic lesions, characterized by an increase in lesional macrophages and skewing of M2 to M1 phenotype.62 In another study, ablation of Arntl specifically in SMCs in an atherosclerotic mouse model led to increased SMC migration, monocyte recruitment, production of reactive oxygen species, and SMC apoptosis, resulting in worsening of atherosclerosis.58 While these observations are generally in line with our findings, we have additionally uncovered a novel genetic basis for ARNTL to be potentially causal to the disease, particularly in SMCs. We also found that ARNTL levels were significantly reduced in patients belonging to Cluster 2 compared with Clusters 1 and 3. However, we were unable to determine the source of this difference among the clusters using the bulk transcriptome, given the broad expression profile of ARNTL in plaques. In agreement with ARNTL transcript correlations in plaques, the expression of contractile SMC markers among the clusters showed an opposite trend to ARNTL, while ARNTL could be positively associated with pathways relevant for modulated SMC function. Namely, adipogenesis and mTOR signalling are involved in regulation of SMC proliferation,63,64 while increased heme absorption in plaques due to leaky neovessels and intraplaque haemorrhage can activate SMCs through increased ROS.45,65 It is also known that Bmal1 controls cell proliferation through PDGF signalling,66 as well as senescence by regulating p16-INK4.49 In line with these reports, we found increased ARNTL expression in SMCs stimulated with PDGF-BB, while silencing of ARNTL led to an increase in gene expression levels of several contractile markers, a decrease in SMC proliferation, and a trend towards increased senescence. Although ARNTL silencing did not lead to significant changes in the corresponding protein levels of its downstream SMC targets, it is known that small changes in expression of transcription factors may result in cumulative functional consequences downstream, as here evidenced by a decrease in SMC proliferation. Overall, repression of ARNTL in SMCs, possibly linked with mesenchymal cell-specific variants in ARNTL gene identified by us, likely aggravates plaque vulnerability67 in symptomatic patients.
Limitations and advantages
The BiKE biobank contains exclusively end-stage atherosclerosis patients with a narrow age-span, which constricts our findings and limits extrapolation of conclusions into processes relevant for earlier stages of atheroprogression. Transcriptomic data from carotid plaques was deconvoluted with scRNAseq data from both coronary and carotid plaques with good agreement in cell fractions and clustering results, validating our findings. We investigated the function of ARNTL in SMCs by modulating gene expression levels, however, further studies are needed to investigate the causal functional effects of the identified SNPs in genome edited SMCs. Finally, it is worth noting that silencing of ARNTL gene expression in our experiments did not translate into a strong BMAL1 protein repression, possibly due to the fact that BMAL1 levels oscillate during the day and can be present concomitantly as phosphorylated and unphosphorylated isoforms in both nucleus and cytosol, based on its complex regulation in association with cellular circadian rhythm.
Conclusions
Taken together, this study shows the potential of combining scRNAseq data with vertically integrated per-patient clinical, genetic and transcriptomic data from a large biobank of human plaques, for refinement of patient vulnerability and risk prediction stratification. Numerous novel mesenchymal cell-specific genetic loci in atherosclerotic plaques were identified, few of them also associated with patient symptoms. ARNTL should be further mechanistically explored as a regulator of SMC proliferation vs. senescence, key processes that lead to profound detrimental changes in plaque morphology and aggravated vulnerability.
Acknowledgements
The Authors acknowledge all staff and surgeons from the Vascular Surgery Unit, Karolinska University Hospital for participating in the collection of patient material.
Supplementary data
Supplementary data are available at European Heart Journal online.
Declarations
Narayanan S et al. Atheroma transcriptomics identifies ARNTL as a smooth muscle cell regulator and with clinical and genetic data improves risk stratification
Disclosure of Interest
C.L.M. has received funding from AstraZeneca for an unrelated project. Bruna Gigante is board member of the European Society of Cardiology, Working Group on Thrombosis. All other authors declare no conflict of interest.
Data Availability
All data generated during this study is available in the paper or Supplementary Files. BiKE microarray dataset is available from Gene Expression Omnibus (GEO accession number GSE21545). Additional data, excluding those that are subject to GDPR regulations, can be provided by the Corresponding Author pending a reasonable request.
Funding
L.M. is the recipient of fellowships and awards from the Swedish Research Council (VR, 2023-02724, 2019-02027), Karolinska Institute Consolidator program, Swedish Heart-Lung Foundation (HLF, 20230357, 20210466, 20200621, 20200520, 20180244, 20180247, and 201602877), Swedish Society for Medical Research (SSMF, P13-0171). Ljubica Matic also acknowledges funding from Mats Kleberg’s, Sven and Ebba-Christina Hagberg’s, Tore Nilsson’s, Magnus Bergvall’s, and Karolinska Institute research (KI Fonder) and doctoral education (KID) foundations. Project funding was also obtained by Ulf Hedin from the Swedish Heart-Lung Foundation (20180036 and 20170584), the Swedish Research Council (2017-01070 and 2019-02027), and King Gustav Vth and Queen Victoria’s Foundation. C.L.M. acknowledges funding from the National Institutes of Health (R01HL148239 and R01HL164577); the Fondation Leducq ‘PlaqOmics’ (18CVD02); and the Single-Cell Data Insights award from the Chan Zuckerberg Initiative, LLC, and Silicon Valley Community Foundation. This work was also funded by a research grant from the European Union’s HORIZON-HLTH-2023-TOOL-05 program with Grant agreement No 101136962 (NextGen).
Ethical Approval
The collection and use of human samples in this study is approved by the regional Ethical Committees and follows the guidelines of the Declaration of Helsinki. All human samples were collected with informed consent from patients or organ donors’ guardians.
Pre-registered Clinical Trial Number
Not applicable.
