Genes I and A exhibit different responses to changes in signal duration
Rate laws were written for each of the reactions in the model and the corresponding chemical master equations were solved using the Gillespie algorithm 
. We explicitly considered stochastic fluctuations in the number of molecules in this system via the Master equations because genes are often present in small copy numbers. The solution of the mean field equations could lead to the identification of signaling thresholds that in reality cannot be met due to fluctuations. For instance, the mean field solution to the Model 1 equations obtained via MATLAB predicts that the number of I molecules never reaches 50 molecules at any point in time for the set of parameters identified in using a signal degradation rate of 0 (). However, the Gillespie algorithm shows that half of the trajectories of I cross 50; therefore it would be inaccurate to conclude that I never reaches the signaling threshold of 50 molecules ().
Table of parameters used in Model 1.
We are more interested in the general qualitative behavior of the genetic network as opposed to fitting the results to actual data, so no attempt was made to use experimentally derived parameters relevant to a specific system. As a result, it is important to use sensitivity analysis to verify that the qualitative results of the model are robust over a range of parameters. Sensitivity analysis was performed according to the procedures outlined in Web S1
. The basic procedure was to divide all the parameters into nine different classes and testing all combinations of low, medium, and high values for each parameter class (Table S1
). We then determined if the response could be switched from pro-inflammatory response to an anti-inflammatory response by changing the signal degradation rate for a given parameter set. In an initial search, 1809 parameter sets were identified that led to the desired qualitative behavior (Web S1
A snapshot of the qualitative behavior of the model can be obtained by considering two values of kdeg,T
with all other parameters fixed at the values given in . A value of 1 for kdeg,T
corresponds to a transient signal, whereas a value of 0 for kdeg,T
corresponds to a sustained signal. Sensitivity analysis revealed that the mRNA degradation rate did not alter the results (Figure S1
), so the rate of degradation for I and A mRNA was set to 0. We consider the case where the rate of transcription of I is the same whether its induced by pSTAT3 or by itself (i.e. ktrx,I,T
) and where the transcription of A is slower than the transcription of I (i.e. ktrx,A,T
). The qualitative results are robust to a range of parameters including the case where the rate of transcription of I and A is the same and the case where the rate of transcription of A is greater than I as determined by parameter sensitivity analysis (Web S1
For sustained signals, the number of I molecules peaks at early times and then degrades, while the number of A molecules increases steadily through time ( and ). There are no negative feedback loops on A, and so, for sustained signals, the increase in its molecular levels is only hindered by its degradation. This is the cause for its steady increase with time. Conversely, I is strongly affected by the negative feedback loop at later times. Protein A is in high concentration at later times and is statistically likely to be bound to both copies of the I genes. In the parameter range chosen, A is able to bind tightly to the I gene and shuts down its production at later times. Once I production is shut down, the level of I gradually decreases due to its degradation/deactivation. The combined effect of the negative feedback loop and the degradation/deactivation reaction causes the peak in the number of I molecules at early times ().
An alternative way of viewing the data is to run 1000 trajectories, determine the molecular levels of I and A at “long times” (time
2000 min), and construct a histogram. For sustained signals, the distributions appear similar in shape; however, I is at much lower levels than A ( and ). Therefore, these cells would appear to be producing mostly protein A if assayed at long times. This implies that the cells have made the decision to activate the anti-inflammatory gene program for sustained signals, an idea which will be explored in more detail in the next section. The same qualitative results were obtained using Model 2 ( and ).
For transient signals, the molecular levels of I increase gradually with time in some of the trajectories, while A remains at or around zero in all of the trajectories ( and ). Protein A is entirely signal-dependent so it is not able to be produced in significant quantities when the signal degrades too quickly ( and ). As a result, the negative feedback loop is not effective for highly transient signals and the production of I is limited only by its ability to initiate the positive feedback loop. For this particular degradation rate, kdeg,T
1, only some of the trajectories of I led to nonzero amounts, which is clearly observed in the long time histogram (). The bistability that results is entirely stochastic in nature, evident by the fact that the mean field time trajectory for I follows the average of the nonzero curves (Figure S2
. Bistability results from the “roll of the dice” in determining whether or not the positive feedback loop is initiated before the signal and resulting I molecules degrade. It appears that bistability develops at approximately kdeg,T
(). As the signal becomes more transient in nature, an increasing number of I trajectories remain at I
0. In many biological systems bistability serves as a tool in making cellular decisions, but in this system it is simply a consequence of the signal degrading too quickly and does not have direct relevance to the decision to be inflammatory versus anti-inflammatory 
. Instead it determines what percentage of cells will be inflammatory or exhibit no response.
Cellular decisions can be made by manipulating signal duration.
In silico knock-out experiments can be performed where either the positive or negative feedback loops are removed. When the positive feedback loop is removed, the results for sustained signals are virtually unchanged since the negative feedback loop dominates (Figure S3
). For the transient case, both I and A would remain at approximately zero since both would be entirely signal dependent. The end result of blocking the positive feedback loop is that the cell would activate the anti-inflammatory gene program for sustained signals or neither program for transient signals. When the negative feedback loop is knocked out, the results for transient signals remain unchanged since the positive feedback loop dominates the behavior (Figure S4
). For sustained signals, both I and A would be produced in significant quantities. There would no longer be a peak in I at short times, but a steady increase as a function of time. The end result of blocking the negative feedback loop is that the cell would choose the pro-inflammatory gene program for transient signals. The model predicts that both programs would be initiated for sustained signals, but this is likely not biologically relevant since the ability of anti-inflammatory genes to act “anti-inflammatory” has been removed by the deletion of the negative feedback loop. Therefore, sustained signals would likely lead to pro-inflammatory cells in this system. If both feedback loops were removed, both I and A would be entirely signal dependent and would exhibit similar qualitative behavior (Figure S5
). The end result of blocking both the positive and negative feedback loops is that the cell would choose both I and A for sustained signals or neither for transient signals. Again, since the negative feedback has been deleted, it would likely be observed that the cells are pro-inflammatory for sustained signals or neither for transient signals.
The same qualitative results of A dominating transcription ( and ) over I ( and ) for sustained signals and I dominating transcription ( and ) over A ( and ) for transient signals is observed for both Model 1 and Model 2. Therefore, whether or not transcription factors dissociate after each transcriptional event is not important in this analysis. All the remaining results are based on the models allowing for either A or I to dominate transcription at different signal durations. We shall only consider results generated using Model 1, although the exact same conclusions would be reached using Model 2.
Tuning signal duration as a tool for cellular decisions
From the two test values of kdeg,T, it appears as if this mechanism could be a viable tool for making cellular decisions based entirely on signal duration modification. To further probe this idea, we constructed the following in silico experiment. Thousands of I and A trajectories similar to those obtained in the last section were generated for kdeg,T values ranging from 0 to 10 with all other parameters fixed at the values given in . It is assumed that the downstream gene programs (inflammatory or anti-inflammatory) activated by I or A are only initiated when the respective molecules cross a threshold value in molecular number at any point in the time trajectory. To measure this computationally, the maximum point in each trajectory was compared to a threshold value and the number of I and A trajectories crossing the threshold was recorded for each value of kdeg,T. The percent of I or A trajectories crossing the threshold directly corresponds to the percent of cells that exhibit pro-inflammatory or anti-inflammatory responses in an equivalent biological experiment since each trajectory corresponds to 1 cell.
Cells with sustained signals (kdeg,T<10−3) clearly exhibit an anti-inflammatory response while cells with transient signals (0.02<kdeg,T<1) clearly exhibit an inflammatory response ( and ). There is a small crossover region for mildly transient signals (10−3<kdeg,T<0.02), the behavior in which is determined by the choice in downstream threshold values for I and A. For the given set of parameters, it is clear that I has the capability of attaining a higher maximal value since its transcription rate is an order of magnitude higher than A and since I has a positive feedback loop ( and ). Therefore, one natural choice is to set the downstream threshold higher for I to account for this difference. With different I and A threshold values, the inflammatory and anti-inflammatory response curves overlap such that a mixture of cells exhibiting either response would be observed in the crossover region (). An alternative choice is to set the thresholds to be the same for both I and A, which leads to no response in the crossover region (). Setting a downstream threshold is an arbitrary process, which merely shifts the inflammatory and anti-inflammatory response curves, but does identify two distinct crossover region behaviors that may both have biological relevance.
One difficulty in performing this analysis is choosing where to set the long time cutoff for when an inflammatory or anti-inflammatory response can occur. For instance, it would not necessarily be biologically relevant if it takes several days of sustained signaling for A to cross its signaling threshold. For the purposes of this study, we arbitrarily chose this cutoff time to be 2000 minutes (~1.5 days). Under most conditions, the decision is made much earlier than 2000 minutes, but we chose a longer time for inspection to ensure that we are identifying the decision for all kdeg,T
values. As an example, consider and . For this kdeg,T
value, the decision to be I and not A has occurred as early as 300 minutes. Given the nature of the curves, the same conclusion can be drawn by looking at 2000 minutes. Different parameter sets could give even earlier decision times or later decision times. Measuring at 2000 minutes allows us to reliably measure what decision has occurred, regardless of when it occurs. Immediate early gene transcription is on the order of minutes 
, while Th17 differentiation is on the order of 5 days 
. Choosing 1.5 days as the cutoff seemed like a reasonable intermediate choice for considering biological relevance. An example of a cellular decision becoming evident after 1 day is in Th17 differentiation, where IL-17 mRNA starts to be transcribed in appreciable quantities at this time (unpublished observations).
Using 3.5 days as the cutoff time revealed an important limitation of our model. For any signal with a degradation rate greater than zero, I will eventually be produced in significant quantities if the simulation is allowed to run for long enough times (Figure S6A
). The cause of this limitation is highlighted by the time course of I for a signal degradation rate of 0 (). Even with a sustained signal, I does not go to zero since its potent positive feedback loop coupled with the sustained signal are still active in the moments when a given A molecule bound to the I gene dissociates or is deactivated, which can counteract the decrease in I due to its own degradation. For all cases where the signal degradation rate is not zero, A will eventually decay to zero since it is signal dependent. If the number of A molecules becomes low enough, it can no longer prevent the few remaining I molecules from initiating the positive feedback loop. Performing the threshold analysis at 3.5 days, I and A thresholds can still be chosen to allow for a cellular decision to be made based on signal duration changes, but the peak in I is broader and the crossover region occurs at smaller values of the signal degradation rate (Figure S6B
). One way to deter this from occurring is to decrease the degradation rate of A, so that the negative feedback loop will remain active even at long times (Figure S6C
). We do not believe that this limitation would be relevant under most biological situations since it would take longer than 3.5 days of steady signaling without any other changes within the cell for it to be an issue. Conceivably, once a pro-inflammatory or anti-inflammatory response was selected, the cell would have other means for shutting down the other set of genes. In this sense, the model is not entirely self-contained.
Amplitude only significantly affects cellular decisions in the crossover region
Intracellular signals are generally described by two variables, signal amplitude and duration. Often in the literature, these two variables are combined into one quantity, the strength of signal. One interpretation of the strength of signal is the integrated signal, which is the area under the curve on a plot of signal strength vs. time. The gene transcriptional network depicted in was formulated to allow for disparate cellular decisions based strictly on differences in signal duration. To declare the cellular decisions made in the transcriptional model as strictly dependent on signal duration, it was also necessary to study the qualitative effect of modifying the signal amplitude.
The effects of amplitude modification were determined by setting the initial value of pSTAT3 molecules to different amounts and using the same procedure that was used to generate . We observed that decreasing the signal amplitude from 100 molecules of pSTAT3 to 10 molecules shifted the I and A curves to the left without significantly altering their shape (). The effects of amplitude modification can be observed within the crossover region from the inflammatory to anti-inflammatory response in that a 10-fold decrease in signal amplitude would shift the crossover from occurring at a degradation rate of 0.01 min−1 to 0.001 min−1. The production of A is entirely signal dependent, so significant decreases in the amount of pSTAT3 molecules translates into significant decreases in the amount of A molecules produced during the lifetime of the signal. As the number of A molecules is reduced, I becomes less regulated and will eventually win out. The effect of reducing pSTAT3 on the number of A molecules produced during the lifetime of the signal can partially be compensated for by increasing the lifetime of the signal, i.e. decreasing kdeg,T, which accounts for the shift in the curves to the left. The crossover region is susceptible to these effects since it is at the edges of the inflammatory and anti-inflammatory response regimes.
Altering the initial signal amplitude shifts the I and A curves, but does not affect their shape.
Amplitude modification also affects significantly transient signals (kdeg,T>0.1 min−1) by shifting the region of bistability to lower kdeg,T values as amplitude is decreased. As an example, consider a kdeg,T value of 0.2 min−1. For a signal amplitude of 100 pSTAT3 molecules, 100% of cells assayed would exhibit an inflammatory response. With a 10-fold decrease in amplitude, 50% of cells would exhibit an inflammatory response, while the other 50% of cells would exhibit no response. Therefore, amplitude only affects the number of cells exhibiting a given response for transient signals, but does not affect the type of response. The shift in the bistable region is also due to the fact that a decrease in signal amplitude can be partially compensated for by increasing the lifetime of the signal. Highly transient, low amplitude signals do not lead to any response since they do not have enough strength to kick off the positive feedback loop. Signals that are moderately sustained (kdeg,T<0.001 min−1) or moderately transient (0.1<kdeg,T<0.01), are unaffected by these changes in amplitude.
Decreasing the signal amplitude to five or less leads to both inflammatory and anti-inflammatory responses for sustained signals in some of the trajectories, unless the threshold for an inflammatory response is increased from that used to generate (data not shown). Since the transcription rate of I is faster than A, I is able to initiate the positive feedback loop before appreciable amounts of A are produced in some of the trajectories. At later times when appreciable amounts of A are produced, the I signal shuts down but had already crossed the inflammatory response threshold leading to an observed inflammatory response starting at early times and also an anti-inflammatory response at later times.
Model predicts that IL-6 will promote an anti-inflammatory response in SOCS−/− cells
Another key difference between the IL-10R and the IL-6R, besides their differing interactions with SOCS3, is that the IL-10R has two motifs for STAT3 activation while the gp130 subunit of the IL-6R has four, leading to the IL-6R inducing a higher amplitude STAT3 signal 
. Given that IL-6 induces a higher amplitude signal than IL-10, one could conclude that strong signals induce an inflammatory response, while weak signals induce an anti-inflammatory response. In this definition of signal strength, amplitude is the defining feature. If an experiment were performed where SOCS3 was no longer able to inhibit IL-6 signaling, the high amplitude signal would be converted from a transient signal into a sustained signal, and thus would have an even higher signal strength than before. Thinking about this problem in terms of signal strength, one would predict that IL-6 would still induce an inflammatory response when SOCS3 inhibition has been removed.
A revealing experiment performed by Yasukawa et. al. in which they exposed SOCS3−/−macrophages to IL-6 and LPS and measured the concentration of two pro-inflammatory cytokines (TNF and IL-12(p40)) in the supernatant showed just the opposite 
. LPS induces an inflammatory response in macrophages, so a decrease in inflammatory cytokine production would be considered anti-inflammatory. As the dosage of IL-6 was increased, the pro-inflammatory cytokine concentration dropped meaning that IL-6 was acting as an anti-inflammatory cytokine on the SOCS3−/− macrophages. This puzzling result cannot be explained by a signal strength model where signal strength is essentially just a measure of the amplitude of the signal.
Our model is capable of making a prediction for how IL-6 will act on SOCS3−/− macrophages by making very few assumptions. We assume that IL-6 induces a moderately transient pSTAT3 signal that can roughly be approximated by a kdeg,T value of 0.1 min−1 in WT macrophages. The second assumption is that IL-10 induces a slowly decaying pSTAT3 signal that can be approximated by a kdeg,T value of 0.003 min−1 in both SOCS3−/− and WT macrophages. IL-6 and IL-10 would conceivably induce similarly decaying pSTAT3 signals in SOC3−/− macrophages, so we can use the same kdeg,T value. In their experiment, Yasukawa et. al. measure concentrations of secreted inflammatory cytokines. In our model this translates to counting the percentage of trajectories (cells) that cross the imposed threshold value of I at any point for a set value of kdeg,T, which is characteristic to the macrophage type and the type of cytokine added along with LPS. By carrying out this analysis for a variety of initial amplitudes, we can make a prediction about how the macrophages will respond to increasing dosages of either IL-6 or IL-10. Since LPS induces an inflammatory response on its own, it is assumed that 100% of the cells would be pro-inflammatory in the absence of IL-6 and IL-10. Therefore, until a large enough signal amplitude is reached where gene A can cause a drop in the number of trajectories which cross the threshold value of I due to the negative feedback loop, it is assumed that all the cells will remain pro-inflammatory.
Using these assumptions, it is clear that the model predicts IL-6 will be anti-inflammatory when acting on SOCS3−/− macrophages (). The model also captures the observed phenomenon that as the pSTAT3 amplitude is increased (i.e. the IL-6 dose is increased), IL-6 has a stronger anti-inflammatory effect. Since there is no distinction between IL-10 acting on SOCS3−/− or WT macrophages, the model predicts that no distinction will be seen experimentally, as was observed 
. For transient signals, the model predicts a pro-inflammatory response, which is why the percentage of inflammatory cells remains at 100% for all amplitudes when WT macrophages are exposed to IL-6. If it were IL-6 alone, there would be a ramping up of the percentage of inflammatory cells at low amplitudes since a basal level of signaling is required to initial the positive feedback loop. However, since LPS induces an inflammatory response on its own, the ramping up period is not observed in the figure. There is no amplitude range for transient signals where A is capable of decreasing the percentage of inflammatory cells.
Model predicts that IL-6 will be anti-inflammatory in SOCS3−/− cells when given some inflammatory stimulus such as LPS for most pSTAT3 signal amplitudes.
has proposed that sustained STAT3 signals are essential for anti-inflammatory responses while transient signals promote inflammation, which is consistent with our model. The model serves as a tool for visualizing how the duration of the signal can lead to two different responses from the same genetic network activated by the same transcription factor from two different pathways. As further evidence of signal duration being a defining feature in determining the inflammatory response, El Kasmi et. al. also developed a system for analyzing this problem 
. They transfected macrophages with an altered form of the EPO receptor (EPOR), which was capable of inducing pSTAT3 and was not susceptible to SOCS3 inhibition. This altered form of the EPOR induced an anti-inflammatory response indistinguishable from IL-10, which also supports the hypothesis that signal duration is a determinant of whether or not a response is inflammatory or anti-inflammatory 
While it is true that many different cytokines activate STAT3 that may or may not induce an inflammatory or anti-inflammatory response, often these cytokines activate multiple STAT molecules, such as STAT1 or STAT5, in addition to STAT3. It is beyond the scope of this analysis to determine how genes activated by other STATs affect the decision to be pro-inflammatory or anti-inflammatory. Since different STATs can be activated by the same receptor, heterodimer formation may be prevalent and it is difficult to predict how this may impact a network with this structure.
One could argue that the proposed transcriptional network model cannot easily account for common gene expression. The key wiring in our model that makes it sensitive to changes in signal duration is the positive and negative feedback loops. Common genes would not necessarily need to have the same wiring as the anti-inflammatory/inflammatory genes. As an example, consider the example of a common gene C that is also induced by pSTAT3. If this gene had a positive feedback loop, significant amounts of C would be produced from both transient and sustained signals. Even in the absence of a positive feedback loop, a transient signal still could produce significant amounts of C for high amplitude signals and low values of the C degradation/deactivation.
We have provided the structure of a gene transcriptional network which is sensitive to changes in signal duration. The essential feature of the model is that both genes of interest are activated by the same transcription factor and that one of the genes has a positive feedback loop to sustain itself in the absence of the main signal, but is also negatively regulated by the opposing gene. We believe that the model may be applicable to IL-6/IL-10 signaling since all essential steps of the model have been observed except for the positive feedback loop, namely that STAT3 induces both inflammatory and anti-inflammatory genes and that the anti-inflammatory genes inhibit the inflammatory genes. Experimentation would be necessary to identify a positive feedback loop in the inflammatory gene network in order for the model to be fully validated. Bioinformatic techniques could be used to determine whether inflammatory genes have the potential to bind to the promoter region of other inflammatory genes. Gene array analysis could then be used to verify the binding. Mutations to these promoter sites could then be used to knockout the positive feedback loop. Our model would then predict that IL-6 would not be able to induce an inflammatory response (Figure S3
Even though the model was presented in the context of IL-6/IL-10 signaling, it may be applicable to a range of biological networks which are sensitive to signal duration. There has been extensive work in the field of systems biology to identify network motifs 
. Our network contains components of previously identified motifs such as an autoregulatory positive feedback loop 
and a negative feedback loop. We are by no means the first to explore the effects of regulatory loops on cellular decision making. A number of studies have investigated the effects of positive feedback loops on decision making 
. Mangan and Alon described the characteristics and features of feed-forward loops 
. Our focus in this work was not to propose a new signaling motif, but instead to show that a very simple genetic network could be used to explain the differences between IL-6 and IL-10 signaling without having to resort to using other more complicated arguments such as receptor-specific STAT3 conformational changes or effects due to unknown species. Further investigation is necessary to see how our model compares to existing models in the literature and to look for other examples in biology where our model may be applicable.