Single-cell ‘omic profiles of human aortic endothelial cells in vitro and human atherosclerotic lesions ex vivo reveal heterogeneity of endothelial subtype and response to activating perturbations

Heterogeneity in endothelial cell (EC) sub-phenotypes is becoming increasingly appreciated in atherosclerosis progression. Still, studies quantifying EC heterogeneity across whole transcriptomes and epigenomes in both in vitro and in vivo models are lacking. Multiomic profiling concurrently measuring transcriptomes and accessible chromatin in the same single cells was performed on six distinct primary cultures of human aortic ECs (HAECs) exposed to activating environments characteristic of the atherosclerotic microenvironment in vitro. Meta-analysis of single-cell transcriptomes across 17 human ex vivo arterial specimens was performed and two computational approaches quantitatively evaluated the similarity in molecular profiles between heterogeneous in vitro and ex vivo cell profiles. HAEC cultures were reproducibly populated by four major clusters with distinct pathway enrichment profiles and modest heterogeneous responses: EC1-angiogenic, EC2-proliferative, EC3-activated/mesenchymal-like, and EC4-mesenchymal. Quantitative comparisons between in vitro and ex vivo transcriptomes confirmed EC1 and EC2 as most canonically EC-like, and EC4 as most mesenchymal with minimal effects elicited by siERG and IL1B. Lastly, accessible chromatin regions unique to EC2 and EC4 were most enriched for coronary artery disease (CAD)-associated single-nucleotide polymorphisms from Genome Wide Association Studies (GWAS), suggesting that these cell phenotypes harbor CAD-modulating mechanisms. Primary EC cultures contain markedly heterogeneous cell subtypes defined by their molecular profiles. Surprisingly, the perturbations used here only modestly shifted cells between subpopulations, suggesting relatively stable molecular phenotypes in culture. Identifying consistently heterogeneous EC subpopulations between in vitro and ex vivo models should pave the way for improving in vitro systems while enabling the mechanisms governing heterogeneous cell state decisions.


Introduction
Endothelial cells (ECs) in the vascular endothelium maintain hemostasis, mediate vasodilation, and regulate the migration of leukocytes into tissues during inflammation.Dysfunctions of the endothelium are a hallmark of the aging process and are also an important feature of diseases, including atherosclerosis.Atherosclerosis is an inflammatory process fueled by cholesterol and leukocyte accumulation in the sub-endothelial layer of arteries.It is the underlying pathobiology of ischemic heart disease and the leading cause of morbidity and mortality worldwide due to heart attack and stroke (Brown et al., 2020;Hajra et al., 2000;Birdsey et al., 2015).Atherosclerosis of the coronary arteries is estimated to be about 50% genetic, with hundreds of genomic loci contributing to genetic risk (Marenberg et al., 1994;Aragam et al., 2022;Tcheandjieu et al., 2022).A major opportunity for better understanding the molecular basis for how disease progresses lies in identifying the genomic and downstream functions impaired by risk variants in disease-relevant cell types.Genetic studies are increasingly suggesting that a significant proportion of genetic risk for atherosclerosis is encoded in perturbed functions of vascular ECs (Aragam et al., 2022;Tcheandjieu et al., 2022;Kessler and Schunkert, 2021).
Single-cell sequencing technologies have begun to characterize the extent of EC molecular diversity in vitro and in vivo (Zhao et al., 2018;Li et al., 2019;Kalluri et al., 2019;Liu et al., 2021;Kalucka et al., 2020;Rohlenova et al., 2020;Zhao et al., 2021a;Xu et al., 2020;Cheng et al., 2021;Khan et al., 2019;Andueza et al., 2020;Tombor et al., 2020).Genetically engineered, lineage-traced mouse models have also been instrumental for identifying which cells in atherosclerotic plaques arose from EC origin.These studies have demonstrated that many cells of EC origin in plaques lack canonical EC marker genes and luminal location (Evrard et al., 2016;Chen et al., 2012).As many as one-third of mesenchymal-like cells in plaques have been reported to be of endothelial origin (Evrard et al., 2016), suggesting that phenotypic transition from endothelial to mesenchymal (EndMT) is a feature of atherosclerosis; however, whether EndMT is a cause or bystander of atherogenesis or plaque rupture is not fully understood.Although lineage tracing is not possible in humans, immunocytochemical techniques suggest that EC heterogeneity is prevalent in atherosclerotic vessels.These studies have described an unexpectedly large number of cells co-expressing pairs of endothelial and mesenchymal proteins, including fibroblast-activating protein/von Willebrand factor (FAP/VWF), fibroblast-specific protein-1/VWF (FSP-1/VWF), FAP/platelet-endothelial cell adhesion molecule-1 (CD31 or PECAM-1), FSP-1/CD31 (Evrard et al., 2016), phosphorylation of TGFB signaling intermediary SMAD2/FGF receptor 1 (p-SMAD2/FGFR1) (Chen et al., 2015), and α-smooth muscle actin (αSMA)/PECAM-1 (Moonen et al., 2015).An important implication of this result is that the use of canonical EC markers to isolate or identify ECs will likely omit certain EC populations.The extent of EC molecular and functional heterogeneity within a tissue during homeostasis and during disease is not well understood.One notable study exemplifying EC heterogeneity demonstrated that the EC-marker gene von Willebrand Factor (VWF) was expressed only in a subset of ECs from the same murine vessel, and the penetrance of VWF expression across ECs was tissue-specific (Yuan et al., 2016).In a related study, expression of the leukocyte adhesion molecule VCAM-1 was found to be upregulated by the proinflammatory cytokine tumor necrosis factor α only in some of the ECs of a monolayer (Turgeon et al., 2020).In both studies, variability in DNA methylation on CpG dinucleotides at the gene promoters negatively correlated with VWF and VCAM-1 expression.These findings raise the question as to how many molecular programs exist within ECs of a same tissue or culture, how this heterogeneity influences response to cellular perturbations, and what factors regulate these cellular states.
There are notable benefits and limitations for studying heterogeneity using in vitro and in vivo approaches in atherosclerosis research.In vitro approaches provide unique opportunities for interrogating consequences of genetic and chemical perturbations in highly controlled environments and are adept at identifying mechanistic relationships on accelerated timelines.In vivo approaches benefit from the complexity of the crosstalk among all cell types and tissues of the organism and are adept for identifying how perturbations manifest in living systems.It reasons that the integration of results from both approaches will best accelerate discovery.However, comprehensive analysis comparing heterogeneity of vascular ECs observed in vivo and in vitro remains unexplored.In the current study, we performed meta-analysis on four human in/ex vivo single-cell transcriptomic datasets (Pan et al., 2020;Alsaigh et al., 2022;Chowdhury et al., 2022;Wirka et al., 2019), containing 17 arterial samples, from mild-to-moderate calcified atherosclerotic plaques to evaluate the ability of the in vitro EC models to recapitulate molecular signatures observed in human atherosclerosis.
Human aortic endothelial cells (HAECs) are among the most appropriate cell type for in vitro modeling of the arterial endothelium in atherosclerosis research insofar as they are human cells, they are more readily available than coronary artery ECs, they are not of venous origin like human umbilical vein ECs, and they can be isolated from explants of healthy donor hearts during transplantation.We set forth in the current study to quantify heterogeneity among HAECs using multimodal sequencing that simultaneously measures transcripts using RNA-seq and accessible chromatin using ATAC-seq from the same barcoded nuclei.To provide estimates for heterogeneity due to genetic background, we molecularly phenotyped HAECs from six genetically distinct human donors.We also quantified single-cell responses to three perturbations known to be important in EC biology and atherosclerosis.The first was activation of transforming growth factor beta (TGFB) signaling, which is a hallmark of phenotypic transition and a regulator of EC heterogeneity (Evrard et al., 2016;van Meeteren and ten Dijke, 2012).The second was stimulation with the pro-inflammatory cytokine interleukin-1 beta (IL1B), which has been shown to model inflammation and EndMT in vitro (Bujak et al., 2008;Bujak and Frangogiannis, 2009;Maleszewska et al., 2013;Chaudhuri et al., 2007;Sánchez-Duffhues et al., 2019), and whose inhibition reduced adverse cardiovascular events in a large clinical trial (Ridker et al., 2017).The third perturbation utilized in our study was knockdown of the ETS-related gene (ERG), which encodes a transcription factor of critical importance for EC fate specification and homeostasis (Sperone et al., 2011;Fish et al., 2017;Lathen et al., 2014;Vijayaraj et al., 2012;Hogan et al., 2017).
Lastly, we examine whether epigenetic landscapes among heterogeneous EC subtypes observed in this study were differentially enriched for coronary artery disease (CAD) genetic risk variants.Taken together, this study provides evidence that EC heterogeneity is prevalent in vivo and in vitro and that not all ECs respond similarly to activating perturbations.

EC single-cell transcriptomic profiles reveal a heterogeneous population
To systematically uncover the heterogeneity of molecular landscapes in ECs at single-cell resolution, we cultured primary HAECs isolated from luminal digests of ascending aortas from six de-identified heart transplant donors at low passage (passages 3-6) (Navab et al., 1988;Figure 1A).Using the 10X Genomics multiome kit (Genomics x, 2022a), single-nucleus mRNA expression (snRNA-seq) and chromatin accessibility (snATAC-seq) data were collected simultaneously for a total of 15,220 nuclei after stringent quality control ('Materials and methods').RNA and ATAC data were integrated separately by treatment condition and then with each other as reported previously ('Materials and methods'; Hao et al., 2021).

EC subtypes exhibit distinct open chromatin profiles and enriched motifs
To investigate how different transcriptional signatures across ECs correspond to distinct chromatin states, we utilized the snATAC-seq portion of the multiome dataset.The snATAC-seq data were sequenced to a median depth of 22,939-126,122 reads with 3480-19,259 peaks called per nucleus (Supplementary file 1b and f).Of the 204,904 total identified peaks, 13,731 were differential across subtypes, with 79-8091 peaks uniquely accessible per EC subtype (Supplementary file 1h).Over 80% of total peaks were intergenic or intronic (Figure 2A and B) and most unique peaks were from EC2 and EC4.
The online version of this article includes the following figure supplement(s) for figure 1:   and 3) and lack of expression in male in vitro donor cells (2, 4, 5, and 6).TEAD3 (adjusted p-value 2.1 × 10 -306 ), and TEAD4 (adjusted p-value 6.9 × 10 -252 ) (Figure 2C and D).Notably, TEAD factors have been found as enriched in vascular smooth muscle cells (VSMCs) (Wirka et al., 2019;Örd et al., 2021), which is consistent with EC4 having the most mesenchymal phenotype of our EC subtypes.Taken together, these data demonstrate that EC1 and EC2 are the subtypes most canonically like 'healthy' or angiogenic ECs insofar as they exhibit ETS motif enrichments.Additionally, we conclude that EC4 is the most mesenchymal EC insofar as it exhibits TEAD factor enrichments.

EC-activating perturbations modestly shift cells into the EC3 subtype
Embedded in the dataset of this study were three experimental conditions known to promote EndMT along with their respective controls.Each experimental condition was administered to between three and five genetically distinct HAEC cultures.The conditions included 7-day exposure to IL1B (10 ng/ ml), 7-day exposure to TGFB2 (10 ng/ml), and 7-day siRNA-mediated knockdown of ERG (siERG).The control for IL1B and TGFB2 treatments was 7-day growth in matched media lacking cytokine and the control for the siERG condition was transfection with scrambled RNA.The UMAP presented in Figure 1 includes all the nuclei profiled across donors and conditions.We hypothesized that EC4, the most mesenchymal cluster, would be enriched for cells exposed to IL1B, TGFB2, and/or siERG relative to the controls thereby consistent with the hypothesis that the EC4 subtype were a consequence of EndMT.Detailed in Figure 3A and B are the relative proportions of cells from each experimental condition and donor by cluster.Contrary to our hypothesis, the EC4 cluster was not enriched for cells that were treated with cytokine or siERG relative to the controls; in fact, there is a nonstatistically significant trend for decreased numbers of EC4 cells from these conditions relative to controls insofar as all the donors with cells in EC4 show diminished proportions upon perturbation (Figure 3).The one cluster exhibiting increased proportions of cells upon perturbations was EC3, with three of four EC IL1B-exposed donors having increased proportions in EC3 (p=0.08 by two-sided paired t-test; Figure 3A), four of five TGFB2-exposed donors having increased proportions (p=0.04 by two-sided paired t-test; Figure 3A), and three of three donors having increased EC3 proportions upon ERG knockdown (Figure 3B).
In addition to heterogeneity across EC clusters, data in Figure 3 underscores that there is heterogeneity among EC cultures.To quantify this effect, we performed principal component analysis (PCA) to evaluate the overall contributions that donor and experimental conditions have on variance in this dataset.We found that pro-EndMT perturbations elicited greater variance in RNA expression (38-56% of variance) than donor (17%-27% variance) (Figure 3-figure supplement 1A-C), supporting that the transcriptional and epigenetic programs elicited by experimental conditions have a greater overall consequence than donor.This finding provides the opportunity to elucidate how different EC clusters respond to pro-EndMT exposures across genetically distinct ECs.

Pro-EndMT perturbations in vitro elicit EC subtype-specific transcriptional responses
We next sought to evaluate the similarities and differences among pro-EndMT perturbations and evaluate the transcriptional response elicited in each EC subtype.Differential gene expression analysis was performed using pseudo-bulked profiles grouped by donor, subcluster, and experimental groupings (Supplementary file 1i).
Overall, we found heterogeneity in transcriptional responses across EC subtypes.While EC1 and EC2 transcripts were predominantly perturbed by siERG, the greatest number of transcripts differentially expressed in EC3 were those responsive to IL1B, though siERG and TGFB2 also regulated tens to hundreds of transcripts in EC3.In contrast, transcripts in EC4 were predominantly responsive to TGFB2 (Figure 4A, Supplementary file 1i).With respect to EC4, we questioned whether transcripts were predominantly responsive to TGFB2 due to differences in expression of TGFB receptors.While we observed increased TGFBR1 expression in EC4, we observed relatively less expression of TGFBR2 and ACVRL1 in EC4 when compared to EC1, EC2, and EC3 (Figure 3-figure supplement 2A).We next questioned whether EC3 transcripts were predominantly responsive to IL1B due to differences in IL1B receptor expression.Notably, we did not observe differences in IL1B receptor expression, suggesting that their transcription is not responsible for divergent EC responses across EC subtypes (Figure 3-figure supplement 2B).Interestingly, we did observe differential expression of IL1RL1 in EC2, which may influence EC2 response to cytokine (Figure 3-figure supplement 2B).

Meta-analysis of ex vivo human atherosclerotic plaque snRNA-seq datasets
To understand the diversity of ECs in human atherosclerotic plaques and evaluate their relationships to our in vitro system, we performed a meta-analysis of data from recent publications that utilized scRNA-seq from human atherosclerotic lesions (Pan et al., 2020;Alsaigh et al., 2022;Chowdhury et al., 2022;Wirka et al., 2019; accessions in Supplementary file 1k).We identified 24 diverse clusters among 58,129 cells after integration of 17 different coronary and carotid samples (Figure 5A and Supplementary file 1l).Clusters were annotated using a combinatorial approach including canonical marker genes, CIPR (Ekiz et al., 2020), and the original publications (Figure 5B).Clusters were annotated as T-lymphocytes, natural killer T-cells, ECs, macrophages, VSMCs, fibroblasts, B-lymphocytes, basophils, neurons, and plasmacytoid dendritic cells (PDCs) (Figure 5A).We find the greatest proportion of cells belonging to each major cell type derive from carotid arteries, except for neurons that derive exclusively from coronary arteries, and PDCs that derive exclusively from carotid arteries (Figure 5-figure supplement 1B and C).Expected pathway enrichments are observed for annotated cell types, including NABA CORE MATRISOME (M5884; p-value 4.8 × 10 -41 ) for fibroblasts, blood vessel development (GO:0001568; p-value 5.6 × 10 -33 ) for ECs, and actin cytoskeleton organization (GO:0030036; p-value 1.3 × 10 -15 ) for VSMCs (Figure 5-figure supplement 1D-G).These observations support the diverse composition of human atherosclerotic lesions.

Ex vivo-derived module score analysis reveals differences among in vitro EC subtypes and EndMT stimuli
To directly evaluate the relationships between the ex vivo and in vitro cell subpopulations, we utilized module scores.These quantitative scores are based on the sum of ex vivo marker genes across each cluster and were used to evaluate similarity to each in vitro cell subcluster.Unexpectedly, the ex vivo cluster that consistently generated the greatest module scores for in vitro ECs is the VSMCs cluster 5 (VSMC5) (Figure 5A, Figure 5-figure supplement 3A).VSMC5 bridges the EC to SMC and fibroblast clusters in the ex vivo analysis (Figure 5A).Marker genes of VSMC5 are expressed across ex vivo and in vitro clusters (Figure 5-figure supplement 4A) and include important regulators of   5-figure supplement 3B-E).This supports that EC3 is a more activated subtype than EC1.Finally, among in vitro cells, those with highest VSMC5 module scores are concentrated in EC4, underscoring that EC4 is a more mesenchymal sub-phenotype in vitro (Figure 5-figure supplement 3B-E).

EC subtype is a major determinant in modeling cell states observed in atherosclerosis
In addition to module score analysis, we applied a complementary approach to quantitatively relate in vitro EC subtypes and pro-EndMT perturbations to ex vivo cell types.We calculate average expression profiles for all major cell populations in both ex vivo and in vitro datasets and examine the comprehensive pairwise relationship among populations with hierarchical clustering using Spearman Correlation (Figure 6A).All in vitro transcripts significantly regulated across all pro-EndMT perturbations at 5% false discovery rate (FDR) (Benjamini and Hochberg, 1995) are used in this analysis, although several additional means to select transcripts showed similar results.This analysis reveals three major observations.First, in vitro EC4 cells are most like mesenchymal ex vivo cell types, including VSMCs and fibroblasts (indicated by the yellow block of correlations in the bottom left of the heatmap in Figure 6A).Second, in vitro EC1, EC2, and EC3 are most like ex vivo Endo1 and Endo2 populations, especially among the siSCR and 7-day control cells.Moreover, cells in the siSCR condition in EC1 are most like ex vivo Endo1, reinforcing that these two populations are the most canonically 'healthy' endothelial populations.Third, while pro-EndMT perturbations did elicit variation in how similar in vitro ECs resembled ex vivo transcriptomic signatures, these effects are secondary to which subtype the cells belonged (Figure 6A).Taken together, these findings underscore that EC subtype, versus perturbation, is a greater determinant of similarity to ex vivo cell types.

CAD-associated genetic variants are enriched across EC subtype epigenomes
Genetic predisposition to CAD is approximately 50% heritable with hundreds to thousands of genetic loci supposed to be involved in shaping an individual's propensity for disease (Drobni et al., 2022; The online version of this article includes the following figure supplement(s) for figure 5:      .Endothelial cell (EC) subtype is a major determinant in the ability to recapitulate 'omic profiles seen in atherosclerosis.(A) Heatmap displaying average expression between in vitro perturbation-subtype combinations and ex vivo cell subtypes using all up-and downregulated genes between IL1B, TGFB2, or siERG versus respective controls.Spearman correlation was used as the distance metric.Rows (in vitro EC subtypes) and columns (ex vivo cell subtypes) are clustered using all significant genes (adjusted p-value<0.05)induced and attenuated across all in vitro EC subtypes for each

Figure 6 continued on next page
McPherson and Tybjaerg-Hansen, 2016).Most CAD-associated variants are not protein coding, suggesting that they perturb cellular function through gene regulatory functions.We therefore asked whether the open chromatin regions in this in vitro dataset coincided with locations of single-nucleotide polymorphisms (SNPs) reported in the latest CAD meta-GWAS analysis from the Millions Veterans Project, which includes datasets from CARDIoGRAMplusC4D 1000G study, UK Biobank CAD study, and Biobank Japan (Tcheandjieu et al., 2022).We found significant enrichment in CAD-associated SNPs for the complete set of accessible regions across all EC subtypes (termed 'panEC'; adjusted p-value 1.5e × 10 -93 ; odds ratio [OR] = 1.8; Figure 6B, Supplementary file 1o and p) when comparing CAD SNPs exceeding the genome-wide significance threshold of p<5 × 10 -8 versus nonsignificant SNPs ('Materials and methods).Among accessible regions unique to EC subtypes, EC4 shows the greatest enrichment (adjusted p-value 7.85 × 10 -6 ; OR = 1.74).Additionally, EC2 is also enriched for CAD SNPs (adjusted p-value 6.3 × 10 -8 ; OR = 2.15), supporting a role for proliferative ECs in CAD.Of all accessible regions influenced by pro-EndMT perturbations, siERG and TGFB2 sets are most enriched for CAD variants (Figure 6B, Supplementary file 1o and p).
The measurement of both gene expression and DNA accessibility in the same cell enables testing for direct correlation, or 'links', between accessibility of noncoding DNA elements and gene expression of their potential regulatory targets (i.e., gene promoters).This is achieved by testing for correlation between DNA accessibility and the expression of a nearby gene across single cells (Stuart et al., 2021;Cao et al., 2018).Focusing on EC4, we search for EC4-specific sites of correlated chromatin accessibility and linked target gene expression.Upon restricting linked peaks overlapping CAD SNPs, we identify 81 significant SNP-peak-gene trios (p<0.05)representing 46 unique genes with specific activity in EC4 (Supplementary file 1q).We submit the 46 unique genes to Metascape (Zhou et al., 2019) and observe enrichment in EndMT-related pathways including blood vessel development (GO:0001568; p-value 2.1 × 10 -10 ), crosslinking of collagen fibrils (R-HSA-2243919; p-value 1.4 × 10 -8 ), and canonical and non-canonical TGFB signaling (WP3874; p-value 2.2 × 10 -6 ) (Figure 6figure supplement 1).Literature review of this gene list further confirms several linked EC4-restricted genes associated with cardiovascular disease, including COL4A1, COL4A2, PECAM1, DSP, and BMP6 (Figure 6C-E; Liu et al., 2019;Yang et al., 2016;Woodfin et al., 2007).
Altogether, these data underscore that common genetic variation influences individual propensities for CAD through ECM-organizing functions evidenced by the EC4 phenotype.

Discussion
The major goals of this study were fourfold: (1) quantitatively assess molecular heterogeneity of cultured HAECs in vitro, (2) evaluate and compare molecular changes elicited by EC-activating perturbations at single-cell resolution, (3) assess similarities between in vitro and ex vivo EC signatures to inform the extent to which in vitro models recapitulate ex vivo biology, and (4) investigate how heterogeneous EC populations are enriched for genetic associations to CAD.The findings for each of these goals are discussed below, along with important implications and questions arising from this work.
The multiomic single-cell profiles of 15,220 cells cultured in vitro from six individuals enabled the discovery of five EC subpopulations, named EC1, EC2, EC3, EC4, and EC5.Except for EC5, EC subpopulations were comprised of cells from multiple donors and perturbations, which lends credence to the reproducibility of these biological states.The loosely defined phenotypes, based on pathway perturbation versus its respective control.(B) Heatmap of coronary artery disease (CAD)-associated single-nucleotide polymorphism (SNP) enrichments across in vitro EC subtypes and perturbation-subtype combinations.Rows (EC subtypes and perturbation-subtype combinations) are clustered using -Log 10 (P) for enrichment in significant CAD-associated SNPs (p-value<5 × 10 -8 ).Note that 'diff' represents peaks common to more than one EC subtype; it is found by subtracting EC1-5 subtype-specific peaks from the entire in vitro peak set (termed 'panEC').(C) Coverage plots displaying links for COL4A1/COL4A2 genes to EC4-specific peaks, including one overlapping with CAD-associated SNP rs9515203.(D) Coverage plot showing links for PECAM1 gene to EC4-specific peaks, including one overlapping with CAD-associated SNP rs1108591.(E) Coverage plot showing links for BMP6 gene to EC4-specific peaks, including one overlapping with CAD-associated SNP rs6597292.
The online version of this article includes the following figure supplement(s) for figure 6:  enrichment analysis, were healthy/angiogenic for EC1, proliferative for EC2, activated for EC3, and mesenchymal for EC4.Angiogenic (Li et al., 2019;Kalluri et al., 2019;Zhao et al., 2021a), proliferative (Tombor et al., 2020;Rodor et al., 2022), and mesenchymal (Tombor et al., 2020) ECs have been previously reported in the literature.The three activating perturbations (TGFB2, IL1B, siERG) had markedly unique effects on different EC subclusters, highlighting the fact that in vitro systems contain populations of discrete cell subtypes, or states, that respond divergently to even reductionistic experimental conditions.Implications of such heterogeneity include both a need to elucidate what factors dictate treatment responsiveness, as well as experimental design and data interpretation that considers heterogeneity of response.The exact origin of EC heterogeneity observed in this study is unknown.We consider it likely that EC1 EC2, EC3, and EC4 subpopulations, which were populated by most donors, date back to the original isolation of ECs from aortic trimmings, implying that different states were preserved across passage in the culture conditions.However, we cannot exclude the possibility that some of the subpopulations have expanded since seeding of the cultures.If that were the case, EC1, EC2, EC3, and EC4 represent reproducible cell states consequent to primary culture of arterial cells.In fact, the limited correlation with ex vivo data supports this interpretation.Future studies will be required to delineate the exact source of heterogeneity in these systems.
In this study, we set out to elucidate whether the mesenchymal phenotype of EC4 was an end-stage result of EndMT and whether TGFB2, IL1B, and/or siERG would increase the proportion of cells in EC4.As shown in Figure 3, this hypothesis was incorrect, and the only cluster with a modest increase in cell proportions upon stimulation was EC3.Moreover, while the percent of cells in EC3 increased with TGFB or IL1B, they decreased in EC4, suggesting trans-differentiation from EC4 into EC3 with these perturbations.We cannot exclude the possibility that EC3 is an EndMT cluster, although we would have expected more significant deviation from clusters EC1 and EC2.It is also possible that the postmortem state experienced by aortic explants prior to EC isolation could induce changes in the ECs, or that the duration and doses of perturbations chosen were not sufficient to elicit complete EndMT.While the duration and doses employed in our study were established based on literature reports reporting EndMT phenotypes (Maleszewska et al., 2013;Nagai et al., 2018;Medici et al., 2011), EndMT was quantified by expression of only a few marker genes rather than complete transcriptomic analysis.This raises an important conclusion of our study, which is that EndMT is not well-defined molecularly and it remains possible that several different molecular profiles may each represent variant flavors of EndMT.
We found that TGFB2, IL1B, and siERG have many distinct effects on EC molecular profiles (Figures 3 and 4).In general, TGFB2 elicits a greater transcriptomic and epigenomic response in the mesenchymal EC subtype, EC4, while siERG and IL1B regulate the greatest numbers of shared transcripts and chromatin regions in more endothelial clusters EC1, EC2, and EC3.One interpretation for this finding is that IL1B treatment and depletion of ERG directly affect rewiring transcription in ECs while TGFB2 may affect other cell types in the vascular wall (or culture plate) that in turn affect ECs through paracrine interactions.Part of the similarities between IL1B and siERG responses may be explained by the fact that ERG depletion increases IL1B production (Hogan et al., 2017).
A major question raised by this work is the origin of cells in the mesenchymal cluster EC4.We originally hypothesized this cluster was the result of EndMT, which led to our investigations as to whether we could leverage EndMT-promoting exposures (IL1B, TGFB2, siERG) in vitro observe an expansion of treated cells in the EC4 population.To our surprise, the EC4 population did not expand.If anything, these exposures reduced the proportion of cells in ECs (Figure 4).Nonetheless, it remains a possibility that EC4 represents cells that had undergone EndMT in vivo prior to culture and that the exposures we presented in vitro were not sufficient to elicit a complete EndMT transition.Another viable hypothesis is that cells in EC4 are of SMC origin and have persisted in culture alongside their EC counterparts.Cells used in this study were isolated by luminal collagenase digestion of explanted aortic segments and were tested at early passage for EC phenotypic markers including VWF expression, cobblestone morphology, and uptake of acetylated LDL.Notably, these rigorous metrics to ensure pure EC isolation occurred prior to our group's studies.In addition, if some of the isolated cells had undergone EndMT in vivo prior to isolation, it would be nearly impossible to distinguish their cell of origin after isolation since their collective molecular phenotypes would appear as an SMC.Without lineage tracing, which is currently not possible in human tissue explants, it would not be possible to distinguish cell origin.Nonetheless, this remains an important issue that is the subject of ongoing investigations.What we can confidently discern from these data is that these distinct cell subpopulations respond differently to the disease-relevant exposures of IL1B, TGFB2, and ERG depletion.
The current study sought to evaluate similarities and differences between in vitro primary cultures of HAECs to ex vivo single-cell signatures of cells from human lesions.First, we leveraged transcriptomic profiles from clusters in the scRNA meta-analysis of human lesions and evaluated each in vitro cluster using a module score (Figure 5, Figure 5-figure supplement 4).The three ex vivo clusters with greatest similarity to in vitro clusters were Endo1, Endo2, and VSMC5.Pathway enrichment analysis suggested that the ex vivo Endo1 cluster is close to the classic 'healthy' EC state relative to Endo2, which returned pathway enrichments consistent with activated endothelium (Figure 5C and  D).Interestingly, Endo2 is depleted in ribosome transcripts as well as transcripts in the Dicer complex (Figure 5C-E), which may serve as hallmarks of dysregulated endothelium in vivo.VSMC5 is an interesting ex vivo cluster insofar as it spans the endothelial, fibroblast, and VSMC clusters (Figure 5A) and is enriched for genes in actin cytoskeleton, extracellular matrix organization, and more (Figure 5figure supplement 4).In vitro EC1, EC2, and EC3 score generally greater in Endo1 and Endo2 relative to the more mesenchymal EC4 (Figure 5-figure supplement 3).Consistent with the intent of the pro-EndMT treatments, they generally decrease Endo1 and Endo2 scores and increase VSMC5 scores.However, these effects are unexceptional in comparison to effects of EC subtype.In addition to module scores, we also utilized unsupervised clustering of Spearman correlation coefficients across ex vivo and in vitro average gene expression profiles, finding again that EC1, EC2, and EC3 are more like Endo1 and Endo2 and EC4 is more like VSMCs (Figure 6A).As expected, the control (siSCR) cells are most correlated to healthy Endo1 transcriptomes; however, the correlation coefficient achieved is modest, at rho = 0.56.We cannot exclude the possibility that the moderate correlation coefficient observed between in vitro and ex vivo ECs may be explained by anatomic differences (i.e., aortic versus coronary and carotid arteries).While reinforcing that in vitro cell cultures best resemble ECs isolated ex vivo, regardless of perturbation, this finding accentuates how different cultured cells are and paves the way for quantitatively evaluating and improving in vitro models.
Finally, GWAS have established that hundreds of independent common genetic variants in human populations affect risk for CAD, yet discovering the causal mechanisms remains a major challenge given that most of the risk is in non-coding regions of the genome.One approach to prioritize causal variants in regulatory elements is through integration of open chromatin regions from the cell type and states of interest followed by expression quantitative trait loci (eQTL) or other linking evidence to target gene (Stolze et al., 2020;Toropainen et al., 2022).In the current study, we find significant enrichment for CAD-risk variants in open chromatin regions across the entire dataset ('panEC') as well as specifically for EC2 and EC4 subpopulations (Figure 6B, Supplementary file 1o-q).While EC3 was found to be more sensitive to perturbations in our in vitro experiments, we did not expect to see CAD-related SNPs enriched in EC3 because plasticity does not necessarily imply a pathological process.Moreover, while EC3 and EC4 both have mesenchymal phenotypes, EC3 may represent a reversible state that is lacking in EC4.This hypothesis would explain the enrichment of EC4, but not EC3, in CAD-related SNPs.
Taken together, these data emphasize the value in multimodal datasets in human samples for prioritizing disease-associated SNPs and mechanisms.

Tissue procurement and cell culture
Primary HAECs were isolated from eight deidentified deceased heart donor aortic trimmings (belonging to three females and five males of Admixed Americans, European, and East Asian ancestries) at the University of California Los Angeles Hospital as described previously (Navab et al., 1988;Supplementary file 1g).The only clinically relevant information collected for each donor was their genotype (see 'Genotyping and multiplexing cell barcodes for donor identification').HAECs were isolated from the luminal surface of the aortic trimmings using collagenase and identified by Navab et al. using their typical cobblestone morphology, presence of Factor VIII-related antigen, and uptake of acetylated LDL labeled with 1,1′-dioctadecyl-1-3,3,3′,3′-tetramethyl-indo-carbocyan-ine perchlorate (Di-acyetl-LD) (Navab et al., 1988).Cells were grown in culture in M-199 (Thermo Fisher Scientific, Waltham, MA, MT-10-060-CV) supplemented with 1.2% sodium pyruvate (Thermo Fisher Scientific, Cat# 11360070), 1% 100× Pen Strep Glutamine (Thermo Fisher Scientific, Cat# 10378016), 20% fetal bovine serum (GE Healthcare, Hyclone, Pittsburgh, PA), 1.6% Endothelial Cell Growth Serum (Corning, Corning, NY, Cat# 356006), 1.6% heparin, and 10 μl/50 ml Amphotericin B (Thermo Fisher Scientific, Cat# 15290018).Mycoplasma testing was not carried out on the cells used.HAECs at low passage (passages 3-6) were treated prior to harvest every 2 d for 7 d with either 10 ng/ml TGFB2 (Thermo Fisher Scientific, Cat# 302B2002CF), IL1B (Thermo Fisher Scientific, Cat# 201LB005CF), or no additional protein, or two doses of small interfering RNA for ERG locus (siERG; Supplementary file 1r), or randomized siRNA (siSCR; Supplementary file 1r).Donors 7 and 8 were treated prior to harvest for 6 hr with either 1 ng/ml IL1B, or no additional protein, and included in the dataset during integration to generate the original UMAP (Figure 1B), but not used for the purposes of downstream analyses in this study (Supplementary file 1g).All HAECs used were authenticated based on morphology, gene expression profiles indicative of ECs, and donor genotypes.No commercially available cell lines were analyzed in this study.
siRNA knockdown, qPCR, and western blotting Knockdown of ERG was performed as previously described (Hogan et al., 2017) using 1 nM siRNA oligonucleotides in OptiMEM (Thermo Fisher Scientific, Cat# 11058021) with Lipofectamine 2000 (Thermo Fisher Scientific, Cat# 11668030).Transfections were performed in serum-free media for 4 hr, then cells were grown in full growth media for 48 hr.All siRNAs and qPCR primers used in this study are listed in Supplementary file 1r.Transfection efficiency for the siRNAs utilized in this study was verified using qPCR 7 d after transfection (Figure 3-figure supplement 3A).Protein knockdown is shown 2 d after transfection using the same siRNAs from a representative experiment (Figure 3figure supplement 3B).Antibodies used included 1:1000 recombinant anti-ERG antibody (ab133264) and 1:5000 anti-histone H3 antibody (ab1791) (Abcam).Western blots were quantified using ImageJ (Schneider et al., 2012).

Genotyping and multiplexing cell barcodes for donor identification
Genotyping of HAEC donors was performed as described previously (Stolze et al., 2020).Briefly, IMPUTE2 (Howie et al., 2009) was used to impute genotypes utilizing all populations from the 1000 Genomes Project reference panel (phase 3) (1000Genomes Project Consortium et al., 2015).Genotypes were called for imputed SNPs with allelic R2 values greater than 0.9.Mapping between genomic coordinates was performed using liftOver (Kuhn et al., 2013).VCF files were subset by genotypes for the donors of interest using VCFtools (Danecek et al., 2011).
To identify donors across the in vitro dataset, snATAC-and snRNA-seq output BAM files from Cell Ranger ARC (10X Genomics, v.2.0.0;Genomics x, 2022a) were concatenated, sorted, and indexed using samtools (Danecek et al., 2021).The concatenated BAM files were input with the genotype VCF file to demuxlet (Kang et al., 2018) to identify best matched donors for each cell barcode, using options '-field GT'.Verification of accurate donor identification was confirmed by visualizing female sex-specific XIST for the known donor sexes (Figure 1-figure supplement 2).

snRNA-seq bioinformatics workflow
A target of 10,000 nuclei were loaded onto each lane.Libraries were sequenced on NovaSeq6000.Reads were aligned to the GRCh38 (hg38) reference genome and quantified using Cell Ranger ARC (10X Genomics, v.2.0.0;Genomics x, 2022a).Datasets were subsequently preprocessed for RNA individually with Seurat version 4.3.0(Hao et al., 2021).Seurat objects were created from each dataset, and cells with <500 counts were removed.This is a quality control step as it is thought that cells with low number of counts are poor data quality.Similarly, for each cell, the percentage of counts that come from mitochondrial genes was determined.Cells with >20% mitochondrial gene percent expression (which are thought to be of low quality, possibly due to membrane rupture) were excluded.Demuxlet (Kang et al., 2018) was next used to remove doublets.The filtered library was subset and merged by pro-EndMT perturbation.Data were normalized with NormalizeData, and cell cycle regression was performed by generating cell cycle phase scores for each cell using CellCycleScoring, followed by regression of these using ScaleData (Luecken and Theis, 2019).Batch effects by treatment were corrected using FindIntegrationAnchors using 10,000 anchors, followed by IntegrateData.

snATAC-seq bioinformatics workflow
A target of 10,000 nuclei were loaded onto each lane.Libraries were sequenced on an NovaSeq 6000 according to manufacturer's specifications at the University of Chicago.Reads were aligned to the GRCh38 (hg38) reference genome and quantified using Cell Ranger ARC (10X Genomics, v.2.0.0;Genomics x, 2022a).Datasets were subsequently preprocessed for ATAC individually with Seurat v4.3.0 (Hao et al., 2021) and Signac v1.6.0 (Heidecker et al., 2010) to remove low-quality nuclei (nucleosome signal >2, transcription start site enrichment <1, ATAC count <500, and % mitochondrial genes >20) (Hao et al., 2021).Next, demuxlet (Kang et al., 2018) was used to remove doublets.A common peak set was quantified across snATAC-seq libraries using FeatureMatrix, prior to merging each lane.A series of two iterative peak calling steps were performed.The first step consisted of calling peaks for every EndMT perturbation, and the second involved calling peaks for every cluster generated from weighted nearest-neighbor analysis (WNN) (see 'Integration and weighted nearestneighbor analyses').Latent semantic indexing (LSI) was computed after each iterative peak calling step using Signac standard workflow (Stuart et al., 2021).Batch effects by treatment were finally corrected using FindIntegrationAnchors using 10,000 anchors, followed by IntegrateData.

Integration and weighted nearest-neighbor analyses
Following snRNA-seq and snATAC-seq quality control filtering, barcodes for each modality were matched, and both datasets were combined by adding the snATAC-seq assay and integrated LSI to the snRNA-seq assay.WNN (Hao et al., 2021) was next calculated on the combined dataset, followed by joint UMAP ( WNN UMAP) visualization using Signac (Stuart et al., 2021) functions FindMul-timodalNeighbors, RunUMAP, and FindClusters, respectively.WNN is an unsupervised framework to learn the relative utility of each data type in each cell, enabling an integrative analysis of multimodal datasets.This process involves learning cell-specific modality 'weights' and constructing a WNN UMAP that integrates the modalities.The subtypes discovered in the first round of WNN were utilized in an additional peak calling step for snATAC-seq, followed by LSI computation, re-integration, and a final round of WNN to achieve optimal peak predictions (see 'Single-Nucleus ATAC sequencing bioinformatics workflow') (Yan et al., 2020).

Differential expression and accessibility region analyses across EC subtypes and EndMT perturbation-subtype combinations
Differential expression between clusters was computed by constructing a logistic regression (LR) model predicting group membership based on the expression of a given gene in the set of cells being compared.The LR model included pro-EndMT perturbation as a latent variable and was compared to a null model using a likelihood ratio test.This was performed using Seurat FindMarkers, with ' test.use= LR' and ' latent.vars' set to perturbation.Differential expression between perturbation and control for each cluster was performed using pseudobulk method with DESeq2 (Love et al., 2014).
Raw RNA counts were extracted for each EndMT perturbation-subtype combination and counts, and metadata were aggregated to the sample level.
Differential accessibility between EC subtypes was performed using FindMarkers, with ' test.use= LR' and latent.vars set to both the number of reads in peaks and perturbation.Finally, differential accessibility between perturbation and control for each cluster was performed using FindMarkers, with ' test.use= LR' and latent.vars set to the number of reads in peaks.
Bonferroni-adjusted p-values were used to determine significance at adjusted p-value<0.05for differential expression, and p-value<0.005for differential accessibility (Benjamini and Hochberg, 1995).

Motif enrichment analysis
A hypergeometric test was used to test for overrepresentation of each DNA motif in the set of differentially accessible peaks compared to a background set of peaks.We tested motifs present in the Jaspar database (2020 release) (Fornes et al., 2020) by first identifying which peaks contained each motif using motifmatchr R package (https://bioconductor.org/packages/motifmatchr).We computed the GC content (percentage of G and C nucleotides) for each differentially accessible peak and sampled a background set of 40,000 peaks matched for GC content (Stuart et al., 2021).Per-cell motif activity scores were computed by running chromVAR (Schep et al., 2017), and visualized using Seurat (Hao et al., 2021) function FeaturePlot.
Human atherosclerosis scRNA-seq public data download, mapping, and integration across samples Count matrices of 17 samples taken from four different published scRNA-seq datasets were downloaded from the NCBI Gene Expression Omnibus (accessions listed in Supplementary file 1k), processed using Cell Ranger (10X Genomics Cell Ranger 6.0.0;Zheng et al., 2017) with reference GRCh38 (version refdata-gex-GRCh38-2020-A, 10X Genomics), and analyzed using Seurat version 4.3.0(Hao et al., 2021).Seurat objects were created from each dataset, and cells with <500 counts and >20% mitochondrial gene percent expression were excluded.Additionally, doublets were removed using DoubletFinder (McGinnis et al., 2019), which predicts doublets according to each real cell's proximity in gene expression space to artificial doublets created by averaging the transcriptional profile of randomly chosen cell pairs.Next, normalization and variance stabilization, followed by PC analysis for 30 PCs, were performed in Seurat (Hao et al., 2021) using default parameters.Batch effects across the 17 samples were corrected using Seurat functions (Hao et al., 2021) FindIn-tegrationAnchors using 10,000 anchors, followed by IntegrateData.During the integration step, cell cycle regression was performed by assigning cell cycle scores with Seurat (Hao et al., 2021) function CellCycleScoring.The ex vivo dataset was first visualized, and canonical markers were identified for annotating cell types using FindAllMarkers.

Module scoring
FindAllMarkers was used to identify the top DEGs between each ex vivo cell subtype.Cells from the in vitro dataset were assigned an ex vivo cell subtype module score using Seurat (Hao et al., 2021) function AddModuleScore.The difference in module score between each in vitro EC subtype was established using Wilcoxon rank sum test with continuity correction and a two-sided alternative hypothesis.
Comparison of ex vivo snRNA-seq data to in vitro snRNA-seq data Meta-analyzed ex vivo human scRNA-seq data and in vitro snRNA-seq data were compared.Gene expression values for each ex vivo cell subtype and in vitro EC subtype-perturbation were produced using the AverageExpression function in Seurat (Hao et al., 2021) (which exponentiates log data, therefore output is depth normalized in non-log space).Figure 6A was generated using hclust function in R (Murtagh and Legendre, 2014).Spearman correlation was used as the distance metric.Sample clustering was performed using all significant genes (adjusted p-value <0.05) induced and attenuated across all in vitro EC subtypes for each pro-EndMT perturbation versus its respective control.Figure 5-figure supplement 4A was made using average expression data for marker genes for each ex vivo cell subtype.Hierarchical clustering across ex vivo cell subtypes was performed using hclust function in R (Murtagh and Legendre, 2014) using average expression as the distance metric for a given gene.

Peak-to-gene linkage
We estimated a linkage score for each peak-gene pair using the LinksPeaks function in Signac (Stuart et al., 2021).For each gene, we computed the Pearson correlation coefficient r between the gene expression and the accessibility of each peak within 500 kb of the gene TSS.For each peak, we then computed a background set of expected correlation coefficients given properties of the peak by randomly sampling 200 peaks located on a different chromosome to the gene, matched for GC content, accessibility, and sequence length (MatchRegionStats function in Signac).We then computed the Pearson correlation between the expression of the gene and the set of background peaks.A z score was computed for each peak as z = (r − μ)/σ, where μ is the background mean correlation coefficient and σ is the SD of the background correlation coefficients for the peak.We computed a p-value for each peak using a one-sided z-test and retained peak-gene links with a p-value<0.05and a Pearson correlation coefficient.The results were restricted to peak regions that overlapped with significant CAD-associated SNPs (see 'GWAS SNP enrichment analysis').
each gene is tested for differential expression between cells of each cluster and the cells outside of that cluster.p_val corresponds to the p-value of the Wilcoxon rank sum test, avg_log2FC is the log2-fold change in expression of the marker for the average cell in the cell type, cell type is the collapsed cluster of cells for which the marker was discovered, gene is the markers, pct.1 is the percentage of cells expressing the marker within the tested cluster, pct.2 is the percentage of cells expressing the marker outside the cluster, and p_val_adj is the Bonferroni adjusted p-value (based on the number of genes tested).(m) Effect of EndMT perturbation on in vitro EC subtypes according to ex vivo cell type module scores: Adjusted p-values<0.05 are generated using Wilcoxon rank sum test with continuity correction setting alternative hypothesis to ' two.sided'.Colored arrows represent significantly (adjusted p-value<0.05)upregulated (green) and downregulated (red) module scores for each EC sub-phenotype and perturbation combination.(n) Results table for effect of EndMT perturbation on in vitro EC subtypes according to ex vivo cell type module scores: Adjusted p-values<0.05 are generated using Wilcoxon rank sum test with continuity correction setting alternative hypothesis to ' two.sided'.(o) SNP enrichment analysis results across clusters and perturbation-cluster combinations: briefly, SNP enrichment analysis was used with traseR with ' test.method' set to 'fisher' and 'alternative' set to 'greater'.(p) Significantly enriched SNPs (adjusted p-value<0.05)across clusters and perturbation-cluster combinations: briefly, SNP enrichment analysis was used with traseR with ' test.method' set to 'fisher' and 'alternative' set to 'greater'.(q) Significant links between genes and EC4-specific peaks which overlap with CAD-associated SNPs (p-value<0.05):Briefly, peak-togene linkage was performed using Signac 'LinkPeaks' function on a dataset consisting of EC4-specific peaks and genes.Results were filtered based on peak regions which overlap with significant (p<5e-8) CAD associated SNPs.(r) siRNAs and qPCR primers used in this study.siERG #1, #2, #4, and #5 were pooled together.Non-targeting siRNA (siSCR) #3 and #4 were pooled together.

Data availability
Sequencing data have been deposited in GEO under accession code GSE228428.This project utilized data deposited previously in GEO accessions GSE155512, GSE159677, and GSE131778.The code used for analysis can be found in GitHub at https://github.com/cromanoski/Adelus_2024_Elife/(copy archived at Romanoski, 2024).
The following dataset was generated: The following previously published datasets were used:

Figure 1 .
Figure 1.Human aortic endothelial cell (HAEC) transcriptomic profiling discover a heterogeneous cell population.(A) Schematic diagram of the experimental design.Endothelial cells (ECs) were isolated from six human heart transplant donor's ascending aortic trimmings and treated with IL1B, TGFB2, or siERG (ERG siRNA) for 7 d.(B) Weighted nearest-neighbor UMAP ( WNN UMAP) of aggregate cells from all perturbations and donors is shown.Each dot represents a cell, and the proximity between each cell approximates their similarity of both transcriptional and epigenetic profiles.Colors Figure 1 continued on next page

Figure supplement 1 .
Figure supplement 1. Marker genes and pathways for endothelial clusters.

Figure supplement 2 .
Figure supplement 2. Violin plot of XIST showing expected expression in female in vitro donor cells (1and 3) and lack of expression in male in vitro donor cells(2, 4, 5, and 6).

Figure 2 .
Figure 2. Endothelial cells (ECs) have epigenetically distinct cell states.(A) Upset plot of differential peaks across EC subtypes.Intersection size represents the number of genes at each intersection, while set size represents the number of genes for each EC subtype.(B) Genomic annotation for the complete peak set.(C) Heatmap of top transcription factors (TFs) from motif enrichment analysis for marker peaks in each EC subtype.Top TFs for each EC subtype are selected based on ascending p-value.Rows (TFs) and columns (EC subtype) are clustered based on enrichment score (ES).(D) Feature plots and position weight matrices (PWMs) for top TF binding motifs for EC1 (TCF12), EC2 (ETV1), EC3 (GATA5), and EC4 (TEAD3).Per-cell motif activity scores are computed with chromVAR, and motif activities per cell are visualized using the Signac function FeaturePlot.(E), PWMs comparing Jaspar 2020 ETV1 motif to ERG motif reported in Hogan et al.

Figure 3 .
Figure 3. Endothelial cell (EC)-activating perturbations modestly shift cells into the EC3 subtype.(A) The proportion of cells in 7-day control and 7-day IL1B treatment are shown per human aortic endothelial cell (HAEC) donor and cluster on the top and for 7-day control and 7-day TGFB2 on the bottom.(B) The proportion of cells in 7-day siSCR control and 7-day siERG knockdown are shown per HAEC donor and cluster.EC1 was omitted in (A) due to lack of cells in both conditions.The online version of this article includes the following source data and figure supplement(s) for figure 3: Source data 1.Annotated western blots for Figure 3-figure supplement 3B where the leftmost six wells are shown from left to right as (1) the protein ladder (labeled in kD), (2) the lipofectamine transfected control, (3) the scrambled siRNA control, (4, 5) two lanes using different siRNAs against the TCF4 gene (not relevant to these studies), and lastly, (6) siRNA against ERG.Source data 2. Original western blots for ERG, annotated and uncropped.Source data 3. Original western blots for H3, annotated and uncropped.Source data 4. Original western blots for ERG, unannotated and uncropped.Source data 5. Original western blots for H3, unannotated and uncropped.

Figure supplement 1 .
Figure supplement 1. Principal component analysis using RNA-seq data.

Figure supplement 2 .
Figure supplement 2. Receptor expression profiles across endothelial clusters.

Figure 4 .
Figure 4. Endothelial cell (EC)-activating perturbations in vitro elicit EC subtype-specific transcriptional responses.(A) Upset plots of up-and downregulated differentially expressed genes (DEGs) across EC subtypes with siERG (gray), IL1B (pink), and TGFB2 (blue).Upset plots visualize intersections between sets in a matrix, where the columns of the matrix correspond to the sets, and the rows correspond to the intersections.Intersection size represents the number of genes at each intersection.(B) Pathway enrichment analysis (PEA) for EC3-4 up-and downregulated DEGs with TGFB2 compared to control media.(C) PEA for EC2-4 up-and downregulated DEGs with IL1B compared to control media.(D) PEA for EC1-4 up-and downregulated DEGs with siERG compared to siSCR.(E) PEA comparing up-and downregulated DEGs that are mutually exclusive and shared between IL1B and siERG in EC3.The online version of this article includes the following figure supplement(s) for figure 4: Figure supplement 1. Profiles of shared peaks and motif enrichments across endothelial subtypes.

Figure 5 .
Figure 5. Endothelial cells (ECs) from ex vivo human atherosclerotic plaques show two major populations.(A) scRNA-seq UMAP of different cell subtypes across 17 samples of ex vivo human atherosclerotic plaques.(B) Dot plot of top markers for each cell type.(C) Heatmap of pathway enrichment analysis (PEA) results generated from submitting 200 differentially expressed genes (DEGs) between endothelial cells 1 (Endo1) and endothelial cells 2 (Endo2).Rows (pathways) and columns (cell subtypes) are clustered based on -Log 10 (P).(E) Heatmap displaying expression of genes belonging to ribosome cytoplasmic pathway for Endo1 and Endo2.

Figure 5
Figure 5 continued on next page

Figure supplement 1 .
Figure supplement 1. Characterization of RNA-seq profiles from ex vivo arterial samples.

Figure supplement 3 .
Figure supplement 3. Cross comparisons between in vitro and ex vivo RNA-based module scores.

Figure supplement 4 .
Figure supplement 4. Heatmap and pathway analysis for marker genes of the VSMC5 ex vivo cluster.

Figure supplement 5 .
Figure supplement 5. Breakdown of ex vivo module scores across in vitro clusters and sample identity.

Figure 5 continued
Figure 5 continued
Figure 6 continued