Fractional response analysis reveals logarithmic cytokine responses in cellular populations

Although we can now measure single-cell signaling responses with multivariate, high-throughput techniques our ability to interpret such measurements is still limited. Even interpretation of dose–response based on single-cell data is not straightforward: signaling responses can differ significantly between cells, encompass multiple signaling effectors, and have dynamic character. Here, we use probabilistic modeling and information-theory to introduce fractional response analysis (FRA), which quantifies changes in fractions of cells with given response levels. FRA can be universally performed for heterogeneous, multivariate, and dynamic measurements and, as we demonstrate, quantifies otherwise hidden patterns in single-cell data. In particular, we show that fractional responses to type I interferon in human peripheral blood mononuclear cells are very similar across different cell types, despite significant differences in mean or median responses and degrees of cell-to-cell heterogeneity. Further, we demonstrate that fractional responses to cytokines scale linearly with the log of the cytokine dose, which uncovers that heterogeneous cellular populations are sensitive to fold-changes in the dose, as opposed to additive changes.

, and has associated Shannon entropy As the above entropy describes uncertainty regarding input value after observing a specific output value, !, and is not representative of all possible values of the output. Averaging ! !|! over all output values, ! , which are distributed as gives the average entropy of the input after observing the output The difference between initial, a priori, entropy of the input, !(!), and average entropy of the input given the output, a posteriori, H(X|Y), quantifies the information gain that defines the Shannon information ! !, ! = ! ! − ! !|! .

6/38
Shannon information depends on the input distribution, !(!). Considering an input distribution that maximizes mutual information, referred to as optimal input distribution, ! * (!), leads to maximal mutual information known as information capacity Information capacity was defined by Shannon 5 to quantify log ! of the maximal number of discrete symbols that can be transferred in a single transmission with a negligible error when messages are encoded in terms long sequences of discrete symbols, through a communication channel described by the probability distribution !(!|!). Such interpretation of Shannon capacity is known as Shannon coding theorem 5 .
Shannon information gained widespread applicability in engineering and has also been adapted to quantify the amount of information transfer in cellular signaling pathways. It can be debated, however, whether it is the most suitable measure for quantification of information flow along cellular signaling relays. For instance, Shannon information depends on the input distribution, !(!), which in most applications is unknown, and a rationale for a choice of a particular distribution is missing. Besides, maximization with respect to the input distribution for the calculation of information capacity raises the question of whether cellular signaling actually operates under optimal conditions. Last but not least, cellular signaling systems do not transfer information with messages encoded in terms of long sequences of discreet symbols, which makes the interpretation of Shannon information difficult to formulate under the circumstances in which cellular signaling operates.

Advantages of Rényi min-information over Shannon information for cellular signaling
The Rényi min-capacity, ! !"# * , does not suffer from the above drawbacks. It does not depend on the input distributions, and it does not involve long sequences of discrete symbols for rigorous interpretation. As explained before, (i) 2 ! !"# * is the sum of fractions of cells simulated with a given dose with responses typical (most likely) for this dose, and (ii) 2 ! !"# * /! is that the probability of a cell generating a response typical for a dose it was exposed to. If an external observer is decoding the input signals using the maximum likelihood principle, than 2 ! !"# * is the sum of fractions of cells that can be decoded correctly, and 2 ! !"# * /m is the overall probability of correct decoding (see Supplementary Table 2 7/38 for the summary). In addition, Rényi min-capacity offers the decomposition into fractions of cells with different response levels, which can be presented along with FRC and as pie-charts.

Quantitative relationship between Rényi min-information capacity and Shannon capacity
The relationship between Shannon capacity and Rényi min-capacity is not well established under general conditions. Nevertheless, it has been shown that 1,4 ! min * ≥ ! * .
The power of 2 of the Shannon information capacity, 2 ! * , ranged from 1.37 and 2.37 being smaller by 0.75 on average from Rényi measure 2 ! !"# * . The question, which quantification more adequately represents signaling performance of the studied systems cannot be debated without reference to the axiomatic definition of Rényi entropies. Nonetheless, in contrast to Shannon, Rényi information is straightforward to interpret, which is one of its main strengths.

Discrimination between binary and graded response
In order to demonstrate how FRA can discriminate between graded and binary responses, we considered a synthetic model in which the same mean response results from different response distributions in cellular populations. Denote the response of a cell in the cellular population as Z, log of the response as Y=log 10 (Z), and stimulation level as x. Assume that the mean log response, µ=E(Y), is described as where ! !"" denotes the mean baseline log-response, i.e., the off-state mean for x=0, and ! !" the maximal log-response, i.e., the on-state mean for saturating x, Supplementary Fig. 15a. The same mean can arise in different modes of stochastic responses, i.e., in graded or binary responses, Supplementary   Fig. 15b,c. In the graded mode, the responses are centered around the mean, which assuming normal distribution of log-responses can be written as where ! ! denotes the log-response variance, Supplementary Fig. 15b. In the binary case, the logresponses can be assumed to center around the basal level, ! !"" , or around the fully induced level, ! !" , which by further assumption of log-response normality gives two possible response distributions If the probability of a cell generating a response following the first distribution is given as ! !!! ! and the second as ! ! !!! ! then the overall probability of responses is given as

12/38
which is shown in Supplementary Fig. 15c and has the mean given by Eq. 6. Having constructed the model, we considered six representative doses, x=0, ½, 1, 1 ½, 2, and performed FRA, Supplementary

Mean responses mask IFN-γ sensitivity in populations of cancer cells
In order to further demonstrate how FRA can uncover patterns masked by population averages, we corresponding to the phosphorylated form of STAT1 was measured. Population distributions are shown in Supplementary Fig. 16a, whereas means of responses (red line) and means of log-responses (black lines) are shown in Supplementary Fig. 16b. Both variants of the mean indicate that A549 cells exhibit stronger responses than CALU1 over the whole range of considered concentrations. Nevertheless, an inspection of response distributions indicates that A549 cells exhibit higher cell-to-cell heterogeneity, compared to CALU1, which leads to stronger overlaps between distributions corresponding to different doses. Therefore, even though the mean responses shift more strongly for A549 than CALU1, it does not imply that the responses distributions are more distinct. It is possible that for increasing doses, a smaller fraction of cells exhibits different responses for A549 than for CALU1 due to higher cell-to-cell heterogeneity. We have therefore performed FRA for both cell lines, Supplementary Fig. 16c,d. As opposed to analysis of means, the FRC increases more rapidly for CALU1 than for A549 cells. The maximal value for CALU1 cells is achieved for lower doses. However, both cell lines reach the same maximal value. The more rapid increase of FRC for CALU1 reflects the fact that even though the increments in means due to dose increase are smaller than for A549 cells, the lower cell-to-cell heterogeneity makes the response distributions more distinct. In terms of fractions of cells that exhibit different responses due to dose change CALU1 is more sensitive to IFN-γ than A549, which is masked by the analysis of population averages. Accounting for cell-to-cell heterogeneity is, therefore, an essential component of sensitivity of cellular populations, which is implemented within FRA.

Sufficient data size is required for accurate estimation
In order to examine how estimation of FRA depends on the data size and data dimensionality, we have considered a simple test model, previously explored in 10 , involving two possible input values, X ∈ {x 1 , x 2 }, and an output Y |X = x i described by a d-dimensional normal distribution with a diagonal covariance matrix, and mean vectors that differ only in the first dimension. Precisely, To account for non-normal distributions, we have also considered the exponent of the above model, which resulted in log-normally distributed data. Further, we have considered three different dimensionalities, d=1, 10, 100, 1000, and performed estimation of the cumulative fraction of cells, , for different sample sizes, N, for the normal and log-normal data, Supplementary Fig. 17a  Nonetheless, FRA depends on the experimental data design, which should be taken into account when interpreting its results. where rows correspond to individual cells. The first column "input" contains stimuli levels, ! ! , whereas subsequent columns represent multivariate responses of individual cells, ! ! !,! , where l varies from 1 to ! ! ,i.e., the number of cells measured for ! ! , and k varies from 1 to d, i.e., the number of response covariates. Once the data-frame is available the model for computing FRA needs to defined, which in the simplest instance is done as follows model=FRA::FRA( data= data, signal="input", response=c("output_1", "output_2", …,"output_d") )

16/38
FRC can be then plotted by calling

Equivalency of Eq. 4 in the main paper and Eq. 11
More formally, the relationship between FRA and ! min * can be derived in the following way. The fraction of cells stimulated with the dose k that can be decoded correctly is the fraction of the probability density, P(y|x k ), with highest values among all doses x 1 ,...,x i , where ! ! is the set of responses y typical for the dose k Therefore, the sum of fractions of cells that can be decoded correctly among doses x 1 ,...,x i , gives the value of the FRC which demonstrates that definitions of FRC given by Eq. 4 and Eq. 11 are equivalent.

Cancer cell lines methods (corresponding to section 3.2 of Supplementary Note 3)
Human  Probability of a cell in the population stimulated with a given does to generate a response that is typical (most frequent) for the encountered dose.
Probability of a cell in the population stimulated with a given does to generate a response that is typical (most frequent) for a given dose other than the encountered one.

Decoding by an external observer
using maximum likelihood principle Probability of correctly discriminating the dose xi from lower doses.
Sum of probabilities of correct decoding of the doses not greater than xi,. This probability divided by i gives the probability of correctly decoding doses not greater than xi.
The overall sum of probabilities of correct decoding. This probability divided by the number of doses, m, gives the overall probability of correct decoding.
Probability of correctly decoding a given dose.
Probability of confusing one dose with another dose.    Lineage markers CD3 and CD19 were firstly used to define CD3+ T cells and CD19+ B cells from the live, single, CD45+ cells.