Two Component Systems: Physiological Effect of a Third Component

Signal transduction systems mediate the response and adaptation of organisms to environmental changes. In prokaryotes, this signal transduction is often done through Two Component Systems (TCS). These TCS are phosphotransfer protein cascades, and in their prototypical form they are composed by a kinase that senses the environmental signals (SK) and by a response regulator (RR) that regulates the cellular response. This basic motif can be modified by the addition of a third protein that interacts either with the SK or the RR in a way that could change the dynamic response of the TCS module. In this work we aim at understanding the effect of such an additional protein (which we call “third component”) on the functional properties of a prototypical TCS. To do so we build mathematical models of TCS with alternative designs for their interaction with that third component. These mathematical models are analyzed in order to identify the differences in dynamic behavior inherent to each design, with respect to functionally relevant properties such as sensitivity to changes in either the parameter values or the molecular concentrations, temporal responsiveness, possibility of multiple steady states, or stochastic fluctuations in the system. The differences are then correlated to the physiological requirements that impinge on the functioning of the TCS. This analysis sheds light on both, the dynamic behavior of synthetically designed TCS, and the conditions under which natural selection might favor each of the designs. We find that a third component that modulates SK activity increases the parameter space where a bistable response of the TCS module to signals is possible, if SK is monofunctional, but decreases it when the SK is bifunctional. The presence of a third component that modulates RR activity decreases the parameter space where a bistable response of the TCS module to signals is possible.


Introduction
Two component systems (TCS) are biochemical signaling modules that are ubiquitous in bacteria and are also present in some eukaryotes. Prototypical TCS are composed of two proteins: a sensor kinase (SK) and a response regulator (RR). The SK phosphorylates a histidine residue and subsequently transfers the phosphate to an aspartate residue in the RR. There are many variations around this prototype, ranging from phosphorelays that can concatenate up to three phosphotransfers (HisRAspRHis-RAsp) between different proteins to hybrid kinases in which the SK and the RR domains are fused in the same protein [1,2]. In prototypical TCS, the SK can be bifunctional if, when unphosphorylated, it increases the dephosphorylation rate of the RR. Otherwise, the SK is monofunctional. The majority of well characterized SKs are bifunctional, with a few, such as the chemotaxis regulating CheA, being monofunctional.
In addition to SKs and RRs, some TCS are also known to interact with specific phosphatases that regulate dephosphorylation of the RR [3]. These core components of TCS and phosphorelays are also complemented by auxiliary proteins that play a regulatory role in the activity of some TCS, transmitting the cognate signal to the SK. For example, the SK CheA is regulated through its interaction with membrane receptors that detect chemical compounds in the medium and direct organisms towards higher concentrations of nutrients [4] and the activity of the SK NRII that regulates nitrogen fixation is modulated through its interaction with the protein PII [5].
In recent years, interactions between the TCS and auxiliary proteins were identified as a strategy to integrate non-cognate signals in the regulation of TCS [6]. For example, the orphan SK RetS interacts with the GacS SK, preventing the response of the latter to its cognate signal [7,8,9,10,11] and the peptide PmrD binds to and protects the phosphorylated form of the RR PmrA from the phosphatase activity of its cognate SK, PmrB [12]. The GacS/GacA TCS regulates virulence in Pseudomonas aeruginosa [13,14], while the PmrB/PmrA TCS is required for resistance of Salmonella to acidic and antibiotic stresses, among others [12,15]. These systems raise the question of understanding the effect of such interactions with the core TCS module in the operating regime of the module and what consequences these effects may have on the influence of the module on the cellular physiology of the organism [16,17,18,19,20,21,22]. Previous studies suggested that a third component that binds to and protects the phosphorylated form of the RR causes delays in the response of autogenous TCS that regulate their own expression [12,17,22]. However, to our knowledge, no studies were made about the effect that binding of a third component to the SK has on the potential dynamic behavior of the TCS module. In addition, the effect of both types of third component proteins were not studied in TCS that do not regulate their own expression.
In previous work we have used mathematical models to characterize the effect of diverse architectures on the signaling response of prototypical TCS. The analysis of such models enables understanding if particular physiological responses are more effectively achieved by one of several alternative designs of the network that executes the biological process of interest [23]. Such studies are difficult, if not impossible to do without the assistance of those mathematical models. In the case of the TCS, we showed that TCS with bifunctional SKs are more effective in buffering the TCS against crosstalk, while monofunctional SK are more effective in integrating different signals [24,25]. We have also identified necessary conditions for the existence of post-translational bistable responses in prototypical TCS [25]. If a system is capable of bistable responses, this means that its output variable can assume one of two possible values as a consequence of the same input. The specific value that the variable assumes depends on the value that the variable had before the stimulus. Posttranslational bistability is only possible in TCS in which the affinity between the phosphorylated SK and unphosphorylated RR is similar to that between the unphosphorylated forms of the proteins. In addition, a large fraction of the dephosphorylation flux of the RR must be independent of any phosphatase activity of the SK [25].
Given these considerations, in this work our goal is to understand the physiological effect of a third protein, such as RetS or PmrD, on the function of canonical TCS in the absence of auto-regulation of gene expression. To achieve this, we built and analyzed mathematical models for the alternative designs of TCS with and without such a third component, and compared the dynamic behavior of the different systems. This analysis identifies specific physiological behaviors that are more effectively executed by each alternative design for the TCS.
Our study reveals that a RR-binding third component (TC RR ) decreases the region in parameter space where a bistable response is possible, while a SK-binding third component (TC SK ) increases the parametric region where a bistable response is possible when the SK is monofunctional and decreases it when the SK is bifunctional.

Results
In order to understand the physiological effect of a third component (TC) on the function of a prototypical TCS, we built models of TCS with and without that TC and compared the dynamical behavior of those models. Figure 1 shows a schematic representation of the three models used in our analysis. These models are mathematically described by using a mass action system of ordinary differential equations (ODE) [26]. The resulting ODE systems for each of the three alternative models can be

Models and Comparisons
The network model that we use to describe the prototypical TCS in our analysis is that defined in Igoshin et al. [25], which is based on earlier work [27]. In Model A, shown in Figure 1, the SK can autophosphorylate and/or autodephosphorylate in response to an external signal. Both phosphorylated and unphosphorylated forms of SK are allowed to bind RR with similar affinities, as reported in [28,29,30]. Binding of unphosphorylated SK and RR is reversible and forms a dead-end complex (SKRR). Phosphorylated SK (SKP) can transfer its phosphate to the RR. The phosphorylated RR (RRP) will modulate the biological levels and activity of relevant proteins.
This network for the prototypical TCS was modified to study the effect of a TC binding to either the SK or the RR. The changes in the network are also shown in Figure 1. Model B represents a TCS where a third component binds to the SK (TC SK ), inactivating it. Model C represents a TCS where a third component binds to the phosphorylated RR (TC RR ) and stabilizes this phosphorylated form. In prototypical TCS modules with bifunctional sensors, the unphosphorylated SK can destabilize the phosphorylated form of the RR and it increases the dephosphorylation rate of RRP (k 8 .0 in Figure 1). In prototypical TCS modules with monofunctional sensors, the unphosphorylated SK has no effect upon the dephosphorylation rate of RRP (k 8 = 0 in Figure 1). The model includes a phosphatase that dephosphorylates RRP independently of the SK. This is done for generality. In the cases where no such phosphatase exists, this set of reactions can be replaced by a single reaction where the unstable RRP phosphate bond hydrolyzes over time. An appropriate choice of parameter values will make the results of the analysis similar to those described for the full model.
In this study we analyze the potential effect of a TC in the physiological behavior of TCS modules with bifunctional and monofunctional sensors independently. If the TC has no effect on the physiological behavior of the TCS, then the presence of TC in particular instances of TCS should be understood as an evolutionary accident. If the TC has an effect on the physiological behavior of the TCS, this could provide a rationale for the selection of a TCS design that includes a TC. To perform the analysis, we compare the dynamical behavior of Model A to that of Models B or C, independently. This comparison is done in two ways.
First, Models A and B (or C) are compared ensuring that the parameter values of all processes that are common are the same in the two models. This guarantees that whatever differences are found are only due to the addition of the TC. This comparison is equivalent to comparing an organism where the TCS interacts with a TC to another where the TC has been deleted from the genome. This situation could occur, for example during the creation of a new biological circuit by genetic manipulation in a biotechnological context. Thus, this type of comparison is relevant for understanding the differences in behavior of biological circuits created using synthetic biology techniques.
Second, we also perform a mathematically controlled comparison between Models A and B (or C). This is a well established method for evaluating the irreducible effect of a change in the design of a biological circuit on the dynamic behavior of the network [31]. In this comparison, in addition to ensuring that Models A and B (or C) have the same values for corresponding parameters of all processes that are common, we use the differences between the designs as degrees of freedom that evolution can use as a substrate to minimize differences between the dynamic behavior of the two systems. If the alternative designs can be made equivalent by using these degrees of freedom, then one may argue that they cannot be distinguished by natural selection. If, after making the systems as equivalent as possible, there are still irreducible differences in the physiological behavior between designs, then one may expect one of them to be preferably selected when its functionality provides better adaptive advantage. In the models under comparison, the difference is the deletion of a protein from the module between Model B (or C) and Model A. In this situation, the protein burden caused by Model A is lower than that caused by its alternative designs. Hence, we allow that the system changes the total concentrations of the remaining proteins (SK and/or RR). The details for this comparison are given in the methods section. This comparison is thus relevant for understanding the differences in the dynamic behavior that are intrinsic to the differences in design between Models A and B (or C), and to those alone, in evolutionary terms.

Effect of a third component on TCS signal amplification and bistability
Signal amplification is an important physiological property of TCS. TCS with appropriate signal amplification can provide evolutionary advantages to organisms harboring them. Thus, understanding how signal amplification is affected by adding a TC to a TCS would help in predicting under which conditions to expect such a design to be selected. Figure 2 shows that all models can achieve the same signal amplification, whether the environmental signal modulates the autophosphorylation (k 1 ) or the autodephosphorylation (k 2 ) of the SK. This can be seen because the difference between the amount of RRP (phosphorylated RR) when k 1 is low (k 2 is high) and when k 1 is high (k 2 is low) can be similar for all models. Nevertheless, Model B responds at higher signal intensities and Model C responds at lower signal intensities than Model A, when the stimulus modulates the SK autophosphorylation reaction rate (compare the curves for k 1 response of Model A to those of Models B and C in Figure 2). When the signal modulates the SK autodephosphorylation reaction rate, Model B responds at lower signal intensity and Model C at higher signal intensity than Model A (compare the curves for k 2 response of Model A to those of Models B and C in Figure 2). However, mostly, the differences in signal intensity at which the systems are turned ON or OFF are small.
In addition, the prototypical TCS shown in Model A can show bistable behavior [25], making it possible that a signal can lead to one of two alternative responses, depending on the history of the system. Such a response may have some evolutionary advantages, for example in situations like sporulation where an irreversible developmental decision is made by cells. Bistable regions in the curves of Figure 2 have three values of RRP for a single value of signal intensity. The two extreme values are the alternative stable steady states, while the middle value is a biologically irrelevant unstable steady state that is mathematically required to exist if two stable steady states are present. In the figure one can see that the signaling ranges where bistability exists are different if the environmental signal modulates the autophosphorylation (k 1 ) or the autodephosphorylation (k 2 ) of the SK.
Necessary, although not sufficient, conditions for the existence of such bistable behavior in the prototypical TCS are i) the formation of a dead-end complex between the dephosphorylated forms of SK and RR and ii) that a sufficiently high fraction of the flux for the dephosphorylation of RRP is independent of SK. To understand how the presence of a TC affects the possibility of a bistable response in the prototypical TCS, we analyzed Models B and C in search of the existence of multiple steady states, followed by a comparison of the physiological behavior between Models A and B, and between Models A and C.
Given that signals can in principle modulate either the autophosphorylation (k 1 ) or the autodephosphorylation (k 2 ) rate of SK, we performed parallel computational experiments independently modulating their intensity. These experiments were done independently for models with monofunctional and bifunctional SK (Figure 2). When the signal modulates SK dephosphorylation (changes in k 2 ), the system behaves symmetrically to when SK phosphorylation (changes in k 1 ) is modulated. In the first case, increases in signal intensity cause the fraction of RRP to decrease, while in the latter, increases in signal intensity cause the fraction of RRP to increase. Our results show that, in an uncontrolled comparison, the range of bistability for the bifunctional prototypical TCS is larger than if a TC binds any of the proteins of the module (compare panel B to panels D and F of Figure 2). Bistability for Model B in panel D is only observed for k 1 signaling, while no bistability is observed for Model C in panel F. On the other hand, the range of bistability for the monofunctional prototypical TCS is larger than if a TC binds the RR of the module (compare panel A to panel E of Figure 2), but smaller than if the TC binds the SK (compare panel A to panel C of Figure 2). Differences among the three systems are more pronounced when the signal induces dephosphorylation of the SK (k 2 ), rather than inducing SK autophosphorylation (k 1 ).
An additional definition is needed before further presenting and discussing the results. Hereafter the system is said to be in an ON state if most of its RR is in the phosphorylated RRP form. If most of the RR is in its dephosphorylated form, the system is said to be in its OFF state. With this in mind, and as one might expect, systems with a TC SK are in an ON state for a smaller signaling range (panels C and D) and systems with a TC RR are in an ON state for a larger signaling range (panels E and F), in comparison with the uncontrolled Model A (panels A and B).
When the comparisons are controlled we see that the response of Model A can become similar to that of Model B or C by adjusting the total amount of available SK. If the response of Model B is to be mimicked, the total amount of SK in Model A is decreased ( Figure 2, panels C and D, see methods for the exact values of the total amount of SK), while mimicking the response of Model C leads to an increase in the concentration of SK (Figure 2, panels E and F, see methods for the exact values of the total amount of SK).
The k 2 -response curves in Figure 2 panels B and C show that the switch from ON to OFF (from high to low levels of RRP) in these models could be irreversible or very difficult to reverse. In other words, modulation of the autodephosphorylation rate of SK by an external signal could generate nearly irreversible biological switches.
Our simulations also show that the necessary conditions for bistability in prototypical TCS remain necessary in the TCS with a TC. If either no independent phosphatase is present in the system (Ph = 0) or no dead end complex is formed (k 10 = 0) all TCS modules analyzed here are monostable (see section ''Effect of changes in SK-independent RRP dephosphorylation and SKRR affinity on bistability'' below).
In summary, a TC RR causes a reduction in the TCS parameter space of bistability and an increase in the signaling range in which the system is in the ON state (responds at lower k 1 -signal intensity and at higher k 2 -signal intensity), whether the SK is monofunctional or bifunctional. This can be more effectively compensated by prototypical TCS through a change (an increase) in the concentration of the SK. In contrast, TC SK increases the signaling range in which the TCS can show a bistable response if and only if the SK is monofunctional and the environment modulates k 2 (SK dephosphorylation rate). The behavior of TCS with a TC SK can be mimicked by prototypical TCS through a change (a decrease) in the concentration of the SK.

Effect of a third component on TCS response time
In addition to signal amplification, the response time to signals is an important physiological property of TCS. In evolutionary terms, a change in response time may have important consequences to the fitness of the system. Therefore, we analyzed the effect of a TC on the response times of the TCS. To do this we performed four independent sets of experiments for each of the models, and independently considering systems with a monofunctional SK and with a bifunctional SK. In experiments 1 and 2 we instantaneously change the signal k 1 and measure how long the system takes to come within 90% of its new steady state. This measures the response time of the system if the physiological signal modulates SK phosphorylation. In experiments 3 and 4, we instantaneously change the signal k 2 and measure how long the system takes to come within 90% of its new steady state. This measures the response time of the system if the physiological signal modulates SK dephosphorylation. The details about how the experiments were run are as follows: 1 -We set each system to its OFF state, with k 1 = 10 25 s 21 . Then, we increased the value of k 1 to a value k 1 higher and measured how long the system took to get to within 90% of its new steady state value. k 1 higher was systematically changed between 10 25 and 10 s 21 . 2 -We set each system to its ON state, with k 1 = 10 s 21 . Then, we decreased the value of k 1 to a value k 1 lower and measured how long the system took to get to within 90% of its new steady state value. k 1 lower was systematically changed between 10 25 and 10 s 21 . 3 -We set each system to its OFF state, with k 2 = 10 s 21 . Then, we decreased the value of k 2 to a value k 2 lower and measured how long the system took to get to within 90% of its new steady state value. k 2 lower was systematically changed between 10 25 and 10 s 21 . 4 -We set each system to its ON state, with k 2 = 10 25 s 21 . Then, we increased the value of k 2 to a value k 2 higher and measured how long the system took to get to within 90% of its new steady state value. k 2 higher was systematically changed between 10 25 and 10 s 21 .
Results are shown in Figure 3. We see that the response times increase by more than two orders of magnitude when the new parameter value k ?lower or k ?higher approaches the threshold value for exiting the bistability region of a system. The peaks of slower response in the curves in Figure 3 are in the region of signal intensity that lies immediately beyond the border of the bistability ranges shown in Figure 2. Given that the peaks of slower response are located at the exit of the bistable region, there is no peak in the signal-response time curve when the response is monostable or when there is an irreversible turning OFF of the system. Model B and Model A|B (A controlled for B) don't have a peak in their OFF to ON k 2 -response times (Panel C of Figure 3) because these models irreversibly turn OFF after an increase in k 2 (as depicted in Figure 2 panel C). Model C also has no peak in the response time (Panels C and D of Figure 3) because this model has a monostable response to changes in k 2 (see Figure 2 panel E). In panels G and H of Figure 3, neither of the three systems shows a peak in their signal-response time curve because of the lack of bistability in their signal-response steady state curve (see Figure 2 panels D and F). When Model A is compared to Model B in an uncontrolled manner, the time response peaks of Model A appear at signal intensities that are always lower than those where the peak appears in the response of Model B. When Model A is compared to Model C in an uncontrolled manner, the time response peaks of Model A appear at signal intensities that are always higher than those where the peak appears in the response of Model C (see Figure S1).
In order to have a proxy of the integral temporal responsiveness of each system, we calculated the area under each of the signalresponse time curves shown in Figure 3. This area is the sum of all the transient response times for each signaling response. The values of these areas are given in Table 1 and show that overall response times are similar between Models A and B. In contrast, Model A has a faster response than Model C. When the comparison is not controlled, differences between integrated response times of the three models are small, when the signal modulates autophosphorylation of SK. However, if SK dephosphorylation is modulated, Model B has the fastest integrated response, followed by Model A. Model C is, again, the slowest responder (Table S1).
In summary, Model B is a faster overall responder than the prototypical TCS when the system is turned ON by modulating the phosphorylation rate of the SK, and it is a slower responder in any other case. In contrast, Model C is always slower to turn ON or turn OFF than the prototypical TCS, under controlled comparison conditions.

Stochastic effects of a third component
Fluctuations in the amount of proteins that participate in biological reactions can lead to stochastic effects in the system's behavior, when the total number of proteins participating in reactions is small. We performed stochastic simulations to understand the role of stochasticity on the effect of the TC on the physiological response of the TCS networks. These simulations take into account that the number of TCS proteins present in the cell are typically in the 10-1000 molecules range.
The simulation experiments performed were similar to those described in experiments 1-4 of the previous section, although with a smaller number of data points. Figures 4 and 5 show the results of these simulations.
The OFFRON plots start with the system at the OFF steadystate (low concentration of active RR) corresponding to a low value of k 1 (k 1 = 10 25 s 21 ) or a high value of k 2 (k 2 = 10 s 21 ), and depict the temporal trajectory of the RRP concentration after an instantaneous increase in k 1 or decrease in k 2 , for three different values of k 1 and k 2 .
The ONROFF plots start with the system at the ON steadystate (high concentration of active response regulator) corresponding to a high value of the signal k 1 (k 1 = 10 s 21 ) or a low value of k 2 (k 2 = 5?10 26 s 21 ), and depict the temporal trajectory of the RRP concentration after an instantaneous decrease in k 1 or increase in k 2 , for three different values of k 1 and k 2 .
The simulation results for three different signal intensities are plotted in Figures 4 and 5. Three independent simulations are shown for each signal intensity. The values of k 1 and k 2 in each trajectory are chosen to be below, next to and above the threshold value at which the system switches from OFF to ON, or from ON to OFF (in the cases in which this threshold exists). Because each system has a different threshold value, the parameter scan is different for each plot.
The results from the analysis of the continuous model are consistent with the stochastic simulations: as discussed in the previous section (Figure 3), in systems with a signal range of bistability the response times increase when the signal intensity is near the threshold value at which the system exits the bistability region. One can see in Figures 4 and 5 that, in many cases, the curves that correspond to a signal that is just outside of the bistability range do not reach steady state during the simulation time. These curves correspond with the peaks in Figure 3.
Furthermore, our simulations predict that the systemic response becomes noisier as the signal intensity approaches the threshold  . Temporal responsiveness curves of Models A, B, and C. The systems are at an initial steady state and, at time zero, the signal, represented in the x axis, changes instantaneously and the time it takes for the system to get to within 90% of the new steady state is measured and plotted in the y axis. A-D: Response times of TCS with monofunctional SK. E-H: Response times of TCS with bifunctional SK. The OFF to ON plots start with the systems at an OFF steady state (low levels of RRP) corresponding to a low value of k 1 (A, C, E, G) or a high value of k 2 (B, D, F, H). The signal is then changed to increase the steady state level of RRP. The ON to OFF plots start with the systems at an ON steady state (high levels of RRP) corresponding to a high value of k 1 or a low value of k 2 . The signal is then changed to decrease the steady state level of RRP. Peaks that indicate slower response times are located immediately outside the range of bistability. The lack of a peak in a curve can be due to monostability or irreversibility. The dashed lines indicate the signal value at which Models B and C exit its bistable range. Absence of a dashed line indicates irreversible turning ON or OFF of the system (Model B in panel C ) or absence of bistability (see the signal-response curves of Figure 2 value for bistability. Just above and just below this value there is an increase in the stochastic fluctuations of the system. This can be seen because the triplicate curves corresponding to these values in Figures 4 and 5 are much more different among themselves than the triplicate curves for the signals away from this threshold. The response in the systems A, B and C is noisier when k 1 is modulated than when k 2 is modulated. The OFF to ON trajectories of Model B after an instantaneous decrease in k 2 confirm that the turn OFF of this system due to an increase in k 2 is irreversible and the system can't return to the ON state (see Figure 2 panel C). The system C does not have a bistability region in its k 2 -response curve (see Figure 2 panels E and F). Therefore, we don't find a range of k 2 values for which the systemic response becomes slower and noisier.

Robustness of the analysis
The analysis thus far was done using the specific set of parameter values reported in Table 2. In order to study the generality of the results we performed sensitivity analyses of the bistability to changes in the different parameter values and concentrations of the systems. The results of the controlled and uncontrolled comparison between Model A and Model B or C with respect to the effect of changing parameter values on a possible bistable response of the TCS are summarized in Table 3. The detailed results are shown in Figure S2, where we show a set of two-dimensional sections of the multidimensional parameter space in which bistability is observed. Overall, a system with a TC SK appears to have a wider parameter range of bistability if the SK is monofunctional, and a lower parameter range of bistability if the SK is bifunctional, while a system with a TC RR appears to have a lower parameter range of bistability, for systems with either a monofunctional or a bifunctional SK, when either system is compared to a prototypical TCS. However, if the comparison between Model A and Model B or C is controlled, then we see that the robustness of the parameter range of bistability is larger in the prototypical TCS (Model A) with only one exception: in systems with a bifunctional SK, Model C has a more robust parameter range of bistability.
Effect of changes in SK-independent RRP dephosphorylation and SKRR affinity on bistability SK-independent RRP dephosphorylation and SKRR complex formation are needed for bistable responses to exist in Models A, B, and C. In order to investigate how quantitatively changing these features affects bistability we performed the following computational experiments (Table 4). We independently and simultaneously changed the values for k 8 (the reaction that regulates dephosphorylation by the SK) and k 9 (changing the rate of dissociation between SK and RR) between 10 26 and 10. Then, we calculated the steady state(s) for each system at different values of the signal represented by the parameters k 1 or k 2 . k 1 and k 2 were independently and systematically scanned between 10 26 and 10 in logarithmic space at intervals of 0.01 units. The results are shown in Table 5 and Figure S3. Table 3 shows that, overall, bistability is possible in Model C in a smaller interval of parameter values than that for Models A and B. However, the picture changes when we analyze only the parameters that directly influence the necessary conditions for bistability (k 8 , k 9 , k 10 ). For these parameters, Model C is the system where overall bistability is possible in a wider range of parameter values, followed by Model B. Model A is the one where bistability is limited to a smaller region of parameter values. Nevertheless, when Model A is controlled to have signal-response curves that are as similar as possible to those of either Model B or Model C, Model A becomes the system where bistable responses can occur in a larger fraction of the space for k 8 , k 9 , and k 10 . For values of k 8 below a threshold that depends on the system and is lower in Model B than in Model A, bistability is present in both models. Within the range of k 8 values that permit bistability, an increase in k 8 causes an increase in the k 2 range of bistability (up to approximately six orders of magnitude for k 2 at the threshold value for k 8 ). This is so, despite the enlargement of the fraction of RRP dephosphorylated by SK, because the increase in k 8 causes an increase in the concentration of the SKRR dead-end complex (see Table 2. Basal values for the parameters and concentrations of the models in Figure 1. The value for this parameter was chosen to be one order of magnitude larger than that representing SK autodephosphorylation, because the TC SK enhances SK autodephosphorylation. c TC SK total is the total amount of the third component in Model B. This third component protein binds the SK of the TCS module. The amount for this protein was chosen taking into account that basal mRNA levels for RetS in GEO micro profiles of Pseudomonas aeruginosa are between 2 and 10 times higher than those of GacS. GacS is an SK and RetS is its cognate TC SK [47]. d TC RR total is the total amount of the third component in Model C. This third component protein binds the phosphorylated RR of the TCS module. The amount for this protein was chosen to be in the same order of magnitude as that of the RR, as is done in reference [43]. e This is the average value for the autophosphorylation catalytic constant between Salmonella typhimurium and Escherichia coli [16]. f It should be noted that, for Model C, this value for the phosphatase rate constant could be as high as 0.14 in Escherichia coli [16]. g Although some measurements have suggested that the affinity between nonphosphorylated forms of the SK and RR is much lower than the affinity between phosphorylated forms of the proteins [48], more recent measurements suggest the opposite [10]. doi:10.1371/journal.pone.0031095.t002 Figure S4). As k 8 decreases, the range of signal k 2 in which the models show bistability decreases steadily for a few orders of magnitude. Then, a lower boundary is reached and bistability is observed for one or less than one order of magnitude of k 2 signal, independently of the value for k 8 . Given that the formation of a dead-end complex between SK and RR is a necessary condition for bistability, we also want to understand the isolated effect of different fractions of RR and SK being sequestered into this complex on bistability. To understand the effect of changing the amount of SKRR dead-end complex on the signaling range in which the systems can be bistable we performed the following numerical experiment. First, we took each model from Figure 1. Then, we systematically scanned the values of the parameters k 9 and k 10 , independently and simultaneously, between 10 26 and 10 in logarithmic space at intervals of 0.01 units. These parameters regulate the amount of SKRR that is formed. Finally, for each pair of values for k 9 and k 10 , we independently calculated the steady state(s) for each system at different values of the signal represented by k 1 or k 2 . Each of these parameters was independently and systematically scanned between 10 26 and 10 in logarithmic space at intervals of 0.01 units. The results are shown in Table 5 and Figure S3.
Bistability can be found only for intermediate steady state concentrations of SKRR. If too little or too much SKRR is formed, then no bistable response is possible. Overall, for bifunctional TCS, Model C has the largest range of SKRR steady state concentrations for which bistability is possible, followed by Model B. In its uncontrolled form Model A has the smallest interval of SKRR steady state concentrations where bistability is permitted. This interval of concentrations decreases further when Model A is controlled to be comparable to Model B. However, when Model A is controlled to be comparable to Model C, the range of SKRR steady state concentrations that enable bistability becomes the largest of the three systems. In monofunctional TCS, Model C has a smaller range of SKRR steady state concentrations for which bistability is possible than Model B.
The notion that Model C is the one in which bistable responses are less sensitive to changes in the steady state concentrations of SKRR (in consequence of changing the affinity between SK and RR) is misleading. Bistability is only found in this model if the affinity between the dephosphorylated forms of SK and RR is much larger than that between SKP and RR or SK and RRP. Given that the affinity between all forms of SK and RR was measured as similar, it is not likely that bistability can be found in vivo in systems that are represented by this model.
A similar experiment was made by changing independently and simultaneously the total amount of SK and RR, followed by independent calculation of the steady state(s) for each system at different values of the signal represented by k 1 or k 2 . Again, each of the parameters was independently and systematically scanned between 10 26 and 10 in logarithmic space at intervals of 0.01 units. The results are shown in Table 5 and Figure S3. They are consistent with the situation described for changes in k 9 and k 10 .

Effect of the SK/TC SK and RR/TC RR concentration ratios on bistability
In order to understand how the relationship between the total amounts of SK (RR) and TC SK (TC RR ) influences the signaling range in which bistable responses are possible, we have performed a number of computational experiments. First, we took Models B and C from Figure 1. Then, we systematically, simultaneously and independently scanned the total amounts of SK (RR) and TC SK (TC RR ) in Model B (Model C), as described in Table 4. Finally, for each total amount of SK (RR) and TC SK (TC RR ), we calculated the steady state(s) for each system at different values of the signal represented by k 2 . This parameter was also systematically scanned between 10 26 and 10 in logarithmic space at intervals of 0.01 Table 4. Experiments to analyze the effect of changes in different parameter values and protein concentrations on the range of bistability for the alternative TCS modules a . Ratio between SK total and RR total . SK total 10 23 -10 3 mM R R total 10 23 -10 3 mM Ratio between SK total and TC SK total . TC SK total 10 23 -10 3 mM S K total 10 23 -10 3 mM Ratio between RR total and TC RR total . RR total 10 23 -10 3 mM T C RR total 10 23 -10 3 mM Formation of the SKRR dead-end complex and rate of RRP dephoshoprylation by SK   Some bidimensional sections of the multidimensional parameter space of bistability are shown in Figure S2. The results show that in TCS with a bifunctional SK, both a TC SK and a TC RR cause a decrease in the size of the parametric region of bistability, with one exception: Model C has a larger parametric region of bistability when the signaling target is SK autophosphorylation (k 1 units. The results are shown in Figure S3. We also performed similar test replacing k 2 by k 1 . The range of signal k 2 for which Model B can show a bistable response is observed to be dependent on the TC. Bistability is observed only within a narrow band of the SK-TC SK concentration space. Outside of this band, a bistable response cannot be observed. The range of total amount of SK in the system that may lead to a bistable response remains approximately constant for low total amounts of TC SK . However, within the band of total SK and TC SK in which bistability is observed, as total TC SK increases, the range of total SK amount that can generate bistable responses also increases. At concentrations of TC SK between approximately 2 and 7 mM, we find bistability for total SK concentrations between 0.2 and 0.001 mM or lower. At higher total TC SK concentrations, only small amounts of SK are available in free form. This prevents formation of the SKRR dead-end complex that is required for bistability.
As is the case in Model B, bistability in Model C can be achieved in a narrow band of the concentration space. However, within the range of values of this simulation, whatever the concentration of TC RR , the system can always show bistability.

Summary of the comparisons
Tables 6 and 7 summarize our findings regarding the different physiological criteria that are relevant for TCS signal transduction and can be asserted from the analysis of our models. In general, if the signaling target is SK autophosphorylation Model C responds at lower signaling intensities, followed by Model A, and finally by Model B. If the signal enhances SK dephosphorylation, Model B is the one that responds at lower signal intensities, followed by Model A, and Model C. This causes Model C to be in an ON state for a wider signaling range, and Model B to be in an ON state for a narrower signaling range, in comparison with Model A.
The system with the largest range of signaling in which it can show a bistable response depends on both, the type of SK in the module and the SK activity (autophosphorylation or autodephosphorylation) that is targeted by the signal. For TCS with monofunctional SK, Model A has the largest signaling range for bistability, as well as the largest fraction of parameter space where such bistability can be observed, if the environment modulates SK phosphorylation. In contrast, Model B has the largest signaling range for bistability, as well as the largest fraction of parameter space where such bistability can be observed, if the environment modulates SK dephosphorylation. For TCS with bifunctional SK, Model B has the largest signaling range for bistability if the environment modulates SK phosphorylation. However, it is Model C that has the largest fraction of parameter space where bistability can be observed. In contrast, Model A has the largest signaling range for bistability, as well as the largest fraction of parameter space where such bistability can be observed, if the environment modulates SK dephosphorylation.
Modulation of SK dephosphorylation leads to responses that have an equally small amount of noise in all Models. However, modulation of SK phosphorylation leads to noisier responses in Model B, followed by Model A and finally Model C.
As is the case with bistability, the model with fastest response times depends on the type of SK in the module and on the SK activity (autophosphorylation or autodephosphorylation) that is targeted by the signal. Both in systems with monofunctional and  bifunctional SK, Model A is the fastest to respond (Model C is the slowest) whether the signaling target is the autophosphorylation or the autodephosphorylation of the SK, with only one exception: Model B turns ON faster if SK autophosphorylation is modulated directly. The response times of Models A and B are similar, but Model C tends to be much slower than Model A.

Biological Relevance
Bacteria often sense and adapt to changes in the environment through TCS and phosphorelays. A question that this work addresses is how variations to the prototypical TCS by means of an accessory third protein that either binds the SK or the RR affect the dynamical behavior of the TCS module.  TCS can, in principle, mediate both gradual and switch like (bistable) responses to environmental stimuli [32,33]. The switchlike response has typically been associated to the positive feedback introduced by genetic regulatory loops in the regulation of autogenous TCS. Nevertheless, such feedback does not necessarily imply the existence of bistability [34]. In fact, genetic positive feedback loops are not strictly necessary for the existence of bistable responses in prototypical TCS. Such responses can also come about through post-translational regulation of bacterial signal transduction networks [25,35]. Namely, bistability is possible in prototypical TCS if a reversible dead-end complex is formed between the dephosphorylated SK and RR and if a sufficient amount of RRP is dephosphorylated independently of the SK phosphatase activity [25].
TC proteins that regulate signal transmission to prototypical TCS have been known for years [36,37]. However, only recently have such interactions been proposed as a way to integrate noncognate signals in the TCS regulated responses. In fact, these interactions have been reported in TCS that are responsible for regulating both, resistance to antibiotics and virulence [6,7,8,9,12,13,14,15].
Biological examples of the first situation can be found in the PmrB/PmrA/PmrD system. The third component PmrD binds and stabilizes the active form of the RR, PmrA. This system regulates antibiotic resistance in Salmonella and other bacteria. Various studies of the PmrA/PmrB/PmrD system suggest that this TC RR could be an intermediate evolutionary step to evolve indirect regulation of the TCS [12,16,17,22,38]. The feedforward connector loop formed by PmrD is presented as a design that speeds up activation and slows deactivation of the gene expression of the proteins in the TCS [17]. Our results suggest that this may not be so in non-autogenous TCS. If the TCS has a TC RR , loss of this protein will make the corresponding prototypical TCS faster to turn ON and OFF (Table S1 and Figure S1). In fact, if the steady state response curve of the prototypical TCS is mathematically controlled to be as similar to that of the TCS with a TC RR as possible, then that prototypical system is always faster. A TC RR appears also to be a feature that decreases the fraction of parameter space in which bistable responses are possible (Tables 3  and 5), except in TCS with a bifunctional SK and when the environment modulates SK autophosphorylation. Thus, a TC RR creates a TCS module that is less likely to show bistable responses and slower in responding to environmental signals, which it can sense at lower intensities than the prototypical TCS without any TC, if SK phosphorylation is modulated.
Antibiotic resistance is arguably a trait whose response should be gradual and proportional to the amount of antibiotic found by the bacteria to increase its survival chances. If this is not so, and a bistable response is possible, bacteria can be made more sensitive to antibiotics [39] and therefore their survival will be hindered. Given that bistability has been observed in the antibiotic resistance of some bacteria [39], a TC that binds the RR would reduce the possibility of such bistable response, potentiating adaptation and tolerance to threatening stress challenges. In addition, having such a TC could enable a response at low antibiotic concentrations, thus increasing the chances of survival for the organism.
The other well studied example of a TC interacting with the TCS is the RetS/GacS/GacA system, where RetS reversibly binds and inactivates the SK GacS. This system regulates virulence in Pseudomonas aeruginosa. Recently, it has been shown that the GacS/ GacA TCS acts exclusively through the regulation of the transcription of two genes, rsmY and rsmZ [40]. The product of these genes are two untranslated small regulatory RNAs (sRNAs), RsmY and RsmZ, that counter translational repression exerted by the RNA-binding protein RsmA on target mRNAs encoding virulence factors. There is an additional SK, LadS, that appears to counter the action of RetS on GacS. However, this effect is indirect, as not direct physical interaction between GacS and LadS was observed [10]. It may be that LadS sequesters RetS, as RetS does with GacS. Our analysis of a TCS with a TC SK reveals that this module will respond at signal intensities that are slightly higher (lower) than those of the prototypical TCS, if SK authophosphorylation (autodephosphorylation) is directly modulated. Furthermore, if one is to synthetically change a TCS module and create an artificial circuit with a TC SK , the engineered circuit will typically respond faster to signals if the environment modulates SK dephosphorylation. However, evolution can eventually equalize response times by changing the SK concentration of the module and making both TCS modules have steady state response curves that are similar. A TC SK can increase the signaling range in which a bistable response is possible (Table 5). Bistability could be advantageous when the system has to choose between two different operational states [35,41], as is often the case for virulent organisms. For example, Mycobacterium tuberculosis is a persistent organism in the lungs of 2/7 of the world population [42]. However, only under certain conditions that are not yet completely clear does this organism causes tuberculosis [42]. Bistability could provide populations with the capacity to sample which type of phenotype is more advantageous at different times and enhance survival of the organisms through bet-hedging strategies [43,44].
Experiments to test the existence of bistability in a TCS with a TC SK could be as follows, taking the RetS/GaS/GacA system as an example. First, determine if the system can show bistable response: incubate two Pseudomonas aeruginosa strains (a wild type strain, with the TC protein RetS, and a RetS mutant strain, without the TC protein) at different environmental conditions of inducing signal intensity, allow the cells to approach a steady state and measure the levels of expression of the sRNAs RsmY and RsmZ. In a TCS module with a monostable gradual response, the level of expression of the output molecules should be proportional to the environmental inducing signal intensity: at intermediate signal intensities there is an intermediate amount of output molecule. However, if the RetS/ GacS/GacA response is bistable, we will find that, for intermediate intensities of inducing signal, the measured levels of RsmY and RsmZ in single cells are distributed in a bimodal manner, with low and high levels (but no intermediate levels) of this sRNAs. If bistability is present and the in vivo effect of RetS is to amplify the signaling range for which a bistable response is possible, we will find that this bimodal distribution of the measured levels of RsmY and RsmZ in single cells of the RetS mutant strain will be observed in a smaller range of signal values. To investigate if the results of our simulations are valid for in vivo conditions and if the RetS/GacS/GacA system could have an irreversible response (as observed in Figure 2 C), we can incubate both strains in a high-stimulus environment, allow the cells to approach a steady state and measure the levels of the sRNAs RsmY and RsmZ . If the in vivo system behaves as its in silico proxy, when the stimulus is removed (transfer the cells to a non-inducing environment), we will find that in RetS mutant cells the levels of RsmY and RsmZ shift from a low value to a high value, but in wild type cells the levels of RsmY and RsmZ remain at a low value.
The arguments discussed thus far explain part of the biological relevance of our work. Another interesting aspect of it regards the modulation of SK autophosphorylation and dephosphorylation. Currently the community is inclined to assume that dephosphorylation is the target of modulation by environmental signals in many cases. However, to our knowledge, conclusive experiments that decide the issue are still lacking in most systems and it is still unclear whether the physiological signal modulates SK autophosphorylation (k 1 ) or SKP dephosphorylation (k 2 ). That is why we have performed our simulations taking as a signal both changes in k 1 and k 2 . An unexpected result of our simulations may shed some light on this issue, and allow us to hypothesize which one of the reaction rates is modulated by the signal in the case of TCS with a TC. We have found that, for TCS with a bifunctional SK, a TC decreases the possibility of a bistable response. For TCS with a monofunctional SK, the same effect is observed if the signal modulates k 1 . However, if the signal modulates k 2 , a TC SK increases the range of signal intensities in which a TCS can show bistability, and a TC RR decreases it. Thus, for TCS with a monofunctional SK, the results suggest that the physiological signal should modulate SK dephosphorylation (k 2 ) both when bistability is an advantageous feature in the function of a TCS with a TC SK component, and when bistability is a disadvantageous feature in the function of a TCS with a TC RR . Conversely, the physiological signal should modulate SK autophosphorylation (k 1 ) when bistability is a disadvantageous feature in the function of a TCS with a TC SK .
The work presented in this paper provides motivation for further analyses of the TCS responsible for regulating virulence and antibiotic resistance, providing clues as to possible mechanisms to both decrease virulence and antibiotic resistance. In the case of virulence, whenever it is regulated by a TCS of the type analyzed here, simultaneously targeting the TC and the SK appropriately could prevent the organism from becoming virulent. In the case of antibiotic resistance, targeting the TC and its interaction with the RR could be used to facilitate locking the bacteria in an antibiotic-sensitive state and facilitate treatment of infections.

Equations
In order to compare the physiological behavior of the three systems in Figure 1, we must create a mathematical representation for each of the networks. The positive and negative terms of each ODE correspond to individual reactions that give rise to the synthesis and degradation of the reactant, respectively. Each reaction is considered to be mass action.
Because the turnover times for protein synthesis and degradation are much higher than those for the phosphorylationdephosphorylation reactions, we consider the total amount of each participating protein to be approximately constant. The parameters for the models are given in Table 1. All these parameters have an experimental basis, clearly presented in Igoshin et al. [25].

Mathematically controlled comparisons
We aim at comparing the physiological behavior of the three models in order to understand if the presence of a TC in a TCS module causes intrinsic differences to the potential physiological responses that the modules can have. To make sure that the differences observed in the behavior of the systems that are being compared are due to the presence of the TC, the comparisons must be made in a controlled way. For this we use the method of mathematically controlled comparisons [31]. This method requires that all components and processes that are common to the alternative models that are to be compared are made numerically equal, making the models internally equivalent. In contrast, the components and processes that are different between the alternative models are degrees of freedom that nature could potentially use to compensate the changes in the physiological responses caused by the differences between systems. In this case, the systems with a TC invest additional resources to synthesize a new protein that binds either the SK or the RR and modulates their phosphorylation state. All new processes of Models B and C with respect to Model A are due to the presence of this TC. In order to control the comparison between TCS with TC and the prototypical TCS, the prototypical system (Model A) should also be allowed to invest additional resources in adjusting the total amount of the SK or the RR. These adjustments will allow the prototypical system to have a physiological response that is as similar as possible to that of the model with a SK-binding or a RRbinding TC (Models B and C, respectively). This control condition ensures maximal external equivalency between the models. Once the maximum equivalency is achieved between the compared models, the remaining behavioral differences can be related to the presence of the TC.
To determine the changes in the total amount of SK or RR that make the physiological responses between Model A and Models B or C as similar as possible, we have used a minimum square differences method. We have calculated the steady state responses of the system in Models B and C to changes in the input phosphorylation or dephosphorylation rate of the modules, by calculating the steady state concentration of RRP in Models B and C, at input signal strengths between 10 26 and 10. These curves were then used individually to fit Model A and calculate the concentration of SK and/or RR that would minimize the differences in the steady state RRP concentration between Model A and Models B or C, independently. All calculations were done using Mathematica. The best fits are achieved by allowing the total amount of SK to change in Model A. The values for the total amount of SK in Model A that minimize the differences between the responses of this model and Model B or Model C are shown in Table 8.

Calculations
All simulations were performed in Mathematica [45] and COPASI [46]. Analyses of regions of bistability were done in Mathematica, using in-house scripts. Figure S1 Temporal responsiveness curves of Models A, B, and C. The systems are at an initial steady state and, at time zero, the signal, represented in the x axis, changes instantaneously and the time it takes for the system to get to within 90% of the new steady state is measured and plotted in the y axis. A-D: Response times of TCS with monofunctional SK. E-H: Response times of TCS with bifunctional SK. The OFF to ON plots start with the systems at an OFF steady state (low levels of RRP) corresponding to a low value of k 1 (A, C, E, G) or a high value of k 2 (B, D, F, H). The signal is then changed to increase the steady state level of RRP. The ON to OFF plots start with the systems at an ON steady state (high levels of RRP) corresponding to a high value of k 1 or a low value of k 2 . The signal is then changed to decrease the steady state level of RRP. Peaks that indicate slower response times are located immediately outside the range of bistability. The lack of a peak in a curve can be due to monostability or irreversibility Absence of a dashed line indicates irreversible turning ON or OFF of the system (Model B in panel C ) or absence of bistability (see the signal-response curves of Figure 2). The difference between this Figure and Table S1. (TIF) Figure S2 Effect of changing the parameter values on the range of bistability in the three TCS modules. In the panels, the x-axis represents values for k 1 (SK autophosphorylation rate constant) or k 2 (SK dephosphorylation rate constant), and the y-axis represents values for each of the other reaction rate constants that are common to the three models (from k 2 to k 13 ). The region where bistability is possible is shaded in blue. The number above each set of plots represents the summation of all areas of bistability in a given model, that is, is a measure of the size of the parametric space of bistability. A, B: Comparison between Models A and B, with a monofunctional SK. C, D: Comparison between Models A and B, with a bifunctional SK. E, F: Comparison between Models A and C, with a monofunctional SK. G, H: Comparison between Models A and C, with a bifunctional SK. (TIF) Figure S3 Percentage of parameter space where a bistable response is possible for Models A, B, and C. Experiments as described in Table 4. The x and y axis represent the values of the scanned parameters, while the z-axis represents the orders of magnitude of signal for which there is a bistable response. The red projection represents the area of parameter space where bistable responses are possible. A -Bifunctional system, signal modulating dephosphorylation of the SK.; B -Bifunctional system, signal modulating the phosphorylation of the SK; C -Monofunctional system, signal modulating dephosphorylation of the SK.; D -Monofunctional system, signal modulating the phosphorylation of the SK. See text for details and discussion. (TIF) Figure S4 Influence of the k 8 value (SK bifunctionality rate constant) on the k 2 range of bistability. Within a k 8 range of values, an increase in k 8 causes an increase in the k 2 range of bistability (panel a and b). This is so, despite an enlargement of the fraction of RRP dephosphorylated by SK (panel c), because of an increase in the SKRR concentration due to a higher value of k 8 (panel d). The simulations were performed using the system represented by Model A.

Supporting Information
(TIF) Author Contributions