PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptNIH Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
Sci Signal. Author manuscript; available in PMC Jan 20, 2014.
Published in final edited form as:
PMCID: PMC3896501
NIHMSID: NIHMS502727
Dissecting common γ chain cytokine family signaling in T cells using cell-to-cell variability analysis
Jesse W. Cotari,#1,2 Guillaume Voisinne,#1,2 Orly Even Dar,1 Volkan Karabacak,1,2 and Grégoire Altan-Bonnet1,2,4
1ImmunoDynamics Group, Programs in Computational Biology and Immunology, Memorial Sloan-Kettering Cancer Center, New York, NY 10065, USA
2Center for Cancer Systems Biology, Memorial Sloan-Kettering Cancer Center, New York, NY 10065, USA
#Contributed equally.
4To whom correspondence should be addressed: altanbog/at/mskcc.org
Natural variability in abundance of signaling regulators can lead to divergence in cell fate, even within genetically identical cells sharing a common differentiation state. To leverage this observation, we introduce cell-to-cell variability analysis (CCVA), an experimental and computational methodology to quantify the correlation between variability in signaling regulator abundance and variation in sensitivity to stimuli. Here, we apply CCVA to investigate the unexpected effects of the interleukin 2 (IL-2) receptor α chain (IL-2Rα) on the sensitivity to common-gamma chain (γc) cytokines in primary T lymphocytes. Our work demonstrates that increased IL-2Rα abundance decreases the concentration of IL-2 but increases the concentrations of IL-7 and IL-15 required for a half-maximal activation (EC50) of downstream signal transducer and activator of transcription 5 (STAT5), without affecting the EC50 of other γc cytokines. To probe the mechanism of IL-2Rα's effect on γc family cytokine EC50s, we introduce a Bayesian-inference computational framework that models the formation of receptor signaling complexes using prior biophysical measurements. Using this framework, we demonstrate that a model in which IL-2Rα drives γc depletion through pre-assembly of complete IL-2 receptors is consistent with both CCVA data and prior measurements. The combination of CCVA and computational modeling yields quantitative understanding of the crosstalk of γc cytokine signaling in T lymphocytes.
Quantifying the impact of protein abundance on cellular function has attracted considerable attention in recent years (14). To do so in bacteria, researchers have changed incrementally the abundance of a chosen protein and measure the functional consequences (59). However, this approach is more cumbersome in primary mammalian cells, such that protein function has been principally studied in an all-or-nothing fashion using genetic mutants or RNAi. As an alternative, we propose that inherent natural variability in protein abundance, as recently observed within populations of genetically identical mammalian cells (1015), can be used to dissect the quantitative regulation of signal transduction.
To assess the phenotypic variability of populations of isogenic cells, researchers can quantify the variability abundance of mRNA or protein with individual cell resolution using single cell RT-qPCR (16) or flow cytometry (12, 13). Of note, studies utilizing these techniques have demonstrated large heterogeneity in the abundance of signaling components (receptors, kinases, adapters, phosphatases and cytokines), with typical coefficients of variation (CV) for the lognormal distribution of mRNA and protein amounts larger than 0.5 within activated T cell clones for example (16, 17). Concretely, in such distributions, 15% of cells will have protein abundances deviating from the median by more than two fold. Even larger variability was uncovered in the case of the interleukin 2 (IL-2) receptor α chain (IL-2Rα), with CVs of up to 3.0 in populations of genetically identical transgenic T cells activated in vitro (17). In these cells, 15% of the population has IL-2Rα abundance that deviates from the median by more than 10 fold in either direction (17). In settings of infection, this variability in T cells' IL-2Rα abundance has been shown to correlate with a split between short-lived effector or memory precursor fates (18). Examples in which a continuous spectrum of surface protein abundance maps onto discrete differentiation paths have been reported in other biological systems as well (12, 19). These observations raise the question of how variability in protein abundance affects signaling thresholds and ultimately cell differentiation decisions.
In this work, we present an experimental methodology to quantitatively correlate such variability in protein abundance with variable regulatory function. As part of this methodology, we present a software package, ScatterSlice, which quantitatively correlates the functional heterogeneity of cellular signaling thresholds with variation in protein abundance. We apply this methodology to probe the signaling responses to common γ chain (γc) cytokines, interleukins (IL-) 2, 4, 7, 9, 15 and 21, whose receptors share the γc chain (CD132) (20). IL-2 and IL-15 also share their β chain (CD122). Each cytokine receptor has an α chain that imparts specificity through enhanced affinity for its respective cytokine. For example, the IL-2 receptor is formed through the assembly of the γc and β chains with a specific α chain, IL-2Rα (CD25) whereas the IL-7 receptor is composed of the γc and its specific α chain, IL-7Rα (CD127). Though it has been established that expression of each cytokine's α chain is required for signaling (21), quantitative understanding of the effects of sharing a receptor chain on the ability of T cells to trigger downstream phosphorylation in response to γc cytokines remains elusive.
Specifically, we investigated whether variation in IL-2Rα abundance influences sensitivity to IL-2 and other γc family cytokines. In agreement with previous results (17), we report that large variation in IL-2Rα correlated with much of the variability in cells' sensitivity to IL-2. More surprisingly, we found that higher abundance in IL-2Rα correlated with inhibited sensitivity to IL-7. Decreases in IL-7Rα abundance at the cell surface or direct inhibition of IL-7Rα signaling were ruled out as mechanisms explaining the dependency of IL-7 sensitivity on IL-2Rα. Instead, a simple thermodynamic model suggested that large abundance of IL-2Rα could inhibit IL-7 signaling indirectly, by sequestering γc chains into preformed heterotrimeric IL-2 receptors in the absence of IL-2. Extending our experimental system to other cytokines of the γc family, we determined that IL-15 sensitivity is slightly dampened by IL-2Rα abundance, whereas sensitivities to IL-4 and IL-21 are unaffected.
Phenotypic variability between T cell populations leads to differences in IL-2 EC50
We sampled isogenic populations of primary T cells at different timepoints following activation and found that a given concentration of IL-2 can trigger greatly varied signaling responses, as reflected in different degrees of STAT5 phosphorylation (Fig. 1A). We then analyzed T cells' pSTAT5 response to a complete dose-response of IL-2 after a 10-minute exposure. In one culture, we estimated the typical concentration (EC50) of IL-2 that yields half-maximal pSTAT5 response to be 10 pM (Fig. 1B, green line); this is consistent with the estimated dissociation constant (KD) for the “high affinity IL-2 receptor” (22), as well as functional readouts such as cell proliferation of CTLL-2 cells (23). However, in other cultures, we observed shifts in IL-2 EC50 by ten-fold— up to 100 pM or down to 1 pM (Fig. 1B). Flow cytometry measurements of IL-2Rα abundance revealed a close anti-correlation between IL-2 EC50 and the median abundance of IL-2Rα (Fig. 1C). Cells from a culture with higher EC50 have lower median abundance of this receptor subunit, while cells from a culture with lower EC50 have greater median abundance. Overall, the median abundance of IL-2Rα per cell correlates well with IL-2 EC50 (R2 = 0.87; Fig. 1C).
Fig. 1
Fig. 1
IL-2 EC50s are anti-correlated with IL-2Rα abundance
Cell-to-cell variation in IL-2Rα abundance leads to intra-population EC50 variation
Cells within a given population have differences in IL-2Rα abundance of the same magnitude as differences in median IL-2Rα abundance between populations (Fig. 1D). Since the variability in IL-2 EC50 between T cell cultures correlates well with variability in median IL-2Rα abundance, we conjectured that the same correlation should exist when examining individual cells within a population.
This cell-by-cell assessment requires that IL-2Rα abundance be measured with precision via flow cytometry. Dual labeling of IL-2Rα with two different antibodies (Fig. 1E) demonstrates broad variation in either channel, but minimal variation in the ratio of two stains (Fig. 1F). Thus, even with the combined noise from staining with two antibodies and measurement across two channels, IL-2Rα abundance can be determined with less than 1.5-fold error across a wide range (up to 105-fold). Cell sorting experiments further confirm that distinct subpopulations can be accurately defined by their IL-2Rα abundance (fig. S1).
To quantify systematically the correlation between cytokine EC50 and receptor abundance, we have developed a three-stage methodology, named cell-to-cell variability analysis (CCVA, Fig. 2A). First, cells are treated with a dose-response of stimulus. Second, protein abundance and downstream responses are measured concomitantly with single-cell resolution using flow cytometry. Finally, the flow cytometry data are analyzed with a software package named ScatterSlice (Fig. 2B–E). Conventional programs for flow cytometry analysis deliver snapshots of cellular response for a given stimulus dose, with the possibility of manually gating for subpopulations based on protein abundance. In contrast, our software was specially designed to automatically parse the heterogeneous population into subpopulations defined by protein abundance (Fig. 2B), quantify each population's downstream phosphorylation (Fig. 2C), then determine stimulus sensitivity using all doses of stimuli to fit an EC50 within each subpopulation (Fig. 2D). As a whole, CCVA delivers a complete map of the relationship between protein abundance and response sensitivity as quantified by the EC50 (Fig. 2E).
Fig. 2
Fig. 2
Cell-to-cell variability analysis (CCVA) methodology
We applied CCVA to quantify the effect of variation in IL-2Rα abundance on IL-2 EC50 within one population of T lymphocytes. We found that, as receptor abundance ranges from 102 to 107, the EC50 of the IL-2 response decreases from 100 pM to 1 pM (Fig. 2E). Extending the same analysis to the different T cell populations shown in Figure 1 reveals that despite large differences in median IL-2Rα abundance and corresponding shifts in the overall IL-2 EC50 from population to population (see Fig. 1C), an identical anti-correlation between IL-2Rα abundance and EC50 is observed (fig. S2).
CCVA uncovers increased IL-7 EC50 correlating with increased IL-2Rα abundance
In effector T cells, IL-2Rα (i.e. CD25) is commonly used as a marker of activation, identifying an active effector phenotype between naïve and memory stages (24). Furthermore, effector T cells with lower abundance of IL-2Rα have been shown to preferentially commit to a memory phenotype, whereas those with higher abundance of IL-2Rα do not (18). As both naïve and memory phenotypes are classically associated with enhanced sensitivity to IL-7 (25) we used CCVA to quantify the relationship between IL-2Rα abundance and IL-7 EC50. This analysis demonstrated a positive correlation between IL-2Rα abundance and IL-7 EC50, with increased IL-2Rα abundance (from 102 to 107 copies) scaling with increased EC50 for IL-7 (from 3 to 60 pM) (Fig. 3A).
Fig. 3
Fig. 3
Positive correlation between IL-7 EC50 and IL-2Rα abundance: CCVA and modeling
Decrease of IL-7 receptor abundance or signaling does not account for observed increase in IL-7 EC50
The impact of IL-2Rα abundance on IL-7 EC50 could be hypothesized to stem from the reported temporal anti-correlation of IL-2Rα and IL-7Rα abundance during the course of differentiation from naïve to effector to memory phenotypes (25, 26). However, in our system, flow cytometric analysis demonstrated that abundances of IL-2Rα and IL-7Rα at the cell surface were not anti-correlated (fig. S3A). We confirmed this result by FACS sorting cells according to IL-2Rα abundance and demonstrating by western blot a lack of correlation with IL-7Rα abundance (fig. S3B). More critically, our investigations of population-level EC50 revealed that although IL-7 EC50 varies from 3 to 60 pM across populations, this variation does not correlate with IL-7Rα abundance with any statistical significance (fig. S3C). This lack of correlation reflects a lack of uniformity in EC50 despite fairly uniform abundance of IL-7Rα across effector cell cultures, suggesting that factors other than IL-7Rα abundance regulate IL-7 EC50 in effector cells.
As an alternate mechanism, we tested whether IL-2Rα abundance correlated with a direct inhibitory modification of the IL-7Rα chain leading to fewer functional IL-7 receptors despite equal abundance. If this were the case, all cytokines signaling through IL-7Rα would have reduced EC50 correlating with higher abundance of IL-2Rα. To test this possibility, we applied CCVA to characterize whether the EC50 of T cells' response to thymic stromal lymphopoietin (TSLP)—a cytokine whose receptor also utilizes IL-7Rα (fig. S4A)— also correlates with IL-2Rα abundance. Our results showed no correlation between TSLP EC50 and IL-2Rα abundance (fig. S4B). Therefore, the increase in IL-7 EC50 cannot be attributed to direct inhibition of IL-7Rα signaling.
An equilibrium binding model of γc chain competition quantitatively accounts for IL-2Rα effects on IL-7 EC50
We then focused on the fact that IL-2 and IL-7 receptors share the γc chain (21, 27) (Fig. 3B). Specifically, we hypothesized that greater abundance of IL-2Rα might drive the pre-formation of complete heterotrimeric IL-2 receptors composed of IL-2Rα, IL-2Rβ, and γc chains. This sequestration of the γc chain in preformed IL-2 receptor complexes could lead to the increase in IL-7Rα EC50 with increased IL-2Rα abundance, as reported in Fig. 3A.
Experimental evidence for existence of complete heterotrimeric IL-2 receptor complexes on the surface of living cells in the absence of cytokine remains equivocal. Indeed, the crystal structure of IL-2 in complex with the soluble extracellular portions of the full receptor indicates no contact between IL-2Rα and either β or γc chains (28). Moreover, isolated extracellular domains of IL-2Rα and γc have no detectable binding affinity in solution (29). On the other hand, the extracellular domains of IL-2Rα and β have a binding affinity of 280 nM, as measured in solution (29). Furthermore, all pairwise interactions between IL-2Rα, β and γc have been demonstrated on the surface of living cells in the absence of IL-2 (20, 30). It is therefore possible that in the absence of IL-2, high abundance of IL-2Rα would drive the preformation of full heterotrimeric IL-2 receptor complexes.
In the following, we set out to investigate whether this hypothesis could realistically and quantitatively account for the correlation between IL-2Rα abundance and IL-7 EC50 observed in our CCVA experiments using an equilibrium binding model. We first compute the EC50 of the response to IL-2 and IL-7 from affinity constants (see Supplementary Material) and quantitatively compare these results to values measured by CCVA. We then identify the simplest equilibrium binding model with realistic affinity constants that can account for the correlation of IL-2 and IL-7 EC50s with the abundance of IL-2Rα uncovered by CCVA.
Estimating the cytokine EC50 using an equilibrium binding model for the receptor chains is appropriate because: 1) the system reaches a steady state within 10 minutes of exposure to cytokines (fig. S5); 2) pSTAT5 correlates linearly with receptor occupancy (fig. S6); and 3) depletion of cytokines is negligible (fig. S7), such that cytokine concentrations can be taken as constant on our experimental timescale. In this framework, a model for cytokine EC50 is determined by the set of possible receptor complexes formed, or equivalently by the set of affinity constants associated with their formation (see Supplementary Material). For this system, all of the receptors chains' abundance at the surface can be measured (fig. S8, Materials and Methods) and all the potentially required affinity constants have been experimentally determined, directly or indirectly (fig. S9) (20, 29, 3136). Thus one should be able to create a model reproducing the data directly from these experimentally determined affinity constants. Based on the potential interactions between subunits of the IL-2 and IL-7 receptors, there exist 553 relevant models that can be enumerated systematically (fig. S10–S11, Supplementary Methods). However, none of the possible models quantitatively reproduce the positive correlation between IL-7 EC50 and IL-2Rα abundance when affinity constants are fixed to their experimentally determined values (Fig. 3D, “Model with measured affinity constants”). The best among these poorly performing models is shown in fig. S12. Alternatively, when the affinity constants are unconstrained and free to vary, very good agreement between data and model output (Fig. 3D, “Model with free affinity constants”) can be obtained for a large number of models, but at the expense of fitted affinity constants that greatly deviate from their measured values (fig. S13).
This discrepancy can be resolved by acknowledging some intrinsic uncertainties in the experimental determination of affinity constants (see fig. S9) and translation of in vitro measurements to living cells. To do so, we introduced a Bayesian strategy to model the CCVA-measured EC50 values while relying on previous measurements of affinity constants as prior estimates of model parameters (see Supplementary Methods). This procedure allows the identification of a set of models that achieve a good agreement with the CCVA-measured EC50 values (Fig. 3D, “Model with constrained affinity constants”), while maintaining the affinity constants in a realistic range around their measured values (fig. S14–18). We relied on Akaike's information criterion (see Supplementary Material) to select the minimal model that accurately reproduces the CCVA data. This results in a model in which external cytokines bind only to the fully assembled receptors IL-7Rα:γc and IL-2Rα:β:γc, and the opposite correlation of IL-2 and IL-7 EC50s with IL-2Rα abundance can be readily interpreted as a consequence of the preformation of IL-2 receptors depriving IL-7 receptors of the γc (Fig. 3E). This prediction of the model could be tested and confirmed experimentally: we sorted effector T cells into highest and lowest 35% of IL-2Rα abundance, and analyzed their steady-state association of γc, and IL-7Rα by immunoprecipitation and western blot analysis (Fig. 3F). This blot demonstrated that IL-7Rα:γc complexes were present in cells with low abundance of IL-2Rα, but almost totally absent in cells with high abundance of IL-2Rα.
An extended equilibrium binding model accounts for IL-2Rα effects on other γc family cytokine EC50s
IL-7 is only one among the γc family of cytokine whose receptors share some components with the IL-2 receptor. Among them, IL-4, -7, -15, and -21 all triggered pSTAT5 response in effector T cells, whereas IL-9 did not (fig. S19). Using CCVA, we found that the EC50 to IL-4, -7, -15, and -21 do not behave identically as a function of IL-2Rα abundance (Fig. 4A). Aside from the case of IL-2, the EC50 of IL-7 response is the most strongly affected by the abundance of IL-2Rα. For IL-15, a slight dependency can be observed whereas the EC50 of IL-4 and IL-21 responses are largely unaffected by IL-2Rα abundance. To account for the variety of these behaviors, we extended the simple model selected above to include these additional receptors (Fig. 4B), using typical protein abundances defined by flow cytometric calibration (fig. S8).
Fig. 4
Fig. 4
Extension of the preassembly model to full γc family
We modeled the full system using our previously estimated affinity constants related to the formation of the IL-2 and IL-7 receptors and adjusted the other affinity constants (fig. S20) to reproduce the effects of IL-2Rα abundance on the EC50 of the responses to these cytokines (Fig. 4C). We find that high affinity constants for the binding of γc to IL-21R and IL-4Rα are required to recapitulate the independence of IL-4 and IL-21 EC50 from IL-2Rα abundance (fig. S20). Consequently, preformation of IL-4 and IL-21 receptors is largely unaffected by IL-2Rα (Fig. 4D). To our knowledge, affinity constants for the binding of γc to IL-21R and IL-4Rα have not been measured and our model predicts high values for these affinity constants on cell surfaces. Average EC50 is tuned through the binding affinity of the preassembled receptor to the cytokine and abundance of IL-21R or IL-4Rα. Similarly, in the case of IL-15, a high affinity binding between IL-15Rα and β:γc (fig. S20) leads to weak correlation of IL-15 EC50 with IL-2Rα abundance. In our model of the γc cytokine family, independence of a given cytokine EC50 from IL-2Rα abundance is achieved when the corresponding α chain has a high affinity for the shared component(s).
In this study, we have presented a comprehensive methodology to leverage the variation of protein abundance within cell populations in order to quantify variability in cytokine sensitivity and investigate its underlying mechanisms. This method relies on the demonstrated ability of flow cytometry to precisely measure the abundance of multiple proteins simultaneously in single cells (Fig. 1). We have used this methodology to study receptor signaling in the γc chain family, uncovering unexpected relationships between related receptors.
During activation of isogenic T cells, the abundance of IL-2Rα varies by up to five orders of magnitude (17). Our co-staining experiments verified that these measurements of IL-2Rα precisely define IL-2Rα abundance in these cells (Fig. 1F, fig. S1). This allowed us to study the effect of this broad range of IL-2Rα abundance on cytokine sensitivity, as defined by the EC50 of STAT5 phosphorylation stimulated by the cytokine. In agreement with previous findings, we found that IL-2Rα abundance strongly correlates with IL-2 sensitivity, with an EC50 decreasing from 100pM to 1pM as IL-2Rα abundance increases from 102 to 107 proteins (Fig. 1, fig. S2). This is consistent with Cantrell & Smith's observation that cells endowed with higher IL-2Rα abundance proliferate more robustly in vitro (22). Strikingly, we observed that this increase of IL-2 sensitivity is concomitant with a 20-fold decrease of IL-7 sensitivity (Fig. 3). The dependency of IL-7 sensitivity on IL-2Rα abundance could not be explained by anti-correlation of IL-7Rα and IL-2Rα abundance or by a potential inhibition of IL-7Rα signaling capacity as a function of IL-2Rα abundance within activated effector cells (fig. S3–S4). Such anti-correlation of receptor abundance has, however, been observed both during T cells' transitions from naïve to effector and effector to memory stages in vivo: naïve cells have abundant IL-7Rα, but no IL-2Rα; memory cell precursors are marked by relatively high abundance of IL-7Rα (25) and low abundance of IL-2Rα (18). Additionally, we found that highly activated effector cells, with large amounts of IL-2Rα, are inhibited from sensing IL-7: this implies that they would rely chiefly on IL-2 for survival. This could reinforce the suppressive effect of IL-2 depletion by T regulatory cells, preventing IL-7 from rescuing effector cells when IL-2 becomes scarce (17, 37). In this context, the inhibitory effect of IL-2Rα abundance on IL-7 sensitivity represents an additional mechanism to enhance separation of cytokine specificities across different T cell subtypes, with naïve and memory cells relying chiefly on IL-7, while effectors rely on IL-2.
Addressing the impact of IL-2Rα abundance on IL-7 sensitivity at the computational level was made possible because this system is well characterized at the biophysical level (20, 27, 29, 3136). Binding interactions between receptor subunits, and between receptor subunits and cytokines have been quantified (29). However, because many of these measurements have been made in different systems and with varying methodologies—on the surface of cells or with soluble fractions— multiple models have been proposed for the function of these receptors (28). Our CCVA methodology (Fig. 2) allows us to add detailed quantitative data, across a large range of protein abundance, thereby placing tight constraints on any proposed model of IL-2 and IL-7 signaling. In the present work, we show that equilibrium binding models incorporating the measured affinity constants as priors quantitatively reproduce aspects of the dependency of IL-7 and IL-2 sensitivities on IL-2Rα (Fig. 3). The switch of sensitivities between IL-2 and IL-7 as a function of IL-2Rα abundance can be readily attributed to the preformation of the full heterotrimeric IL-2 receptor, which sequesters the γc chain away from the IL-7 receptor. By comparing all possible equilibrium binding models of formation of signaling cytokine-receptor complexes, we have identified a minimal equilibrium binding model with affinity constants constrained around the measured values that accounts quantitatively for the CCVA data. This model was selected using a statistical criterion comprising a cost for deviation from both biochemically measured parameters and CCVA data, and penalizing for inclusion of a larger set of parameters. As such, this selection procedure does not exclude other models with more parameters, representing other paths of receptor formation (35, 38, 39), but more data would be needed to justify their inclusion. The selected model presents the practical advantage of providing a straightforward interpretation of the scaling of IL-7 sensitivity as a function of IL-2Rα abundance. The increased EC50 at high IL-2Rα abundance is a direct consequence of the sequestration of the γc chain in the preformed IL-2 receptor. As a consequence, a cytokine receptor subunit with a high affinity for a shared component (i.e. γc or β) should make it difficult for the IL-2 receptor to sequester this component, and the effect of IL-2Rα abundance on the sensitivity to the corresponding cytokine should be weak.
Extending this minimal model to other cytokines of the γc family, we recapitulated different measured effects of IL-2Rα on cytokine sensitivities (Fig. 4). Specifically, we found that IL-2Rα abundance negatively impacts sensitivity to IL-7 and to a lesser extent to IL-15, but not at all to IL-4 or IL-21. In our model, non-dependence of sensitivity to IL-2Rα for these two cytokines implies high affinity between the corresponding α chains and the shared component(s). The ability of each α chain to tune not only the sensitivity to its respective cytokine, but also cross-antagonism of other signals provides a functional separation between γc-utilizing cytokines that matches their respective functions in the context of T cell differentiation (40). IL-7 is dedicated to the homeostasis of naïve and memory cells, but not to the expansion of effectors, whereas IL-15 affects both effectors and memory cells (41). Therefore, segregating the homeostatic functions of IL-7 from the expansion function of IL-2 may be advantageous (24). In contrast, IL-4 and IL-21 are both involved in the differentiation of effector cells into helper subtypes, and are experienced during expansion (21). For these two cytokines, cross-antagonism from large amounts of IL-2Rα would thus be functionally counter-productive.
Our CCVA method can be applied to virtually any cellular signaling investigation wherein cellular heterogeneity and signaling can be measured, such as other cytokines, chemokines, or growth factors (42, 43). The number of signaling regulators that can be simultaneously measured is limited only by technical constraints of flow cytometry. As the limitations of flow cytometry recede with the introduction of machines capable of recording dozens of parameters simultaneously (10, 11), our method is poised to identify and quantify regulatory effects in signaling networks. This technical framework, accompanied by quantitative modeling efforts, will improve mechanistic understanding of signal transduction cascades in living cells (2, 3).
Reagents and antibodies
Recombinant mouse IL-2, IL-4, IL-7, IL-12, IL-21, interferon gamma, and thymic stromal lymphopoietin (TSLP) were from eBioscience (San Diego, CA); IL-15 was from R&D systems (Minneapolis, MN). See table below for details of antibodies used.

EpitopeCloneConjugateUsageSupplier

CD4RM4-5Alexa-700FCeBio

IL-2Rα (CD25)PC61.5PE-Cy7FCeBio

PC61.5Per-CP-Cy5.5FCeBio

PC61.5PEFACSMiltenyi

7D4FITCFCeBio

IL-7Rα (CD127)eBioSB/199PEFCeBio

P14noneWBscbt

γc (CD132)4G3PEFCeBio

H-300noneIPscbt

IL-2Rβ (CD122)TM-βlPEFCeBio

IL-4Rα (CD124)mIL4R-M1PEFCBD

IL-21R4A9PEFCBD

IL-15-Ra (CD125)DNT15RaAPCFCeBio

pY694-STAT5C11C5noneFCCell Signal

rabbit IgGNAAPCFCJackson IR

goat IgGNAHRPWBscbt

FC - Flow Cytometry, FACS - Fluorescence-Activated Cell Sorting, WB - Western Blot, IP - Immunoprecipitation
PE – Phycoerythrin, APC – Allophycocyanin, FITC – Fluorescein Isothiocyanate, HRP – Horseradish Peroxidase, Per-CP – peridinin chlorophyll protein, Cy - Cyanine

eBio - eBioscience, San Diego, CA
Miltenyi - Miltenyi Biosciences, Cambridge, MA
Cell Signal - Cell Signaling Technology, Danvers, MA
BD - BD Biosciences, San Jose, CA
scbt - Santa Cruz Biotechnology, Santa Cruz, CA
Jackson IR - Jackson Immunoresearch, West Grove, PA

Media
All in vitro experiments were performed in complemented RPMI medium (prepared by the Media facility at MSKCC), which consists of RPMI 1640 supplemented with 10% heat-inactivated fetal bovine serum, 2 mM L-glutamine, 10 mM HEPES (pH 7.4), 0.1 mM non-essential amino acids, 1 mM sodium pyruvate, 100 μg/ml of penicillin and 100 μg/ml of streptomycin and 50 μM β-mercaptoethanol. All cell cultures were maintained in an incubator at 37 °C with 5% CO2.
Cell culture
5C.C7 TCR transgenic, Rag2−/− IL-2−/− CD4 T cells were isolated and expanded as described previously (13). Briefly, cells were expanded with irradiated B10A CD3−/− splenocytes pulsed with 1 μM MCC K5 peptide: ANERADLIAYFKAATKF (GenScript, Piscataway, NJ). IL-2 (10 nM) was given at one and four days post stimulation, then every three days for maintenance of cultures. Cells were used between day 5 and 20 of culture, unless otherwise noted. Live cells were purified prior to analysis by separation on a Ficoll gradient (GE, Piscataway, NJ).
Cytokine titrations
Before measuring cytokine responses, bound IL-2 was stripped from cell surface by a 2 minute incubation in 0.1 M glycine buffer equilibrated at pH 4.0. After two washes in RPMI, cells were rested for 30 minutes at 37 °C (44). Cytokine titrations were prepared in 96-well V-bottom plates on ice. Cells (105 per well) were added to cytokines and incubated for ten minutes at 37 °C before processing for flow cytometry.
Flow Cytometry
Cells were prepared for single cell staining following protocols optimized by the Nolan group (45). Briefly, cells were fixed for 10 minutes on ice in 1.6% paraformaldehyde followed by permeabilization on ice in 90% methanol. Cells were then washed twice in fluorescence-activated cell sorting (FACS) buffer (PBS with 4% FCS and 0.1% sodium azide) and stained intracellularly for phospho-STAT5, followed by incubation with labeled secondary antibody and conjugated antibodies to receptors (CD4 and IL-2Rα). All staining was performed at room temperature in FACS buffer. Before flow cytometric acquisition, cells were stained with 5 μM 4',6-diamidino-2-phenylindole (DAPI; Sigma-Aldrich, St. Louis, MO) for exclusion of dividing (G2) and apoptotic (sub-G1) cells. Flow cytometry acquisition was performed on a LSR-II cytometer (BD Biosciences, San Jose, CA).
Receptor abundance calibration
A standard curve correlating fluorescence as measured by flow cytometry to PE molecules was calculated using Bangs Laboratories (Fishers, IN) Quantum R-PE MESF beads. “Corrected median fluorescence intensity” (cMedFI) values were calculated by subtracting the median fluorescence of unstained beads from the median fluorescence value of each PE standard. A linear regression of log cMedFI to log PE bead MESF values was used to establish a standard curve. To calibrate typical abundances of receptor chains, cells were stained with saturating concentrations (as established by titration) of the respective PE-labeled antibodies or a nonspecific PE-labeled rat IgG control. Corrected median fluorescence was calculated by subtracting the median fluorescence of the isotype control from the median fluorescence of each antibody staining. The calibration curve was used directly to calculate numbers of molecules from the corrected median fluorescence (fig. S8), using a ratio of 1 PE molecule per antibody, as specified by the manufacturer (eBioscience or BD Biosciences). For calibration of IL-2Rα abundance in ScatterSlice data, the calibration curve was applied to each bin's fluorescence value.
Population sensitivity fitting
Population EC50 values were determined from flow cytometry data, by fitting the dose response for the geometric mean of pSTAT5 as a function of the cytokine concentration with the following equation:
equation M1
where pSTAT5low and pSTAT5high are the baseline and the plateau of pSTAT5 signals. Fits were performed with Prism (GraphPad Software; La Jolla, CA). The fitting algorithm uses the Hessian matrix (46) to estimate errors on the fitted parameters.
Fluorescence-Activated Cell Sorting (FACS)
Cells were Ficoll purified, then stained with PE- anti-CD25 (Clone PC61.5, Miltenyi Biosciences, Cambridge, MA) and APC-anti-CD4 on ice for 30 minutes. Prior to sorting, 5 μM DAPI was added for exclusion of dead cells. Sorting was performed on a FACSAria from BD Biosciences (San Jose, CA). Sorted cells were kept on ice until further processing.
Western Blotting
Equal numbers of sorted cells were rinsed in PBS, then lysed in buffer containing Tris HCl 25 mM, EDTA 1 mM, NaCl 130 mM, NP-40 1%; supplemented with protease inhibitor cocktail from Sigma Aldrich (St. Louis, MO) containing to final concentrations of 1.04 mM AEBSF, 800 μM Aprotinin, 40 μM Bestatin, 14 μM E-64, 20 μM Leupeptin, and15 μM Pepstatin A. Lysates were mixed with NuPAGE LDS loading buffer with 2.5% beta-mercaptoethanol, incubated at 95 °C for 2 minutes before loading and running on NuPAGE Novex 4–12% bis-tris gels (all products from Life Technologies, Grand Island, NY). Blots were transferred to Immobilon membrane (Millipore, Billerica, MA), prior to blotting with primary (IL-7Rα) followed by secondary (anti-rabbit-IgG-HRP) antibodies in TBS with 0.05% Tween-20 and 2% BSA.
Immunoprecipitation of γc from sorted cells
T cells from (5 × 108) were stripped in 5 ml of 0.1 M glycine buffer equilibrated at pH 4.0. After two washes in RPMI, cells were rested for 30 minutes at 37 °C. Cells were then pelleted and resuspended in 10 ml fresh RPMI for staining, and processed for FACS as described above. 35% highest and lowest abundance IL-2Rα fractions were collected. Following sorting, immunoprecipitation was performed using the Dynabeads co-immunoprecipitation kit (Life Technologies), following the manufacturer's protocol, with the following modifications. Cells were pelleted, weighed, and suspended in a 13:1 ratio of Dynabeads lysis buffer (“Buffering Salts” pH 7.4, 110 mM potassium acetate, 0.5 % Triton X-100, 100 μM sodium chloride). After incubation on ice for 20 min, cell lysates were passed through a 27-gauge syringe 10 times on ice. Immunoprecipitation was carried out according to the manufacturer's instructions, using Dynabeads freshly conjugated to anti-γc antibody (clone H-300). Western blotting was carried out as described above.
ScatterSlice analysis
Initial analysis of flow cytometry data was performed using FlowJo (Treestar, Ashland, OR). Single, CD4 positive cells with 2n DNA content were identified and exported as text files for analysis by our processing R program ScatterSlice.
ScatterSlice was written in R (47). The software is freely available upon request. Details of the algorithm are as follows:
Text files are imported into R, and then divided into user-specified bins. Within each bin, a three parameter Hill equation—fitting base, amplitude and EC50—is fit by minimizing the sum across all doses in a dose-response (from i = 1 to Ndoses) of the squared difference between the log of each dose's geometric mean fluorescence intensity of the response channel (pSTAT5 in our studies) (ziobs) minus the log of the three parameter hill equation (fitting base, amplitude, and EC50), normalized by multiplying by the number of cells in that dose (Ni) divided by the phospho-response variance (σi2). Errors on fitted parameters are estimated using the Hessian matrix (46). Tables of fitting results (values for base, amplitude, EC50 and their errors for different amounts of receptors) are presented as color maps and can be further analyzed with R or exported as text files for analysis with other software.
Supplementary Material
Modeling & Supplementary Materials
Acknowledgments
The authors would like to thank James Sethna, K. Christopher Garcia, Karen Tkach, Nihal Altan-Bonnet, Massimo Vergassola, and Thierry Rose for useful comments and discussion. This work was supported by NSF 0848030, NIH R01-AI083408, NIH U54-CA143798, NIH T32-AIO7621 (JWC).
1. Loewer A, Lahav G. We are all individuals: causes and consequences of non-genetic heterogeneity in mammalian cells. Current opinion in genetics & development. 2011;21:753–758. [PMC free article] [PubMed]
2. Munsky B, Neuert G, van Oudenaarden A. Using Gene Expression Noise to Understand Gene Regulation. Science. 2012;336:183–187. [PMC free article] [PubMed]
3. Pelkmans L. Cell Biology. Using cell-to-cell variability--a new era in molecular biology. Science. 2012;336:425–426. [PubMed]
4. Gaudet S, Spencer SL, Chen WW, Sorger PK. Exploring the contextual sensitivity of factors that determine cell-to-cell variability in receptor-mediated apoptosis. PLoS Comput Biol. 2012;8:e1002482. [PMC free article] [PubMed]
5. Kollmann M, Lovdok L, Bartholome K, Timmer J, Sourjik V. Design principles of a bacterial signalling network. Nature. 2005;438:504–507. [PubMed]
6. Rosenfeld N, Young JW, Alon U, Swain PS, Elowitz MB. Gene regulation at the single-cell level. Science. 2005;307:1962–1965. [PubMed]
7. Suel GM, Kulkarni RP, Dworkin J, Garcia-Ojalvo J, Elowitz MB. Tunability and noise dependence in differentiation dynamics. 2007;315:1716–1719. [PubMed]
8. Swain PS, Elowitz MB, Siggia ED. Intrinsic and extrinsic contributions to stochasticity in gene expression. Proc Natl Acad Sci USA. 2002;99:12795–12800. [PubMed]
9. Raj A, van Oudenaarden A. Nature, nurture, or chance: stochastic gene expression and its consequences. 2008;135:216–226. [PMC free article] [PubMed]
10. Bendall SC, Simonds EF, Qiu P, Amir el AD, Krutzik PO, Finck R, Bruggner RV, Melamed R, Trejo A, Ornatsky OI, Balderas RS, Plevritis SK, Sachs K, Pe'er D, Tanner SD, Nolan GP. Single-cell mass cytometry of differential immune and drug responses across a human hematopoietic continuum. Science. 2011;332:687–696. [PMC free article] [PubMed]
11. Newell EW, Sigal N, Bendall SC, Nolan GP, Davis MM. Cytometry by time-of-flight shows combinatorial cytokine expression and virus-specific cell niches within a continuum of CD8+ T cell phenotypes. Immunity. 2012;36:142–152. [PMC free article] [PubMed]
12. Chang HH, Hemberg M, Barahona M, Ingber DE, Huang S. Transcriptome-wide noise controls lineage choice in mammalian progenitor cells. Nature. 2008;453:544–547. [PubMed]
13. Feinerman O, Veiga J, Dorfman JR, Germain RN, Altan-Bonnet G. Variability and robustness in T cell activation from regulated heterogeneity in protein levels. Science. 2008;321:1081–1084. [PMC free article] [PubMed]
14. Spencer SL, Gaudet S, Albeck JG, Burke JM, Sorger PK. Non-genetic origins of cell-to-cell variability in TRAIL-induced apoptosis. Nature. 2009;459:428–432. [PMC free article] [PubMed]
15. Lee YK, Mukasa R, Hatton RD, Weaver CT. Developmental plasticity of Th17 and Treg cells. Curr Opin Immunol. 2009;21:274–280. [PubMed]
16. Peixoto A, Evaristo C, Munitic I, Monteiro M, Charbit A, Rocha B, Veiga-Fernandes H. CD8 single-cell gene coexpression reveals three different effector types present at distinct phases of the immune response. J Exp Med. 2007;204:1193–1205. [PMC free article] [PubMed]
17. Feinerman O, Jentsch G, Tkach KE, Coward JW, Hathorn MM, Sneddon MW, Emonet T, Smith KA, Altan-Bonnet G. Single-cell quantification of IL-2 response by effector and regulatory T cells reveals critical plasticity in immune response. Molecular systems biology. 2010;6:437. [PMC free article] [PubMed]
18. Kalia V, Sarkar S, Subramaniam S, Haining WN, Smith KA, Ahmed R. Prolonged interleukin-2Ralpha expression on virus-specific CD8+ T cells favors terminal-effector differentiation in vivo. Immunity. 2010;32:91–103. [PubMed]
19. Chang JT, Palanivel VR, Kinjyo I, Schambach F, Intlekofer AM, Banerjee A, Longworth SA, Vinup KE, Mrass P, Oliaro J, Killeen N, Orange JS, Russell SM, Weninger W, Reiner SL. Asymmetric T lymphocyte division in the initiation of adaptive immune responses. Science. 2007;315:1687–1691. [PubMed]
20. Pillet A-H, Lavergne V, Pasquier V, Gesbert F, Thèze J, Rose T. IL-2 Induces Conformational Changes in Its Preassembled Receptor Core, Which Then Migrates in Lipid Raft and Binds to the Cytoskeleton Meshwork. Journal of Molecular Biology. 2010;403:671–692. [PubMed]
21. Rochman Y, Spolski R, Leonard WJ. New insights into the regulation of T cells by gamma(c) family cytokines. Nature reviews. Immunology. 2009;9:480–490. [PMC free article] [PubMed]
22. Cantrell DA, Smith KA. The interleukin-2 T-cell system: a new cell growth model. 1984;224:1312–1316. [PubMed]
23. Smith KA. The interleukin 2 receptor. 1988;42:165–179. [PubMed]
24. Mazzucchelli R, Durum SK. Interleukin-7 receptor expression: intelligent design. Nat Rev Immunol. 2007;7:144–154. [PubMed]
25. Kaech SM, Tan JT, Wherry EJ, Konieczny BT, Surh CD, Ahmed R. Selective expression of the interleukin 7 receptor identifies effector CD8 T cells that give rise to long-lived memory cells. Nature immunology. 2003;4:1191–1198. [PubMed]
26. Palmer MJ, Mahajan VS, Chen J, Irvine DJ, Lauffenburger DA. Signaling thresholds govern heterogeneity in IL-7-receptor-mediated responses of naive CD8(+) T cells. Immunology and cell biology. 2011;89:581–594. [PMC free article] [PubMed]
27. Wang X, Lupardus P, Laporte SL, Garcia KC. Structural biology of shared cytokine receptors. Annual review of immunology. 2009;27:29–60. [PubMed]
28. Wang X, Rickert M, Garcia KC. Structure of the quaternary complex of interleukin-2 with its alpha, beta, and gammac receptors. Science. 2005;310:1159–1163. [PubMed]
29. Rickert M, Boulanger MJ, Goriatcheva N, Garcia KC. Compensatory energetic mechanisms mediating the assembly of signaling complexes between interleukin-2 and its alpha, beta, and gamma(c) receptors. Journal of molecular biology. 2004;339:1115–1128. [PubMed]
30. Damjanovich S, Bene L, Matkó J, Alileche A, Goldman CK, Sharrow S, Waldmann TA. Preassembly of interleukin 2 (IL-2) receptor subunits on resting Kit 225 K6 T cells and their modulation by IL-2, IL-7, and IL-15: a fluorescence resonance energy transfer study. Proceedings of the National Academy of Sciences of the United States of America. 1997;94:13134–13139. [PubMed]
31. Wang H-MM, Smith KA. The interleukin 2 receptor. Functional consequences of its bimolecular structure. Journal of Experimental Medicine. 1987;166:1055–1069. [PMC free article] [PubMed]
32. Myszka DG, Arulanantham PR, Sana T, Wu Z, Morton TA, Ciardelli TL. Kinetic analysis of ligand binding to interleukin-2 receptor complexes created on an optical biosensor surface. Protein science : a publication of the Protein Society. 1996;5:2468–2478. [PubMed]
33. Balasubramanian S, Chernov-Rogan T, Davis AM, Whitehorn E, Tate E, Bell MP, Zurawski G, Barrett RW. Ligand binding kinetics of IL-2 and IL-15 to heteromers formed by extracellular domains of the three IL-2 receptor subunits. International Immunology. 1995;7:1839–1849. [PubMed]
34. McElroy CA, Dohm JA, Walsh STR. Structural and Biophysical Studies of the Human IL-7/IL-7R[alpha] Complex. Structure. 2009;17:54–65. [PMC free article] [PubMed]
35. Rose T, Pillet A-H, Lavergne V, Tamarit B, Lenormand P, Rousselle J-C, Namane A, Thèze J. Interleukin-7 Compartmentalizes Its Receptor Signaling Complex to Initiate CD4 T Lymphocyte Response. Journal of Biological Chemistry. 2010;285:14898–14908. [PubMed]
36. Noguchi M, Nakamura Y, Russell SM, Ziegler SF, Tsang M, Cao X, Leonard WJ. Interleukin-2 receptor gamma chain: a functional component of the interleukin-7 receptor. Science (New York, N.Y.) 1993;262:1877–1880. [PubMed]
37. Pandiyan P, Zheng L, Ishihara S, Reed J, Lenardo MJ. CD4+CD25+Foxp3+ regulatory T cells induce cytokine deprivation-mediated apoptosis of effector CD4+ T cells. Nat Immunol. 2007;8:1353–1362. [PubMed]
38. de Bakker BI, Bodnár A, van Dijk EMHP, Vámosi G, Damjanovich S, Waldmann TA, van Hulst NF, Jenei A, Garcia-Parajo MF. Nanometer-scale organization of the alpha subunits of the receptors for IL2 and IL15 in human T lymphoma cells. Journal of cell science. 2008;121:627–633. [PubMed]
39. Hand TW, Morre M, Kaech SM. Expression of IL-7 receptor alpha is necessary but not sufficient for the formation of memory CD8 T cells during viral infection. Proceedings of the National Academy of Sciences of the United States of America. 2007;104:11730–11735. [PubMed]
40. Junttila IS, Mizukami K, Dickensheets H, Meier-Schellersheim M, Yamane H, Donnelly RP, Paul WE. Tuning sensitivity to IL-4 and IL-13: differential expression of IL-4Ralpha, IL-13Ralpha1, and gammac regulates relative cytokine sensitivity. J Exp Med. 2008;205:2595–2608. [PMC free article] [PubMed]
41. Romano E, Cotari JW, Barreira da Silva R, Betts BC, Chung DJ, Avogadri F, Fink MJ, St Angelo ET, Mehrara B, Heller G, Munz C, Altan-Bonnet G, Young JW. Human Langerhans cells use an IL15Ralpha/IL15/pSTAT5-dependent mechanism to break T-cell tolerance against the self-differentiation tumor antigen, WT1. Blood. 2012 [PubMed]
42. Zikherman J, Jenne C, Watson S, Doan K, Raschke W, Goodnow CC, Weiss A. CD45-Csk phosphatase-kinase titration uncouples basal and inducible T cell receptor signaling during thymic development. Immunity. 2010;32:342–354. [PMC free article] [PubMed]
43. Saini M, Sinclair C, Marshall D, Tolaini M, Sakaguchi S, Seddon B. Regulation of Zap70 expression during thymocyte development enables temporal separation of CD4 and CD8 repertoire selection at different signaling thresholds. Sci Signal. 2010;3:ra23. [PubMed]
44. Duprez V, Cornet V, Dautry-Varsat A. Down-regulation of high affinity interleukin 2 receptors in a human tumor T cell line. Interleukin 2 increases the rate of surface receptor decay. J Biol Chem. 1988;263:12860–12865. [PubMed]
45. Krutzik PO, Clutter MR, Nolan GP. Coordinate analysis of murine immune cell surface markers and intracellular phosphoproteins by flow cytometry. J Immunol. 2005;175:2357–2365. [PubMed]
46. Gutenkunst RN, Waterfall JJ, Casey FP, Brown KS, Myers CR, Sethna JP. Universally sloppy parameter sensitivities in systems biology models. PLoS computational biology. 2007;3:1871–1878. [PubMed]
47. R Development Core Team R: A Language and Environment for Statistical Computing. Vienna, Austria: 2011.