Behaviors and Strategies of Bacterial Navigation in Chemical and Nonchemical Gradients

Navigation of cells to the optimal environmental condition is critical for their survival and growth. Escherichia coli cells, for example, can detect various chemicals and move up or down those chemical gradients (i.e., chemotaxis). Using the same signaling machinery, they can also sense other external factors such as pH and temperature and navigate from both sides toward some intermediate levels of those stimuli. This mode of precision sensing is more sophisticated than the (unidirectional) chemotaxis strategy and requires distinctive molecular mechanisms to encode and track the preferred external conditions. To systematically study these different bacterial taxis behaviors, we develop a continuum model that incorporates microscopic signaling events in single cells into macroscopic population dynamics. A simple theoretical result is obtained for the steady state cell distribution in general. In particular, we find the cell distribution is controlled by the intracellular sensory dynamics as well as the dependence of the cells' speed on external factors. The model is verified by available experimental data in various taxis behaviors (including bacterial chemotaxis, pH taxis, and thermotaxis), and it also leads to predictions that can be tested by future experiments. Our analysis help reveal the key conditions/mechanisms for bacterial precision-sensing behaviors and directly connects the cellular taxis performances with the underlying molecular parameters. It provides a unified framework to study bacterial navigation in complex environments with chemical and non-chemical stimuli.


The Unified Model for Bacterial Taxis
In this supplementary section, we present the analytical details about our unified model for bacterial taxis. As mentioned in the maintext, the internal state of an individual cell, represented by the receptor-kinase activity (a), can be described by the MWC model: . (S1) where S(x) = {S 1 , S 2 , ...} represents the external condition and m denotes the average (aggregate) methylation level of the receptor cluster. The receptor methylation-demethylation process takes place in a slow time-scale, τ a , and restores the total activity, a, to the adapted level, a 0 . Such receptormodification reactions serves as the cell's memory about the environment and can be described by the following differential equation: The receptor activity, a, regulates the level of the intracellular response regulator (CheY-P) which promotes the clockwise (CW) rotation of motors. A swimming cell tumbles when its motors rotate clockwise or when it experiences rotational diffusion (denoted by z 0 ). The total tumbling rate is z(a) = z 0 + z 1 (a) = z 0 + τ −1 (a/K 1/2 ) H .
Here, 1/τ sets the duration time of the tumbling state, K 1/2 represents the activity level at which the CW bias is 0.5, and H denote the Hill coefficient of the motor response function, respectively. In the one-dimensional setup, a cell moves toward either the right (+) or the left (−) direction with speed v(S) and changes its direction randomly with the tumbling frequency z(a).
In response to the environmental stimuli, an ensemble of bacterial cells will migrate in the physical space and also distribute in the internal state space (due to memory). We define P ± (x, a, t) as the probability density of cells which is in state a and move in the "±" direction at (x, t). The master equation governing P ± (x, a, t) is where G ± (a, S) are the rate functions for the internal state variable a as the cell moves in the "±" direction. Since a is a function of m and S, its change over time is contributed by the local methylation changes and the cell migration in space, i.e., We define the density (ρ) and the flux (J) of cells in the spatial coordinate: Summing Eqs. (S4) and (S5) and integrating over a leads to Similarly, subtracting these two equations and integrating over a yields The system over a finite region [x 0 , x 1 ] with reflecting boundaries should satisfy the zero flux condition (i.e., J = 0) in steady state. Thus we should have P + da = P − da = ρ(x)/2 in steady state and Eq. (S10) becomes where z ± (x) = zP ± da/ P ± da represents the average tumbling rate for the right or left moving cells at position x. From Eq. (S11), one can easily find the equilibrium cell density: For the simple case where both v and ∆z are (nonzero) constant, the steadystate cell density takes an exponential profile in space. Here, the chemotactic drift arises purely from the tumbling rate difference between the left and right moving populations. Intuitively, if z + (x) < z − (x) for example, then on average cells tend to move in the right (+) direction because it is more difficult for cells to enter a region where they tend to tumble more frequently.
The origin of ∆z is that the average activity over the left-moving cells (a + = aP + da P + da ) differs from the average over the right-moving population (a − = aP − da P − da ). The average receptor activity in the steady state is We can estimate the activity difference (a + − a − ) ≡ ∆a by considering the scenario in Fig. S1: If the representative cell is moving to the right at (x, t), it may have swum all the way along the path (x − l + , t − δt + ) → (x, t) with the initial activity a + (x − l + , t − δt + ) = a(x, t); if the cell is moving to the left at (x, t), it may have come from the path ( Then, we have the following approximations along the two paths: Since the direction of motion is randomized during each tumbling event, the time scale of δt + and δt − (durations of swimming without tumbling) is set by the short run time 1/z(a). Then one can make Taylor expansions around S(x) and derive the difference between a + and a − at (x, t): It is worth remarking that the function G m has a slow time scale set by the intrinsic methylation time scale τ a which is much longer than the average run time z −1 ∼ 1s. This suggests that in the last line of Eq. (S15) we have which becomes negligible given τ a z −1 . Therefore, we arrive at the following key approximation for ∆a: Then we can use the approximation ∆z ≈ ∂z ∂a a=a ∆a and Eq. (S17) to rewrite Eq. (S12) as follows: Eq. (S18) shows that the equilibrium distribution is ultimately shaped by the two motility behaviors, tumbling (z) and swimming (v), both of which directly or indirectly depend on the external conditions. To complete this model, we still need the dynamic equation for a. First, we notice the following Then by using Eqs. (S4) and (S5), we find that Using Eq. (S17), one can see that a(x) in steady state (J = 0) should satisfy: where, to the first order approximation in the shallow gradients, we have Given ∂a ∂m = 0, Eq. (S23) indicates that a(x) ≈ a 0 in shallow gradients. If one keeps the second-order term in Eq. (S22), one will see that the average activity a(x) tends to increase with the gradient steepness.

Application to Bacterial Chemotaxis under Dual Chemical Gradients
In the natural environment, cells are often exposed to multiple chemical stimuli. Here we extend our analysis to the case where bacterial cells respond to two opposing chemoattractant gradients (aspartate: d[L] 1 dx > 0 and serine: dx < 0 ) that are sensed by the two most abundant E. coli receptors, Tar and Tsr, respectively. The free energy energy for the receptor-kinase activity in the MWC model is given by: The steady-state cell density is found to be Clearly, the density profile ρ depends on the Tar/Tsr ratio (r 1 /r 2 ). When Tar (or Tsr) dominates in the receptor complex, cells are expected to respond mainly to aspartate (or serine) gradient. For intermediate values of r 1 /r 2 , the competition between Tar and Tsr may lead to the accumulation of cells around some position x * . This position x * can be determined by the first- Then, in the logarithmic sensing regime, the first order condition suggests that . If the solution of x * exists within the interval [x 0 , x 1 ], then it is given by

Application to Bacterial pH Taxis
In this supplementary section, we consider the case of bacterial pH taxis. The free energy for the total receptor-kinase activity is The steady-state cell distribution in a pH gradient: with the effective potential V eff = −η(r 1 f 1 + r 2 f 2 ). We can use the firstorder condition, V eff (pH * ) = −η[r 1 f 1 (pH * ) + r 2 f 2 (pH * )] = 0, to determine the preferred pH. Let w * ≡ 10 pH * , w A q ≡ 10 K A q and w I q ≡ 10 K I q for q = 1, 2, we find that which is simply a quadratic equation of w * . The observed opposite responses to pH changes indicate that K A 1 < K I 1 for Tar and K A 2 > K I 2 for Tsr. Thus, we expect that w I 1 w A 1 , w A 2 w I 2 , and w I 1 w * w I 2 , which together simplify Eq. (S28) as Therefore, for a given Tar/Tsr ratio (r 1 /r 2 ), the preferred pH is mostly determined by the values of w A 1 and w A 2 (i.e., K A 1 and K A 2 ). The exact solution to Eq. (S29) is given by For the special case that r 1 = r 2 , Eq. (S29) has a simple solution, w * = w A 1 w A 2 , which means pH * = (K A 1 + K A 2 )/2 at r 1 /r 2 = 1. Considering the simple scenario where K I , we can obtain from Eq. (S29) that w * ≈ w A 2 r 2 /r 1 and correspondingly the preferred pH level is pH * = K A 2 − log 10 (r 1 /r 2 ).
Thus, we can conjecture an empirical relationship between pH * and r 1 /r 2 : where λ represents the sensitivity of pH * to the change of the Tar/Tsr ratio. The parameter λ depends on the relative values of K A 1 and K A 2 , and can be calculated numerically. The empirical Eq. (S32) has been tested through extensive numerical experiments and turns to be useful for our data analysis.

Application to Bacterial Thermotaxis
In this supplementary section, we apply our unified model to bacterial thermotaxis. The free energy for the Tar receptor activity is assumed to be Here, g(T ) satisfies g (T ) > 0 and g(T 0 ) = 0 at a reference temperature Thus, consistent with the experiment, the Tar receptor in our model acts as a warm sensor ( ∂a ∂T < 0) when m < m c and as a cold sensor ( ∂a ∂T > 0) sensor when m > m c . The (steady-state) adapted activity depends on temperature and can be modeled as: It is observed that a 0 ≈ 1/3 at room temperature T = 22 • C and a 0 ≈ 1/2 at T 0 = 32 • C, from which we estimate that β ≈ 0.07/ • C. The critical temperature T c at which the Tar receptors invert response is set by ∂a ∂T | T =Tc = 0 which is equivalent to m(T c ) = m c by Eq. (S34). The steady state methylation level m(T ) can be determined by dm/dt = 0 or a(m, T ) = a 0 (T ). It amounts to solve N f a (m, T ) = −β(T − T 0 ) which yields Thus, using Eq. (S36) in the condition m(T c ) = m c , one can find which shows how the critical temperature is encoded in the signaling pathway. Combining Eq. (S36) and Eq. (S37), one can show that which changes sign at T c given that g(T ) + E m > 0 holds for the range of temperature under consideration. By Eq. (S18), the steady-state cell distribution in a linear temperature gradient over the physiological range [T − , T + ] can be calculated as follows In deriving Eq. (S39), we have used Eq. (S34) and Eq. (S38). The expression of ρ(T ) is composed of two parts: the temperature-dependent swim speed v(T ) and the temperature-dependent effect from the tumbling behavior. From Eq. (S39), we can see that cells are able to accumulate near the critical temperature T c as long as β > 0, g (T ) > 0, and g(T ) + E m > 0. The first condition β > 0 reflects the temperature-dependent imperfect adaptation kinetics for E. coli and is required for the existence of the critical temperature; see Eq. (S37). The second condition g (T ) > 0 ensures that Tar acts as a warm sensor when its methylation level m is below m c and switches to a cold sensor when m > m c ; see Eq. (S34). The third condition g(T ) + E m > 0 means that an increase of the methylation level will always reduce the free energy in Eq. (S33) and ensures that the methylation level m(T ) increases with the temperature near T c ; see Eq. (S38). For convenience, we can define a general temperature-dependent function: Since g (T 0 ) = 0, the above equation determines the following transformation which allows us to rewrite the free energy Eq. (S33) as follows Thus, θ(T ) determines g(T ) and contains all the information about how the receptor activity depends on temperature. We further introduce which represents the accumulative effect of the tumbling rate difference ∆z on the cell migration. Using Eqs. (S40) and (S43), we rewrite Eq. (S39) as Thus, the effective potential function for bacterial thermotaxis is From the condition V eff (T * ) = 0, we can find the preferred temperature Thus, if v is constant, the first-order condition gives that T * = T c and the second-order condition, V eff (T * ) = β N η(T * )θ(T * ), indicates that T * = T c is the preferred temperature only if θ(T ) > 0. If, however, the swimming speed v increases with T , then Eq. (S46) suggests that T * < T c . Last, for mutant cells lacking all receptors (θ(T ) = 0), the density becomes ρ(T ) = Ω/v(T ) and T * corresponds to the temperature where v(T * ) reaches its minimum.

Bacterial Taxis under Both Chemical and Non-Chemical Stimuli
Our unified model can be applied to study bacterial taxis in more complex environments. As a final example, we investigate the bacterial response to an integration of both chemical and nonchemical stimuli. Specifically, we consider the decision-making process of the Tar-only mutant cells under two opposing linear gradients over the interval [x 0 , x 1 ]: a temperature gradient ( dT dx = G T > 0) and a chemoattractant gradient ( d[L] dx = −G L < 0). Both signals are sensed by Tar whose receptor-kinase activity is given by a = 1/ 1 + e N fa(m,[L],T ) , with the free energy difference between two states, In the presence of a chemical gradient, the critical temperature at which Tar inverts response is no longer defined intrinsically because it depends on the external chemical condition. However, we can still use Eq. (S37) which defines T c as an intrinsic parameter (solely encoded by the signal transduction pathway). To avoid confusion, T c in this section always refers to Eq. (S37) and is not necessarily the critical temperature for Tar to change its sensing mode. Then we can use Eq. (S37) to rewrite Eq. (S48) as follows Since f L decreases with x, it is possible that the steady-state methylation level m also decreases with x if the chemical gradient is strong enough. For appropriate gradient profiles of [L](x) and T (x), there may exist a particular position x c satisfying m(x c ) = m c . Then, by Eq. (S49) we find that from which x c can be determined, although the solution may be out of the interval [x 0 , x 1 ] where cells are constrained. Again, we use Eq. (S18) to derive the steady-state cell distribution: Clearly, without any chemical interference (i.e., f L = 0), Eq. (S51) recovers our earlier result Eq. (S44) for thermotaxis. For uniform chemical background (f L > 0 and df L /dx = 0), Eq. (S51) becomes Compared to Eq. (S44), the extra exponential term in Eq. (S52) is a decreasing function of T if θ(T ) > 0. This will suppress the accumulation of cells at high temperatures. For constant v, it is easy to verify that the presence of a uniform chemoattractant background can shift the preferred temperature from T * = T c to a lower temperature T * = T c − N f L /β, in agreement with Eq. (S50). Now, if we increase the chemical gradient (df L /dx < 0) against the temperature gradient, more and more cells will be lured away to chase the chemoattractant.