Migration and accumulation of bacteria with chemotaxis and chemokinesis

TL;DR

This paper extends a mathematical model to include bacterial chemokinesis alongside chemotaxis, revealing how combined behaviors influence bacterial migration and improve responses to nutrient sources in complex environments.

Contribution

It introduces a Keller-Segel-type model incorporating chemokinesis and analyzes its effects on bacterial migration dynamics and response to chemoeffector fields.

Findings

Chemokinesis enhances bacterial population response when combined with chemotaxis.

Chemokinesis modifies the migration pattern, potentially reducing bias in nutrient gradient measurements.

The model predicts bacterial behavior around dynamic nutrient sources like marine particles or roots.

Abstract

Bacteria can chemotactically migrate up attractant gradients by controlling run-and-tumble motility patterns. In addition to this well-known chemotactic behaviour, several soil and marine bacterial species perform chemokinesis: they adjust their swimming speed according to the local concentration of chemoeffector, with higher speed at higher concentration. A field of attractant then induces a spatially varying swimming speed, which results in a drift towards lower attractant concentrations -- contrary to the drift created by chemotaxis. Here, to explore the biological benefits of chemokinesis and investigate its impact on the chemotactic response, we extend a Keller-Segel-type model to include chemokinesis. We apply the model to predict the dynamics of bacterial populations capable of chemokinesis and chemotaxis in chemoeffector fields inspired by microfluidic and agar plate migration assays. We find that chemokinesis combined with chemotaxis not only may enhance the population response with respect to pure chemotaxis, but also modifies it qualitatively. We conclude presenting predictions for bacteria around dynamic finite-size nutrient sources, simulating, e.g., a marine particle or a root. We show that chemokinesis can reduce the measuring bias that is created by a decaying attractant gradient.

Figures

Tables

Table 1. Table 1 : Parameters used in simulations. Please see main text for parameter definitions.

Parameter	Value	Figure	Reference
$H$	$0.0$	3, 4	-
	$3.5$	S2, S3, 8	Croze2011
$N$	$0.0$	3, 4	-
	$0.5$	S2, S3, 8	Croze2011
${\delta }_0$	$50.0$	3, 4, S3	-
	$105.0$	S2	Croze2011
$\omega$	$0.2$	3, 4, S3	-
	$0.5$	S2	-
$\eta$	$2.0$	3, 4, S3	-
	$0.5$	S2	Garren2014
$n$	$1$	4, S3	-
	$5$	3, 4,S2	-
	$10$	4	-
	$40$	4	-
$\zeta$	$0.0$	3, 4,S2	-
	$8.164$$\cdot {10}^{-3}$	S3	Croze2011
$K_S$	$1.0$	3, 4,S2, S3	Croze2011
$K_{\chi }$	$0.53$	3, 4,S2, S3	Croze2011
$S$	$0.5$	S3 8

Equations

$$J(x)=-D_b(x)\frac{\partial b}{\partial x}+V(x)b$$

$$\text{centering}\text{with }{D}_b(x)=\frac{v^2}{2\alpha },V(x)=V_k+V_{\chi },\text{@add@centering}$$

$$V_k=-\frac{v}{\alpha }\frac{\partial v}{\partial x}$$

$$V_{\chi }=\frac{v}{2\alpha }({\alpha }_L-{\alpha }_R),$$

$$V_{\chi }=\chi \left(\nabla f_{\chi }+\frac{1}{v}{\partial }_t{f}_{\chi }\right),$$

$$f_{\chi }=\frac{c(x,t)}{c(x,t)+k_{\chi }}.$$

$$\chi =\frac{v^2}{\alpha }\beta ,$$

$$\chi (x)={\chi }_0\frac{v(x{)}^2}{v_0^2},$$

$$v(x)=v_0+v_c\frac{c(x{)}^n}{c(x{)}^n+k_c^n}$$

$$\frac{\partial c}{\partial t}$$

$$\frac{\partial b}{\partial t}$$

$$J(x)$$

$$\frac{\partial C}{\partial T}$$

$$\frac{\partial B}{\partial T}$$

$$J$$

$${\mathsf{V}}_k$$

$${\mathsf{V}}_{\chi }$$

$$\mathcal{V}(X)=1+\eta \frac{C^n}{C^n+{\omega }^n},$$

$$\nabla B$$

$${\nabla }^{2}B$$

$$\frac{n\eta {\omega }^n}{(C^n+{\omega }^n{)}^2}{\mathcal{V}}^{-1}{C}^{n-1}>\frac{{\delta }_0{K}_{\chi }}{(C+K_{\chi }{)}^2}\forall X\in \Omega ,$$

$$n>\frac{4{\delta }_0\omega K_{\chi }}{(\omega +K_{\chi }{)}^2}\left(\frac{1}{\eta }+\frac{1}{2}\right)\forall X\in \Omega :C(X)=\omega .$$

$$\frac{\partial B}{\partial X}+B\left[\frac{1}{2D}\frac{\partial D}{\partial X}-{\delta }_0\frac{d}{dX}\left(\frac{C}{C+K_{\chi }}\right)\right]=0,$$

$$\frac{B}{B^{\ast }}=\frac{{\mathcal{V}}^{\ast }}{\mathcal{V}}\exp \left\{{\delta }_0\left(\frac{C}{C+K_{\chi }}-\frac{C^{\ast }}{C^{\ast }+K_{\chi }}\right)\right\},$$

$$\Delta (t)={\int }_{-\infty }^{t}d{t}^{\prime }K(t-t^{\prime })f_{\chi }(x(t^{\prime })).$$

$$\alpha (t)\approx {\alpha }_0(1-\Delta (t)),$$

$$\tau ={⟨{\int }_0^{\infty }\alpha (t)\exp \left(-{\int }_0^{t}d{t}^{\prime }\alpha (t^{\prime })\right)tdt⟩}_{paths},$$

$$\tau ={⟨{\int }_0^{\infty }{e}^{-{\alpha }_{0}t}\exp \left({\int }_0^{t}d{t}^{\prime }\Delta (t^{\prime })\right)dt⟩}_{paths}.$$

$$\tau =\frac{1}{{\alpha }_0}+{\alpha }_0{\int }_0^{\infty }dt{e}^{-{\alpha }_{0}t}⟨{\int }_0^{t}d{t}^{\prime }\Delta (t^{\prime })⟩,$$

$$\begin{array}{rl}& \tau =\frac{1}{{\alpha }_0}+{\alpha }_0{\int }_0^{\infty }dt{e}^{-{\alpha }_{0}t}\times \\ & ⟨{\int }_0^{t}d{t}^{\prime }{\int }_0^{\infty }duK(u)f_{\chi }(x(t^{\prime }-u))⟩.\end{array}$$

$$f_{\chi }(x(t-u))\approx \text{const}+x(t-u)\nabla f_{\chi }{\mid }_t+(t-u){\partial }_t{f}_{\chi }{\mid }_{x(t)}.$$

$$K(u)=A\delta (u-\theta ),$$

$$\tau \approx \frac{1}{{\alpha }_0}+A{\alpha }_0{\int }_0^{\infty }{e}^{-{\alpha }_{0}t}[I_1(t)+I_2(t)]dt$$

$$\begin{array}{rl}{I}_1(t)& =⟨{\int }_0^t\nabla f_{\chi }x(t^{\prime }-\theta )d{t}^{\prime }⟩\\ & =\nabla f_{\chi }{\int }_0^t⟨x(t^{\prime }-\theta )⟩d{t}^{\prime }.\end{array}$$

Full text

I Introduction

Many bacteria are able to swim by rotating helical flagella distributed on their cell body, and control their swimming pattern by modulating the speed or direction of rotation of their flagellar motors. For example, in the model run-and-tumble motion of Escherichia coli Berg1993 , a bacterium swims approximately straight in a ‘run’ by rotating its flagella in a bundle. When some flagella change their rotation direction, the bundle comes apart and the cell randomly changes direction in a ‘tumble’. In the absence of external bias, this microscopic pattern resembles a random walk, and leads to macroscopic diffusion of a bacterial population. In the presence of a chemical gradient, the random walk is biased, a response known as chemotaxis. As illustrated in Fig. 1(a), a bacterium achieves the biased motion up a gradient of attractant by varying the frequency of tumbles in its random walk: if the bacterium moves up the gradient, the tumbling rate $\alpha$ decreases and, thus, the run length increases, while the swimming speed remains constant.

The run-and-tumble model was originally introduced for enteric bacteria such as E. coli and Salmonella typhimurium Berg2004 ; Tindall2008 , which commonly live in nutrient-rich environments, such as the gut. Marine and soil bacteria, however, often experience heterogeneous and nutrient-scarce environments, and have been found to display different motility patterns. For example, several species living in these harsher environments respond to higher concentrations of attractant by increasing their speed Barbara2003 ; Garren2014 ; Armitage1997 . This response, known as ‘chemokinesis’, modifies the swimming speed in response to the local chemical concentration without affecting the tumbling rate, as illustrated in Fig 1(b). A positive chemokinetic response leads to a higher swimming speed at higher attractant concentrations, whereas a negative response lowers the speed at those concentrations. The strength of positive chemokinetic response can be defined as relative increase in swimming speed over the speed in the absence of chemokinetic effector. A wide range of the chemokinetic response strength has been reported, even for a single species. The responses have been found to vary for the symbiotic soil bacteria Sinorhizobium meliloti and Azospirillum basilense from $7.5$ to $35\%$ Attmannspacher2005 ; Meier2007 ; Sourjik1996 and $40$ to $77\%$ Zhulin1993 , respectively; $7.5$ to $73\%$ for the soil and freshwater purple bacterium Rhodobacter sphaeroides Brown1993 ; Packer1994 ; $26$ to $53\%$ for the enterobacterium E.coli Deepika2015 ; $48\%$ Garren2014 or $6$ to $64\%$ Garren for the marine pathogen Vibrio coralliilyticus. The marine bacterium V.alginolyticus showed an increase of up to $80\%$ upon stimulation with glucose in Son2016 . However, to the best of our knowledge, the pure chemokinetic speed increase as a function of attractant concentration has not been systematically measured for any of these species.

While the role of chemokinesis has been studied extensively in Paramecium spp. and other protozoa vanHouten1978 ; vanHouten1981 ; vanHouten1982 ; Giuffre2011 , the biological significance of the chemokinetic response of marine and soil bacteria has yet to be fully elucidated. Based on the environment that chemokinetic bacteria have been found in, (positive) chemokinesis might be beneficial in heterogeneous environments with scarce sources of nutrients (attractants). For example, alga-sized microbeads coated with various amino acids were used to study the response of marine bacteria to point-like sources of attractants Barbara2003 . All marine bacteria studied were observed to accumulate in bands around the point-like sources while displaying a chemokinetic response. Furthermore, chemokinesis could allow marine bacteria to track algae, helping to foster symbioses with these microorganisms, as well as permitting to respond quickly to short bursts of nutrients, such as those generated from lysing algae Barbara2003a . Another example of a chemokinetic marine bacterium is the coral pathogen Vibrio coralliilyticus. Microfluidic experiments on this pathogenic bacterium in combination with mathematical modelling have suggested that the maximum accumulation in response to chemical cues produced by heat-stressed coral hosts is larger and is reached faster than in the absence of chemokinesis Garren2014 ; Garren . As heat-stressed corals are more susceptible to pathogens, chemokinesis could be a crucial evolutionary advantage in oceans heating up due to climate change. In fact, the chemokinetic response was shown to be even stronger at elevated temperatures; increasing from $6\%$ at $20^{\circ}$C to $64\%$ at $30^{\circ}$C Garren .

Recent interest in chemokinesis has also been sparked by synthetic microswimmers, such as Janus particles. Janus particles are synthetic colloids in a bath of fuel (e.g. $\rm{H}_{2}{O}_{2}$) that propel due to an asymmetric chemical reaction on their surface Walther . These particles show a positive chemokinetic response since their swimming speed increases with increasing fuel concentration Howse2007 , and therefore accumulate in areas of lower fuel concentration.

To date, theoretical work on the combination of positive chemokinesis and chemotaxis has focussed on single-cell level using agent-based models Garren2014 ; Son2016 ; Jackson1987 . The chemokinetic response of the marine pathogen V. coralliilyticus has been modelled as a step increase in swimming speed beyond a threshold attractant concentration Garren2014 . This model was used to analyse the chemotactic response to a transient attractant gradient in a microfluidic device after an initial release of attractant. The chemotactic index (i.e. the enhancement in cell concentration over a control region) suggests that chemokinesis enables a stronger and faster response. This model has been further adapted to include speed dependent changes in the probability of flicking and reorientation frequency of V. alginolyticus Son2016 ; Son2013 . In this particular case, the speed induced changes in motility pattern are responsible for a significant part of the chemotaxis improvement as shown by the agent-based model in Son2016 .

In this work, we use a continuum model to study the spatio-temporal dynamics of bacterial populations with chemokinetic and chemotactic responses. We incorporate chemokinesis into the standard Keller-Segel model for chemotaxis by deriving the model from microscopic run and tumble dynamics. The model is then used to obtain analytical conditions for chemokinetic drift, and solved numerically for three different example attractant distributions that are inspired by existing experimental systems.

II Model

We derive a model for chemotaxis in combination with positive chemokinesis by considering a one-dimensional system in which a cell can move either to the right or left with speed $v$. In this system, a right (left) moving particle changes direction with rate $\alpha_{R}$ ($\alpha_{L}$). Following previous approaches Schnitzer1993 ; Othmer2002 , the one-particle probability density for bacteria evolves according to $\partial b/\partial t=-\partial J/\partial x$, and the bacterial flux $J$ can be derived as

[TABLE]

[TABLE]

where $\alpha_{R}+\alpha_{L}=2\alpha$, and where we have defined the chemokinetic

[TABLE]

and chemotactic

[TABLE]

drift speeds, respectively.

We will now connect the bacterial flux (1) to the commonly used Keller-Segel model of chemotaxis, adapted to the description of chemokinetic populations in dynamic environments.

II.1 Chemotaxis

For chemotactic populations, the drift-diffusion flux (1) is coupled to a chemoattractant density field $c$ via the chemotactic drift speed, given by Eq. (4). In the standard Keller-Segel model it is phenomenologically asserted that this chemotactic drift speed is proportional to the change in the chemical attractant in space, $V_{\chi}=\chi\nabla f_{\chi}$, where $\chi$ is the chemotactic sensitivity parameter and $f_{\chi}$ is a function of $c$ that ensures that the chemotactic drift is biased towards higher attractant concentrations Cates2012 . This definition of the chemotactic drift speed assumes that the chemical attractant profile changes in space but not in time. However, as mentioned, known chemokinetic bacteria inhabit dynamic environments, such as in the ocean or soil. As pointed out by Hein et al Hein2016 , the effective gradient perceived by a bacterium changes in a temporally varying attractant profile depending on the direction of its run. Consider a source of attractant at one point of space, and the associated negative gradients of concentration as one gets away from it. If these gradients are steady in time, for example as in a microfluidic setting, a bacterium exploiting the attractant landscape will correctly detect the gradient and swim towards the source, as in standard chemotaxis. However, if the source corresponds to a single point-like release of attractant, and the concentration at the source position also decays due to diffusion, a bacterium travelling toward the source perceives a smaller increase (or even a decrease) in concentration compared to the steady case. On the other hand, a bacterium moving away from the source perceives a decrease, which is reduced in magnitude compared to the steady case. We show in Appendix A how this influences the mean run length of a bacterium, and, in light of this, modify the chemotactic drift speed to

[TABLE]

where we see that the presence of a temporally increasing (decreasing) gradient increases (decreases) the chemotactic drift speed. Further, we note that the perturbation to the chemotactic drift speed is smaller the larger the value of the swimming speed. For the chemotactic function $f_{\chi}$, we may choose, for example Croze2011 ,

[TABLE]

Following DeGennes2004 , using relationship (4) for the chemotactic speed, the microscopic swimming parameters $v$ and $\alpha$ can be related to the macroscopic parameter $\chi$ (see Appendix A for details) as

[TABLE]

where $\beta$ is constant dependent on the chemotactic response (again, see Appendix A for details).

II.2 Chemokinesis

In this section we consider how chemokinesis can be modelled. We can modify Eq. (7) to include chemokinesis (a spatially varying swimming speed) to give

[TABLE]

in which the subscript refers to the parameters in the absence of chemokinesis (the interested reader can find a detailed derivation in Appendix A). Here we assumed that the swimming speed is constant during a run.

In addition to modifying the chemotactic drift $\chi$, the spatially varying swimming speed of chemokinetic bacteria causes an additional chemokinetic drift to arise with a speed given by Eq (3). When acting alone, this drift drives the cells towards regions of lower speed Schnitzer1993 .

Next we consider how to quantify the chemokinetic coupling between speed and local concentration of attractant. We need to make assumptions about the relationship between speed and attractant concentration, as experimental studies have not been carried out to provide this. Firstly, we assume that cells swim at a base level speed, $v_{0}$. Secondly, we reasonably posit that chemokinesis monotonically increases the swimming speed up to a maximum speed denoted as $v_{0}+v_{c}$. The dynamics between these limits are given by an unknown function characterising the chemokinetic response. In the following, we choose a Hill-type equation

[TABLE]

to approximate the chemokinetic response. The Hill parameter $n$ allows us to introduce an inflection point and change the gradient $\partial v/\partial c$. As can be seen from Fig 2, the speed increases monotonically with increasing attractant concentrations $c$ for all $n$ and $v_{c}>0$, where the half-maximum speed is reached at the attractant concentration $k_{c}$. For any $n$, $v=v_{0}+k_{c}/2$ at $c=k_{c}$ as all functions covered by Eq. (9) have the same half-saturation constant. Note that upon setting $n=1$, one recovers a Michaelis-Menten type response; for $n\rightarrow\infty$, Eq. (9) approaches a step function, which has been used previously to approximate the chemokinetic response Garren2014 ; Son2016 .

Finally, in situations where bacterial growth cannot be ignored (e.g. migration across agar plates considered in section III.2), bacteria are also assumed to undergo logistic growth, which comprises growth and death terms Croze2011 . We note, however, that the chemoeffectors eliciting chemotaxis and chemokinesis need not in general be metabolisable nutrients that induce growth. The full model equations for the chemical attractant field and the chemotactic and chemokinetic bacterial population field, therefore, are

[TABLE]

where $k_{g}$ is the maximum growth rate, $Y$ is the bacterial yield, and $k_{b}$ is the carrying capacity. Furthermore, $g(c)$ is chosen as a Monod-type growth function $g(c)=c/(c+k_{s})$ with the half-saturation constant $k_{s}$. Note that $\chi(x)$ is given by Eq. (8), $f_{\chi}$ by Eq. (6) and $v(x)$ by Eq. (9).

Let us summarise the effects of chemokinesis in our model. A spatially varying speed affects all three terms of the bacterial flux $J$ in Eq. (10c) as: (i) there are regions with a higher diffusivity since $D_{b}\sim v^{2}$ (first term); (ii) it introduces a drift where $\partial v/\partial x\neq 0$ (second term); and (iii) there are regions with a larger chemotactic drift as $\chi\sim v^{2}$ (third term).

II.3 Non-dimensionalisation

The system of partial differential equations (10) is non-dimensionalised using the characteristic time and length scales $t_{0}=k_{g}^{-1}$ and $x_{0}=\sqrt{t_{0}D^{0}_{b}}$, where the bacterial diffusivity is $D^{0}_{b}=v_{0}^{2}\alpha^{-1}$. We rescale the attractant and bacterial densities by their respective initial densities, $c_{0}$ and $b_{0}$. The system of PDEs in dimensionless form thus reads

[TABLE]

with non-dimensional parameters $B=b/b_{0}$, $C=c/c_{0}$, $N=D_{c}/D^{0}_{b}$, $H=b_{0}/(Yc_{0})$, $K_{S}=k_{s}/c_{0}$, $\eta=v_{c}/v_{0}$, $\omega=k_{c}/c_{0}$, $\delta_{0}=\chi_{0}/D_{b}^{0}$, $\zeta=v_{0}/(\alpha x_{0})$ and $K_{\chi}=k_{\chi}/c_{0}$, and non-dimensional time $T=t/t_{0}$ and space $X=x/x_{0}$. The non-dimensional speed function in Eq. (11) is

[TABLE]

where $\eta$ is the maximum increase in swimming speed, $\omega$ is the attractant concentration at which the half maximum speed increase is reached, and $n$ is the Hill parameter. Note that $\omega$ corresponds to the inflection point of Eq. (12). The chemokinetic response is positive for $\eta>0$.

Finally, the model is extended to the 2D axisymmetric case by introducing polar coordinates as

[TABLE]

and the equivalent equations for the chemical field, $C$. The details of the numerical solution and simulations of the model are described in appendix B, including a summary of parameters used. The parameters for the results presented in the main text were chosen to best illustrate the chemokinetic effect (see Appendix B for a discussion of parameter values and the Supplementary Information for further simulations with different parameter sets, including $n$). In section III we will present the simulation results for three different types of attractant gradient.

II.4 Condition for dominant chemokinetic drift

Before solving the extended Keller-Segel model numerically, we can use Eq. (11) to analytically derive a condition on the relative importance of chemokinetic and chemotactic contributions to the drift of the bacterial population. The drift due to a spatially varying swimming speed causes cells to accumulate in regions where they have low speeds. For $\eta>0$, by construction of the velocity function (12), the speed is low at low attractant concentrations. The chemotactic drift, on the other hand, is directed towards higher attractant concentrations by virtue of Eq. (6). Hence, the bacterial density is governed by two competing drifts, as can be seen from the opposing signs in Eq. (11d) and (11e). If the chemokinetic drift is larger than the chemotactic drift for a large part of the spatial domain, this could lead to accumulation at low attractant concentrations, instead of the biologically-desirable accumulation at high concentrations. Assuming a stationary and linear attractant profile (i.e. $\partial C/\partial T=0$ and $\partial C/\partial X=\mathrm{const}$), we have from Eq. (11d) and (11e) that the chemokinetic drift is larger than the chemotactic drift if

[TABLE]

where $\Omega$ is the spatial domain. In the case of a linear attractant profile, we know from Eq. (12) (and Fig. 2) that the gradient $\partial V/\partial C$ is maximum close to the half-saturation constant $\omega$. Thus, we evaluate condition (14) with $C=\omega$, which yields the condition

[TABLE]

If the Hill parameter $n$ exceeds this threshold, the chemokinetic drift is predicted to be larger than the chemotactic drift at the attractant concentration $C^{*}=\omega$ for $\eta>0$. Note that, for a step function, condition (15) is always met at the threshold concentration $C^{*}=\omega$ since $n\rightarrow\infty$. Conversely, if $\eta<0$ (i.e. modelling a negative chemokinetic response), condition (15) will never be met as the chemotactic and chemokinetic drift have the same direction (chemokinesis in this case is stabilising).

III Results

III.1 Steady linear attractant profile

The full model Eq. (11) includes the effect of growth and consumption as well as chemotaxis and chemokinesis. Thus, any change in the spatial distribution of the bacterial population due to chemotaxis and chemokinesis will feed back onto the attractant distribution due to consumption by the bacteria. In order to identify the influence of chemokinesis without such additional complications, we first solve the model for a steady attractant gradient, i.e. $\partial C/\partial T=0$. We furthermore assume a linear profile such that $\partial C/\partial X=\mathrm{const}$ and ignore consumption, and thus population growth. This situation might be achieved experimentally in microfluidic devices Rusconi2014 , where the gradient may be fixed and bacterial growth can be neglected on experimental timescales short compared to growth timescales $\sim k_{g}^{-1}$.

Fig. 3 compares the response of a purely chemotactic population to the response of a chemotactic-chemokinetic population. Chemokinesis leads to a stronger and faster accumulation than in the purely chemotactic case. At the critical concentration, $C=\omega$, however, the chemokinetic drift holds back a subset of the population because it is directed towards lower attractant concentrations as described in section II.4. This can be seen in the form of an accumulation of cells at low attractant concentrations. As the chemotactic sensitivity parameter, $\delta_{0}$, is large in this simulation, the population subset overcomes the drift and accumulates at high attractant concentrations at long times. However, if condition (15) is fulfilled, the chemokinetic drift is larger than the chemotactic drift. Thus, there is a subpopulation driven to small attractant concentrations by chemokinesis. The effect of varying the Hill parameter $n$ in Eq. (12) is illustrated in Fig. 4. For the parameters chosen in Figure 3, condition (15) is met for $n>39.78$. In experiments, this would require observing the transient bacterial concentration profiles in addition to the commonly reported steady-state profiles.

As the attractant concentration is fixed, we can determine a steady-state for the bacterial population, i.e. by setting $\partial B/\partial T=\partial J/\partial X=0$. Due to the homogeneous Neumann boundary conditions of the problem, we have $J=0$ in Eq. (11), which yields the differential equation

[TABLE]

where we used $D(X)=\mathcal{V}(X)^{2}$. This equation can be integrated to give the steady-state

[TABLE]

where $B^{*},C^{*},V^{*}$ are reference values at a chosen reference point $X^{*}$. It is thus clear that, in addition to the influence on the dynamics, chemokinesis affects the steady-state solution via the term $\mathcal{V}^{*}/\mathcal{V}$, where $\mathcal{V}$ varies in space due to chemokinesis. In the case of chemokinesis but no chemotaxis (i.e. $\delta_{0}=0$), the steady-state is determined by the inverse of the speed distribution, i.e. the bacteria accumulate at low speed, as expected and shown previously (Schnitzer1993, ; Othmer2013, ; Arlt, ). For non-zero $\delta_{0}$, if the speed is uniform in space, $\mathcal{V}^{*}/\mathcal{V}=1$ and Eq. (17) reduces to the chemotactic steady-state solution. The exponential term in Eq. (17) represents the chemotactic contribution to the steady state, which does not depend on the swimming speed. Thus, the increase in chemotactic sensitivity $\chi$ (see Eq. (8)) must be balanced by the increase in diffusivity at steady state in a fixed chemical gradient. However, chemokinesis still affects the steady state via the term $\mathcal{V}^{*}/\mathcal{V}$. This chemokinetic effect in the steady-state may only be detectable in experiments with small $\delta_{0}=\chi_{0}/D_{b}^{0}$, since the chemotactic exponential term will dominate $\mathcal{V}^{*}/\mathcal{V}$ for large $\delta_{0}$.

III.2 Self-generated gradient: agar plate

We now consider an evolving attractant field with consumption and growth of bacteria. The attractant is initially uniformly distributed in a 2D axisymmetric setting. This set-up is reminiscent of the classical agar plate experiments, in which bacteria are inoculated in the centre of a nutrient agar plate, see e.g. Adler66 ; Croze2011 . While growing and consuming nutrient, the population creates a gradient of attractant, which it then follows outwards in a chemotactic wave. The attractant profile is a travelling wave itself, and we assume here that the profile relative to the bacterial travelling wave is stationary, i.e. $\zeta=0$ in the chemotactic drift (11e).

In Fig. S2, we compare two chemotactic populations to a chemotactic-chemokinetic population. The chemotactic populations travel at a constant speed, either $\mathcal{V}=1$ (blue curve) or $\mathcal{V}=1+\eta$ (purple curve). Both populations develop a sharp travelling wave, with a larger wave speed for the population at speed $\mathcal{V}=1+\eta$. The chemotactic-chemokinetic bacterial population, on the other hand, develops a broader wave profile. The peak of the wave front is smaller and is followed by a plateau. This effect is more pronounced at late times, as can be seen in the third lower panel of Fig. S2. The reduced pulse also travels slower than the pulse of the chemotactic population at elevated speed $\mathcal{V}=1+\eta$ because the front speed scales with the number of bacteria in the pulse Seyrich2019 . This observation might explain why in agar plate experiments testing for chemotaxis, chemokinetic species such as Sinorhizobium meliloti lack the sharp bands Platzer1997 ; Meier2007 , which are known to be a hallmark of chemotaxis for other species, e.g. E. coli Adler66 ; Croze2011 .

Chemokinesis confers an additional biological advantage in the form of increased population growth as can be seen in Fig. 6, which shows the integrated number of cells over time. At any point in time, the chemotactic population is smaller than the chemotactic-chemokinetic population. However, the chemotactic population at elevated speed, $\mathcal{V}=1+\eta$, (purple curve) shows a stronger population growth than the chemotactic-chemokinetic population due to the faster travelling pulse, which is caused by the increased swimming speed. This increase in swimming speed is associated with metabolic cost Martinez2014 . Thus, permanently swimming faster independent of the attractant concentration could be a beneficial strategy, if the metabolically available energy is not constrained by nutrient supply. When nutrient concentrations are low, on the other hand, increasing swimming speed provides no benefit to bacteria and metabolism is a limiting factor. Chemokinesis could provide an advantageous speed enhancement when it is both metabolically affordable and beneficial Barbara2003 ; Barbara2003a . While the situation considered in this section assumed an abundant supply of chemoeffectors (and in this case nutrients), the next section will consider the response to a transient burst of chemoeffectors.

III.3 Transient source

A localised burst of chemoeffector may, e.g., occur in the sea if algae/phytoplankton lyse and release their content, as has been recently studied in the laboratory Smriga2016 , or when marine particles exude plumes of chemoeffector Stocker4209 . In soil, plant roots exude sugars and other potential nutrients, which locally create a high concentration of chemical attractants. In the following, we consider a single strong pulse of chemoeffector originating from a finite-size axisymmetric source that dissipates via diffusion, modelling a potential dynamic environment around roots or marine particles. The attractant profile that develops is $C(R,T)=S(4\pi NT)^{-1}\exp(-R^{2}/4NT)$, with $S$ representing the amount of chemoeffector contained in the pulse in non-dimensional units. The attractant profile is illustrated in Fig. S3 (top row), while the bacteria are initially uniformly distributed in the domain at concentration $0.2$. To model the response to such a transient attractant profile, we need to include the chemotactic drift velocity, Eq. (11e) with $\zeta>0$, modified to account for the effective gradient perceived by bacteria as they traverse the temporally varying pulse (see Appendix A).

As can be seen from Fig. S3, bacteria with chemotaxis and chemokinesis display a faster and stronger response to a chemoeffector pulse than those with chemotaxis alone. This strong accumulation occurs in spite of the fact that, because of chemokinesis, diffusivity close to the source is higher for these bacteria. For the transient pulse under consideration, temporal variations in the chemoeffector concentration ($\sim\partial_{t}f_{\chi}$) need to be considered when modelling the chemotactic response. Indeed, the run duration for a bacterium travelling up/down the gradient is $\tau^{L,R}\propto\pm(v\nabla f_{\chi}+\partial_{t}f_{\chi})$ (see Appendix A). Thus, when $v$ is higher, as for chemokinetic bacteria, the effect of temporal variation on the chemotactic bias of tumbles is reduced, increasing the accuracy of the chemotactic response. This effect contributes to accounting for the stronger accumulation of chemotactic-chemokinetic bacteria. To illustrate this further, in Fig. 8 we plot the number of bacteria accumulated at the attractant source, $B_{S}$ (the maximum of the bacterial profiles shown in Fig. S3), as a function of time, for chemotactic and chemotactic-chemokinetic populations. The plot displays bacterial accumulation when the temporal perturbation to the chemotactic response is included in the model ($\zeta>0$ in (11e)) and when it is not ($\zeta=0$). In the case of a purely chemotactic population, it can be seen that the predicted amount of bacteria accumulated at the source is lower for a model that ignores the temporal perturbation than for one that includes it. Chemokinesis, on the other hand, reduces the relative effect of temporal perturbation so much that there is very little difference between model predictions with $\zeta=0$ and $\zeta>0$.

These results suggest that a chemokinetic population might be able to overtake purely chemotactic competitors in response to a sudden nutrient release. While the difference observed in Fig.8 of at most $\sim 12\%$ may seem small, a corresponding boost to the growth rate can be sufficient to outcompete a purely chemotactic strain within a few generations. In such transient nutrient landscapes, we have further shown that chemokinesis can reduce the adverse effect that a temporal change in attractant profile can have on the chemotactic response.

IV Discussion

Chemokinesis is a known response for many environmentally relevant bacteria, yet its consequences for bacterial population dynamics have been little explored. In this work, using a modified Keller-Segel continuum model we have shown how chemokinesis significantly affects both the dynamics and steady-states of bacterial populations capable of chemotactic and chemokinetic behaviour. Our model incorporates the effects of a concentration dependent speed, including an increase in the chemotactic sensitivity, and a recently suggested modification to the chemotactic response in dynamic environments Hein2016 , which was derived adapting the microscopic model first suggested by de Gennes DeGennes2004 .

We have solved our model numerically to explore the effect of chemokinesis on migration and accumulation in experimentally realistic gradients. In a fixed attractant gradient, our results show that chemokinesis can lead to two subpopulations travelling at different speeds, with the slower one being held back by the chemokinetic drift. In the case of agar plate migration, where bacteria inoculated onto the plate generate their own gradient by consuming nutrients, we find that populations with chemokinesis migrate out from the inoculation point in waves that are faster, but broader than purely chemotactic migrating populations. While the increase in front speed could be explained by a population at a uniformly increased swimming speed, the broadening is only observed in the chemokinetic population. It is a new feature not predicted by previous studies using agent based simulations Son2016 . It is possible that this broadening might explain why the classic chemotactic Adler bands observed for E. coli Adler66 are not observed for chemokinetic soil bacteria such as Sinorhizobium meliloti Platzer1997 ; Meier2007 . Furthermore, chemokinesis increases the population growth significantly in comparison to purely chemotactic migration.

Our simulations also considered the case of a transient source of nutrients, e.g., a lysed algal cell. In this case, our results show how bacteria with chemokinesis and chemotaxis accumulate faster and more strongly around the source, while concentrations of nutrients are high, with respect to purely chemotactic bacteria. This chemokinetic advantage is both due to the enhanced migration discussed in the previous examples, but also to the fact that chemokinesis mitigates the perturbation to the chemotactic response due to the transient nature of the attractant profile. We note that, while our model includes the effect of transient chemical fields on tumbles, it does not include a recently modelled effect of chemokinesis on the precision of chemosensing Brumley10792 . It will be interesting to include this additional effect, which could lead to further enhancements in chemokinetic accumulations, in future elaborations of our model.

The role of metabolism is an important consideration for chemokinesis. For example, positive chemokinesis might be caused purely by an increase in nutrient availability. A resulting increase in the energy level of the cell may lead to more energy being available for flagella rotation, which allows the cell to swim faster. As the chemotactic sensitivity scales with the swimming speed as $\chi\propto v^{2}$, a faster swimming population will always show a stronger chemotactic drift. However, swimming is associated with a considerable metabolic cost Martinez2014 ; Mitchell2002 ; Ni2020 . Indeed, the energetic cost of swimming increases quadratically with the swimming speed Taylor2012 . Therefore, it might be beneficial to swim faster (and, thereby, improve chemotaxis) only if favourable nutrient conditions are available. While the energy level of a cell may influence the extent of the chemokinetic response, there are several indications of metabolism-independent chemokinesis for different bacterial species Deepika2015 ; Garren2014 ; Packer1994 . In a dynamic environment such as considered in section III.3, chemokinesis can then help to reduce temporal bias and improve the chemotactic response. As pointed out by Hein et al, the threshold for detecting absolute concentration is smaller than for gradient detection Hein2016 . Thus, chemokinesis can take place at lower background concentrations to improve chemotaxis only when needed.

The predictions of our model include interesting qualitative effects, which have not previously been observed in agent-based models: the slower subpopulation in the fixed attractant profile, and a broadening of the travelling wave in a self-generated gradient. To test these predictions experimentally, chemotaxis and chemokinesis should be addressed independently. For example, a recently developed E.coli system with a swimming speed that is controlled via light Arlt could be used to engineer populations with a swimming speed that can be controlled independently of chemotaxis.

Future theoretical investigations would benefit significantly from the experimental measurement of the chemokinetic response function relating swimming speed and local attractant concentration. In this work we assumed a smooth change from a reference speed to an increased speed, where the degree of change in swimming speed changed with a single parameter. As we have shown, a very steep change in swimming speed (e.g. in form of a step change as assumed previously in agent-based models Garren ; Son2016 ), could actually inhibit chemotaxis rather than promote it. Experimental work so far has been restricted to measurements at very few different attractant concentrations, which makes it difficult to deduce a functional relationship between speed and attractant concentration. Thus, further work is required to determine the function $v(c)$ for chemokinetic bacterial species. Such measurements would also allow elucidation of the rate of adaptation, i.e. how quickly the swimming speed adapts to its new value both for an increase and decrease in attractant concentration, and which we have here assumed instantaneous. The increase in the population-averaged swimming speed in response to a uniform addition of effector occurred on the time-scale of $100-200$s for V. alginolyticus Son2016 . In Packer1994 , on the other hand, the chemokinetic response of R. sphaeroides was measured within $10$s upon uniform addition of chemoeffector, at which point the swimming speed had already adapted to its increased level. Furthermore, a desensitisation to a sustained higher level of attractant, as observed for chemotaxis, does not seem to occur as the swimming speed remained at elevated levels for hours in R. sphaeoroides and A. brasilense Packer1994 ; Zhulin1993 . The discrepancy in the order of magnitude in the response times might be caused by the experimental set-ups, e.g. in Son2016 , the effector first needs to diffuse from two sides in a wide microfluidic channel before a uniform population response can be measured. To conclude, further experiments could shed some light on the chemokinetic response function and adaptation time, which would benefit the further development of the model and its assumptions to understand how bacteria make use of chemokinesis in dynamic environments.

Data accessibility

Code and data are available from the Zenodo repository: https://doi.org/10.5281/zenodo.4323421.

Acknowledgements

We acknowledge financial support from EPSRC EP/L504920/1 and EP/N509620/1 (T.J.) and the Winton Programme for the Physics of Sustainability (T.J., O.A.C.), and funding from the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 798411 (F.J.P.).

Authors’ contributions

TJ carried out analytical calculations and simulations. JL carried out initial calculations and simulations. FP assisted with initial simulations. TJ, FP and OC interpreted the simulations. TJ and OC wrote the manuscript, which was critically revised by FP. All authors gave final approval for publication and agree to be held accountable for the work performed therein.

Appendix A Chemotactic drift and sensitivity

Here we derive new relations for the chemotactic drift speed $V_{\chi}$ and sensitivity parameter $\chi$ for bacteria undergoing chemotaxis and chemokinesis in dynamic environments, connecting the macroscopic parameter $\chi$ with the microscopic swimming speed $v$ and tumbling rate $\alpha$. We follow the approach of de Gennes DeGennes2004 ; Locsei07 ; Croze2011 , considering run-and-tumble bacteria with a simplified E. coli-like chemotactic response. Our contribution for the purpose of this work is to include: i) a swimming speed that is a function of position; ii) a temporally varying contribution to the gradient perceived by the bacteria. As in DeGennes2004 ; Croze2011 , we consider a one-dimensional model of run-and-tumble bacteria propagating in gradient of attractant, to which bacteria are exhibiting a ‘small response’, that is the tumbling response remains close to the adapted value Xue2015 . Chemotactic memory is modelled through kernel integral, reflecting the fact that bacteria ‘remember’ their chemical environment over a characteristic delay time of a few seconds. We neglect directional persistence and rotational diffusion, which can be addressed in the same framework Locsei07 , but are here ignored for simplicity. The modified de Gennes model we employ could be applicable to chemokinetic species like rhizobia, which display run and tumble dynamics Armitage1997 . The model neglects features of the motility pattern of chemokinetic marine bacteria, such as a run-reverse-flick mechanism or the influence of a change in swimming speed on reorientation frequency, as observed for V. alginolyticus Son2016 , but instead focus on the pure effect of speed change.

As the tumbling frequency follows an exponential distribution, tumbling events can be treated as independent events of a Poisson process with rate $\alpha(t)=\alpha_{0}e^{-\Delta(t)}$, where $\alpha_{0}$ is the tumbling rate in the absence of a gradient and where $\Delta(t)$ is the chemotactic bias given by the memory integral

[TABLE]

The concentration function was chosen as in the main text as $f_{\chi}=c(x,t)/(c(x,t)+k_{\chi})$. The memory kernel $K(t)$ was first measured for E. coli Segall and more recently also for V. alginolyticus Xie2015 . Based on these experiments we assume that the kernel obeys $\int_{0}^{\infty}K(t)dt=0$, i.e. the tumbling rate perfectly adapts. As mentioned, we shall assume a small response ($\Delta(t)\ll 1$), so that the bacterial tumble rate can be linearised to

[TABLE]

Considering a run that starts at $t=0$, the probability density of a tumble event in the interval $[t,t+dt]$ is $\alpha(t)\exp(-\int_{0}^{t}dt^{\prime}\alpha(t^{\prime}))$. The mean run duration is given by

[TABLE]

where angled brackets denote averaging over all possible bacterial swimming paths, as the tumble rate is path dependent (the ‘paths’ subscript will be omitted henceforth). Integrating (20) by parts and recalling (19) we can write

[TABLE]

Next, since $\Delta(t)\ll 1$, we can linearise the exponential integral to obtain

[TABLE]

where we have brought the angled brackets inside the time integral to surround only path dependent quantities. Substituting the expression for the bias (18) and performing a change of variables in the memory integral (19) by defining $u=t-t^{\prime}$, we then obtain

[TABLE]

Next, as the time interval of interest is small compared to the gradient variations, we can Taylor expand the concentration function about a reference position and time:

[TABLE]

We note that the constant term does not influence the integral in (23), since the response function $K(u)$ integrates to zero. Further analysis is simplified by considering the special response function

[TABLE]

where $\theta$ is a single delay time DeGennes2004 . We note that this function does not represent a physical response and is considered solely for the purpose of facilitating the calculation. Then, substituting equations (25) and (24) into (23) and integrating over the delta function we obtain

[TABLE]

where we have defined

[TABLE]

The average over paths can be taken inside the time integral since tumbles for $t>0$ are treated as having no effect on cell motion, so that paths are independent of time on the interval considered. Proceeding similarly we define

[TABLE]

To evaluate (27), we recall that the bacteria we are considering do not undergo rotational diffusion and do not possess directional persistence. Thus, following DeGennes2004 , we see that for times preceding a run (when a bacterium is tumbling) the position of a bacterium is, when averaging over paths, not correlated to the velocity. On the other hand, during a run, the position correlates with velocity. In this way, in (27) we have

[TABLE]

So that integral (27) becomes:

[TABLE]

We proceed analogously to evaluate integral (28). In this case, the integral over paths for times preceding a run does not average to zero, but sums up the contributions from the temporal variation of the gradient during tumbles. During a run, on the other hand, the temporal variation is evaluated for when bacteria are travelling up (down) the gradient, providing positive (negative) weights to the path integration, which provides

[TABLE]

where $F(t^{\prime},\theta)$ is a function of time and delay. In this way, we can carry out integral (28)

[TABLE]

where we have defined $G(\theta,t)=\int_{0}^{\theta}F(t^{\prime},\theta)(t^{\prime}-\theta)\,dt^{\prime}$. Because it does not depend on direction up/down the gradient, this function (and its integrals) will cancel out in the evaluation of the run time, and will henceforth be omitted from the derivation. Inserting expressions (30) and (A) back into (26), we can thus write the mean run durations up ($R$) and down ($L$) the gradient as

[TABLE]

which can be integrated to obtain

[TABLE]

Extending to general response function $K(\theta)$ with a distribution of delays DeGennes2004 , equation (34)

[TABLE]

Since $\alpha_{L,R}=1/\tau_{L,R}$ and (see main text)

[TABLE]

we then arrive at an expression for the chemotactic drift

[TABLE]

where we have defined the constant $\beta=\int_{0}^{\infty}d\theta K(\theta)e^{-\alpha_{0}\theta}$. The standard definition of the chemotactic sensitivity is obtained from the empirical drift $V_{\chi}=\chi\nabla f_{\chi}$, where as previously $\chi$ is the chemotactic sensitivity. We extend here the definition to include temporal gradients, and define

[TABLE]

so that chemotactic drift is given by Eq. (5) in the main text. In the absence of chemokinetic alterations to the swimming speed, the sensitivity is simply

[TABLE]

where subscripts denote a constant swimming speed. Dividing (38) by (39) provides relationship (8) in the main text.

Appendix B Numerical solution and parameter values

A finite difference scheme was chosen to compute the numerical solution of the system of PDEs (11). The domain $\Omega$ is represented by a vector of $M$ equally spaced grid points, $X_{1},...,X_{M}$. A fourth order scheme in space and a first order scheme in time are used to approximate the derivatives. Homogeneous Neumann boundary conditions are imposed for both substrate and bacterial concentration, i.e. $\frac{\partial C}{\partial X}|_{X=1,M}=0$ and $\frac{\partial B}{\partial X}|_{X=1,M}=0$.

In the case of steady attractant gradient (Fig. 3 and Fig. 4) the initial attractant distribution is $C(X,0)=0.01X$, while the bacteria are uniformly distributed $B(X,0)=1.0$. For the self-generated gradient (Fig. S2 and Fig. 6) we assume a Gaussian bacterial inoculum, $B(R,0)=\exp(-R^{2}/\sigma^{2})$, with a corresponding attractant distribution, $C(R,0)=1-\exp(-R^{2}/\sigma^{2})$ (Croze2011, ). Finally, for Fig. S3 and Fig. 8 the bacteria are initially uniformly distributed $B(R,0)=0.2$, while the attractant distribution is $C(R,0)=S(4\pi NT)^{-1}\exp(-R^{2}/4TN)$, where $S=0.5$ and $T=0.02$.

The parameters used in the non-dimensional model are summarised in Table 1. The parameters are mostly based on literature values for E. coli. Due to the lack of experimental data required for a full estimate of the functional form of the chemokinetic response of, e.g., a single species, the values for the chemokinetic response ($\eta$, $\omega$, $n$) are motivated by several studies and chosen to illustrate the response. For example, $\eta=0.5$ corresponds to a $50\%$ maximum increase in the swimming speed, which is on the order of magnitude that has been reported for several species Packer1994 ; Deepika2015 ; Garren2014 . The largest increase reported, to the best of our knowledge, is $80\%$ for V. alginolyticus Son2016 . Thus, $\eta=2$ is chosen as an extreme value to illustrate the effect of chemokinesis more clearly. The half-saturation constant $\omega$, which is required for the chemokinetic function given by Eq. (12), has not been reported directly. Upon (visual) inspection of results, half-saturation constants may be about $0.1$mM of glucose and $0.1$mM of acetate for E.coli Deepika2015 and R. sphaeroides Packer1994 , respectively. In Garren2014 , the half maximum speed seems to be reached at about $30\%$ relative mucus concentration, while the maximum is reached at about $60\%$. As pointed out in Hein2016 , absolute concentration above a threshold can be detected before gradients can be accurately measured. Therefore, we assumed that the half-saturation constant of chemokinetic response ($\omega$) should be below or on the same order as the chemotactic half-saturation constant ($K_{\chi}$). As there has not been a functional fit to chemokinesis measurements, an estimate for the Hill factor $n$ is difficult. From visual inspection of results in Garren2014 and Deepika2015 we may assume a factor between $1$ and $3$. Simulations with larger Hill factors were performed to compare to results based on the assumption of a step-change in swimming speed (i.e. $n\rightarrow\infty$) as done in agent-based models Garren2014 ; Son2016 . The parameters used to produce the figures in the main text were chosen to illustrate the effect of chemokinesis best. As a comparison, other parameter combinations are given in the Supplementary Material.

SUPPLEMENTARY MATERIALS

Additional simulations

We present further simulation results for the scenarios described in the manuscript by varying crucial and unknown parameters determining the bacterial response. As the chemokinetic response function is unknown, these include the strength of the chemokinetic speed increase $\eta$, the half-saturation concentration of the chemokinetic response $\omega$ and the Hill factor $n$ as well as the strength of the chemotactic response $\delta_{0}$.

Steady linear gradient

The effect of modifying the chemotactic sensitivity $\delta_{0}$ and the chemokinetic parameters $\eta$ and $\omega$ is shown in Figure S1.

Self-generated gradient

The effect of varying the Hill parameter $n$ of the chemokinetic function $V(C)=1+\dfrac{\eta C^{n}}{C^{n}+\omega^{n}}$ is shown in Figure S2(a). Figure S2(b) illustrates the effect of the chemotactic sensitivity $\delta_{0}$ and the chemokinetic parameters $\omega$ and $\eta$.

Transient source

The parameters used in the main text are the same as those used for the steady linear attractant gradient. The results of a simulation that uses the parameters of the self-generated attractant gradient instead are shown in Figure S3.