Investigating Functional Roles for Positive Feedback and Cellular Heterogeneity in the Type I Interferon Response to Viral Infection

Secretion of type I interferons (IFN) by infected cells mediates protection against many viruses, but prolonged or excessive type I IFN secretion can lead to immune pathology. A proper type I IFN response must therefore maintain a balance between protection and excessive IFN secretion. It has been widely noted that the type I IFN response is driven by positive feedback and is heterogeneous, with only a fraction of infected cells upregulating IFN expression even in clonal cell lines, but the functional roles of feedback and heterogeneity in balancing protection and excessive IFN secretion are not clear. To investigate the functional roles for feedback and heterogeneity, we constructed a mathematical model coupling IFN and viral dynamics that extends existing mathematical models by accounting for feedback and heterogeneity. We fit our model to five existing datasets, reflecting different experimental systems. Fitting across datasets allowed us to compare the IFN response across the systems and suggested different signatures of feedback and heterogeneity in the different systems. Further, through numerical experiments, we generated hypotheses of functional roles for IFN feedback and heterogeneity consistent with our mathematical model. We hypothesize an inherent tradeoff in the IFN response: a positive feedback loop prevents excessive IFN secretion, but also makes the IFN response vulnerable to viral antagonism. We hypothesize that cellular heterogeneity of the IFN response functions to protect the feedback loop from viral antagonism. Verification of our hypotheses will require further experimental studies. Our work provides a basis for analyzing the type I IFN response across systems.


Introduction
Type I interferons (IFN) play a key role in the early immune response to viral infection [1]. The type I IFN response can be divided into two components: IFN induction, involving sensing of virus and consequent transcription and secretion of IFNs, and IFN signaling, involving the binding of extracellular IFN to type I IFN receptors and consequent upregulation of antiviral genes. Type I IFN induction begins when viral molecules are sensed by membrane-bound or cytosolic proteins known as pattern recognition receptors (PRRs). Binding of PRRs to their cognate viral molecules activates cellular pathways, leading to transcription of type I IFN genes and secretion of type I IFNs [2][3][4][5]. Type I IFN signaling begins when extracellular IFN binds to type I IFN receptors, which are found on the surface of almost all cells, leading to activation of the JAK/STAT signaling pathway and upregulation of interferon-stimulated genes (ISGs), a collection of hundreds of genes, many with antiviral properties [6,7]. IFN signaling can occur in an autocrine or paracrine manner, depending on whether a particular IFN protein binds to receptors of the cell that secreted it or other cells, respectively. Type I IFNs are potent and crucial viral restrictors. IFN-mediated ISGs restrict many viruses [8,9], and IFN deficiencies in mice and humans significantly compromise immune response to viral infections [10]. However, in counterpoint to its key role in restricting viral infection, the type I IFN response is associated with pathology, including auto-immune disease, neuropathology and tissue damage [11][12][13]. The potential harmful effects of type I IFNs are underscored by the existence of extensive negative regulation of IFN pathways and the pathologies that arise when negative regulation is constitutively absent or experimentally blocked [14,15]. A proper IFN response must therefore achieve a balance between protection and excessive IFN secretion. Further, this balance must be achieved in the face of different viruses, which exhibit a range of viral dynamics and a variety of IFN antagonism strategies. For example, HIV-1 has an estimated R 0 of eight during the first days of infection [16], while New Castle disease virus (NDV) in human cell lines can infect cells, but replicates poorly [17]. Similarly, many viruses encode proteins that antagonize IFN induction and IFN signaling, but the degree and method of IFN antagonism vary greatly across viruses [1,18,19].
The type I IFN response has components that are virus and cellular specific; for example, many ISGs have effector function against particular viruses [9], and cells vary in the pathways that induce IFN [20]. Nevertheless, the IFN response to viral infection has a general form: viruses infect cells; infected cells secrete IFN; and extracellular IFN leads to upregulation of ISGs and protection of cells. At the intercellular level, the functionality of the IFN response is determined by the number of infected cells, the amount of IFN secreted by infected cells and the sensitivity of cellular ISG upregulation to extracellular IFN concentrations. How these factors come together to balance protection and excessive IFN secretion across different viral infections is not well understood.
In particular, two features of the type I IFN response have been noted across many experimental studies: IFN secretion rates are driven by a positive feedback loop and IFN induction is heterogeneous across infected cells. The positive feedback loop is driven by an interplay between IFN induction and signaling. Many of the PRRs and signaling proteins involved in the induction and signaling pathways are themselves ISGs. For example, the proteins IRF7, RIG-I and STAT1 play central roles in the type I IFN pathway and are all ISGs. In fibroblasts and epithelial cells, the type I IFN response typically begins with induction and secretion of IFN-beta, leading to IFN signaling through JAK/STAT, which raises the expression of proteins such as RIG-I, IRF7 and STAT1, which in turn leads to higher levels of IFN-beta induction, as well as IFN-alpha induction [21]. A similar positive feedback loop that is IRF7 dependent exists in dendritic cells [20].
Heterogeneity in cellular IFN expression has been noted in multiple experimental studies, with only a small fraction of infected cells upregulating IFN expression [22][23][24][25][26][27][28][29][30]. In [29], Zhao et al. considered in vitro Sendai virus infection of human and murine fibroblasts. Although almost all fibroblasts were infected, IFN-beta expression was raised above constitutive levels in less than 20% of infected cells. Zhao et al. showed that raising the expression level of RIG-I, TRIM25 and IRF3/7, proteins in the IFN induction pathway, significantly increased the percentage of infected cells that upregulated IFN-beta expression. In [25], Hu et al. demonstrated a similar result using dendritic cells exposed to New Castle disease virus (NDV). Hu et al. found that RIG-I levels were constitutively raised in the subset of dendritic cells that expressed IFN-beta early in infection. More recently, Shalek et al. [28] considered murine dendritic cells exposed to double-stranded RNA and LPS, both of which can stimulate IFN expression. While all cells sensed and altered their gene expression due to RNA or LPS stimulation, only 8% of cells initially showed raised IFN-beta expression.
While positive feedback and heterogeneity in the IFN response have been broadly demonstrated, what functionality these two features have is not clear. A positive feedback loop is thought to restrict a strong IFN response to only significant infections, but what implications does that have on IFN-mediated protection across all infections? What tradeoffs are implicit in the presence of a feedback loop? Previous groups have suggested molecular mechanisms potentially mediating IFN heterogeneity [25,27,[31][32][33], but the functional effect of heterogeneity on the type I IFN response has not been studied.
Understanding the functional roles of feedback and heterogeneity benefits from quantitative models of the IFN response. Many authors have constructed intercellular mathematical models of IFN-viral dynamics, allowing for quantitative investigation of the IFN response [34][35][36][37][38][39][40][41][42][43][44][45][46]. All these models, in some way, compartmentalize cells into target cells, infected cells and refractory cells, with refractory cells representing cells that resist viral infection due to ISG expression. Infected and refractory cells are sometimes further split into different compartments, representing different stages of viral infection and innate response, with the exact compartmentalization varying across models. For example, Saenz et al. [46] split infected cells into three compartments-infected cells that have yet to become viral productive, infected cells that have halted infection through ISG upregulation and infected cells that are virally productive-and they split refractory cells into two compartments-cells that are in an initial state of ISG upregulation and can still be infected and cells that are completely resistant to infection (see Figure 2 in Saenz et al.). In contrast, Rand et al. [40] split infected and refractory cells into eight compartmentA (see Figure 6A in Rand et al.). Intercellular models typically also include compartments representing extracellular viral load and IFN concentrations, with interactions between the cellular, viral and IFN compartments described through differential equations.
Notably In order to investigate the functional roles of positive feedback and heterogeneity in the IFN response, we constructed an intercellular mathematical model of IFN-viral dynamics, similar in structure to existing models, but with the addition of parameters to describe IFN feedback and heterogeneity. We fit the model across five datasets, comparing the fit of the model with and without IFN feedback and estimating the heterogeneity of the IFN response.
Our fits across datasets provided us with a biologically-relevant parameter range over which to investigate IFN-viral dynamics. The parameter range modeled a range of viral infections and IFN secretion rates, allowing us to investigate possible functional roles for a feedback loop and heterogeneity. This portion of our work was based on numerical simulation, and experimental data only entered in defining a proper range for parameter values. Using numerical simulation, we were able to generate hypotheses of functional roles for IFN feedback and heterogeneity. Verification of our hypotheses will require further experimental studies.

Model
We constructed a differential equation model of viral-IFN dynamics. The model is presented in detail in the Results section, see model 2 below. Parameters of the model are given in Table 2.

Fitting to Data
We fit model 2 to datasets presented in Rand et al. [40], Patil et al. [47], Saenz et al. [46] and Schmid et al. [41]. For each dataset, we fixed some parameters to dataset authors' values and fit all other parameters. Fitted parameters are shown in Tables 3-6 below. The fixed parameters are given in the Supplementary Information in Tables S1-S4. All parameters that were fixed were set to the equivalent parameters specified by the dataset authors.
To fit datasets, we used a scaled sum of squares as the negative log likelihood.
where D represents the number of data types we considered from the given dataset (e.g., from the Rand et al. dataset, we considered data describing (1) the fraction of cells upregulating IFN-beta mRNA, (2) the fraction of cells upregulating IRF7 mRNA and (3) the extracellular IFN-beta level (U/mL), so D = 3), θ is a vector of the parameters we are fitting (e.g., θ = (q, β, φ, h, r)),d j,i is sample i from the j-th data type, d j.i (θ) is our predicted value of d j,i based on numerical simulations of the extended model and max i (d j,i ) is the largest sample value in the j-th data type. The 0.1 factor in the denominator reflects an assumption that the standard error is 10% on the scale of the data. We follow Schmid et al., who made the same assumption in constructing their likelihood. We minimized the negative log likelihood over the fitted parameters using the least_squares function from the Python scipy.optimize package under the default settings. The least_squares function requires initial values for the fitted parameters, off of which the optimization starts. Whenever possible, we used author dataset values as starting values, but when no such values were available, we selected values that gave a reasonable initial fit. We then called least_squares 30 times, perturbing the initial value of each fitted parameter by 10% in an effort to find a global minimum.

Construction of Confidence Intervals and Identifiability
We constructed confidence intervals using likelihood profiles as described in [48]. The likelihood profile for a given parameter is the minimum of the negative log likelihood as a function of the given parameter with all other parameters being allowed to vary: where LP(θ i ) is the likelihood profile and θ i is the i-th coordinate in the θ vector, i.e., the i-th fitted parameter. In practice, for each parameter, we selected a grid of values centered around its fitted value and computed the likelihood curve at those grid points. As an example, likelihood curves for each fitted parameter in the Rand  Following [48], we assumed: where θ * i is the fitted (optimal) value of θ i and χ 2 (1) is the chi-squared distribution with 1 degree of freedom. We then calculated 95% confidence intervals as: where w 0.95 is the 0.95 quantile of χ 2 (1), which is roughly 3.85. In [49], Raue et al. distinguish between three types of parameter identifiability: structural non-identifiability, practical non-identifiability and identifiability. In terms of our computations, structural non-identifiability corresponds to likelihood profiles with no global minimum and represents parameters that cannot be identified even with less noise in the sampling values. Practical non-identifiability corresponds to likelihood profiles with a global minimum, but for which the confidence interval is unbounded because the likelihood curve does not exceed 3.85 and represents parameters that can be precisely identified with less noise in the sampling values. Finally, identifiable parameters have global minimums and finite confidence intervals. The distinction between practical non-identifiable and identifiable parameters depends on the 0.1 standard error that we built into our log likelihood above. Lowering that value would convert the practical non-identifiable parameter to identifiable parameters. Table 1 shows the identifiability of each parameter for each dataset.

Numerical Simulations
Numerical simulations of Model 2 were performed using the odeint function from the Python package PyDDE. For the dataset fits, we set T(0) = 1 and all other initial values to 0, except for V(0), which was either fit or set according to dataset author values, as specified above. For the numerical experiments, we set T(0) = 0.99 and E 1 (0) = 0.01 with all other initial conditions starting at 0, reflecting an infection at MOI of 0.01.

Results
As a base model of IFN-viral dynamics, we considered a mathematical model similar to the model of Saenz et al. [46]. Figure 1a shows our baseline model. We modeled uninfected cells through two compartments: target cells, denoted by T, which can be infected, and refractory cells, denoted by R, which are resistant to viral replication. In many viral infections, when target cells become infected, an initial struggle occurs between the intracellular immune response and the virus for control of the cell [40,[50][51][52][53]. For simplicity, we assumed that either the immune system wins and the cell restricts viral replication or the virus wins and the cell becomes virally productive. Accordingly, we modeled infected cells through three compartments: eclipse cells, denoted by E 1 , which are infected cells whose outcome has yet to be determined; effector cells, denoted by E 2 , which are cells that were infected, but have entered an antiviral state that restricts viral replication; and viral cells, denoted by I, which are virally productive. Viral cells release virions into the extracellular space, and virions can infect target cells, transforming them into eclipse cells (E 1 ). We let V be the virion concentration level.
We let F be the extracellular IFN concentration. Uninfected cells express and secrete low levels of constitutive IFN, but only infected cells upregulate IFN expression through the IFN induction pathway [2]. We ignored constitutive IFN expression in our model and assumed that IFN is secreted only by infected cells. We assumed that infected cells secrete IFN at a fixed rate, but we allowed the secretion rate to vary depending on the infected cell compartment, i.e., E 1 , E 2 and I cells. IFN has a fast diffusion rate relative to other infection dynamics [54], so many authors have assumed that all secreted IFN immediately diffuses across the extracellular region, e.g., [36,40,41]. Importantly, we modified this assumption by assuming that some secreted IFN is reabsorbed by the secreting cell, a dynamic seen in experiments and theoretical analysis [55][56][57][58][59]. To model this effect, we assumed that for every unit of IFN that is secreted and diffuses into the extracellular region, a units of IFN act in an autocrine fashion. Given these assumptions, we constructed Model 1, which serves as a base model and is defined by the following system of equations.
Viruses infect target cells (T) and transform them into eclipse cells (E 1 ) at rate β. Extracellular IFN transforms target cells into refractory cells (R) at rate φ. Eclipse cells are transformed into viral cells (I) at rate k and are transformed into effector cells (E 2 ) at rate φF, due to paracrine signaling, and at rate φarq, due to autocrine signaling. Effector, eclipse and viral cells secrete IFN at rates q, rq and (1 − σ)q, respectively. The parameter r models lower secretion rates in eclipse cells relative to effector cells due to either delayed IFN induction or viral antagonism. σ models the level of viral antagonism in viral cells on a scale of 0-1, with σ = 1 representing complete antagonism. Effector cells die or stop producing IFN, and viral cells die or stop producing virions at rate δ. Viral cells produce virions at rate p. Virions and IFN are cleared from the extracellular space at rates c and d, respectively. Parameter definitions and units are shown in Table 2. Most of our notation is taken from Saenz et al. (compare their Figure 2 to our Figure 1a). Importantly, our base model does not account for IFN feedback or heterogeneity.
We extended our base model to account for a positive IFN feedback loop and heterogeneity. To model a feedback loop, we replaced the parameter q, which gives the IFN secretion rate of E 2 cells, by a function q(F) that specifies IFN secretion rates as a function of the extracellular IFN concentration. Intuitively, a feedback loop leads to higher IFN secretion rates as the extracellular IFN concentration F rises. We investigated three forms for q(F): a constant rate model parameterized by q(F) = q 0 , an increasing rate model parameterized by q(F) = q 0 * F/(F + F) + q 1 * F/(F + F) where q 0 is less than q 1 and a pulsed rate model parameterized by q(F) = q 0 * F/(F + F) The parameterF is the extracellular IFN concentration at which q(F) = (q 0 + q 1 )/2. For the increasing rate model, the secretion rate starts at q 0 for F = 0 and rises to q 1 for F significantly greater thanF, implicitly modeling a positive IFN feedback loop in which cells secrete more IFN as the extracellular IFN levels rise. The pulsed rate model starts at q 0 for F = 0, rises to (q 0 + q 1 )/2 for F =F and then drops and stays at q 0 as F increases further, implicitly modeling a positive feedback loop for 0 ≤ F ≤F and a negative feedback loop for F ≥F. Figure 2 shows the profile of each secretion rate model for particular choices of q 0 , q 1 andF. Figure 2. Secretion rate models. Profiles for the constant (q 0 = 100), increasing (q 0 = 500, q 1 = 2000, F = 100) and the pulsed (q 0 = 200, q 1 = 2000,F = 100) secretion rate models. The models specify the cellular secretion rate (q(F)) as a function of the extracellular IFN level (F). For the increasing and pulsed model, q(F) = (q 0 + q 1 )/2. For the pulsed model, the maximum secretion rate is achieved atF.
To model a heterogeneous IFN response, we introduced the parameter h, which is the fraction of cells that can induce IFN in response to infection. We split E 1 (eclipse) cells into two compartments, E 1r and E 1n , reflecting cells that can respond to infection by inducing IFN and those that cannot, respectively, and we similarly split the I (viral) cells. Putting our new parameterization of IFN secretion and heterogeneity, we arrive at the following model.
Figure 1b shows the extended model graphically, and Table 2 includes definitions and units for parameters specific to the extended model.

Fitting the Extended Model to Data
We fit our extended model, Model (2), to datasets presented in Rand et al. [40], Patil et al. [47], Saenz et al. [46] and Schmid et al. [41]. For each dataset, we fit our model assuming each of our three secretion models: constant, increasing and pulsed. Results for the increasing rate model were roughly intermediate to the other two, and here, we only present results for the constant and pulsed secretion models. Except for Patil et al., the authors fit their own mathematical intracellular model of IFN-viral dynamics to their dataset. Many of the authors' models include parameters that connected directly to our own, and in those cases, we used the authors' inferred values; however, we always obtained our own fit to the parameters describing IFN secretion rates (q 0 , q 1 ,F), IFN induction and signaling rates (r, φ), viral infectivity rates (β) and the heterogeneity level of IFN induction (h). We also fit p and σ, except for the Rand et al. and Patil et al. datasets. Those datasets involved New Castle Disease virus infection, for which p = 0 and σ = 0; see below for details.
Given the parameter values that best fit the data, we constructed a confidence interval for each fitted parameter; see the Methods for details. For Saenz et al., we could not construct meaningful confidence intervals (i.e., confidence intervals other than [0, ∞)) for almost all of the parameters, so we set β and φ to the values found by Saenz et al., after which we were able to generate meaningful confidence intervals for most of the parameters. Overall, we inferred 5-9 parameters depending on the dataset and secretion rate model.
Below, we discuss each dataset in detail. Overall, the pulsed model was a statistically better fit for the Rand Figure 3a), while the pulsed model infers a more sudden increase (see the right panel in Figure 3b), but both models fit the data well. Importantly, NDV in murine cells does not antagonize the IFN response and cannot produce infection particles, although NDV can enter murine cells and initiate viral protein synthesis [40]. Over the first 48 h of infection, Rand et al. determined the fraction of cells upregulating IFN-beta mRNA, the fraction of cells upregulating IRF7 mRNA, the extracellular IFN-beta level (U/mL) and the total viral expression level of hemagglutinin (HN) across all cells. To connect our model to the data, we assumed that the IFN-beta mRNA level was proportional to cellular IFN secretion rates and that IRF7 was upregulated in our E 2 and R cellular compartments. Since our model does not track viral protein expression and since NDV virions could not initiate second round infections, we did not fit the HN data. Figure 4 shows the fit of our model to the Rand et al. data, and Table 3 provides the fitted parameter values.  For the constant secretion model, we inferred a secretion rate of 1190 U/mL(hour), while for the pulsed secretion rate, we inferred a baseline secretion rate of 100 U/mL(hour), which can rise to 2760 U/mL(hour) when extracellular IFN levels reach 1120 U/mL, which occurred between 18 h and 24 h of infection (Rand et al. estimated a secretion rate of 1300 U/mL(hour), roughly comparable to our values). The pulsed model had a better fit that was statistically significant (F-test, p-value = 5.6 × 10 −6 ). Notably, the pulsed model is able to capture the rapid upregulation in IFN concentrations between 12 h and 24 h, while the constant model cannot; compare the far right panels of Figure 4a,b.
We inferred h = 0.45 and h = 0.66 under constant and pulsed secretion rate models, respectively, suggesting that roughly 45-66% of target cells were capable of inducing IFN upon infection. The confidence intervals for h extended to one, but were bounded below by 0.38 and 0.50, so the presence of heterogeneity was not statistically significant, but the absence of extreme heterogeneity was significant. The dataset shows that roughly 30% of cells induce IFN (see fraction IFN+ cells panel in Figure 4). The difference is explained by refractory cells. Our model predicts that roughly half of target cells become refractory. Given 60% of target cells capable of inducing IFN, half become refractory and half become infected, leading to the sampled 30% levels of IFN induction.

Patil et al. studied in vitro NDV infection of human, monocyte-derived, dendritic cells (DCs).
As is the case for the murine cells of Rand et al., NDV does not antagonize IFN expression in human DCs, and NDV infection can lead to viral protein synthesis, but not infectious virions. Patil et al. measured the fraction of cells upregulating IFN-beta mRNA, the fraction of cells expressing viral HN, the fraction of cells upregulating both IFN-beta mRNA and expressing viral HN and the total IFN mRNA expression level across all cells. To connect our model to the data, we assumed that cells upregulating IFN-beta expression were secreting IFN, that only our viral cells (I compartment) were expressing viral protein at levels high enough to be detected and that the total IFN mRNA level across all cells was proportional to IFN extracellular concentrations. Figure 5 shows the fit of our models to the Patil et al. data, and Table 4 provides the fitted parameter values. For both the constant and pulsed models, we could not fit the β parameter. In both cases, fits improved as β increased, suggesting a model in which all cells were immediately infected. Patil et al. reported infection MOIs of roughly one and, based on the presence of viral HN protein, suggested that roughly 75% of cells were infected. Our modeling suggests that all cells were infected, but roughly 25% of cells were either still in the eclipse phase (E 1 cells) or had become effector cells (E 2 cells). Under the constant secretion model, we inferred a secretion rate of 7.3 mRNA/mL(hour), while under the pulsed secretion rate model, we inferred a baseline secretion rate of roughly 9.6 mRNA/mL(hour). The best fits for the parameters q 1 andF, 0 and 4.5 mRNA/mL(hour), respectively, suggest downregulation of secretion rates from the baseline secretion rate of 9.6 mRNA/mL(hour). However, our confidence intervals for q 1 andF are wide (indeed, forF, the confidence interval is unbounded). When we fixed q 1 to values significantly exceeding q 0 , the best fits pushedF to high values, so that IFN secretion rates were never upregulated in the model. Overall, the data is not consistent with a positive feedback loop for IFN secretion rates. We speculate that the short duration of the experiment, 10 h, may not have been long enough to capture positive feedback.
For both constant and pulsed secretion models, we inferred heterogeneity values of roughly 50% with confidence intervals restricted to roughly 45%-55%, matching the results of Patil et al. The pulsed model had a better fit, but the fit was not statistically significant (F-test, p-value = 0.16). In contrast to the other datasets, the Saenz et al. dataset is more limited, including only extracellular IFN and virus concentrations and lacking cellular compartment frequencies. As a result, we were not able to identify the parameters β and φ for either the constant or pulsed models, meaning that different values of β and φ achieved equivalent fits (see the Methods for details). In order to fit our extended model, we therefore used the β and φ values inferred by Saenz et al. (see their Table 2). With β and φ fixed, we were able to identify the constant rate model parameters, but we still could not identify several of the pulsed rate model parameters. Below, we present the best fit of the data under the pulsed rate model; however, this fit serves simply as a hypothesis consistent with the data, and we cannot assign any statistical certainty to it. Figure 3 shows the fit of our models to the data, and Table 5 provides the fitted parameter values and their confidence intervals (note that several of the pulsed rate model parameters have unbounded confidence intervals, reflecting our inability to identify the parameters). For the constant secretion rate model, we inferred a secretion rate of 4.7 mRNA/mL(hour). For the pulsed secretion rate model, we inferred a low baseline secretion rate of 0.05, but the secretion rate could be upregulated to roughly 12.5 mRNA/mL(hour) when extracellular IFN levels reached 3.6 mRNA/mL, which occurred at Day 3 of infection. We inferred a heterogeneity level of h = 1 for both the constant and pulsed models, meaning that all cells could induce IFN.
Notably, under both the constant and pulsed secretion rate model, we inferred a viral IFN antagonism level of σ = 1, meaning complete antagonism of IFN secretion in viral cells (I cells). We speculate that our estimates of σ reflect the known antagonism properties of the NS1 protein of influenza [61]. The pulsed model had a better fit, but the better fit was not statistically significant (F-test, p-value = 0.71).
The parameter r, which specifies the IFN secretion rate in eclipse cells relative to effector cells, differed significantly between the two models: 0.00467 and 0.44 for the constant and pulsed models, respectively. Although the difference in r estimates is not statistically significant, as seen by the overlap in confidence intervals in Table 5, we do have statistical support that r is near zero for the constant model, while we cannot reject r = 0.46 for the pulsed model. The difference in these two r estimates provides an instructive example of how the two models of secretion rates influence parameter estimates. We find that under the constant model, r must be nearly zero to keep IFN concentrations low during the first two days of infection, when many cells are in the eclipse phase. In contrast, under the pulsed model, r can be relatively large, i.e., 0.46, because the baseline secretion rate q 0 is small, allowing IFN concentrations to remain low during the first days of infection, but still rise quickly at later times due to upregulation of IFN secretion rates.  Table 6 gives the fitted parameter values and confidence intervals. Figure 6 shows the fit of our models to a portion of the   For both the WT and MT datasets, our inferred baseline IFN secretion rate for the constant rate model was roughly double the pulsed rate model baseline secretion rate (WT: 95 pg/mL(hour) for constant and 58 for pulsed, MT: 470 pg/mL(hour) for constant and 250 for pulsed). Although these differences did not rise to statistical significance, the resultant differences in the predicted IFN dynamics are instructive. Figure 6 shows the predicted IFN concentration and secretion rate dynamics under the constant model (Subplot (a)) and the pulsed model (Subplot (b)) for the MT dataset. The sampled IFN concentrations are also shown. The pulsed model suggests the presence of strong IFN feedback between 48 h and 72 h, when the secretion rates rise from the baseline of 250 pg/mL(hour) to 1500, followed by downregulation back to baseline between hours 72 and 96. The sampled IFN concentrations show an increase between hours 48 and 72 followed by a leveling between hours 72 and 96, dynamics which are well fit by the pulsed rate model. In contrast, the constant rate model cannot capture the leveling between hours 72 and 96. The fit of the pulsed rate model to the full WT and MT datasets was significantly better than the constant rate model (F-test, p-value: 0.03 (WT) and 7 × 10 −9 (MT)).
Inferred heterogeneity levels suggested a roughly homogeneous response across both WT and MT datasets and both secretion rate models, although our confidence intervals were too broad to provide statistical significance. In contrast, the parameter r, which sets the IFN secretion in eclipse cells (E 1 ) relative to effector cells (E 2 ), differed across the WT and MT datasets. For the WT dataset, r was essentially zero for both constant and pulsed rate models, while for the MT dataset, r was roughly one for both rate models, with confidence intervals that did not include zero. Our r estimates are in line with the results of Schmid et al., suggesting that the MT strain leads to higher levels of IFN transcription at early stages of intracellular infection. To circumvent this difficulty, we expressed IFN in units defined so that a single unit of IFN transforms target cells into refractory cells at the rate of one per hour. Put another way, one unit of IFN reduces the frequency of target cells by one log (natural) every hour. Similarly, we expressed viral load in units defined so that a single unit of virions reduces target cell frequency at a rate of one per hour. Biologically, measuring in these units means that we are measuring IFN and virions in terms of their effector function on the given cell types. In contrast, measuring IFN using, for example, international units (IU) corresponds to measuring IFN in terms of its effector function on a standardized cell type and virus [63]. Measuring IFN by density, e.g., pg/mL, standardizes the amount of IFN secreted, but not its effect. With this in mind, we refer to our units as IFN effective units (iEU) and viral effective units (vEU).

Cross Dataset Comparisons
Mathematically, measuring in effective units leads to β = 1 and φ = 1, and the only other parameters affected are p, which is replaced by p × β, and q 0 , q 1 andF, which are replaced by q 0 × φ, q 1 × φ andF × φ. We can then use p as a measure of viral infectivity and the q parameters as measures of IFN response strength. Schmid et al. and Rand et al. measured the rate target cells become refractory when exposed to IFN, so for those datasets, we can connect our iEU to standard units. For Rand et al., we estimated φ as 2.8 × 10 −4 (constant model) and 0.63 × 10 −4 (pulsed model). Since 1 iEU is 1/φ U/mL., 1 iEU is roughly 4000-16,000 U/mL. Rand et al. measured that 100 U/mL of IFN-beta leads to one log drop in target cell frequency in roughly 24 h (see their Figure 5A), so that using their measurements, we would expect 1 iEU to be 2400 U/mL, at the same magnitude, but slightly below our 4000-16,000 U/mL estimate. Similarly, for the Schmid et al. dataset, we estimate φ as roughly 10 −5 , meaning that 1 iEU would equal 10 5 pg/mL, or in nanograms 100 ng/mL. Schmid et al. found that 10 ng/mL drops target cell frequency by one log after 6 h (see their Figure 7C); using their measurements, 1 EU would equal 60 ng/mL, roughly in line with our estimate of 100 ng/mL. Below, we use EU derived through our fits, but we emphasize that EU can be measured experimentally and can serve to standardize the comparison of IFN response across systems.
In effective units, comparison between datasets is possible. Table 7

Functional Roles for IFN Feedback and Heterogeneity Based on Model Simulations
Using numerical simulations, we investigated potential functional roles for IFN feedback and heterogeneity by considering a range of parameterizations of the IFN response in our extended model. As the IFN response must cope with a range of viral infections, we investigated the relationship between protection and IFN secretion of each particular parameterization of the IFN response across a range of viral parameterizations. In fitting datasets, we found that IFN secretion rates varied roughly on the order of 0.01-1.0 iEU/hour, while the infectivity rate p varied between 0.01 and 0.3 vEU/hour. We therefore considered parameterizations roughly within those ranges, although we extended the viral range slightly to 0.5 to account for more virulent viruses than IAV. Below, to place extracellular IFN levels in a more concrete context, we converted our IFN EU to U/mL based on the relation: 1 iEU equals 7500 U/mL (recall, for the Rand et al. dataset we estimated 1 iEU as 4000-16,000 U/mL; here, we pick 7500 within this range). As a touchstone against which to assess whether IFN levels in our simulations below are biologically reasonable, we note that Rand et al. measured maximum extracellular IFN concentration at roughly 2000 U/mL in their in vitro experiment.

Our Model Feedback Loop Mediates Protection without Excessive IFN Secretion
To investigate the relationship between protection and IFN secretion for the constant secretion rate model, i.e., q(F) = q 0 , we considered a range of secretion rates 0 ≤ q 0 ≤ 1 over a range of viral infectivity 0 ≤ p ≤ 0.5. Figure (2) numerically for each combination of q 0 , p. Each point within the panels represents a result from a single simulation. Protection, as quantified by the fraction of cells that become virally productive (with a lower frequency reflecting higher protection), improves across all infectivity rates p as the secretion rate q 0 rises. Functionally, raising secretion rates leads to an increase in the fraction of refractory cells, which leads to the drop in virally-productive cell frequency. Notice that maximum IFN concentrations are on the order of 1000 U/mL for the higher end infectivity rates, in line with the IFN concentrations seen in Rand et al. The relationship between protection and secretion rates changes markedly when p is very small, reflecting a weak infection. Figure 7b shows cellular frequencies and extracellular IFN levels when p = 0, an infectivity rate reflecting abortive infection in which the virus can infect cells, but not replicate, e.g., NDV in human cell lines. For all q 0 values, the fraction of virally infected cells is less than 0.02, but as q 0 rises, the frequency of refractory cells and extracellular IFN levels rise significantly, with roughly 75% of cells becoming refractory under the highest secretion rate, q 0 = 1.0. Since raising secretion rates has little functional effect on restricting infection, we view the IFN secreted and the antiviral cell frequency as measures of excessive IFN secretion. In particular, for q 0 = 1, maximum IFN concentrations rose to roughly 300 U/mL in our simulations. In their in vitro study, Rand et al. found that cells exposed to 250 U/mL upregulated their ISGs within 24 h, suggesting that 300 U/mL is a significant IFN level. Intuitively, in low infectivity rate settings, infected cells are largely restricted to cells infected due to an exposure event, and few target cells are subsequently infected. Raising q 0 raises the amount of IFN produced by these infected cells roughly linearly, leading to a roughly linear increase in the extracellular IFN level and the frequency of antiviral cells, but with little effect on the frequency of infected cells. To investigate the relationship between protection and IFN secretion rates under an increasing rate model, i.e., q(F) = q 0 * F/(F + F) + q 1 * F/(F + F), we considered q 1 = 0.01, 0.05, 0.2, 1 and fixed q 0 = 0.01. As shown in Figure 8a, raising q 1 increases protection. As with constant secretion rates, the rise in protection is mediated by an increase in extracellular IFN levels and refractory cell frequency. However, in contrast to the constant secretion rate model, when q 1 = 1, extracellular IFN levels reach levels of 15,000 U/mL, far greater than the 2000 U/mL seen in Rand et al.'s in vitro experiments and possibly biologically unrealistic. Figure 8b,c shows the results for the case p = 0. Unlike the constant secretion rate, the fraction of cells that become refractory and the extracellular IFN concentration stay low across the different q 1 . We repeated the same numerical experiments assuming a pulsed secretion rate. Infection outcomes were essentially identical as under the increasing secretion rate, except that extracellular IFN concentrations were much lower for a pulsed secretion rate. Figure 9 compares the fraction of virally-infected cells and extracellular IFN levels under increasing and pulsed secretion rate models with q 0 = 0.01 and q 1 = 1. Strikingly, the fraction of virally-infected cells is roughly identical for both models, but the pulsed secretion rate model produces much less extracellular IFN, on the order of 3000 for the highest infectivity rate p = 0.5, in line with the concentration seen by Rand et al. Figure 9. A pulsed IFN secretion rate provides the same protection as an increasing IFN secretion rate, but with lower levels of extracellular IFN. See Figure 8 and the text for details.
We sought to better understand how the pulsed rate model mediates protection while avoiding excessive IFN secretion for low infectivity rates p. Figure 10 shows the dynamics of extracellular IFN levels and the frequency of refractory and target cells across a range of infectivity rates: p = 0.02, 0.05, 0, 1, 0.5. For p = 0.01, extracellular IFN levels stay low throughout infection, reflecting non-activation of the feedback loop. In contrast, for the higher infectivity rates, extracellular IFN levels experience transitions from low levels to significant levels, reflecting activation of the feedback loop.  In summary, a constant secretion rate model is unable to provide protection from high infectivity viruses without excessive IFN secretion against low infectivity viruses. In contrast, a positive feedback loop, modeled by our pulsed secretion rate model, is able to provide protection without excessive IFN secretion. Protection is mediated by activation of the feedback loop, which leads to high levels of extracellular IFN when infectivity rates are high. Excessive IFN secretion is avoided by non-activation when infectivity rates are low.

Our Model Feedback Loop Is Vulnerable to Viral Antagonism
We modeled viral antagonism of IFN induction using the parameters r and σ. σ quantifies antagonism in I (viral) cells, with zero and one signifying no antagonism and complete antagonism, respectively. The parameter r models reduced IFN secretion by E 1 (eclipse) cells relative to E 2 (effector) cells, either due to cellular mediated delays in IFN induction or early antagonism of IFN induction by viral proteins.
Just as we investigated the relationship between protection and IFN secretion across a range of infectivity rates, we investigated the relationship between protection and IFN secretion across a range of antagonism strengths. As an example of strong antagonism, we set σ = 0.99 and r = 0.06, reflecting 99% and 80% antagonism of virally-productive and eclipse cells, respectively, as compared to our baseline values of σ = 0 and r = 0.3, which we took as modeling no antagonism. We also considered interpolations between strong and no antagonism: σ = 0.95 and r = 0.1, corresponding to intermediate antagonism, and σ = 0.5 and r = 0.2, corresponding to weak antagonism. Figure 11 shows infection outcome for the four cases assuming a pulsed secretion model with q 0 = 0.01 and q 1 = 1. Strong antagonism blocks the IFN response, with IFN levels only reaching roughly 50 U/mL for the highest infectivity of p = 0.5, leading to complete loss of protection. Intermediate antagonism leads to significant maximal IFN concentrations, but the response does not mediate significant protection. For weak and no antagonism, once p > 0.1, extracellular IFN levels are nearly maximal, reflecting full activation of the feedback loop, leading to significant levels of protection. To understand the mechanism by which antagonism mediates loss of protection under our model more precisely, we considered extracellular IFN dynamics for the specific case of p = 0.3 and a pulsed secretion rate with q 0 = 0.01 and q 1 = 1. Figure 12 shows the dynamics of extracellular IFN levels. Under strong antagonism, extracellular IFN levels stay low. In contrast, for all other levels of antagonism, the feedback loop activates. However, under intermediate antagonism, the IFN feedback loop activates relatively late and does not reach full activation levels, leading to the loss of protection discussed above. Using 50 U/mL as a mark of IFN activation (as we did above), IFN activation occurs at 20, 23   For all levels, we set p = 0.3, and we assumed a pulsed secretion rate model with q 0 = 0.1 and q 1 = 1. As the level of antagonism is increased, activation of the feedback loop is delayed, and in the case of strong antagonism, the feedback loop is blocked. The delay or blockage of activation leads to the loss of protection, as shown in Figure 11.
In summary, viral antagonism of IFN induction can delay or block feedback loop activation, leading to a loss of protection. Our increasing and pulsed secretion rate models avoid excessive extracellular IFN levels through low secretion rates in early infection, but viral antagonism can exploit the low secretion rates to prevent feedback loop activation.

Cellular Heterogeneity Can Protect Our Model Feedback Loop from Viral Antagonism
Restricting IFN induction to a fraction of cells reduces the average secretion rate per cell, thereby convolving the heterogeneity level with the average secretion rate. To modulate heterogeneity with average secretion rates fixed, we considered a scaled heterogeneity that modulates both the fraction of cells that can induce IFN and the IFN secretion rate of those cells. We assumed that a fraction h of cells can induce IFN and that those cells secrete IFN at rate q(F)/h, thereby keeping average secretion per cell fixed at q(F).
To investigate the effect of heterogeneity on protection and IFN secretion levels, we considered a pulsed secretion rate q 0 = 0.1 and q 1 = 1. Figure  Functionally, under our model, heterogeneity serves to raise the autocrine-mediated IFN signal in eclipse cells while keeping paracrine-mediated signaling in target cells roughly fixed. Eclipse cells are transformed into viral and effector cells at rate k and φ(F + arq(F)), respectively. The rate of transformation into effector cells can be further decomposed into a paracrine-mediated rate, φF, and an autocrine mediated rate φarq(F). In early infection, F ≈ 0, and transformation of eclipse cells into effector cells is autocrine mediated and occurs at rate φarq 0 . Putting all this together, in early infection, the probability that an eclipse cell becomes an effector cell rather than a viral cell is roughly φarq 0 /(k + φarq 0 ). The presence of antagonism in eclipse cells, parameterized by small values for r, lowers this probability to near zero, leading to delay or blockage of feedback loop activation. Our heterogeneity scaling replaces a q 0 secretion rate by q 0 /h, leading to a higher probability that eclipse cells transform into effector cells. However, solely scaling the secretion rate by q 0 /h leads to excessive IFN secretion. By restricting IFN induction to a fraction h of eclipse cells, heterogeneity functionally creates a small fraction of eclipse cells that are protected by the autocrine response, while maintaining low secretion rates in early infection that keep paracrine signaling in target cells low, thereby avoiding excessive IFN secretion for low infectivity viruses. Figure 14 shows infection dynamics for three responses: a scaled heterogeneous response with h = 0.1 and secretion rates q(F)/h, a homogeneous response (h = 1) and a response in which we raised secretion rates to q(F)/h with h = 0.1, but assumed that all eclipse cells can induce IFN. The raised secretion rate response serves to demonstrate the importance of scaling both the IFN secretion rate by q(F)/h and restricting IFN secretion to a fraction h of eclipse cells. We first investigated the dynamics of these three responses assuming strong antagonism and p = 0.1 (Panel a). The heterogeneous and raised responses activate the feedback loop, while the homogeneous feedback loop is blocked. Underlying this difference is autocrine-mediated protection. For the heterogeneous and raised responses, eclipse cells have a 0.20 probability of transforming into effector cells when infection begins, while for the homogeneous response, the probability is ≈0.01. As the feedback loop activates in the heterogeneous and raised responses, autocrine-mediated protection of eclipse cells rises, and eclipse cells have an ≈0.60 probability of transforming into effector cells, while autocrine-mediated protection stays low for the homogeneous response. The result is significantly increased protection relative to the homogeneous response. The raised response activates the feedback loop earlier than the heterogeneous response, leading to higher levels of extracellular IFN and greater protection. We then considered the same three responses assuming no antagonism and p = 0.01 (Panel b). In this case, the heterogeneous and homogeneous feedback loops do not activate, while the raised response feedback loop does. There is little difference in protection between the three response, but the raised response transforms almost all target cells into refractory cells, reflecting excessive IFN secretion.
(a) strong antagonism, p = 0.1 (b) no antagonism, p = 0.01 Figure 14. A scaled heterogeneous response protects eclipse cells through autocrine-mediated signaling, while not excessively raising paracrine-mediated signaling to target cells. We considered a scaled heterogeneous response, a homogeneous response and a response in which we raised secretion rates, but did not restrict the fraction of cells that could induce IFN. (a) In the presence of strong antagonism and an infectivity rate p = 0.1, the feedback loop is activated, and extracellular IFN concentration rises for the heterogeneous and raised responses, but not for the homogeneous response. Rising extracellular IFN levels are mediated by strong autocrine signaling, which leads to a high probability that eclipse cells transform into effector cells (autocrine prob). Activation of the feedback loop leads to a lower frequency of viral cells and a higher frequency of refractory cells, particularly for the raised response. (b) In the absence of antagonism and an infectivity rate p = 0.01, the heterogeneous and homogeneous feedback loops do not activate, but the feedback loop activates for the raised response. For all three responses, viral cells have low frequency, but the raised response transforms nearly all target cells into refractory cells and raises extracellular IFN levels to 2000 U/mL, reflecting excessive IFN secretion.
In summary, a scaled heterogeneity raises autocrine-mediated signaling in a fraction of cells that secrete IFN at relatively high rates. Autocrine-mediated signaling protects these cells from relatively high levels of viral antagonism, allowing the cells to transform into effector cells, which can activate the feedback loop. At the same time, the average secretion rate per cell is kept fixed, meaning that extracellular IFN levels remain low in early infection, and excessive IFN secretion does not occur for low infectivity viruses.

Discussion
To investigate potential functional roles for feedback and heterogeneity, we constructed, fit and analyzed a mathematical model coupling type I IFN and viral dynamics. Our model extends current intercellular-based models by accounting for IFN feedback and heterogeneity. We model IFN feedback by increasing cellular IFN secretion rates as extracellular IFN concentrations rise. We fit our model to several datasets involving different viruses and host cells: in vitro infection of human, epithelial carcinoma cells by mutant and wild-type DENV (Schmid et al. [41]), in vitro infection of murine fibroblasts by NDV (Rand et al. [40]), in vitro infection of human, monocyte-derived dendritic cells by NDV (Patil et al. [47]) and an in vivo infection of ponies by equine IAV (Saenz et al. [46] Our current understanding of the IFN response suffers from a lack of cross system analysis, so we sought to compare the IFN response across the fitted datasets. Extracellular IFN concentrations are typically reported as standardized densities, either using weights, e.g., nanograms per mL, or using effector function on some standard assay, e.g., IU/mL [62,63]. We were interested in comparing the strength of the IFN response across datasets. In that context, using standardized units to compare IFN secretion rates assesses which system produced the most IFN. In some contexts, this may be the appropriate measure of strength, but we wanted to measure the ability of the IFN response to mediate protection. Since different cell types respond differently to the same concentrations of IFN, suggested that peak upregulation of the feedback loop led to IFN secretion rates of 0.15 effective units per hour, which corresponds to protection of roughly 15% of target cells every hour. In contrast, our fit of DENV infection of epithelial cells suggested peak extracellular IFN secretion rates of 0.005, which corresponds to protection of roughly 1% of target cells per hour.
The notion of effective units also allowed us to combine our fits to define a biologically meaningful parameter range for IFN secretion and viral infectivity rates. Using this range, we investigated potential functional roles for feedback and heterogeneity through numerical simulations. We found that a relatively high constant secretion rate leads to excessive IFN secretion in abortive viral infections, while a relatively low constant secretion rate cannot protect against high infectivity viruses. In contrast, a feedback loop was able to avoid this tradeoff by setting a low baseline secretion rate that rises only under relatively high viral infectivity rates. Further numerical simulations showed that modulation of secretion rates by a feedback loop creates a vulnerability that can be exploited by viral antagonism. To avoid excessive IFN secretion, a feedback loop must keep secretion rates low in early infection, but this essential feature is vulnerable to viral antagonism. When secretion rates are low, viral antagonism of IFN induction can keep extracellular IFN levels from rising to levels that activate the feedback loop, leading to a delay or blockage in activation of the feedback loop. Putting all this together, based on our numerical simulations, we hypothesize that the IFN feedback loop balances protection and excessive IFN secretion, but with the price of vulnerability to viral antagonism.
We next used numerical simulations to consider a functional role for heterogeneity. We modeled a scaled heterogeneity, in which a fraction h of infected cells can induce IFN with secretion rates that are raised by a factor 1/h. The higher secretion rates in this subset of infected cells increases autocrine-mediated IFN signaling, thereby making it more likely that these cells will upregulate their ISGs and become effector cells, which secrete IFN free from viral antagonism. On the other hand, restricting IFN induction to a fraction h of cells keeps the average secretion rate per infected cell fixed, which preserves the extracellular IFN levels mediated by the feedback loop. Putting these two features together, based on our numerical simulations, we hypothesize that heterogeneity protects the feedback loop from viral antagonism through autocrine-mediated IFN signaling in the fraction of cells capable of inducing IFN.
Work by Hu et al. [25] and Shalek et al. [28] suggests that a small fraction of cells, termed initiator cells, are responsible for activating the IFN response. In particular, Shalek et al. experimentally demonstrated, in an in vitro system, that stimulation of dendritic cells by PRR agonists led to high levels of IFN secretion by a small fraction of cells. When autocrine-mediated IFN signaling was blocked, this response was lost. Our modeling suggests that the functional role of limiting IFN induction to initiator cells is to provide a protective autocrine response without raising extracellular IFN levels to excessive levels.
While we are not aware of previous modeling analysis of the effects of heterogeneity on infection outcome, several authors have suggested possible benefits associated with a heterogeneous response. Rand et al. demonstrated heterogeneity in the IFN-beta response time of infected cells and suggested that this heterogeneity may allow for a sustained IFN response [40]. Our modeling suggests that IFN mediated protection is dependent on activation of the feedback loop, rather than sustaining an IFN response, but it is possible that protection against viruses with extended intracellular times may require a sustained response. Zhao et al. suggested that heterogeneity may allow cells to achieve an optimal level of extracellular type I IFN by modulating the frequency of expressing cells [29]. Our modeling supports exactly such a view. Ivashkiv and Donlin and, more recently, Czerkies et al. suggested that heterogeneity may allow for a small fraction of infected cells to mass produce IFN and then enter apoptosis [30,64]. Such mass producers would protect other cells through paracrine signaling and prevent further viral replication through their own apoptosis. While our modeling also suggests that mass producers of IFN protect other cells through paracrine signaling, we suggest that mass producers are protected by autocrine signaling, and we do not model apoptosis. Czerkies et al. noted apoptosis in the case of fibroblasts exposed to RIG-I and MDA5 agonists, but whether this would also be the case under viral infection for mass producers is unclear. Further, the extent to which cells exposed to ongoing viral replication are able to secrete IFN is dependent on virus and cell types [65]. It may be that in some cases, mass producers are protected, while in other cases, apoptosis results.
Given the complexity of the type I IFN response and viral infection, our modeling has many limitations. We have not distinguished between the different type I IFNs, in particular IFN-β and IFN-α, although both are crucial to the IFN response and are induced with different dynamics and by different cell types [20,66]. We have not considered spatial effects. In particular, we modeled the IFN-response as either autocrine or paracrine, but in reality, IFN secreted by a cell is more likely to bind to IFN receptors of nearby cells, a phenomenon observed in experimental studies [26,56,67]. We have also not considered the interaction between IFN and other early cytokines such as IL-1, e.g., [68,69] On the other hand, our results have some general features that transcend our model. Under our model, low secretion rates in early infection leave the feedback loop vulnerable to antagonism, and heterogeneity makes the feedback loop more robust by increasing autocrine-mediated protection. Regardless of the model, if the IFN response is driven by a positive feedback loop, then there must be cases when the feedback loop is not activated in response to infection, otherwise excessive extracellular IFN levels would arise in response to abortive infections. If activation of the feedback loop is largely driven by extracellular IFN levels, then our results would likely still hold since viral antagonism can modulate extracellular IFN levels and therefore modulate feedback loop activation. However, the extent to which IFN induction and signaling within a specific cell is dependent on extracellular IFN levels is unclear. IFN independent pathways that upregulate ISGs exist [70,71]. Individual cells can alter their IFN secretion rates through autocrine-mediated IFN signaling, an effect that is not present in our model, which assumes that secretion rates are determined by a single global extracellular IFN level. Cells are also sensitive to signals of ongoing infection, for example the presence of defective interfering viral particles [72,73], which serve to activate the IFN response independent of extracellular IFN levels. Further work is needed to investigate such effects and their impact on the IFN response.