Chemotactic behaviour of Escherichia coli at high cell density

At high cell density, swimming bacteria exhibit collective motility patterns, self-organized through physical interactions of a however still debated nature. Although high-density behaviours are frequent in natural situations, it remained unknown how collective motion affects chemotaxis, the main physiological function of motility, which enables bacteria to follow environmental gradients in their habitats. Here, we systematically investigate this question in the model organism Escherichia coli, varying cell density, cell length, and suspension confinement. The characteristics of the collective motion indicate that hydrodynamic interactions between swimmers made the primary contribution to its emergence. We observe that the chemotactic drift is moderately enhanced at intermediate cell densities, peaks, and is then strongly suppressed at higher densities. Numerical simulations reveal that this suppression occurs because the collective motion disturbs the choreography necessary for chemotactic sensing. We suggest that this physical hindrance imposes a fundamental constraint on high-density behaviours of motile bacteria, including swarming and the formation of multicellular aggregates and biofilms.

The spatial power spectral density E(q) measured in the 50 µm high channels exhibited a peak at π/q str = 50 µm, as in the ∆f lu mutant, growing with cell density. (b) The maximum value of E(q) grew in exactly the same manner for both strains. (c-d) The normalized chemotactic drift had the same dependence on the cell body volume fraction (c) and on the amplitude of the collective motion E(q str ) (d) in both strains. (d) E(q str ), when corrected for its low density value E min (q str (h)) and normalized to the value at which it saturates at large density E max (h), appears to be a single function of the cell body volume fraction times typical vortex size times swimming speed (Φ c π/q str v 0 ).

Supplementary Note 1. Peclet number
The Peclet number compares diffusion and advection of the chemical attractant in the fluid. It is defined as P e = hv f /D, since the typical size of the flow scales with the height of the channel. Here v f is the velocity of the fluid. It can be estimated from the difference between the high density collective velocity and the free swimming speed assumed by the cell at low density. It never exceeded 10 µm/s (Supplementary Fig. 1).
Taking the maximal h = 50 µm, and considering that D = 500 µm 2 /s, leads to P e ≤ 1 for all experiments, no matter how strong the collective motion is. P e = 1 indicates that distortions of the gradient are to be expected, as we observe in the case of the strongest collective motion.

Supplementary Note 2. Chemotactic drift in shallow gradients
Model of the pathway -The chemotaxis pathway of E. coli is composed of a large superstructure (array) of coupled receptor dimers, embedded in the inner membrane of the cell, and coassembled with the cytoplasmic kinase CheA and adaptor protein CheW. The most aboundant chemoreceptors of E. coli are Tar (the main sensor of MeAsp) and Tsr. Chemoeffectors bind the receptors on their periplasmic domain, inducing conformational changes which cooperatively modify the autophosphorylation activity of the kinase CheA, attractants reducing kinase activity. The cooperative array organization enables to integrate and amplify signals from different receptors. The kinase then transmits its phosphates to the small diffusible molecule CheY. Phosphorylated CheY binds to the motor to induce tumbles and is dephosphorylated by CheZ. The previous reactions are subsecond so that a sudden increase in attractant concentration induces a fast drop in the probability of tumbling. Two enzymes, CheR and CheB, then respectively add and remove methyl groups to specific amino-acids of the respectively inactive and active chemoreceptors in few seconds. This slowly offsets chemoeffectors action, adapting the average kinase activity -and thus the tumbling rate -back to the intermediate value it assumes in homogeneous environments, so that subsequent stimulations can be sensed. The current methylation level of the chemoreceptors then acts as a physical memory which represents the environment experienced by the cell a few seconds before, with which the current situation is compared.
We recall here the model of the chemotaxis pathway constructed in [4]. The chemoreceptor dimers were modeled as two-state variables interacting following a Monod-Wyman-Changeux allosteric model. The probability P on of a signaling team of N a Tar and N s Tsr receptor dimers (and associated kinases) to be active is given by the free energy difference F between the active and inactive states as: with the free energy difference: where K off i (resp. K on i ) is the binding affinity of the chemoattractant, present at concentration c, to the receptor i in its inactive (resp. active) state. The methylation dependent free energy difference ϵ(m) is linear by part, and defined as: The methylation enzyme CheR methylates only inactive receptors with average rate k R , and CheB demethylates only active ones with average rate k B , so that the methylation level evolves according to: The previous set of equations makes the fraction of active teams P on evolve in time, responding and adapting to the history of concentrations experienced by the cell c(t). The balance between autophosphorylation of active kinases and phosphotransfer to CheY sets the fraction of phosphorylated CheA dimers A p as: Finally, the fast phosphorylation -dephosphorylation cycle of CheY sets the fraction of phosphorylated CheY as: If the cells are in the run state, the probability to tumble during the time step δt is given by: If the cell is in a tumbling state, the probability to start running again is however independent of the CheY phosphorylation level: All parameter values are given in Supplementary Table 2 and are as in [4] except for N a and N s which were chosen to match the values expected at the OD 600 to which our cells are grown [5].
Prediction for the chemotactic drift in shallow gradients -We write here in our notations the results of Dufour et al. [6]. From Eq. 3 of this paper, because our simulations are 2D and the adapted value of Y p is 1, we have: were (1 − T B) = τ r /(τ t + τ r ) is the fraction of running cells at any given time in absence of a gradient, τ m = 0.5(N s + N a )P on (1 − P on )(k R + k B ) is the relaxation time of F according to Eq. 2 and 4, 1/τ dec = 1/τ R + 1/τ T , and the derivative is taken at the adapted value F = ln 2. We have also defined the Brownian reorientation time τ R = 1/D r , and the tumbling reorientation time τ T = τ r (Y p ) −H /(1 − exp(−D T τ t )). The latter accounts for the incomplete randomization due to finite tumbling time [7], which was not accounted for in Dufour et al. [6] because complete randomization was assumed. We now considere that: We also have from Eq. 1, 5 and 6 that: since the other factors are very close to 1, when one considers that they have to be evaluated for F = ln 2.
We then have: which corresponds to Eq. 2 (and equivalently 4) of the main text, with the coefficients given in Supplementary Table 2.

Supplementary Note 3. Chemotactic behavior in simulations in absence of hydrodynamic interactions
Supplementary Fig. 7 shows the main characteristics of the collective motility and the drift in the simulations where hydrodynamic interactions are neglected. Note that the rotational Brownian motion was also neglected. The characteristics of the collective motion, as measured by the flow structure factor E(q), were set by the cell length L ( Supplementary Fig. 7b-d). The chemotactic drift first decreased with area fraction, before increasing, peaking and decreasing again for L = 3 µm and 4 µm ( Supplementary Fig. 7e). Contrary to the experiments, there was no scaling with E(q str ) ( Supplementary Fig. 7f). Interestingly, the velocity decorrelation time τ dec extracted from the time autocorrelation of the single cell velocity (Supplementary Fig. 7g) decreased monotonously in all conditions as a function of area fraction ( Supplementary Fig. 10a), as it did when hydrodynamics was accounted for. However, in the dry case, the chemotactic coefficient, considering its dependence in τ dec , did not follow Eq. 4 of the main text ( Supplementary Fig. 7h), contrary to the full simulations (main Fig. 4h). Clearly, at the densities were the drift started to reincrease, the system started to depart from the framework of Eq. 4 of the main text, which assumes a Brownian motion like reorientation process and shallow gradients, i.e. that the average pathway activity P on departs only slightly from its adapted value. Because it occurs at a fairly high cell density, it is however not clear if the peak observed here occurs for the same reasons as the one observed in the experiments.

Supplementary Note 4. Dual effect of physical interactions on the chemotactic drift
An extension of the model we used so far to explain chemotactic drift reduction (i.e. [6] leading to Eq. 2-4 of the main text) was proposed in Long et al. [1]. There, this extension was used to explain non-linear reinforcement of the chemotactic drift when the gradient is sharp enough. It therefore does not assume shallow gradients anymore, and is based on a Schmolukovski equation for the space-averaged probability P (t, F, s) to be at time t with a chemoreceptor free energy F and an orientation s = n · ∇c/||∇c|| [1]: Here time is normalized to the adaptation time τ m and r(F ) is the probability of the cell being running, given the receptor free energy F . The first right hand side term describes the chemoreceptor free energy actuation according to adaptation (term −(F − F 0 ) with F 0 the adapted, unstimulated, value of F ; F 0 = ln 2 in our case) and stimulation due to swimming in the gradient (with 1/Hτ E being the normalized gradient, defined in Supplementary Table 2 for our case). The second term describes reorientations due to Brownian rotational diffusion and tumbles, withL s the rotational diffusion operator and τ D = 1/ (τ m (rD r + (1 − r)D T )) comparing reorientation and adaptation times. Assuming that the detailed balance in the angular fluxes of cell orientation hold, we can include in this second term the effect of collective reorientations in the form of an effective enhanced rotational diffusion D r → D r + D eff coll . Integrating this equation over orientations, in Long et al. [1] it is shown that at steady state in the gradient: Here the chemotactic ratio v ch /v 0 is the product of the normalized gradient strength and of the average over all positions, orientations and times of the shift in chemoreceptor free energy from its adapted value (CFES, ⟨F − F 0 ⟩). Key assumptions for derivating this equation are steady-state and detailed balance holding. In this equation, the effect of reorientations (tumbles, Brownian motion and possibly collective reorientations) on the drift is accounted for by a reduction of the CFES. In the simulations (contrary to experiments), the CFES is readily accessible, and it evolves similarly to the drift as a function of volume fraction ( Supplementary Figs. 9b and 10b). The chemotactic ratio v ch /v 0 follows Eq. 14 as a function of ⟨F − F 0 ⟩ fairly well in the simulations accounting for hydrodynamics ( Supplementary Fig. 9c), but less so for the dry simulations ( Supplementary Fig. 10c). Upon closer inspection using the ratio we find that the ratio decreases in both cases below 1 when collective reorientations increase, as measured by the decrease of the decay time of the autocorrelation of the cell velocity, τ dec (Supplementary Figs. S9d and S10d). This deviation from Eq. 14 is weaker when hydrodynamics is included. Since we ensured that steady state is reached (in 1000 frames, 10 s of real time equivalent) before measuring the averages, we deduce that this discrepancy must come from detailed balance in the angular fluxes not being respected. Indeed, detailed balance assumes that: P (s)ϕ(s → s + ds) = P (s + ds)ϕ(s + ds → s) where ϕ(s → s + ds) is the flux of cells changing their orientation from s to s + ds over a small amount of time, and it is not necessarily satisfied when hydrodynamics-and collision-induced reorientations are concerned. Indeed, vortices are a net flux of orientation which thus cannot satisfy Eq. 16. Therefore, when hydrodynamics is included and even more so in the dry simulations, the effect of cell-cell interactions is not fully comparable to Brownian motion (for which detailed balance holds), and we prefer in this sense to talk about active enhancement of rotational diffusion, where Hτ E ⟨F − F 0 ⟩ represents the drift if detailed balance was holding, and the ratio (15) is a measure of the effect of detailed balance breakdown.