{'fraction_non_alphanumeric': 0.020080751777611032, 'fraction_numerical': 0.008968449637331618, 'mean_word_length': 5.32497175141243, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Highly multiplexed simultaneous detection of RNAs and proteins in single cells

Precise gene expression measurement has been fundamental to developing an advanced understanding of the roles of biological networks in health and disease. To enable detection of expression signatures specific to individual cells we developed PLAYR (Proximity Ligation Assay for RNA). PLAYR enables highly multiplexed quantification of transcripts in single cells by flow- and mass-cytometry and is compatible with standard antibody staining of proteins. With mass cytometry, this currently enables simultaneous quantification of more than 40 different mRNAs and proteins. The technology was demonstrated in primary cells to be capable of quantifying multiple gene expression transcripts while the identity and the functional state of each analyzed cell was defined based on the expression of other transcripts or proteins. PLAYR now enables high throughput deep phenotyping of cells to readily expand beyond protein epitopes to include RNA expression, thereby opening a new venue on the characterization of cellular metabolism.


Introduction
Biological systems operate through the functional interaction and coordination of multiple cell types. Whether one is trying to delineate the complexity of an immune response, or characterize the intrinsic cellular diversity of cancer, the ability to perform single-cell measurements of gene expression within such complex samples can lead to a better understanding of system-wide interactions and overall function.
A current method of choice for study of transcript expression in individual cells is single-cell RNA-seq. This approach involves physical separation of cells, followed by lysis and library preparation with protocols that have been optimized for small amounts of input RNA 1 -11 .
Barcoding of physically separated cells before sequence analysis makes possible the analysis of thousands of individual cells in a single experiment 12 . However, sample handling (such as separation of live cells before lysis) has been shown to induce significant alterations in the transcriptome 13 . Moreover RNA-seq requires cDNA synthesis and does not enable simultaneous detection of protein epitopes and transcripts. The complexity of protocols and the associated costs further limit the applicability of this technology in studies where sample throughput is essential. Finally, the number of cells that can be analyzed is limited by the overall sequencing depth available. These limitations notwithstanding, the possibility of taking a genome-wide approach to the study of gene expression in single cells, coupled with precise quantification through the use of Unique Molecular Identifiers, make single-cell RNA-seq an exceptionally promising technology 14 .
A complementary approach is to quantify a smaller number of transcripts while increasing the number of cells that can be analyzed. Flow cytometry allows multiple parameters to be measured in hundreds to thousands of cells per second. For such a purpose, fluorescence in situ hybridization (FISH) protocols have been adapted to quantify gene expression on cytometry platforms 15 -20 . In such experiments bright FISH signals with excellent signal-tonoise ratios are necessary since flow cytometry does not provide the subcellular imaging resolution necessary to distinguish individual RNA signals from diffuse background. Different techniques have been adapted for the generation and amplification of specific hybridization signals including DNA padlock probes in combination with rolling circle amplification (RCA) 21 , 22 or branched DNA technology 23 . Recently the branched DNA approach has been successfully applied to flow cytometry 24 but the availability of only three non-interfering branched DNA amplification systems and the spectral overlap of fluorescent reporters complicates multiplexing. What was missing for higher parameter purposes was a technology that allowed full access to the parameterization enabled by mass cytometry 25 and also allowed for protein epitopes to be simultaneously measured.
The Proximity Ligation Assay for RNA (PLAYR) system as described here addresses these limitations by enabling routine analyses of thousands of cells per second by flow cytometric approaches and simultaneous detection of protein epitopes and multiple RNA targets. The method preserves the native state of cells in the first step of the protocol, detects transcripts in intact cells without the need for cDNA synthesis, and is compatible with flow cytometry, mass cytometry, as well as microscope-based imaging systems. Making use of the different measurement channels of mass cytometry, this enables the simultaneous quantitative acquisition of more than 40 different proteins and RNAs. Thus, PLAYR adds a unique and flexible capability to the growing list of technologies that merge 'omics datasets (transcript, protein, and signaling levels) in single cells. We expect that PLAYR will lead to a better understanding of stochastic processes in gene expression 26 -28 and allow for deeper insights into complex cell populations.

Overview of the technology and PLAYR probe design
PLAYR uses the concept of proximity ligation 29 , 30 to detect individual transcripts in single cells, as shown schematically in Fig. 1a, and is compatible with immunostaining. Pairs of DNA oligonucleotide probes (probe pairs) are designed to hybridize to two adjacent regions of target transcripts in fixed and permeabilized cells. Each probe in a pair is composed of two regions with distinct function. The role of the first region is to selectively hybridize to its cognate target RNA sequence. The second region, separated from the first by a short spacer, acts as template for the binding and circularization of two additional oligonucleotides (termed backbone and insert). When hybridized to an adjacent probe pair, the backbone and insert oligonucleotides form a single-stranded DNA circle that can be ligated. The ligated, closed circle is then amplified through rolling circle amplification by phi29 polymerase initiated by the 3′ OH of one member of the probe pair. As phi29 continues to polymerize, it creates a linear molecule that contains hundreds of concatenated complementary copies of the original circle 31 . Then, using a labeled oligonucleotide that is complementary to the insert region of the amplicon, one can detect any given probe pair through binding to the amplified product. For analysis by flow cytometry fluorescently labeled oligonucleotides are used for detection. Alternatively, metal-conjugated oligonucleotides enable mass cytometric analyses using a CyTOF instrument 25 .
Lowering of background binding events and increased specificity result from the fact that both PLAYR probes must hybridize independently to adjacent locations of a target RNA in order for the two independent ligation events and subsequent RCA to take place. Nonspecific, off-target binding of single probes does not result in a signal. This is in contrast to FISH approaches, where many different labeled oligonucleotides bind to transcripts directly and any unspecific binding of individual probes can lead to background signals that cannot be distinguished from specific RNA signals using flow cytometry platforms. PLAYR can be multiplexed by designing oligonucleotides with different insert regions that act as cognate barcodes for given transcripts 32 . Different insert sequences are designed to have the same melting temperatures and base compositions and have been found to lead to the formation of RCA products with similar efficiency ( Supplementary Fig. 1). To ensure that the resulting RCA products uniquely barcode a particular transcript the insert sequences do not have common substrings longer than 4 bases.
An open-source R software package has been developed for rapid design of PLAYR probes (Supplementary Fig. 2 and Supplementary Software). The package comes with point-bypoint installation instructions and includes a graphical user interface front end for userfriendly probe design. Candidate probe pairs with similar thermodynamic properties are first produced using the Primer3 software 33 . The application then displays the location of the probes along the target transcript sequence and other characteristics including BLAST matches to other transcripts or to repetitive sequences of the genome 34 and the position of non-constitutively spliced exons (see Supplementary Fig. 2 for additional comments regarding probe design). These features are used to guide the selection of specific probe pairs. For each gene, the user can then manually select the best probe pairs in combination with one of the PLAYR insert systems for multiplexing. Based on these selection criteria the software outputs the complete sequences of probes that can be used to detect transcripts of interest. The sequences of all probes and backbone-insert systems used in this manuscript can be found in Supplementary Table 1.

Simultaneous quantification of proteins and RNAs in single cells
In a first experiment probe pairs specific for beta-actin (ACTB) were designed. In Jurkat T cells that had been fixed and permeabilized, the PLAYR protocol led to a signal that was detected well above background by flow cytometry (Supplementary Fig. 3). No signal was observed if any component of the signal generation or amplification cascade was omitted. Similarly, no signal was observed when sense probes, probes with the same half of the insert-targeted sequence, or combinations of probes targeting different genes were used ( Supplementary Fig. 3). To further demonstrate the specificity of the approach, the protocol was used with one or several probe pairs designed to detect CD10 and CD3E transcripts, which are known to be expressed in pre-B cells and T cells, respectively (Supplementary Fig. 4). As expected, signal intensities for these transcripts increased when multiple PLAYR probe pairs were used simultaneously. Interestingly, the resulting signal increase was in certain cases more than additive. This may be due to formation of RCA products generated from probes of two different pairs on the same transcript even though the target regions of the two probes are not adjacent as they would be for a cognate pair. For instance, bound probes may be brought into proximity in unexpected manners by the structure of RNA molecules in three dimensions. Supporting this there was an inverse relationship between the strength of the signal and the distance in the target hybridization regions of two probes ( Supplementary Fig. 4). Using more than one probe pair per transcript simultaneously as a mixture leads to an increase in signal and improves the sensitivity of the method. Using multiple probe pairs that target different regions of a transcript also makes results for individual genes more reproducible as it limits variability due to differences in probe accessibility to target sites, secondary RNA structure, and alternative splicing.
In general we found that a mixture of four or five probe pairs per gene led to reliable and robust detection of both rare and highly abundant transcripts and we recommend using four or five probe pairs whenever possible. Fewer probes can be used if a transcript is known to be highly abundant. Since PLAYR is primarily used to compare the abundance of a given transcript across different biological conditions rather than to compare expression levels of different transcripts, varying numbers of probe pairs can be used to simultaneously detect different transcripts with different expression levels. Compared to single molecule FISH and branched-DNA approaches, the small number of probes used in PLAYR experiments allows for selection of probes with ideal properties and for targeting of short transcripts. Though we did not apply the principle in this report, we note that careful design of probe pairs could be used to delineate splice variant complexities and genomic translocations in genes of interest.
Using 5 probe pairs per gene we detected the three housekeeping genes HMBS, PPIB, and GAPDH in U937 cells by mass cytometry (Fig. 1b). This application demonstrates a dynamic range of PLAYR that enables the detection of highly abundant transcripts (GAPDH) as well as low abundant transcripts (HMBS) that have been detected at only about 10 copies per U937 cell using branched DNA technology 24 . To further investigate to what extent PLAYR signals correlate with the underlying abundance of a transcript, results obtained with PLAYR and with RT-qPCR were compared for the induction of the cytokines interferon-γ (IFNG) and chemokine ligand 4 (CCL4) in NKL cells at different time points after stimulation with PMA-ionomycin. As shown in Fig. 1c, PLAYR and qPCR measurements were correlated (R 2 values of 0.93 (CCL4) and 0.72 (IFNG)), indicating that PLAYR can be used to assess relative changes in transcript abundance across different biological conditions.
An additional requirement in the PLAYR protocol was the simultaneous detection of transcripts and proteins. The protocol was therefore optimized using conditions that preserve binding of antibodies. Best results were obtained when antibody staining was performed immediately after cell fixation (i.e. at the beginning of the protocol). After antibody staining, amine-to-amine crosslinking using the BS 3 crosslinker was used to prevent antibodies from being washed away during the procedure. This crosslinking step made PLAYR compatible with the majority of antibodies tested. We found that transient permeabilization of cells by the addition of 0.2% saponin in the presence of RNase inhibitors during antibody staining greatly enhanced the preservation of RNA integrity ( Supplementary Fig. 5). Furthermore, this transient permeabilization can be leveraged to stain intracellular proteins with antibodies. Using this protocol NKL cells were stimulated with PMA-ionomycin, in presence of protein-secretion inhibitors, and changes in IFNG protein and transcript levels were determined as a function of time (Fig. 1d). Thus, PLAYR allows studies of the dynamic nature of transcription and translation at the single-cell level. Both RNA and protein expression consistently showed bimodal distributions, which suggested that not all cells in this supposedly homogeneous cell line responded equally to stimulation and highlighted the potential of PLAYR to assess the functional capacity of individual cells within complex populations.

Highly multiplexed detection of specific proteins and transcripts
Using the insert-based multiplexing strategy illustrated in Fig. 1a, multiple targets can be detected simultaneously within individual cells. We designed probes to target 14 different transcripts and first evaluated them individually and then simultaneously in Jurkat T cells by mass cytometry (Fig. 2a). For this experiment cells were incubated either with probes against individual transcripts or with a mixture of all probes. Appropriate control combinations of non-cognate probe pairs were included to demonstrate probe pair specificity. Critically, the presence of insert-backbone oligonucleotides did not lead to observable signals if corresponding cognate probes were not also present in the reaction. Furthermore, the signal amplitude for any given target in the multiplexed sample was not affected by the presence of oligonucleotides against non-cognate targets and corresponding amplification products. This suggests that the number of transcripts that can be quantified within the same cell is only limited by the number of reporters that can be conjugated to detection oligonucleotides and analyzed simultaneously with a given platform.
We made use of the multiplexing capability of PLAYR to simultaneously detect the transcripts of 11 different cytokines and other effector molecules in NKL cells that had been primed with three cytokines (IL2, IL12 and IL18) and stimulated with PMA-ionomycin. Instead of a uniform cellular response, simultaneous transcript quantification revealed complex combinatorial RNA expression patterns in this supposedly homogenous clonal cell line (Fig. 2b). Based on such multiplexed measurements, high-dimensional analysis methods can be leveraged to identify functional NKL subpopulations based on transcript expression profiles. To that end we clustered cells based on the expression of induced effector transcripts, which revealed a remarkable complexity of cellular responses (Fig. 2c). A number of studies 35 , 36 have shown that supposedly homogenous cell populations in primary samples also express such a remarkable diversity of cytokine combinations. While the functional implications of this observation are still poorly understood, the study of any such combinatorial phenomenon clearly benefits from the increased parameterization enabled by PLAYR.
The increased multiplexing capabilities of PLAYR also enable RNA-only experiments, where transcript expression is used to define cell types in which expression patterns of other transcripts can then be studied. Such experiments can be set up at a fraction of the costs typically associated with antibody-based experiments and are not limited by the availability of antibodies for genes of interest. We analyzed an artificial mixture of cells that contained mouse embryonic fibroblasts (MEFs), mouse embryonic stem cells (mESCs), and differentiating mESCs based on the expression of 15 different transcripts. We then visualized the data using viSNE 37 , an algorithm that maps high-dimensional cytometry data onto two dimensions in a manner that best separates cell populations from the original highdimensional space. This type of analysis clearly defined the three different populations of cells in the mixture based on RNA expression (Fig. 3a). Subsequently, different markers of pluripotency (e.g. Nanog), differentiation (e.g. Thy1), proliferation (Mki67), as well as pluripotency-associated long intergenic non-coding RNAs (Lincenc1) could be studied in the context of this cellular system (Fig. 3b). In an independent experiment with mESCs, 22 different transcripts including negative controls were measured in 16 technical replicates, which demonstrated the very high reproducibility of PLAYR for low as well as highly abundant transcripts ( Supplementary Fig. 6).
We further validated this approach by making use of the protein co-detection and multiplexing capabilities of PLAYR. For this experiment we analyzed primary human peripheral blood mononuclear cells (PBMCs) for 10 cell surface proteins and corresponding transcripts. In contrast to the previous experiment, antibody stained protein markers were used to create a viSNE 37 analysis, which enabled the visualization of the major cell types in human peripheral blood (Fig 3c). Subsequent addition of the data on expression of corresponding transcripts demonstrated remarkable cell-type specificity in mRNA expression patterns (Fig 3d). In general, transcripts showed more gradual differences in expression levels than the corresponding proteins and the analysis revealed a clear discrepancy in the case of ITGAX, for which the protein was detected whereas essentially no transcript was observed in a distinct subpopulation of cells (Fig. 3e). Such differences can be explained by longer half-lives of proteins compared to transcripts, post-transcriptional regulation, and the fact that transcription has been observed to occur in bursts 26 -28 .
Accordingly, this type of analysis demonstrates some of the challenges when defining cell types based on individual transcripts and the potential of PLAYR to study the relationship of transcripts and proteins in individual cells within complex primary samples.

Profiling of cytokine transcript induction in primary samples
We next used PLAYR to monitor cytokine transcript induction in PBMCs upon stimulation with lipopolysaccharide (LPS) to correlate protein marker expression with the functional capacity of individual cells. Cytokine expression in single cells is traditionally evaluated on the protein-level by flow cytometry after treatment with secretion inhibitors that lead to accumulation of cytokines in the cells. This approach precludes the study of (and is complicated by) paracrine effects, such as intercellular communication and feedback loops. We used antibodies against surface markers to distinguish different cell populations within human PBMCs while monitoring the expression of a panel of cytokine genes at the transcript level with PLAYR. Similar experiments were performed using fluorescence-based flow cytometry and mass cytometry. The fluorescence experiment involved the detection of four transcripts and four surface markers, whereas mass cytometry allowed for the simultaneous quantification of 8 transcripts and 18 protein epitopes.
In both experiments, antibody staining enabled gating of distinct cell populations (gating for mass cytometry shown in Fig 4a, see Supplementary Fig. 7 for flow cytometry). As expected cytokine production was restricted to the monocyte compartment and therein confined mostly to individual cells that expressed the LPS co-receptor CD14 (Fig 4b). Moreover, different cytokines consistently exhibited distinct expression dynamics. For example, tumor necrosis factor alpha (TNF) and interleukin 8 (CXCL8) were induced early and the former peaked between 2 and 4 hours, while the latter continued to increase during the entire time course. Conversely, expression of interleukin 6 (IL6) was only strongly induced after 4 hours (Fig 4c and 4d). These results recapitulated previous individual observations 38 , 39 and confirmed that PLAYR effectively detects RNA expression in specific cellular subpopulations. Interestingly, while CXCL8 at its peak is expressed in the entire CD14 + monocyte compartment, there was a distinct population of CD14 + cells that did not express TNF (Fig 4e). This observation underscores the utility of protein and RNA co-detection in identifying functional differences in cellular populations. viSNE 37 analysis using the CyTOF data for the cytokine induction experiment demonstrated that all major PBMC populations clustered in unique areas of the viSNE map (Fig. 5a) and could be identified by observing the restricted expression of canonical markers (Fig. 5b).
Similarly, MAP kinase signaling as measured by p38 MAP kinase phosphorylation could be monitored and was restricted to the myeloid compartment. When cytokine transcript expression was overlaid on the map, cells that responded to LPS were mostly restricted to the CD14 + monocytes region (Fig. 5c). This analysis provides a single-cell resolution map of cytokine induction and MAP kinase signaling in PBMCs, highlighting the potential of PLAYR in combination with mass cytometry for system-wide analyses of transcriptional networks in complex samples.

Discussion
PLAYR enables highly multiplexed measurement of gene expression in hundreds to thousands of intact cells per second. On the protein level, single cell measurements have been shown to have prognostic and diagnostic value in multiple clinical settings 40 -43 .
PLAYR extends such analyses to include measurements on the transcript level and can supplement the use of antibodies especially where exon-specific expression is concerned and no relevant antibody reagents exist. Immediate measurement of mRNA as enabled by PLAYR can overcome issues introduced with ex vivo processing of live cells and experimental artifacts are further minimized since PLAYR assays for RNA molecules directly without the need for cDNA synthesis 44 .
PLAYR can simultaneously measure transcripts and their encoded proteins, thus enabling the characterization of the interplay between transcription and translation at the single-cell level. Other applications will include clustering of complex cellular populations purely on the basis of transcript abundance, which could be particularly useful when the availability or quality of antibodies is limiting. We believe that such an approach will help in the definition of yet undetermined cellular populations, which share specific patterns of temporal or spatial regulation of RNA expression. Of relevance to this last point, PLAYR can be deployed for index sorting 45 , 46 and imaging approaches such as fluorescence microscopy and multiplexed ion beam imaging 47 , making it a flexible tool to study gene expression in single intact cells on a variety of platforms.

Cell fixation, permeabilization, and antibody staining
Cells at a density of ~1×10 6 /mL were fixed in RPMI medium without serum in 1.6% paraformaldehyde (Electronic Microscopy Sciences) for 10 min at room temperature under gentle agitation as described previously 49   oligonucleotides were used as detection probes. Also used were a biotinylated oligonucleotide and an oligonucleotide labeled with a single BrdU nucleotide at the 5′ end; cells were then incubated with PE-Cy7-streptavidin or an anti-BrdU-PE antibody conjugate as appropriate. For mass cytometry, cells were incubated with metal-conjugated detection oligonucleotides at a concentration of 10 nM for 30 minutes at 37 °C in PBS, 5 mg/mL BSA, 0.02 % sodium azide. After washing, cells were processed immediately on a fluorescence-based flow cytometer or further processed for CyTOF acquisition as described elsewhere 50 .

Preparation of metal-conjugated detection oligonucleotides
Maleimide-activated Maxpar metal chelating X8 polymers (Fluidigm, Maxpar labeling kit) were loaded with metals and purified using centrifugal filters as per the manufacturer's instructions. Detection oligonucleotides carrying a 5′ Thiol-Modifier C6 S-S (Glen Research) were synthesized at the Stanford Protein and Nucleic Acid Facility. Oligonucleotides were resuspended in DEPC-treated water at 250 μM, and the thiol was reduced by treatment with 50 mM TCEP (Pierce) for 30 min at room temperature. After ethanol precipitation, oligonucleotides were resuspended in C buffer (Fluidigm, Maxpar labeling kit) and conjugation reactions were performed with 2 nmol of oligonucleotide per reaction with X8 polymer. After 2 h at room temperature, TCEP was added to a final concentration of 5 mM, and samples were incubated for 30 min to reduce unconjugated oligonucleotides. Conjugates were filtered through 30-kDa centrifugal filter units (EMD Millipore) in a total of 500 μl water, spun at 14000 g for 12 min, and washed twice with DEPC-treated water (Life Technologies). Purified detection oligonucleotide conjugates were resuspended in DEPC-treated water at a concentration of 1 μM and stored at 4 °C.

RT-qPCR
RNA was extracted using RNeasy Plus Micro Kit (Qiagen), following the manufacturer's instructions. Reverse transcription was performed using SuperScript III First-Strand Synthesis System (Life Technologies), following the manufacturer's instructions. PCR was carried out in a LightCycler 480 II (Roche) using SYBRGreen I Master (Roche).  Each row represents a sample to which probe pairs for one gene only or all genes simultaneously (bottom row) were added. Each column represents a mass cytometry acquisition channel that monitors a metal reporter used to detect transcripts of a given gene. Non-cognate probes that are using the same insert system but bind to different target transcripts were included as an additional control (CTL). b) NKL cells were primed with IL2, IL12 and IL18 and stimulated with PMA-ionomycin for 3 hours. Contour plots display co-expression of NKL effector transcripts as measured by mass cytometry. c) 10000 cells were randomly sampled from the data in (b) and transcript expression was represented in heat map format. Each column corresponds to a single cell and rows denote different effector transcripts. Rows and columns of the heat map were clustered for visual clarity (dendrograms not shown). Experiments were run multiple times and representative examples are shown.