New Research In
Physical Sciences
Social Sciences
Featured Portals
Articles by Topic
Biological Sciences
Featured Portals
Articles by Topic
 Agricultural Sciences
 Anthropology
 Applied Biological Sciences
 Biochemistry
 Biophysics and Computational Biology
 Cell Biology
 Developmental Biology
 Ecology
 Environmental Sciences
 Evolution
 Genetics
 Immunology and Inflammation
 Medical Sciences
 Microbiology
 Neuroscience
 Pharmacology
 Physiology
 Plant Biology
 Population Biology
 Psychological and Cognitive Sciences
 Sustainability Science
 Systems Biology
Relationship between cellular response and behavioral variability in bacterial chemotaxis

Edited by Curtis G. Callan, Jr., Princeton University, Princeton, NJ, and approved December 27, 2007 (received for review June 11, 2007)
Abstract
Over the last decades, bacterial chemotaxis in Escherichia coli has emerged as a canonical system for the study of signal transduction. A remarkable feature of this system is the coexistence of a robust adaptive behavior observed at the population level with a large fluctuating behavior in single cells [Korobkova E, Emonet T, Vilar JMG, Shimizu TS, Cluzel P (2004) Nature 428:574–578]. Using a unified stochastic model, we demonstrate that this coexistence is not fortuitous but a direct consequence of the architecture of this adaptive system. The methylation and demethylation cycles that regulate the activity of receptorkinase complexes are ultrasensitive because they operate outside the region of firstorder kinetics. As a result, the receptorkinase that governs cellular behavior exhibits a sigmoidal activation curve. We propose that the steepness of this kinase activation curve simultaneously controls the behavioral variability in nonstimulated individual bacteria and the duration of the adaptive response to small stimuli. We predict that the fluctuating behavior and the chemotactic response of individual cells both peak within the transition region of this sigmoidal curve. Largescale simulations of digital bacteria suggest that the chemotaxis network is tuned to simultaneously maximize both the random spread of cells in the absence of nutrients and the cellular response to gradients of attractant. This study highlights a fundamental relation from which the behavioral variability of nonstimulated cells is used to infer the timing of the cellular response to small stimuli.
Molecular noise (i.e., stochastic fluctuations) has been largely reported as one important source of phenotypic variability in a number of biological systems as diverse as gene expression and signal transduction in prokaryotes and eukaryotes (1–3). A standard method to characterize noise in biological systems is to analyze the distribution of behaviors across a population of cells at steady state. This approach has been powerful in identifying specific molecular mechanisms responsible for controlling phenotypic variability. Alternatively, the temporal evolution of noise within a single cell also contains key information about underlying cellular dynamics, but this approach is seldom used to characterize cellular behavior outside the steadystate regime (4). It is conceivable, however, that biological systems that are sensitive to intracellular spontaneous noise are also sensitive to small extracellular perturbations such as a sudden change of environmental conditions (5). We wish to investigate whether there exists in bacterial chemotaxis a general relationship between the fluctuation of cellular behavior in single cells and the timing of the cellular response to a small external stimulus.
Bacterial chemotaxis, a cellular locomotion system, has become a classic model for signal transduction. Although the chemotaxis network consists of just a few molecular species, it can perform complex cellular functions such as adaptation in response to environmental changes. Escherichia coli bacteria swim using flagella activated by rotary motors. When motors spin counterclockwise (CCW), the flagella form a corkscrew bundle whose rotation propels the cell in a smooth trajectory called a run. Clockwise (CW) spinning of motors favors the disruption of the bundle and causes the cell to tumble. The trajectory of a swimming cell resembles that of a random walk, which consists of a succession of runs and tumbles. When exposed to a sudden increase of chemical attractant, the cell lengthens the time interval between consecutive tumbles. This modulation of tumbling rate allows bacteria to bias their random walk toward a source of attractants and to perform chemotaxis.
A combination of experiments and models has demonstrated that, at the population level, E. coli bacteria display robust exact adaptation: After an initial response to a chemical stimulus, the average tumbling rate returns precisely to its steady prestimulus level (6–8). Remarkably, these studies showed that the property of exact adaptation is robust to variations of biochemical parameters such as the concentration and reaction rates of the chemotaxis proteins (7, 8). In contrast with population measurements, recent singlecell experiments have revealed that the switching behavior of a single motor, from a cell adapted to a homogeneous environment, exhibits large temporal fluctuations (9, 10). It was found that the characteristic time scale of these fluctuations was so large that a steady tumbling rate could be defined only when taking time averages longer than several hundred seconds. The amplitude of these fluctuations, hereafter called noise, exceeded that expected from Poisson statistics. Very long time series of switching events from individual flagellar motors exhibited distributions of run (CCW) intervals with long tails. These distributions are in general not exponentially distributed as it was previously believed (9, 10).
Here, we develop a unified stochastic model of adaptation in bacterial chemotaxis and use the fluctuation–dissipation theorem (11, 12) to interpret how the noise in nonstimulated cells is related to the cellular response to a small chemical stimulus. We show that the design of the activitydependent (de)methylation cycles that provide robust adaptation is equivalent to the architecture of a covalent modification cycle (13). We demonstrate how the large behavioral variability observed experimentally in single nonstimulated bacteria arises from the amplification of molecular noise in the ultrasensitive covalent modification cycle. As a result, the steadystate activity of the kinase as a function of the ratio of the concentrations of the modifying enzymes CheR and CheB is sigmoidal (Hill coefficient >1). We therefore analyze how the slope of the kinase activation curve simultaneously controls both the amplitude of spontaneous fluctuations and the relaxation time of the system in response to a small stimulus. The prediction emerging from this analysis is that the noisiest cells should also exhibit the largest chemotactic response.
Results
Kinetics of the Adaptation Module.
The adaptation system in bacterial chemotaxis is governed by a well defined module of specific interacting proteins (Fig. 1). This adaptation module consists of two antagonistic enzymes, CheR and phosphorylated CheB (CheBP), which respectively add and remove methyl groups at multiple receptor residues of the receptorkinase complex. This series of methylation and demethylation cycles controls the activity of the histidine kinase CheA in the complex. Methylation and demethylation are complex biochemical processes that involve the enhanced recruitment of CheR and CheB to the receptors by a tethering site that is distinct from methylation sites (14, 15). As in previous models of chemotaxis (7, 16–23), we coarsegrain the underlying biochemical details to capture the behavior of the full chemotaxis system obtained from experiments on living cells. In particular, we assume that receptorkinase complexes can be either active or inactive (23), and we use the Barkai and Leibler activitydependent feedback in the receptor modification system that yields robust exact adaptation (7, 8) (Fig. 1 B). We use Michaelis–Menten kinetics for the methylation–demethylation steps of the complexes. Here, CheR (CheBP) binds the methylation sites only of inactive (active) receptors (16–23). But in Sect. 4.2 of supporting information (SI) Appendix , we also consider two alternative hypotheses where CheR methylates the receptors in both active and inactive conformations, or where only the catalytic step of the methylation reaction is activity dependent. Only the latter assumption is incompatible with singlecell experiments. We report separately the concentration of free receptor complexes and the concentrations of receptor complexes bound to CheR and/or CheBP. Under these conditions, the temporal evolution of the average methylation level, M, of free receptors throughout the methylation–demethylation cycles (19) obeys Eq. 1 (see Materials and Methods): Here, A* and A are the concentrations of free active and inactive receptor complexes, and r and b are the rates of methylation and demethylation of inactive and active receptors complexes. The parameters ε _{r} and ε _{bp} are the concentrations of CheR and CheBP. K_{r} and K_{b} are the effective Michaelis–Menten constants for the methylation–demethylation of a receptor complex, and k_{r} and k_{b} are the corresponding catalytic rates. The total concentration of receptor complexes is constant: The first (second) term in Eq. 2 represents the total concentration of inactive (active) receptor complexes including those bound to CheR and CheB. Eq. 2 is equal to 1 because we normalize all of the concentrations (A, ε _{r} , K_{r} , etc.) with the concentration, N, of receptor complexes present in the system. The concentration of active and inactive receptors at steady state, A̅_{tot}* and A̅_{tot}, is independent from the ligand concentration (see Materials and Methods), which is consistent with the property of perfect adaptation (7, 18). We model the phosphorylation cascade using mechanisms and parameter values similar to those proposed by Sourjik and others (19, 20, 24): The rate of autophosphorylation of the kinase CheA is proportional to the total kinase activity A _{tot}* and the phosphate transfer from the kinase to CheY and CheB is proportional to the concentration of phosphorylated kinases (Eqs. S6–S8 of SI Appendix ).
Importantly, we find that the equation (Eq. S15 of SI Appendix ) defining the steady state of the kinase activity is identical to the equation obtained by Goldbeter and Koshland (13) for covalent modification of a common substrate by two antagonistic enzymes (compare figure 3 in ref. 13 with Eq. S15 of SI Appendix ). Thus, the adaptation module (Fig. 1) should share some of the properties of a covalent modification cycle. As in ref. 13, the normalized constants, K_{r} and K_{b} , control the steepness of the activation curve of the kinase. When K_{r} and/or K_{b} are <1, the activation curve is steeper than a hyperbolic function (Hill coefficient >1). The ratio α = k_{r} ε _{r} /k_{b} ε _{bp} of the maximal enzymatic velocities determines the steadystate activity of the kinase (13). Because the ratio α is proportional to the ratio of [CheR] to [CheBP], it can be easily adjusted by changing the relative concentration of these two proteins.
Using recent biochemical parameters (Table S1 of SI Appendix ), we plotted the steadystate activity of the kinase, A̅_{tot}*, as a function of [CheR] and for [CheB] fixed at wildtype level (Fig. 2 A). In the absence of the CheBP feedback loop, the kinase activation curve exhibits an effective Hill coefficient of 3.5. In the presence of the CheBP feedback loop, the curve is less steep (Hill coefficient ≈2.5) (Figs. 2 B). This feedback loop also helps maintaining the kinase activity within the narrow functioning range of the rotary motor (20). The kinase activity as a function of 1/[CheB], and [CheR] fixed at wildtype level, exhibits a Hill coefficient >1 as well (data not shown). These results suggest that the adaptation module in bacterial chemotaxis is more sensitive, namely ultrasensitive (13), to small variations of [CheR] and [CheBP] than a standard hyperbolic activation curve. Outside of the transition region of the sigmoid (Fig. 2 B), the kinase would be either fully active or fully inactive. Thus, we expect the ratio α of the methylation and demethylation velocities from wildtype cells to lie within the transition region of the sigmoid where the activity of the kinase can vary (see alternative models in Sect. 4 of SI Appendix ) (16).
Stochastic Model of the Adaptation Module.
Recent experiments showed that individual E. coli bacteria adapted to a homogeneous environment exhibit large temporal variations in their behavior (10). It was hypothesized that the observed behavioral variability was due to the fluctuations in kinase activity. More specifically, the slow methylation–demethylation processes in the adaptation module govern the fluctuations in kinase activity. Three complementary experiments in which cells did not display large fluctuations at long time scales support this hypothesis: (i) when [CheYP] was not regulated by the chemotaxis network but substituted by the active mutant CheYD13K stably expressed from an inducible plasmid; (ii) in mutant cells whose receptors have a fixed methylation level; (iii) behavioral variability was found tunable in ΔcheR mutant cells complemented with various levels of [CheR]. Behavioral variability decreased from its maximal to minimal value when [CheR] increased from 1 to 4fold wildtype level (10).
In this section, we wish to discuss how the steep kinase activation curve (Hill coefficient >1) (Fig. 2 B) could amplify the spontaneous stochastic fluctuations associated with the methylation and demethylation reactions. We show that the amplification of fluctuations in kinase activity is large enough to produce the large behavioral variability observed experimentally in nonstimulated cells (10). Chemotaxis assays are usually performed in media that do not support growth but that provide in excess the essential elements to produce energy. Consequently, we consider bacteria that are adapted to a homogenous environment to have reached the equilibrium (25). Under this condition, thermal fluctuations cause spontaneous fluctuations of the rate constants in the chemotaxis network (11, 26).
We use a linear noise approximation (see Sect. 2.5 of SI Appendix ) (11, 12, 27, 28) to derive the spontaneous noise in the kinase activity (and in [CheYP]) around the steady state (Fig. 2 A and B). This approximation allows us to compute the relaxation time, τ _{a} , and the power spectrum of the spontaneous fluctuations in the output signal of the chemotaxis system. The relaxation time is the typical time scale that defines the steady state of a nonstimulated cell. At time scales shorter than τ _{a} , the spontaneous fluctuations of the system are strongly correlated and the behavior is unsteady.
For the sake of simplicity, we first ignore the CheB phosphorylation step. The spontaneous fluctuations in kinase activity, δA*, include “fast” fluctuations in receptor activity due to the binding and unbinding of ligand and also slow fluctuations βδM associated with the methylation–demethylation reactions (see Materials and Methods). At long time scales, only the slow fluctuations are relevant and δA* ≅ βδM. Insertion of the latter relation into Eq. 1 and the use of the linear noise approximation yield the noise, δA _{tot}*, associated with the total kinase activity in nonstimulated cells: Here, δη _{a} is a source of stochastic white noise and D_{a} is the strength of the spontaneous fluctuations associated with (de)methylation reactions of the receptors. Within the linear noise approximation, the individual methylation and demethylation steps follow simple Poisson statistics. Summing up the contributions from all of the methylation and demethylation steps, we find that D_{a} is proportional to the sum of the methylation and demethylation velocities (see Sect. 2.5 of SI Appendix ): Interestingly, D_{a} is independent of the slope (i.e., the Hill coefficient) of the kinase activation curve [b̅A̅* is approximately k_{b} ε̅ _{bp} /(1 + 2K_{b} ) ≅ k_{b} ε̅ _{bp} in the middle of the transition]. By contrast, we find that the relaxation time, τ _{a} , is proportional to the slope of the kinase activation curve (see Materials and Methods): Eq. 3 is analogous to the Langevin equation γẋ = −κx + f(t) that describes the fluctuations of a massspring system in a viscous fluid, with spring constant κ = 2k_{B}T/(τ _{a}D_{a} ), damping constant γ = 2k_{B}T/D_{a} , and random force f(t). Fluctuations and dissipation are related through the fluctuation–dissipation theorem 〈f(t)f(t′)〉 = 2k_{B}Tγδ(t − t′), where T is the temperature and k_{B} is the Boltzmann constant (11). As a consequence, the variance of the spontaneous fluctuations in kinase activity δA _{tot}* is σ_{a} ^{2} = k_{B}T/κ = τ _{a}D_{a} /2. Thus, increasing the Hill coefficient of the kinase activation curve is equivalent to decreasing the spring constant κ in the massspring system without changing the damping constant γ. The higher the Hill coefficient, the larger the fluctuations of the kinase activity, σ_{a} ^{2}, and the longer the relaxation time, τ _{a} . Importantly, we find that τ _{a} and 1/σ_{a} ^{2} calculated for the full adaptation module are proportional to the relaxation time and the inverse of the noise for one covalent modification cycle [e.g., the Goldbeter and Koshland system (13)] with the proportionality factor 1/β.
A practical way to study the behavioral variability of a nonstimulated bacterium is to plot the power spectrum associated with the spontaneous fluctuations in the output of the full chemotaxis system. The Fourier transform of Eq. 3 gives the power spectrum of the spontaneous fluctuations of the kinase activity. Similarly, we also calculate the power spectrum of the fluctuations of [CheYP] from the full pathway that includes the CheBP feedback loop (see Sect. 2.6 of SI Appendix ). In Fig. 3, we plot the power spectra of the fluctuations of [CheYP] for various values of [CheR] and the Michaelis–Menten constants K_{r} and K_{b} . Using recently published biochemical parameters (29), we find that this model best reproduces singlecell experiments (power spectra) (10) when the following two conditions are fulfilled: (i) K_{r} and K_{b} are ≈10^{−1} and (ii) the methylation and demethylation catalytic rates k_{b} and k_{r} are adjusted such that wildtype cells, like in experiments, exhibit the noisiest behavior (see Table S1 of SI Appendix ). The values of these parameters are also in agreement with biochemical data and parameters used in earlier models (see Sects. 3.2 and 4 of SI Appendix ). By contrast, when K_{r} , K_{b} ≈ 1, we find that the behavioral variability is not as sensitive to variations in [CheR] as in experiments (10). Similarly, when K_{r} , K_{b} ≈ 10^{−2}, the difference in behavioral variability between wildtype and mutants cells that express 4fold wildtype level of [CheR] is too large in comparison with experiments (10).
The power spectra in Fig. 3 exhibit two characteristic (knee) frequencies. At low frequency, the fluctuations of kinase activity govern the fluctuating behavior of the full system. The position of the knee at long time scales reflects the magnitude of σ_{a} ^{2} and τ _{a} , which depend on [CheR] and [CheBP]. At short time scales, the second knee frequency corresponds to the spontaneous fluctuations occurring within the phosphorylation cascade, and its position depends only on the time scale associated with the phosphorylation of CheY (Fig. 1 A). Next, we plotted the noise in kinase activity, σ_{a} ^{2}, and the relaxation time, τ _{a} , as a function of [CheR], in the absence or presence of the CheBP feedback loop (Fig. 2 C and D). As expected for a covalent modification cycle with a Hill coefficient >1 (13, 30–32), σ_{a} ^{2} and τ _{a} peak within the transition region of the kinase activation curve (Fig. 2 C). This effect is similar to the one observed in refs. 30 and 32. In the presence of the CheBP feedback loop, the overall profiles of the variance and relaxation time as a function of [CheR] are conserved, with the exception that τ _{a} decreases at very low levels of CheR (Fig. 2 D). We also find that the CheBP feedback loop slightly reduces the noise and the relaxation time (see detailed analysis in Sect. 2.4 of SI Appendix ). We obtain similar results when we compute σ_{a} ^{2} and τ _{a} as a function of 1/[CheB]. In Fig. 2 E, we plot σ_{a} ^{2} and τ _{a} as a function of both [CheR] and 1/[CheB] for the whole system including the CheBP feedback loop. The presence of the crest over the hyperbolic surface is the signature of the deviation from firstorder kinetics. As expected, large fluctuations in kinase activity coincide with long relaxation times.
Figs. 2 and 3 suggest that in wildtype cells, the ratio α of the methylation and demethylation velocities is tuned to the transition region of the kinase activation curve (Fig. 2 B) to exhibit fluctuations (Fig. 3) similar to those observed experimentally (10). This model is also compatible with the experimental observation in single cells where behavioral variability is maximal when CheR is expressed at wildtype level and decreases to the variability expected from Poisson statistics when [CheR] is expressed at 4fold wildtype level (10).
Relationship Between Behavioral Variability in Nonstimulated Cells and the Timing of the Chemotactic Response to Small Stimuli.
The critical hypothesis in this section is that the behavioral variability observed in nonstimulated single cells is fundamentally related to the timing of the response of stimulated cells. Because stimuli encountered by the bacterium in its natural habitat may be small, we calculate the response to stimulus using a linear perturbation analysis of the kinetic system. Although we chose to keep our study independent from the actual chemical stimulus present in the environment, the chemotactic response in the real system also depends on the initial amplification of the input stimulus mediated by complex allosteric mechanisms taking place in the receptorkinase complex (17, 22, 33–37). In this model, we consider that a small external perturbation, such as a sudden exposure to attractant, causes an “instantaneous” change of the receptor activity, ΔA _{input}*. During the subsequent adaptation, the changes ΔM in the methylation level of the receptors govern the changes in kinase activity. We linearize Eq. 1 and eliminate ΔM (see Materials and Methods) to obtain the input–output relationship of the adaptation module solely in terms of changes in total kinase activity ΔA*_{input}: Within this linear approximation, the pathway relaxes with the same time scale τ _{a} in response to either an external stimulus or internal spontaneous noise (Eqs. 3 and 6 ). Thus, one can infer the sensitivity of the chemotactic response ΔA _{tot}*(t) to a small external perturbation ΔA _{input}* by characterizing the time correlation of the noise, δA _{tot}*, in nonstimulated cells. In this picture, a large τ _{a} will cause on average long runs following a small step stimulus ΔA _{input}*. Consequently, we interpret the relaxation time τ _{a} as a relative measure of the sensitivity of the chemotaxis system in response to a small external perturbation ΔA _{input}*. We expect the chemotactic response to peak with the relaxation time τ _{a} at about wildtype level of [CheR] and 1/[CheB] (Fig. 2 D and E). For small or large values of [CheR] and/or 1/[CheB], it is conceivable that [CheYP] becomes too small or too large to fall within the narrow functioning range of the motor. Under this extreme condition, the motor does not switch and the system is not chemotactic. In population measurements (8), it is difficult to observe the peaking of the relaxation time because of the inherent celltocell variability of cheR and cheB expression levels (Fig. 2 D Inset and Sect. 5 of SI Appendix ), but this prediction should be testable by the means of singlecell experiments.
To highlight the significance of the relaxation time for chemotaxis, we hypothesize that the average duration of runs varies like the time τ _{a} when cells respond to small input stimuli (Fig. 4 A). We anticipate that the chemotactic drift associated with cells swimming upgradient of nutrients will be larger for wildtype cells with a larger τ _{a} than that of CheR mutants (Fig. 4 B). We tested this hypothesis by performing largescale computer simulations with the agentbased simulator AgentCell (38) (www.agentcell.org). We find that the chemotaxis response decreases with small variations in [CheR] relative to wildtype level (Fig. 4 C). Small changes in the concentration of CheB from the wildtype level produced similar results (Fig. S8 of SI Appendix ). In these simulations, we adjusted the narrow functioning range of the motor so that the CW bias would remain the same in all populations expressing various levels of [CheR] (CW bias = 0.23). In this way, we ensured that the changes in chemotactic response were solely caused by the changes of the relaxation time and not by the differences in tumbling rate.
Because wildtype cells display the largest behavioral variability, they spread farther than mutant cells in the direction perpendicular to the gradient of nutrient (Fig. 4 D). It is then conceivable that the nonlinearity (i.e., Hill coefficient >1) of the adaptation module (10) confers two advantages. First, it enhances the chemotactic drift by producing a long relaxation time. Second, it produces a large behavioral variability that allows a population of bacteria to explore a wider area. These two predictions should be testable experimentally. In real bacterial populations, celltocell variations of [CheR] and [CheB] arise naturally. Given the finetuning of the chemotactic response on [CheR] and [CheBP], we anticipate the existence of molecular mechanisms to limit the variation of the ratio of the methylation–demethylation velocities. For example, cheR and cheB genes are adjacent on the multicistronic meche operon, a strategy known to reduce independent variations between the expression of the two genes (Fig. S6 of SI Appendix ). Another mechanism that helps maintaining the [CheR]/[CheBP] ratio is the existence of the negativefeedback loop on CheB as illustrated in ref. 20.
This study reconciles the presence of large behavioral variability observed at the singlecell level with the chemotactic response of a cellular population. From the design analysis of the chemotaxis system, we showed that the nonlinearity of the adaptation mechanism simultaneously amplifies both the noise and the relaxation time in the chemotaxis system. This approach highlights a key property of the adaptive system in chemotaxis, which consists of a relationship between cellular response and behavioral variability. Largescale simulations of digital bacteria suggest that the chemotaxis network is tuned to simultaneously maximize both the random spread of cells in the absence of attractants and the cellular response to gradients of attractant. The ability to infer the timing of cellular response to small stimuli from the noise in nonstimulated cells should be applicable to a wide range of biological systems.
Materials and Methods
Kinetic System.
Derivation of the full kinetic and stochastic models is available in Sect. 2 of SI Appendix . Ligand binding and conformational changes are much faster than the (de)methylation reactions. Therefore, we describe the activity of receptor complexes using equilibrium probabilities a_{m} (L) (Fig. 1 B) that are functions of the concentration of ligand in the external medium. The subscript m = 0, …, m _{max} indicates the methylation level of the receptor complex. To ensure perfect adaptation, we assume that a _{0}(L) = 0 and a _{mmax}(L) = 1 (18). The methylation–demethylation steps and the time evolution of the concentrations X_{m} of free (not bound to enzyme) receptor complexes with m methyl groups (Fig. 1 B) are governed by (19) Multiplying Eq. 7 with m and summing over all m's gives Eq. 1 , where A* = Σ _{m}a_{m}X_{m} , A = Σ _{m} (1 − a_{m} )X_{m} , and M = Σ _{m}mX_{m} .
Perturbation Analysis of the Adaptation Module.
The perturbation of the kinase activity ΔA* around steadystate Ā* is To keep our analysis independent from various models of receptors, we use ΔA _{input}* directly as the input of the adaptation module without detailing the relationship of Δa_{m} to changes in ligand concentration (Fig. 1 B). The term ΔA _{adapt}* represents the change in kinase activity due to small changes in the methylation levels of the receptors at constant external concentration of ligand, L̅. In previous models of bacterial chemotaxis (16, 19, 20), the probability of activation of receptors complexes at steady state increases approximately linearly with the number of methyl groups when the external concentration of ligand is small (Fig. S7 of SI Appendix ). In our linear perturbation analysis, we take ā _{m} ≈ m/m _{max}. The second term in Eq. 8 becomes ΔA _{adapt}* = βΔM with β = 1/m _{max}. When considering the slow stochastic fluctuations in nonstimulated cell we have δA* = δA _{input}* + βδM ≅ βδM. This linear approximation simplifies the analytical treatment while capturing the basic dependence of the kinase activity on methylation level as established by biochemical data (20, 34, 39). We validate our analytical results with stochastic simulations that include nonlinear activation probabilities a_{m} (L) (see Sects. 4 and 7 of SI Appendix ). In the cases where the activation probability deviates significantly from ā _{m} ≈ m/m _{max}, we can obtain an estimate for the parameter β by expanding the activation probability around the mean steadystate methylation, ā _{m} = ā(M̅) + ā′(M̅)(m − M̅) + … Inserting the expansion in Eq. 8 and neglecting the contributions from small changes in the second and higher moments of the distribution of methyl groups, we obtain β = ā′(M̅) = (∂a/∂m)_{L=L̄, m=M̄} (see Sect. 2.3 of SI Appendix ).
Relaxation Time.
For the case without the CheBP feedback loop, the linear expansion of Eqs. 1 and 2 yields dδM/dt = −δA*/τ_{GK}, where τ_{GK} ^{−1} = k_{r} θ _{r} (1 + θ _{b} )/(1 + θ _{r} ) + k_{b} θ _{b} with θ _{r} = K_{r} ε _{r} /(K_{r} + A̅)^{2}, θ _{b} = K_{b} ε _{b} /(K_{b} + A̅*)^{2} is the inverse relaxation time for only one modification cycle [e.g., the Goldbeter and Koshland system (13)]. Using δA* ≅ βδM calculated above, we then obtain the relaxation time τ_{a} ^{−1} = βτ_{GK} ^{−1}, as well as Eqs. 3 and 6 . The relationship between the relaxation time τ _{a} and the slope of the kinase activation curve (Eq. 5 ) follows from considering the equilibrium relation r̅A̅ = b̅A̅* and its total derivative τ_{GK} ^{−1} dA̅* ≅ b̅A̅*(d ln ε _{r} − d ln ε̄ _{bp} ) (Eq. S30 of SI Appendix ). Without the CheBP feedback loop and at fixed [CheB], we have d ln ε̅ _{bp} = 0 and therefore τ _{GK} ≅ (∂A̅*/∂ ln ε _{r} )_{εbp}/(b̅A̅*) from which Eq. 5 follows. In the system with the CheBP feedback loop, Eq. 5 remains valid, but now τ_{a} ^{−1} ≅ β(τ_{GK} ^{−1} + b̅μ _{a} ), where μ _{a} represents the strength of the CheBP feedback loop (Eqs. S26 and S33 of SI Appendix ). As in Eq. 5 , we can calculate how the relaxation time depends on [CheB] for fixed [CheR]: τ _{a} (ε_{b} ^{−1}) ≅ (∂A̅*/∂ ln ε_{b} ^{−1})_{εr}/(βb̅A̅*μ _{b} ) (Eq. S34 of SI Appendix ). We neglected the receptor complexes bound to CheR and CheBP. The full derivations are in Sect. 2.4 of SI Appendix .
Filtering Properties of the Pathway.
Eq. 6 represents a negative integral feedback system (40). The Fourier transform of Eq. 6 indicates that the frequency response of the adaptation module resembles that of a highpass filter with a cutoff frequency τ_{a} ^{−1}. Similarly, the response regulator module (Fig. 1 A) behaves like a lowpass filter with cutoff frequency τ_{y} ^{−1}, where τ _{y} is the relaxation time of the response regulator module (see Sect. 2.6 of SI Appendix ). Thus, the chemotaxis system resembles a bandpass filter. The value of the low cutoff frequency τ_{a} ^{−1} predicted by our model for wildtype level of [CheR] is in agreement with measurements from tethered cells (25). A recent analysis of the noisefiltering properties of the chemotaxis system suggested that a lower cutoff frequency τ_{a} ^{−1} corresponds to longer time integration of the input signal and therefore better removal of highfrequency fluctuations (21). However, that study did not take into account the slow fluctuations in kinase activity due to stochastic fluctuations in the methylation–demethylation process.
LargeScale Simulations.
Within each digital cell, the chemotaxis network was simulated by using stochastic methods (16, 38) with rates as in Table S1 of SI Appendix . Digital cells for each population were first adapted to a constant external concentration of aspartate of 1 μM for >1,000 s. Then, at time point t = 0, cells were placed at z = 0 in a linear gradient of attractant. This gradient (10^{−8} M/μm) increases along the z direction (L = 1 μM at z = 0) in a 3D infinite medium. The trajectories of individual cells were subject to rotational diffusion with diffusion coefficient D _{r} = 0.062 rad^{2}/s (41).
Acknowledgments
We thank H. Park, S. Needham, C. Guet, P. Oikonomou, T. Shimizu, and R. Alexander for insightful discussions, C. Macal and M. North for computing time on the ANL Jazz Linux cluster; and Eric Lyons for providing Stochastirator. This work was supported by the National Institutes of Health and joint research funding under H.28 of U.S. Department of Energy Contract W31109ENG38.
Footnotes
 ^{†}To whom correspondence should be sent at the present address: Department of Molecular, Cellular, and Developmental Biology, Yale University, P.O. Box 208103, New Haven, CT 06520. Email: thierry.emonet{at}yale.edu

Author contributions: T.E. and P.C. designed research, performed research, and wrote the paper.

The authors declare no conflict of interest.

This article is a PNAS Direct Submission.

This article contains supporting information online at www.pnas.org/cgi/content/full/0705463105/DC1.
 © 2008 by The National Academy of Sciences of the USA
References

↵
 Pedraza JM ,
 van Oudenaarden A

↵
 Raser JM ,
 O'Shea EK

↵
 Raser JM ,
 O'Shea EK
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Fahrner KA
 ↵

↵
 Bialek W ,
 Setayeshgar S

↵
 Kampen NG v.

↵
 Goldbeter A ,
 Koshland DE
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Endres RG ,
 Wingreen NS
 ↵

↵
 Sourjik V ,
 Berg HC

↵
 Block SM ,
 Segall JE ,
 Berg HC
 ↵
 ↵

↵
 Elf J ,
 Ehrenberg M

↵
 Li MS ,
 Hazelbauer GL
 ↵
 ↵
 ↵

↵
 Duke TAJ ,
 Bray D

↵
 Sourjik V ,
 Berg HC

↵
 Mello BA ,
 Tu Y

↵
 Keymer JE ,
 Endres RG ,
 Skoge M ,
 Meir Y ,
 Wingreen NS
 ↵

↵
 Emonet T ,
 Macal CM ,
 North MJ ,
 Wickersham CE ,
 Cluzel P

↵
 Borkovich KA ,
 Alex LA ,
 Simon MI

↵
 Yi TM ,
 Huang Y ,
 Simon MI ,
 Doyle J
 ↵

↵
 Bourret RB ,
 Stock AM

↵
 Shapiro MJ ,
 Panomitros D ,
 Koshland DE, Jr

↵
 Sanders DA ,
 Koshland DE, Jr

↵
 Cluzel P ,
 Surette M ,
 Leibler S