Balancing speed and accuracy of polyclonal T cell activation: a role for extracellular feedback

Background Extracellular feedback is an abundant module of intercellular communication networks, yet a detailed understanding of its role is still lacking. Here, we study interactions between polyclonal activated T cells that are mediated by IL-2 extracellular feedback as a model system. Results Using mathematical modeling we show that extracellular feedback can give rise to opposite outcomes: competition or cooperation between interacting T cells, depending on their relative levels of activation. Furthermore, the outcome of the interaction also depends on the relative timing of activation of the cells. A critical time window exists after which a cell that has been more strongly activated nevertheless cannot exclude an inferior competitor. Conclusions In a number of experimental studies of polyclonal T-cell systems, outcomes ranging from cooperation to competition as well as time dependent competition were observed. Our model suggests that extracellular feedback can contribute to these observed behaviors as it translates quantitative differences in T cells’ activation strength and in their relative activation time into qualitatively different outcomes. We propose extracellular feedback as a general mechanism that can balance speed and accuracy – choosing the most suitable responders out of a polyclonal population under the clock of an escalating threat.

I is the amount of free IL-2, U and B are the amounts of unbound and bound IL-2R, respectively. Table S1 summarizes the parameters of the model. These parameters were used for the model presented in the main text.

A constant induced synthesis rate of receptors
Here, we consider the case in which receptor synthesis rate depends only on m, and not on B, hence is independent of IL2. In this case, we can use Eqs. S1, with K f = 0. In general, the fixed points of (S1) for only one T-cell obey: With K f = 0, Eqs. S2 can be reduced into a quadratic equation for U, which can be solved as: In the case of strong induction (mg f >> g c , U >> dα), the number of unbound receptors in steady state is U ≈ m(g f -s)/a U ≈ 10000 and the number of bound IL2R is B ≈ ms/a B ≈ 2500. The recognition curve (i.e. bound IL2R a function of m) is linear as long as U >> dα that is when m >> m c1 = (a U dαg c )/(g fs). For any reasonable parameters, m c1 << 1 and thus the recognition curve is linear for almost all range of m values ( Fig. 2A in the main text). In this case, the T-cells always benefit from the interaction (Fig. S1) and the interaction does not depend on time.
Figure S1: Phase space diagrams for the case of a constant induced receptor synthesis rate. The normalized interaction index, C (as defined in Eq. 5 in the main text) is presented for both interacting cells (color bar). The interaction does not depend on time. The cells always benefits, in a negligible way, from the interaction.

A linear induced synthesis rate of receptors
Here, we consider the case in which receptor synthesis rate depends on both m and B. However, the IL2 induced feedback is linear, unlike the non-linear feedback that is discussed in the main text. Thus, the term m·g f ·(B 2 / (B 2 + K f 2 ) is replaced with m·g f ·B. The fixed points are given by: In a similar way to the case of constant synthesis, B is almost always a linear function of m (

A constant secretion rate
In our model, the secretion rate of IL-2 depends on the strength of the TCR activation, m. The overall secretion rate is s•m where s is the maximal secretion rate. However, it is insightful to consider also two variants of the model: 1) The case of constant secretion rate, that is, assuming that all activated cells secrete IL-2 with the same rate ( Fig. S2A) and 2) Assuming that in addition to the IL2 secretion that depends on TCR activation there is a constant secretion rate that can, for example, be the result of other unrelated cells that create a 'background' IL2 pool (Fig, S2B).
In the first case, the term s•m is in Eq. S1 is replaced with the term, s, while in the second case it is replaced with the term, s•(m + constant) and the initial value of IL2 is (constant)/d.
The phase-space diagrams for these two cases are shown in figure S2 and are similar to the result obtained with varying secretion rate (Figs. 3, 4 in the main text). This result implies that the main determinant of interaction is the interplay between the TCR strength and the receptor feedback, and not TCR mediated changes in IL2 secretion levels. This is related the well mixed environment assumption. In this environment, the cells share the same IL-2 pool and are not limited by their own secretion rate. However, cooperation and competition still occurs, due to the dependence of the receptor production rate on m. Figure S2: Phase space diagrams for the case of a constant IL-2 secretion rate (A) and the case of adding of a constant secretion rate which is 20% of the maximal secretion rate (B). The normalized interaction index, C (as defined in Eq. 5 in the main text) is presented for both interacting cells (color bar). The black lines are the critical TCR signal strength, m c , for a single cell. Cell1 activates Δt minutes before Cell2. As Δt is larger, the area of Cell2 that can be excluded (blue) increases. In (A) all cells secrete IL-2 with the same ratethe maximal secretion rate, s (see table S1); while in (B) the IL-2 secretion rate in Eq. S1, s•m, is replaced with s•m + 0.2•s.

Sensitivity Analysis
We varied the parameters of the model (table S1). The variation is at least ± 20%. The main effect of changing model parameters is to vary the critical m threshold value. However, the basic features of the system, as indicated by the phase space interaction diagrams, are more stable ( Figure S3-S9). As the scale for m is arbitrarily chosen to be 0-1, the actual range of parameters in which the system shows a similar phase-space diagram, with areas of time-dependent cooperation, co-existence and competition, is larger.
Changing IL-2 secretion and degradation rates have mild effect. As expected, when the level of available IL-2 is higher (higher secretion / lower degradation), the critical threshold value is slightly reduced (Figs. S3,4). In addition, in the case of increased IL-2 secretion rate the exclusion regions for cells above activation threshold are smaller (Fig. S3, top).  We also checked the effect of changing the parameters of constitutive IL-2R dynamics, g c and a U while keeping their ratio, which is the initial number of unbound IL-2R, U, fixed (Fig. S5). Reducing / increasing the constitutive synthesis rate results in a higher / lower threshold, respectively. The overall pattern of cooperation / competition remains similar. Figure S5: Variation in the constitutive receptor synthesis rate, g c and in the unbound receptor internalization rate, a U . The parameters were changed while keeping the number of initial unbound receptors fixed.
The parameters that strongly affect the IL-2R threshold are related to the induced IL-2R synthesis rate, g f , a B and K f (Fig. S6-S8). Reducing the induced synthesis rate, g f, or increasing the induced internalization rate, a B, strongly increases the threshold. Increasing the feedback saturation constant, K f , has a similar effect to increasing a B . As above, the overall pattern of cooperation / competition remains similar around the shifted threshold.  Finally, we check the effect of changing the Hill parameter in the term describing the positive feedback that upregulates the expression of IL-2R. In the main text we used a value oh n = 2. Figure S9 shows that other values of the Hill coefficient, as long as n > 1, provide similar behavior of the system of two interacting cells (see also figure 3B in the MS).

Mapping activation time to initial conditions
An intriguing phenomenon that emerges from our model is the dependence on time of activation. If a strong cell is activated a long time before a mediocre competitor, the weaker cell is excluded. On the other hand, if the mediocre cell is activated a long time before the strong one, they will co-exist. However, in some cases -when the difference between times of activation is small, there could be a scenario in which increasing the time advance of the mediocre cell leads to exclusion (Fig. S10). Figure S10: The competition between m 1 = 1 and m 2 = 0.7. Negative Δt represents that the stronger cell is activated first. If the strong cell is activated a long time before the weaker one (Δt << 0), exclusion of the weaker cell occurs. If the weaker cell is activated a long time before the strong one (Δt >>0), there is co-existence. When Δt is around zero, the behavior is more sensitive to initial conditions.
To further understand this behavior, we studied the propagation of cell induction in a three-dimensional phase space that fully describes the dynamics of the system. At each time point, a cell is characterized by its levels of U and B, and the level of IL-2 in the environment. When the second cell is activated, the fate of the system, exclusion or co-existence, depends on the current values of I, U and B. We plot for the above case (m 1 =1, m 2 =0.7), the trajectory of a single cell in the 3D concentrations phase-space (black line). The second cell joins at some time, at which the first cell is at a particular point along its trajectory. This point determines the outcome of the interaction. Hence, the phase-space (Fig. S11) is composed of points which results in exclusion (red), and points which results in co-existence (green).
In the case that the strong cell, m =1, is activated first (Fig. S11A), the single cell orbit crosses the separator only once. Note that in the region of small concentrations, the orbit is very close to the separator and thus maybe sensitive to changes in initial conditions or to noise. In the case that m = 0.7 activates first, the orbit is so close to the separator that it crosses it twice leading to the "spike" which appears in Figure S10. In these cases, where the concentrations are low, that is the difference in time of activation is small, the system is sensitive to changes in initial conditions. The stronger cell, m = 1, activates first, and its dynamics evolves along the orbit (black line) determined by Eq. (S1). The fate of the system, exclusion of m = 0.7 (red) or co-existence (green), depends on the position of m = 1 along its orbit at the time of arrival of m = 0.7. (B) The weaker cell, m = 0.7, activates first. The cross of the phase separator by the orbit is manifested by the 'spike' in figure S10. The lower plate is a "zoomin" on the separator crossing region. For clarity, we present only a small fraction of the phase space, around the single cell trajectory.

The role of IL-2R subunits
In Eq. S1, we consider the dynamics of an effective IL-2 receptor. The actual IL-2 receptor is composed out of three subunits, α, β and γ. While activation signal is mediated by a complex of these three subunits, only the α subunit synthesis rate is subject to positive feedback.
The dynamics of the full system, assuming that β and γ vary slowly in time, is given by: The corresponding parameters appear in table S2.

Parameter
It is easy to see that these equations are the same as the steady state equations that result from Eq. S1, if we use: U = α; B = Iαβγ; a 1 = a B ; a 2 = a U ; k on = k on,1 and k off = k off,1 /C (Figs. S12, S13).
Figure S12: One cell fixed points results for full and effective model. The full model (Eq. S6) and the effective one (Eq. S1) give the same fixed points if k off = k off,1 /C (see text). In particular, the number of bound signaling receptors as a function of TCR strength m of the two models coincide when k off = k off,1 /C (black and red lines). If k off = k off,2 , the effective model (green line) yields a slightly different TCR threshold.
Moreover, since k on2 βγ >> k off2, a B , Iαβγ follows rapidly the dynamics of Iα, Iαβγ(t) ≈ Iα(t)·C. Thus, as long as C >> 1, Eq. S1 describes, to a good approximation, the interaction between two cells as given by the full receptor model (Fig. S12).

Interaction between cells with the same m value
When the cells are activated at the same time, interaction between two cells that are slightly below commitment threshold leads to their commitment, akin to quorum sensing (Fig. S14A). When the two cells are activated at different times and the m value of both cells is slightly below threshold, the later activated T-cell become committed while using the IL-2 secreted by the earlier activated T-cell which remains uncommitted (Fig. S14B).
If both are above threshold, there is a competition between the T-cells that could lead to exclusion of the earlier cell, if the m values are slightly above threshold, or the later cell, for intermediate m values (Fig. S14B).

Interaction between regulatory T-cells (Tregs) and effector T cells
We model Treg dynamics with equations similar to equations S1. In our model we assume two differences between an effector T cell and a Treg: 1) the latter does not secrete IL-2 and 2) has a higher constitutive IL-2R synthesis rate. In our calculation we assumed that the constitutive IL-2R synthesis rate of Tregs synthesis rate is 5 times higher than this of effector cells. We further assume that the non-linear positive feedback is the same in both cell types.
The phase-space map describing the interaction between an effector cell and a regulatory cell is shown in Figure  S15. The presence of the Treg strongly reduces the number of bound receptors on the effector cell in a large range of m values, when the two cells are activated at the same time (Fig. S15A). However, this effect diminishes if the Treg is activated after the effector cell. In this case, full inhibition occurs only when the effector cell is not strongly activated (Fig. S15B). It is insightful in this context to examine the result of an interaction between an effector cell and another effector cell that does not secretes IL-2 ( Fig S16). As in the case of Tregs, the number bound IL-2R is reduced due to the interaction and this reduction is less significant if the effector cell arrives before its competitor.

Stochasticity in the number of initial IL-2R
Since the number of initial IL-2R is small (around 100, see table S1) and thus prone to fluctuations, in this section we investigate the effect of stochasticity in the initial number of IL-2R on the interaction between the cells. For each iteration in the simulation, the number of initial receptors is drawn out of a normal distribution such that the mean is g c /a U = 100 and the standard deviation is 50. If this number is negative it is replaced with zero. Figure  S18 shows the result for the three representative points from figure 3 in the MS. In the case of cooperation, points I (Fig. S18A) and II (Fig. S18B), the effect of the fluctuations is negligible. However, in the case of exclusion, point III (Fig. S18C), the result is a bimodal destruction where exclusion happens only in 70% of the interactions. Figure S18: The distribution of B when the initial number of IL-2R is noisy. The interaction shown here are the same as the 3 cases shown in figure 3A in the main text, but with random initial levels of U for the two cells. Levels of U were drawn out of a Gaussian distribution with a mean of 100 and SD of 50. Blue denotes the result when the cell is alone and red when the cells are interacting. In the case of cooperation, (A) and (B) the outcome of the interaction is not sensitive to fluctuations in the initial IL-2R. In the case of exclusion (C) a bi-modal distribution emerges and only on about 70% of the interactions exclusion happens.
A B C 13. The effect of interaction on the total number of IL-2R Figure S19: The normalized number of total IL-2R (B + U), which is defined in a similar way to the interaction index C (as defined in Eq. 5 in the main text). In the case of strong cooperation (point II in Figure 3 in the MS) and exclusion (point III in figure 3 in the MS), the relative change in the total number of IL-2R is similar to the relative change in the number of bound IL-2R. In the case of mildcooperation (point I in figure 3A in the MS) the number of bound receptor increases while the number of total receptors decreases.