Single-cell profiling reveals heterogeneity and functional patterning of GPCR expression in the vascular system

G-protein-coupled receptor (GPCR) expression is extensively studied in bulk cDNA, but heterogeneity and functional patterning of GPCR expression in individual vascular cells is poorly understood. Here, we perform a microfluidic-based single-cell GPCR expression analysis in primary smooth muscle cells (SMC) and endothelial cells (EC). GPCR expression is highly heterogeneous in all cell types, which is confirmed in reporter mice, on the protein level and in human cells. Inflammatory activation in murine models of sepsis or atherosclerosis results in characteristic changes in the GPCR repertoire, and we identify functionally relevant subgroups of cells that are characterized by specific GPCR patterns. We further show that dedifferentiating SMC upregulate GPCRs such as Gpr39, Gprc5b, Gprc5c or Gpr124, and that selective targeting of Gprc5b modulates their differentiation state. Taken together, single-cell profiling identifies receptors expressed on pathologically relevant subpopulations and provides a basis for the development of new therapeutic strategies in vascular diseases.

G -protein-coupled receptors (GPCRs) are the largest family of transmembrane receptors in eukaryotes; they transduce signals of numerous physicochemical stimuli including neurotransmitters, hormones, local mediators, metabolic or olfactory cues, and light 1 . The human genome encodes B800 GPCRs, the majority of them being olfactory receptors. Among the 367 non-olfactory GPCRs are still B150 orphan receptors, that is, GPCRs for which ligand and function are still unknown [2][3][4] . GPCRs have regulatory functions in almost all organ systems, and dysregulation of GPCR signalling has been implicated in the pathogenesis of many diseases [5][6][7][8][9] . The unique combination of diversity and specificity within the GPCR family, together with the fact that they are readily targetable by exogenous agonists and antagonists, has made GPCRs a most successful group of drug targets 4 . In the vascular system, GPCRs modulate critical parameters such as vessel tone or endothelial permeability [10][11][12] ; pharmacological modulation of receptors for angiotensin II, catecholamines or histamine are therefore crucial in the therapy of arterial hypertension or allergic fluid extravasation, respectively. However, although the market share of GPCR-targeting drugs is with 30-40% high, the overall number of targeted GPCR is with B30 rather low, suggesting that the potential of GPCRs as drug targets is not yet fully exploited 4 .
The search for new GPCR-based therapies encompasses various strategies, among them the identification of new biased or allosteric ligands at GPCRs with known function and ligand or the deorphanization of GPCRs for which ligands are not yet known 13,14 . A third approach is the identification of new pathophysiologically relevant functions for specific GPCRs, in particular for orphan receptors. One useful strategy in the identification of new GPCR functions is the ever more detailed expression analysis in functionally relevant cell subpopulations, for example, in subgroups of vascular cell types. While GPCR expression has been studied in bulk cDNA of whole vessels or cultured vascular cells 15,16 , our knowledge about GPCR expression in freshly isolated vascular cell types, in particular on the single-cell level, is insufficient. This is a major limitation, since current interpretations of expression data rely on the assumption that all cells of a given population are equal, or at least comparable, with respect to their GPCR repertoire. In contrast to this, studies in other fields indicate that gene expression is rather heterogeneous in vivo 17,18 , and it is crucial for further development of GPCR-based therapeutic strategies to understand whether this is also true for GPCR expression, and if so, whether certain GPCRs are associated with functional subgroups within a given cell population 19 . Activated EC, for example, upregulate adhesion molecules and proinflammatory mediators and are crucial for the pathogenesis of inflammatory vascular disease 20 ; whether this inflamed subpopulation of ECs is endowed with a specific GPCR repertoire, is unclear. Also vascular SMC undergo phenotypical changes in response to noxious stimuli, resulting in a dedifferentiated state characterized by downregulation of contractile proteins and upregulation of pro-proliferative, promigratory and proinflammatory markers as well as enhanced matrix synthesis 21,22 . Understanding in what way the GPCR repertoire differs between healthy and diseased vascular cells will be crucial to tailor new therapeutic strategies selectively targeting the latter.
To address this issue, we used here a microfluidic-based system for single-cell reverse transcription polymerase chain reaction (RT-PCR) to determine GPCR expression in individual freshly isolated EC and SMC. We found GPCR expression to be highly heterogeneous in all cell types and describe characteristic changes in single-cell GPCR repertoires in response to inflammatory activation. Furthermore, we identify functionally relevant subgroups of cells that are characterized by specific GPCR patterns and show examples of how this information can be used to modulate cell differentiation.

Results
Array design and quality control. To narrow down the spectrum of GPCRs to be tested on the single-cell level, we determined GPCR expression in bulk RNA of vascular EC and SMC. Of 424 tested GPCRs, 122 GPCRs were expressed in at least one of these cell types ( Supplementary Fig. 1A). For these 122 GPCRs, as well as for 32 additional GPCRs, a primer array for single-cell expression analysis was designed (Table 1). In addition, 13 genes Table 1 | Overview over the genes included in the array.
identifying individual cell types and 36 function-defining genes, including the three reference genes b-actin (Actb), glyceraldehyde 3-phosphate dehydrogenase (Gapdh) and hypoxanthine phosphoribosyltransferase (Hprt), were included in the array (Table 1). Whenever possible intron-spanning primer design was used, but for 22 single exon GPCRs this was not possible. These 22 non-intron-spanning primer pairs gave positive results in all tested cell types ( Supplementary Fig. 1B,C), indicative of contamination with genomic DNA. GPCRs not allowing intronspanning primer design were therefore excluded from the analysis (strikethrough in Table 1). Individual EC or SMC were pre-sorted by flow cytometry from enzymatically digested tissues of mice expressing enhanced green fluorescent protein (EGFP) in EC or SMC and subjected to single-cell isolation, cDNA synthesis and RT-PCR amplification using a microfluidic-based system. Various vascular beds were analysed, including aorta (ao), skeletal muscle vasculature (sk), lung (lu) or brain (br). All primer pairs included in the array produced amplification products with the predicted melting behaviour in at least one of three cell types tested (EC, SMC and leukocytes) (Supplementary Table 3). Despite stringent sorting criteria, single-cell analysis of individual EC or SMC showed 10-18% contamination with other cell types ( Table 2). Only cells with clear identity ('Cdh5 only', 'Myh11 only') and positive expression of reference genes Gapdh and Hprt as quality control were included in further analyses (Fig. 1a). A systematic comparison of expression data obtained by bulk RNA analysis or single-cell RT-PCR in SMC from skeletal muscle vasculature (SMsk) showed that of 74 GPCRs undetectable in bulk RNA, 26 showed expression in individual SMsk cells (Fig. 1b, left). Of 11 GPCRs with uncertain result in bulk cDNA, all showed amplification in individual cells (Fig. 1b, middle). Of 37 GPCRs clearly expressed in bulk RNA, two were completely absent in individual cell analysis (Fig. 1b, right). These two GPCRs, Eltd1 and Gpr116, are known to be highly specific for EC 23 , indicating that these discrepancies are due to contamination of SMsk with EC.
Reference genes Actb and Gapdh as well as SMC marker gene Myh11 were homogenously expressed (Fig. 2a, top). On average, individual cells expressed 20.3 ± 0.9 of the tested 132 GPCRs, though the individual values varied between 3 and 38 GPCRs per cell (Fig. 2b). mRNA sequencing of individual SMao showed even lower frequencies of GPCR expression (selected data in Supplementary Fig. 2, whole data set in Supplementary Data 1), which led us to verify single-cell expression data in GPCR reporter mice that are genetically engineered to express b-galactosidase (bgal) under control of different GPCR promoters. Flow cytometric analysis of bgal expression in freshly isolated SMao from Mrgprf-, Gabbr1-or Gprc5b-reporter mice closely matched the results of the single-cell RT-PCR expression analysis, while mRNA sequencing data underestimated GPCR expression (Fig. 2c,d). We also sequenced single-cell RT-PCR amplicons to exclude off-target amplification or amplification of highly homologous GPCRs and found that the amplified sequences were in all cases specific for the targeted receptor (Supplementary Fig. 3 and Supplementary Data 2). To investigate whether the same degree of heterogeneity was present on the protein level, we analysed expression of select GPCRs in individual SMao by flow cytometry. We found that also on the protein level GPCR expression was heterogeneous, and that the percentages roughly matched the values obtained by single-cell RT-PCR ( Supplementary Fig. 4A,B). Ex vivo culture of primary SMao (passage 1) resulted not only in an upregulation of genes indicative of a dedifferentiated SMC phenotype, for example, intercellular adhesion molecule 1 (Icam1), vascular cell adhesion molecule 1 (Vcam1), marker of proliferation Ki-67 (Mki67) or interleukin-6 (Il6), but also in strongly increased expression frequency for the majority of GPCRs (Fig. 2e). Some GPCRs also showed decreased expression, such as Lgr6, Npy1r, Crhr2 or Bdkrb2 (Fig. 2e), but the average number of GPCRs per cell was still strongly increased (Fig. 2f). GPCR expression frequencies in cultured human aortal SMC (passage 1) matched in most cases the murine data ( Fig. 2g middle), though some GPCRs showed notable difference between the species or between individual human donors ( Fig. 2g right).
GPCR repertoire in different types of SMC. We next investigated whether the strong heterogeneity of GPCR expression was a special feature of conductance arteries such as the aorta, or whether it was also present in resistance arteries, for example, from SMsk. Also in SMsk GPCR expression was very heterogeneous (Fig. 3a), although the overall number of GPCRs per cell was higher than in SMao (Fig. 3a,b). The repertoire of GPCRs expressed in the two SMC types differed strongly: GPCRs Lphn2, Cmklr1, Lpar1, Gpr133, P2ry6, Lgr6 and F2r were mainly present in SMao, while a large number of other GPCRs were predominantly expressed in SMsk (Fig. 3a). Among those receptors preferentially or selectively expressed in SMsk, the largest group were peptide hormone receptors (Fig. 3c), for example, receptors for endothelin (Ednra, Ednrb), angiotensin II (Agtr1a), vasopressin (Avpr1a), neuropeptide Y (Npy1r), pituitary adenylate cyclase-activating polypeptide (PACAP, acting on Adyap1r1), vasoactive intestinal polypeptide (VIP, acting on Vipr1 and Vipr2), corticotropin releasing hormone (Crhr2), parathyroid hormone (Pth1r) or calcitonin gene-related peptide (CGRP, acting on Calcrl) (Fig. 3d). Also a number of chemokine and orphan receptors were preferentially expressed in SMsk (Fig. 3e). We also analysed GPCR expression in other SMC types, for example, from the mesenteric vasculature (SMmes) or urinary bladder (SMub). K-means cluster analysis revealed that SMC from mesenteric or skeletal muscle vasculature, two regions rich in resistance arteries/arterioles, were quite similar in their GPCR profile, but differed strongly from SMao or SMub (Fig. 3f,g). An analysis of differential gene expression for all SMC types is shown in Supplementary Fig. 5A-C; a comparison of GPCR frequencies for all SMC types can be found in Supplementary Fig. 5D. An overview of the total number of cells, mice and independent experiments is shown in Supplementary Table 4.
Endothelial GPCR pattern after acute inflammatory activation. We next investigated how the GPCR repertoire changes in individual EC in response to acute inflammatory activation. To do so, we used a murine sepsis model induced by intraperitoneal (i.p.) injection of bacterial endotoxin lipopolysaccharide (LPS), which results in massive direct and indirect activation of EC 24 . ECbr showed 12 h after LPS injection a clear upregulation of markers of inflammatory activation such as Icam1 or E-selectin (Sele), whereas expression of vascular endothelial growth factor receptor 2 (VEGFR2, encoded by Kdr) or endothelin-1 (Edn1) was reduced (Fig. 5a). Also the GPCR pattern changed significantly: numerous orphan receptors (Gprc5a, Gpr97, Gpr111, Gpr153, Lphn2, Gpr107, Gprc5b and Gpr56) were upregulated, as well as select chemokine receptors (Darc, Ccrl2 and Cxcr7) and purinergic receptors (P2ry2, P2ry6 and Adora2a).
To also assess the endothelial GPCR repertoire under conditions of chronic vascular inflammation, we investigated aortal EC from ApoE-deficient mice kept for 16 weeks on a Western-style high-fat diet 25 . Surprisingly, all Cdh5-positive EC obtained from aortae were also positive for SM marker Myh11, although Myh11 levels were very low compared to normal SMC  ( Fig. 6a). Cdh5 pos ; Myh11 low aortic cells expressed EC-specific GPCRs such as Gpr116 or Eltd1, whereas GPCRs normally present in SMao, such as Lgr6 or Mrgprf, were low or absent (Fig. 6b). We therefore concluded that these Cdh5-positive cells are despite low levels of Myh11 expression aortic EC (ECao). Compared to cells from healthy aortae, ECao from atherosclerotic aortae showed clear changes in the expression of both functionrelated genes and GPCRs (Fig. 6c). In detail, markers of inflammatory activation such as Icam1, Vcam1 or Il6 were upregulated, whereas expression of Kdr and Edn1 was reduced (Fig. 6c,d). Most GPCRs showed diminished expression (Fig. 6c,d), resulting in a reduction of the GPCR number expressed per cell (Fig. 6e). However, a trend to increased gene expression was observed for Ccrl2, S1pr3, Darc and Gpr153 (Fig. 6c,d). A comparison of expression changes in acutely and chronically activated EC showed that adhesion molecules and selected GPCRs such as Ccrl2, Darc, F2r, Gpr153 or Lphn2 were upregulated both in acute and chronic inflammation, whereas Kdr, Edn1, Pdgfb and various GPCRs were downregulated in both conditions (Fig. 6f). Interestingly, a number of GPCRs were upregulated in acute inflammatory activation, but downregulated in chronic inflammation, for example, Cxcr7, Gpr56, Gprc5a and others (Fig. 6f).
GPCR repertoire in SMC from atherosclerotic mice. We next studied GPCR patterns in SMC from atherosclerotic aortae. After 16 weeks of high-fat diet, SMao showed clear signs of inflammatory activation and dedifferentiation 26,27 : Icam1, Vcam1 and Il6 were upregulated, whereas contractile proteins such as Acta2 and Myh11 were reduced in expression strength, though not in frequency (Fig. 7a, upper part). In addition, numerous GPCRs showed increased expression frequency (Fig. 7a, lower part), resulting in a significantly increased number of GPCRs per individual cell (Fig. 7b). K-means cluster analysis assigned the majority of SMao from atherosclerotic mice (squares in Fig. 7c) to a cluster characterized by increased expression of GPCRs such as Olfr78, Ednrb, Adcyap1r1, Ptger2, Avpra1a, S1pr3, Ptgir, Vipr2, Ccrl2, Cxcr7, Gprc5b or Gprc5c (Fig. 7d). ApoE-deficient mice without high-fat diet showed an intermediate phenotype, but no clear upregulation of inflammatory genes such as Icam1 or Vcam1; comparable changes were observed in aged but otherwise healthy C57BL/6 mice (Fig. 7e). To understand the mechanisms regulating GPCR expression in healthy and dedifferentiating SMao, we analysed transcription factor (TF) binding sites in promoters of GPCRs that were either up-or downregulated in dedifferentiating SMC. GPCRs that were upregulated in atherosclerotic SMao were more likely to contain binding sites Gpr116 S1pr1 Eltd1 for TFs such as heat shock factors proteins 1,2,4 (HSF1,2,4), retinoic acid receptor a/b (RARA, RARB), NF-kB2, AP2 (TFAP2A,B,C), KLF5 and others (Fig. 7f, left side). GPCRs that were downregulated in atherosclerotic SMao were more likely to contain bindings sites for, among others, estrogen-related nuclear receptors ESRRA and ESRRB, sterol regulatory element-binding transcription factor 1 (SREBF1), or T box TFs TBX1 and TBX19 (Fig. 7f, right side; complete list in Supplementary Table 5).  Dedifferentiating SMC of the healthy aorta. We noticed that some SMao from healthy mice were assigned to the cluster of dedifferentiating SMC in Fig. 7c (marked by arrows), indicating that a spontaneously dedifferentiating population exists in the healthy aorta. Indeed, k-means analysis of SMao from healthy mice identified a small subpopulation of cells as being clearly different from the rest (cluster 3 in Fig. 8a). These cells were characterized by a reduced expression of contractile proteins Myh11 or Acta2, but increased expression of Icam1, Vcam1, Col1a2 and Col3a1 (Fig. 8b). This dedifferentiated phenotype was associated with reduced expression of GPCRs such as Lgr6 and Adra1d, and increased expression of Ptgir, Vipr2, Gpr39, Lpar1, Gprc5b, Gpr124, Cxcr7, Gpr137b, Ednra and others (Fig. 8b). Spearman's rank correlation analysis confirmed the positive correlation between Myh11, Acta2, Lgr6, Adra1d on the one hand and Icam1, Vcam1, CD44, Ptgir, Gpr39, Vipr2, Cxcr7, Gprc5b and so on on the other hand (Fig. 8c). A direct comparison of changes in dedifferentiating SMC from healthy mice and atherosclerotic mice showed that both cell types behaved similarly with respect to upregulation of adhesion molecules and GPCRs such as Ptgir, Vipr2, Gprc5b, Gprc5c, Cxcr7 and so on, as well as downregulation of Acta2, Myh11, Adra1d and Lgr6 (Fig. 8d, left side). However, a number of GPCRs were upregulated only in the context of atherosclerotic dedifferentiation, and not in spontaneous dedifferentiation, such as Olfr78, Ednrb, Adcyap1r, Ptgfr and others (Fig. 8d, right side).
Since vascular SMC dedifferentiation is believed to occur mainly in regions of disturbed flow, for example, the inner curvature of the aortic arch 28 , we used Gprc5b-bgal reporter mice to investigate the localization of spontaneously dedifferentiating SMC in vivo. Immunohistochemical analysis of the expression of bgal and SMC marker a-smooth muscle actin (aSMA) in transverse section of the aortic arch (Fig. 8e) showed only sparse bgal expression in SMC of the outer curvature, but an enrichment of bgal/aSMA-positive cells in the inner curvature (Fig. 8f,g). These data confirm the notion that GPCRs such as Gprc5b are selectively expressed in SMC in atheroprone regions. We finally investigated whether selective targeting of GPCRs that are specifically expressed in dedifferentiating SMao may be used to modulate their functional state. Knockdown of Gprc5b resulted in freshly isolated SMao from ApoE-deficient mice in a significant increase in expression of Icam1 and Il6, indicating that this receptor plays a modulatory role in inflammatory gene expression in dedifferentiating SMC (Fig. 8h).

Discussion
We report here the analysis of GPCR expression in primary vascular cells on the single-cell level. Compared to conventional expression analysis in pooled cell RNA/cDNA, single-cell expression analysis has three crucial advantages: First, it allows to rigorously exclude contaminating cells and thereby precludes misinterpretation of data. Second, single-cell expression analysis is able to detect transcripts in rare cell populations that might be below threshold in bulk RNA/cDNA analyses: 37 GPCRs that were judged negative or uncertain in SMsk by NanoString analysis were detected in individual SMsk, though in some cases only in 1.5% of cells. However, due to the limited amount of RNA obtained from an individual cell also single-cell expression may overlook very low abundance transcripts, in particular if the RNA quality is compromised. Since GPCRs are in many cases expressed at low levels 29 , we directly compared the two major readout systems for single-cell expression analysis, RT-PCR and mRNA sequencing, with respect to their sensitivity for GPCR detection in SMao. RT-PCR detected most GPCRs with higher frequency than mRNA sequencing, and sequencing of the single-cell amplificates as well as analyses in GPCR reporter mice and on the protein level largely confirmed the RT-PCR data. A possible explanation for this difference in sensitivity lies in the fact that single-cell RT-PCR uses target-specific pre-amplification of mRNAs, while pre-amplification for mRNA sequencing is unbiased.
The third and probably most relevant advantage of single-cell expression analysis is that it allows to estimate the degree of GPCR heterogeneity within a cell population and, consecutively, to identify correlations between GPCR profile and functional state of a given cell. We found that all types of SMC and EC showed a surprisingly high heterogeneity of GPCR expression, and reporter analysis, as well as studies on the protein level confirmed these findings. Studies in other fields, mainly developmental biology, suggest that this is not a specialty of GPCRs, but is also observed in other gene families 17,18 . These findings have major implications for pharmacotherapy, since current interpretations of GPCR expression data rely on the assumption that all cells of a given population are equal, or at least comparable 15,16,30,31 . Our data not only clearly disprove this assumption for the majority of GPCRs, they also open up the possibility to selectively target pathologically altered cells based on their specific GPCR repertoire.
We show, for example, that dedifferentiated SMao differ from normal SMao not only in the expression of typical markers indicating inflammatory activation and dedifferentiation, but also in their GPCR repertoire. Among those receptors that are preferentially expressed in dedifferentiating SMao are a number of G s -coupled receptors with known anti-inflammatory and anti-proliferative properties, such as the prostacyclin receptor Ptgir or the VIP receptor Vipr2 (refs 32,33). It is tempting to speculate that also other G s -coupled receptors upregulated in dedifferentiating SMao will exert anti-proliferative effects, for example, the PACAP receptor Adcyap1r1, or the adenosine  Actb Gapdh Hprt receptor Adora2a. In line with this notion, PACAP was shown to inhibit SMC proliferation 34 , and enhanced SMC proliferation was observed in Adora2a-deficient mice 35 . Furthermore, it will be particularly interesting to investigate whether those orphan receptors that are upregulated in dedifferentiating SMC, for example, Gpr39, Gpr124, Gpr153 or Gprc5b have the potential to positively or negatively regulate SMC differentiation. In support of this idea, we found that knockdown of Gprc5b enhanced proinflammatory gene expression in freshly isolated SMao, suggesting that this receptor negatively modulates inflammatory gene expression in SMC. How the GPCR repertoire of an individual vascular cell is shaped and how stable it is over time, is unclear. While numerous studies analysed the posttranslational regulation of GPCRs by phosphorylation, internalization or dimerization, their transcriptional control is little understood. To address the mechanisms regulating GPCR expression in healthy and dedifferentiating SMao, we analysed TF binding sites in promoters of GPCRs upregulated in dedifferentiating SMao. The promoters of these GPCRs were, among others, enriched in binding sites for TFs AP2, KLF5, RARA/RARB, HSF1/2/4 and NF-kB2. Some of these TFs have been implicated in SMC dedifferentiation: KLF5 has been shown to promote proliferation of vascular SMC 36 , whereas activation of retinoic acid receptor a (RARA) increased migration and tissue-type plasminogen activator activity in SMC 37 . Both HSF1 and NF-kB show increased activity in Avpr1a S1pr3 S1pr3 Adora2a S1pr1 Ccrl2 atherosclerotic vessels 22,38 , and NF-kB activation enhances expression of inflammatory genes such as Vcam1 and reduces expression of SM-specific genes 39 . These data suggest that TFs implicated in the regulation of SMC dedifferentiation also control changes in the GPCR repertoire, which in turn might modulate the dedifferentiation process.
Our data also show that both SMC and EC express distinct GPCR repertoires depending on the vascular bed they originate from. SMC from vascular beds rich in resistance arteries, here SMsk and SMmes, express significantly more GPCRs than aortic SMC, in particular more peptide hormone receptors and orphan receptors. Most of the peptide hormones in question have been shown to affect vascular tone: vasorelaxation is, for example, mediated by the gastrointestinal hormones VIP and PACAP (acting on receptors encoded by Vipr1, Vipr2 and Adcyap1r1), by urocortin and corticotropin releasing hormone (acting on the gene product of Crhr2), by parathyroid hormone (acting on the receptor encoded by Pth1r), or by calcitonin gene-related peptide or adrenomedullin (acting on the gene product of Calcrl) [40][41][42][43][44] ; all these receptors are G s coupled 45 . Vasoconstriction, in contrast, is induced by peptide mediators such as endothelin-1, angiotensin II, vasopressin or neuropeptide Y acting on the G i/o -, G q/11 -or G 12/13 -coupled receptors encoded by Ednra, Ednrb, Agtr1a, Avpra1a and Npy1r 46-48 . The second-largest group of GPCRs with preferential expression in SMsk/SMmes are the orphan receptors, for example, Gpr19, Gpr21, Mrgprf, Gprc5b, Gprc5c, Gpr124, Gpr126 or Lphn1. This suggests that also orphan GPCRs expressed in SMsk are potential targets for modulation of blood pressure, but due to the lack of known agonists/antagonists so far no data are available with respect to their function in regulation of vascular tone or other SMC functions.
Also EC show, depending on their anatomical location, remarkable differences in the GPCR repertoire, and also their responses to acute inflammatory activation differ. Brain EC, for example, were characterized by an upregulation of chemokine receptor Darc and orphan receptor Gpr153, a response that was absent in lung EC. In line with this finding it was recently shown that DARC is upregulated in brain EC during neuroinflammation, where it shuttles inflammatory chemokines across the blood-brain barrier 49,50 . The role of GPR153 in activated brain EC, in contrast, is completely unclear. Interestingly, we also detected a small subgroup of brain EC that were characterized by a myeloid-like GPCR expression pattern, but were devoid of any myeloid lineage markers. The function of this subpopulation is currently unclear, but the fact that expression of these GPCRs was negatively correlated with expression of Icam1 or chemokine receptors involved in leukocyte trafficking, such as Darc or Cxcr7(ref. 50), suggests that this population is less supportive of leukocyte extravasation than other LPS-activated ECbr. It is also interesting to note that chronic inflammatory activation during atherosclerosis results in other changes than acute activation by LPS, for example, with respect to atypical chemokine receptor Cxcr7, or orphan GPCRs Gpr111, Gpr107, Gprc5a and Gprc5b. It is also noteworthy that ECao from atherosclerotic aortae showed reduced expression of endothelial genes such as Cdh5, Pecam1, Gpr116 or Gpr56, suggesting that EC, very much like SMC, undergo dedifferentiation in atherosclerotic aortae in vivo.
Taken together, our results show that expression profiling on the single-cell level allows the identification of receptors that are predominantly or selectively expressed on pathologically relevant vascular subpopulations. Understanding the GPCR profile of these cells will not only significantly enhance our knowledge about the pathomechanisms of disease, it will also allow a more selective pharmacological targeting of these cells.
Single-cell quantitative real-time PCR. For single-cell expression analysis, B1,500 cells were loaded in a volume of 5 ml onto the microfluidic C1 Single-Cell Auto Prep System (5-10 mm mRNA Arrays, Fluidigm, San Francisco, CA, USA), followed by RNA isolation and cDNA synthesis according to the manufacturer's protocol. cDNAs derived from visually empty chambers or chambers containing more than one cell were excluded from further analysis; also cDNAs that showed poor expression of reference gene Hprt in conventional RT-PCR were excluded. High-throughput quantitative PCR on single-cell cDNAs was performed on 96.96 Dynamic Array IFCs on the BioMark system (Fluidigm) using Sso-Fast EvaGreen Supermix low ROX (BioRad, Hercules, CA, USA) and Delta Gene primer assays (Fluidigm) as listed in Supplementary Table 1. Primers used for the single-cell analysis of human genes are shown in Supplementary Table 2. Only single-cell cDNAs negative for lineage markers Ptprc (CD45R, all leukocytes), Itgam (CD11b, monocytes/macrophages), Cd4 (CD4 T cells), Cd8 (CD8 T cells), Cd19 (B cells), Ly6g (neutrophils), Cdh1 (E-Cadherin, epithelial cells), Tnni2 (skeletal muscle), but positive for Gapdh, Myh11 (for SMC), or Cdh5 (for EC) were included into analyses (exception: Cdh5 pos ; Myh11 low ECao). Because C1 cDNA samples are potentially contaminated with genomic DNA, intron-spanning primer design was used for all genes. The limit of detection for the BioMark HD System has been estimated to be at a Ct value of 24 cycles (Limit of detection (LoD) Ct); all sample Ct values were therefore subtracted from the (LoD) Ct using the formula: gene expression ¼ 2 (LoD Ct-sample Ct) 54 .
Tissue digestion and cell sorting. Animals were killed, perfused with PBS and different tissues (lung, skeletal muscle of lower limbs, brain, mesentery, whole aorta including common carotid arteries) were dissected, minced and enzymatically digested for 90-120 min while shaking at 37°C in a digestion mix containing collagenase II (2 mg ml -1 ; Worthington), elastase-I (0.04 mg ml -1 ; Sigma), and DNase1 (5 U ml -1 ; New England Biolabs). For lung digestion, also dispase II (1.2 units per ml; Sigma) was added to the mix. Cell suspensions were serially filtered through 70 and 40 mm cell strainers followed by washing twice with PBS. Cells from Cdh5-Cre pos ; dTom/EGFP-reporter pos mice or Myh11-CreERT2 pos ; dTom/EGFP-reporter pos mice were sorted on a JSAN cell sorter (Bay Biosciences, Japan) based on their EGFP expression; in other cases antibodies directed against CD31 (Serotec MCA2388PE) or CD144-PE (ebioscience #12144180) were used to identify EC.
Bioinformatic analysis of single-cell RT-PCR data. RT-PCR-based expression values for individual cells were analysed by customized functions from the RaceID R package 18 . For k-means clustering, we utilized the clustexp and clustheatmap functions; cluster specific genes (Po0.05, according to expected transcript count probability from binomial testing, see Methods section from 18 were exported by a customized clustdiffgenes function. Graphs show fold changes (41.5 or o0.7) from cluster specific genes. T-SNE plots were generated by the plottsne function. Heat maps based on RT-PCR expression data were generated in Perseus software 55 . The correlation between expression of individual genes was determined using Spearman's rank correlation coefficient; graphical representations of Spearman's rank correlation coefficients were generated using Cytoscape 56  To investigate the potential role of transcription factors (TFs) in cell type-specific GPCR sets, we utilized the pscan tool 57 . We investigated the promoter region of GPCRs within a given set at a size of À 950 to þ 50 bp around the respective transcription start site and tested for significantly (Po0.05) overrepresented/underrepresented TF binding sites listed in the current Jaspar 2016 release (http://jaspar.genereg.net/) 58 compared to a global promoter set of the mouse organism. Resulting P values were z transformed and illustrated in a multi column heat map generated by the heatmap2 R function (package gplots, version 3.01).
Single-cell mRNA sequencing and expression analysis. Single-cell transcriptome analysis was performed on a C1 Autoprep station (Fluidigm) using 5-10 mm mRNA arrays and standard protocols. In total, 15,000 events sorted on a JSAN swift sorter were suspended in 10 ml resuspension buffer (PBS þ 0.5% bovine serum albumin þ 2 nM EDTA, sterile filtered) and further diluted (6:4) in suspension buffer. A volume of 3 ml was loaded to the C1 Array while remaining cells (B10,000) were used for total RNA isolation using mRNeasy micro kit (Qiagen) combined with on-column DNase digest (Qiagen). RNA was quantified by Qubit HS RNA Assay (Thermo Fisher) and used for tube control. Pre-amplification and cDNA synthesis were done with Smarter ultra-low RNA Kit (Clontech) without changes to the standard protocol. cDNA obtained from C1 system was quantified using PicoGreen assay (Thermo Fisher) and diluted to 0.1-0.3 ng ml -1 . Picking of positive cDNA samples and library preparation were done with NGS Express pipetting robot (Perkin Elmer) using custom made protocols. Barcoded and amplified libraries were pooled in groups of 16 libraries for final SRRI-bead based cleanup and quantified by Qubit HS dsDNA Assay (Thermo Fisher) and BioAnalyzer 2100 HS DNA Assay (Agilent). Sequencing was performed on Nextseq500 Sequencer (Illumina) using v2 chemistry and High Output Flow Cell with 75 bp single end protocol. Raw data files were trimmed using Reaper with a quality cutoff of 53 (length 20 bp, prefix 50 bp) and a minimum clean read length of 50 bp. Reads were then mapped to the genome (mm10) using STAR. FeatureCounts was used to assign reads to genes from the Gencode vM6 annotation. Genes with less than 50 counts across all cells were discarded before normalization. The knn.error.models function from SCDE 59 was used to construct cell-specific error models by evaluating the gene counts for a observed cell versus the expected gene counts from the k most similar cells. The parameter k was chosen to resemble roughly three subpopulations. Estimated fragments per million mapped reads (FPM, in log scale) were then obtained using the scde.expression.magnitude function and are given in Supplementary Data 1.
NanoString analysis of GPCR expression in bulk RNA. NanoString analyses were performed as described previously 60,61 . In brief, 250-500 ng RNA from sorted cells was applied in a total volume of 30 ml in the assay. Barcodes were counted for B1,150 fields of view per sample. Counts were first normalized to the geometric mean of the positive control spike count, then a background correction was done by subtracting the mean þ two s.d. of the eight negative control counts for each lane. Data were not normalized to reference genes because none of the reference genes showed sufficiently stable expression in all cell types according to the geNorm algorithm. Values that were o20 were fixed to background level.
Cell culture and transfection. Isolated human aortal SMC (passage 1) were obtained from Cell Applications (Cat. No. 354-05a; San Diego, CA, USA) and Innoprot (Cat. No. P10456; Derio, Spain). Cells were thawed and immediately processed for expression analysis without further culture. Primary murine aortic SMC were isolated as described above and cultured in collagen I-coated (100 mg ml À 1 ; Corning) 24-well plates for 7 days. Cells were transfected on day 0 with siRNAs directed against murine Gprc5b (Mm_Gprc5b_6; target sequence: 5 0 -TCGGGCCTACATGGAGAACAA-3 0 ) or AllStars negative control siRNA using Lipofectamine RNAiMAX. All siRNAs were obtained from Qiagen.
Flow cytometric analysis. b-galactosidase expression in murine aortic SMC was analysed using the FlouReporter lacZ flow cytometry kit (ThermoFischer) with propidium iodide to identify viable cells and anti-aSMA antibody (ebioscience, #53976080) to identify SMC. For antibody-mediated detection of GPCR expression, aortae of wild-type mice were digested, fixed in 0.01% formaldehyde and permeabilized using Tween20 (0.5% v/v in PBS), followed by incubation of resulting single-cell suspensions with PE-or FITC-labelled antibodies directed against aSMA (ebiosciences) and APC-labelled antibodies directed against receptors Cmklr1 (Miltenyi), Ccrl2 (R&D), Cxcr7 (R&D), Celsr2 (R&D) (in all cases 10 ml antibody per 10 6 cells per 200 ml). For each GPCR-specific antibody the corresponding isotype control was used. Analyses were performed on a FACS Canto II (Becton Dickinson) and data were analysed by FlowJo software (Tree Star).
General statistical analyses. Data are presented as means±s.e.m. if not otherwise indicated. Comparisons between two groups were performed using two-sample t-test; normalized data (control group set to 1) were analysed by one-sample t-test. P values are indicated as follows: *Po0.05; **Po0.01; ***Po0.001.
Data availability. Single-cell RNA-Seq data have been deposited in GEO database (https://www.ncbi.nlm.nih.gov/geo, under the accession code GSE97955. Full data sets generated during the current study are available from the corresponding author on reasonable request.