Learning run-and-tumble chemotaxis with support vector machines

Motile biological cells can steer their swimming paths upwards a concentration gradient of chemoattractant molecules in a process termed chemotaxis [1–6]. At the microscopic scale of cells, sensing and motility noise severely limit the precision of chemotaxis [7–14]. Biological cells measure chemoattractant concentration by registering binding events of diffusing chemoattractant molecules to receptors on the cell surface. The number of molecules registered in a short time interval is thus a random variable, resulting in inevitable sensing noise. Brownian motion as well as active fluctuations of cell motility randomize the orientation of the cell, resulting in motility noise [15,16]. In a seminal paper, Berg and Purcell argued that the relative error in sensing a concentration c due to molecular shot noise scales as $\delta c/c \sim (r\,T)^{-1/2}$, where $r\sim D c L$ is the detection rate of chemoattractant molecules, which depends on the diffusivity D of these molecules and the size L of the cell, while T is measurement time [7]. The choice of optimal measurement time implies a trade-off between precision and accuracy: longer measurement times increase the precision of concentration measurements, as measurement noise can be effectively averaged out, but reduce the accuracy as the true concentration dynamically changes in time [17]. For a concentration that fluctuates on a characteristic time-scale τ, the optimal choice of measurement time $T_{\textrm{eff}}\sim \sqrt {\tau /r}$ results in the modified scaling $\delta c/c \sim (r\tau )^{-1/4}$ [18]. An analogous scaling was later found for optimal gradient sensing by spatial comparison, where a stationary cell determines the direction of a spatial concentration gradient by comparing local concentration across its diameter [19].

To determine not only a concentration, but a spatial concentration gradient, biological cells use different strategies [3,20]. Possibly the most natural strategy, spatial comparison, compares concentrations across the diameter of the cell, and is employed, e.g., by crawling immune cells [4]. Intriguingly, many smaller cells such as the bacterium E. coli do not use spatial comparison, but instead determine concentration gradients by temporal comparison, i.e., by estimating the time derivative of the local concentration signal $c(\textbf{r}(t))$, while the cell actively moves along a path $\textbf{r}(t)$ in a concentration gradient [1]. Generally, small and comparatively fast cells $(L\lesssim 10\,\mathrm{\mu } \textrm {m},\,v\gtrsim 10\,\mathrm{\mu } \textrm {m}/\textrm {s})$ employ temporal comparison for gradient sensing, while larger and slower cells often employ spatial comparison [20]. In fact, the signal-to-noise ratio of spatial comparison $(\textrm{SNR} = D L^{3} T\,|\nabla c|^{2}/c)$ depends strongly on agent size L. Moreover, the measurement time T is ultimately limited by the persistence time $1/D_{\textrm{rot}}$, which is set by the rotational diffusion coefficient $D_{\textrm{rot}} \sim L^{-3}$ of Brownian motion for small cells [15]. For micrometer-sized bacteria, gradient sensing by spatial comparison would be simply too noisy. Instead, these cells employ a chemotaxis strategy known as run-and-tumble, where short episodes of approximately straight swimming (runs) are interrupted by random reorientation events (tumbles), during which the cell picks a new swimming direction at random [1]. By estimating the time derivative of the chemoattractant concentration signal during a run, and extending the duration of runs if this temporal change is positive, these cells realize a biased random walk with net motion up the concentration gradient. Several works investigated optimal gradient-sensing strategies in the presence of noise [21,22]. A key question concerns the optimal measurement time T_eff for taking the binary decision whether to tumble or not.

Here, we employ supervised machine learning to determine kernels for optimal decision filters for run-and-tumble chemotaxis using support vector machines. Time-continuous decision filters are effectively represented by finite weight vectors that can be learned efficiently. The learned navigation strategy enables search agents to find a target with a chemotactic index close to that of an ideal run-and-tumble agent, demonstrating that machine learning can identify near-optimal functional behavior. We empirically show that the optimal measurement time follows a power law $T_{\textrm{eff}} \sim D_{\textrm{rot}}^{-\alpha }$ with exponent $0.2,\ldots ,0.3$.

To describe a minimal chemotactic agent, we extend the popular model of an Active Brownian Particle (ABP) [23] by introducing tumbles as random re-orientation events, which are triggered by a decision rule that takes noisy concentration measurements of the ABP along its trajectory as input. Specifically, we consider a point-like agent moving along a trajectory $\textbf{r}(t)$ in two-dimensional space in a spatial concentration field $c(\textbf{x})$ of chemoattractant molecules released by a hidden target at the origin, see fig. 1(a). The diffusion of molecules to an absorbing surface, here the circumference of the agent, is well described by a Poisson point process of arriving molecule with rate $r(t) = \lambda \, c(\textbf{r}(t))$ (see footnote ¹ ). Here, λ is a binding rate constant (that depends as $\lambda \sim D L$ on the diffusion coefficient D of the chemoattractant molecules, and the size L of the agent and its surface coverage with chemoattractant receptors [15]). We write $s(t) = \sum _{k} \delta (t-t_{k})$ for the train of stochastic binding events occurring at times t_k ; s(t) represents an (inhomogeneous) Poisson point process with rate r(t): $\langle s(t) \rangle = r(t)$.

Zoom In Zoom Out Reset image size

For the sake of explicitness, we assume a radially symmetric concentration profile with exponential decay ² $c(\textbf{x}) \sim \exp (-k |\textbf{x}|)$. The agent performs runs interrupted by tumbling events that randomize its swimming direction; specifically, the swimming path $\textbf{r}(t)$ obeys $\dot {\textbf{r}} = v_{0}\,\textbf{t}$, where v₀ is a constant swimming speed and $\textbf{t}(t)$ the local tangent vector. The agent is subject to motility noise modeled as rotational diffusion with rotational diffusion coefficient D_rot, i.e., the direction angle ψ characterizing the direction of the tangent $\textbf{t} = \cos \psi \,\textbf{e}_{x} + \sin \psi \,\textbf{e}_{y}$ undergoes Brownian motion as $\dot {\psi } = \sqrt {2D_{\textrm{rot}}}\,\xi (t)$, where ξ (t) is Gaussian white noise with zero mean and unit variance $\langle \xi (t)\xi (t^{\prime }) \rangle = \delta (t-t^{\prime })$. Thus, the agent moves like an Active Brownian Particle (ABP) [23] during runs. Tumbling events are assumed to be instantaneous and cause a randomization of swimming direction with uniformly distributed new swimming direction angle $\psi \in [0,2\pi )$. After tumbling, agents exhibit a refractory period of T_max before they may tumble again.

Run-and-tumble chemotaxis relies on dynamically adjusting the length of runs, by triggering tumbling events earlier or later, depending on whether the agent is moving down-gradient or up-gradient, respectively. For this, the Machine-Learning (ML) agent needs to decide whether the concentration signal measured during its current run is decreasing or increasing with time. This is a non-trivial decision task, as concentration measurements are corrupted by both sensing and motility noise. As a result, run-and-tumble chemotaxis is highly stochastic on time-scales of a single run, but, if tuned correctly, can give robust net motion up-gradient on longer time-scales.

Deciding whether a train of binding events s(t) represents a time-dependent concentration $c(\textbf{r}(t))$ that decreases or increases with time t, can be efficiently addressed by a linear decision filter with memory kernel $w(\tau )$ and bias b

$$\begin{equation} \int_0^T \!\textrm{d}\tau\, w(\tau) s(t-\tau) + b \stackrel{?}{\gtrless} 0, \end{equation} \tag{ 1 }$$

where the integral represents a weighted sum of the recent past of the input s(t) up to a maximal measurement time T. We use supervised machine learning to learn optimal memory kernels $w_{j}=w(t_{j})$ using an appropriate discretization of time, and systematically explore how these memory kernels depend on noise levels.

For efficient learning, we use training datasets consisting of a large number of short trajectories that were evenly distributed over the search domain to ensure unbiased learning, see fig. 1(b). Specifically, we constructed training datasets consisting of short trajectories $\textbf{r}(t)$ of duration T_max = 1 s representing uninterrupted runs with rotational diffusion coefficient D_rot. Start points of training trajectories were uniformly distributed in a disk-like search domain with isotropically distributed initial tangent directions. For each training trajectory, we simulated time-discrete realizations of the train of stochastic binding events s(t), see fig. 1(c). Specifically, we represent s(t) as a vector with components $s_{j}=\int _{t_{j-1}}^{t_{j}}\!\textrm{d}t\,s(t)$, $j=1,\ldots ,m$, which are determined as Poissonian random variables with mean and variance $\lambda \,c(\textbf{r}(t_{j}))$, while $t_{j}=j\Delta t$ with $\Delta t=T_{\textrm{max}}/m$, where m is the number of temporal bins (here, m = 100).

Next, we implemented supervised learning by first assigning each realization $\{s_{j}\}$ an m-dimensional input vector u with components $u_{m+1-j} = s_{j}/\langle s_{j} \rangle$, where we normalized by $\langle s_{j} \rangle = \sum _{j=1}^{m} s_{j}/m$, see fig. 1(c). Second, we assign a class label $y=\pm 1$, where y = 1 if the tangent vector $\textbf{t}(T_{\textrm{max}})$ at the end of run is pointing up-gradient, and $y=-1$ otherwise. The rational behind this is that it would be favorable for a ML agent to continue its run beyond $t>T_{\textrm{max}}$ along its currently favorable direction $\textbf{t}(T_{\textrm{max}})$ if the class label is $y=+1$, but it would be more favorable to tumble if $y=-1$. Our goal is now to find out an m-dimensional weight vector w and scalar bias b, such that the predicted class label given by the sign of $\textbf{w}^{T}\textbf{u} + b$ matches the true label y as often as possible.

We used a training dataset consisting of $n=10^{4}$ short training trajectories, with corresponding input vectors $\textbf{u}_{i}$ and class labels y_i for $i=1, \ldots , n$.

To find the optimal weight vector w and bias b, we used LinearSVC to solve the optimization problem [24]

$$\begin{equation} \min_{\textbf{w}, b} \frac{1}{2} \textbf{w}^T \textbf{w} + C \sum_{i=1}^n \textrm{max}(0, \, 1- y_i(\textbf{w}^T \textbf{u}_i + b)), \end{equation} \tag{ 2 }$$

where C acts as an inverse regularization parameter (here C = 10). The learned weight vector w represents a time-discrete representation of a memory kernel $w(\tau )$ with $w(j\Delta t) = w_{j}$, see fig. 1(d) for a typical example. Memory kernels for additional values of λ and D_rot are shown in fig. S1 in the Supplementary Material Supplementarymaterial.pdf (SM).

On a technical note, we expect $w(\tau )$ to vary smoothly with τ. This could be achieved by using very large training datasets, or, more efficiently, by averaging learned $w(\tau )$ for 10³ learning sessions each using 10⁴ input vectors (economically constructed from 10 sets of $n=10^{4}$ short trajectories, for which 100 different realizations of noisy concentration measurements were simulated). Additional technical details on the learning of memory kernels is provided in the SM, sect. 1.

We tested the accuracy of the trained machine-learning (ML) agent on a different, analogously defined test dataset comprising $n=10^{4}$ short trajectories, see fig. S2 in the SM.

Although the memory kernel $w(\tau )$ was learned on short trajectory segments only, agents can use this kernel perpetually in time to find a hidden target, see fig. 1(e). Specifically, these agents record local concentrations $c(\textbf{r}(t))$ while performing runs as ABP with rotational diffusion coefficient D_rot, and execute a tumbling reorientation event, whenever $\int _{0}^{T_{\textrm{max}}}\!\textrm{d}\tau \, w(\tau ) s(t-\tau ) < 0$, where the memory kernel $w(\tau )$ is applied to the most recent past of noisy concentration measurements up to a maximum delay of T_max = 1 s.

We can characterize the performance of chemotactic agents in terms of the distribution p(t) of first-passage times (FPT) to find a target at the origin, see fig. 2(a). As expected, increasing motility noise increases typical FPTs, and reduces the probability to find the target within the simulation time t_sim = 200 s. We imposed absorbing boundary conditions both at the target boundary at radius R_target = 10 μm and at the domain boundary R_max = 500 μm. For comparison, we consider an ideal agent that has perfect knowledge whether it swims up- or down-gradient at any point in time. Remarkably, the chemotactic performance of the ML agent subject to sensing noise (with λ = 10⁵ μM⁻¹ s⁻¹), albeit smaller than that of the ideal agent, is still comparatively high.

Zoom In Zoom Out Reset image size

The chemotactic index $\textrm{CI}(R) = \langle \textbf{t}\cdot \nabla c / | \nabla c | \rangle _{|\textbf{x}|=R}$ characterizes the ensemble-averaged net motion up-gradient at different radial distances R, see fig. 2(b). A value of 1 would correspond to agents that head straight to the target, whereas unbiased random walks with isotropically distributed swimming direction correspond to $\textrm{CI} = 0$. In fig. 2(b), the CI decreases and even becomes negative close to the domain boundary, reflecting the net flux of agents to this absorbing boundary. The CI also decreases close to the target, presumably due to geometric nonlinearities of the radial concentration field, which are not sufficiently represented in the training. Results are comparable for ML agents and ideal agents, reaching maximal values of $\textrm{CI} \approx 0.6$. We find $\textrm{CI}<1$ because run-and-tumble agents cannot head directly to the target, but execute a biased random walk at best. The observed values are typical for E. coli bacteria at optimal conditions [25]. Consistent with the mostly positive CI, we find that the mean residence time t_res of agents as a function of the distance $|\textbf{x}|$ from the target is maximal in the vicinity of the target (and drops near the absorbing boundaries right at the target and the domain boundary, respectively), see fig. 2(c).

Binary decision making to tumble or not implies decision errors in the presence of noise: false-negative (FN) decisions to wrongly continue a down-gradient run, and false-positive (FP) decisions to wrongly tumble during an up-gradient run. For the training data, due to its symmetry, FN and FP errors have equal probability. But for the long trajectories of chemotactic agents, the positivity of their CI implies a class imbalance between up-gradient and down-gradient swimming. The result are strikingly different dependencies of FN and FP errors as a function of target distance, see fig. 2(d). The relative frequency of FP errors among all tumbling events is maximal far from the target, because chemoattractant concentrations are low there, implying a strong impact of sensing noise. The relative frequency of FN errors is high both close to the target and close to the domain boundary, where it even exceeds the chance level of 0.5. This can be understood from the fact that $\textrm{CI}<0$ close to the outer absorbing boundary, which implies that only few trajectories point up-gradient, i.e., the prevalence of true negative decisions (TN) is low. Similarly, the relative frequency of FN decisions increases again close to the target, because agents spend more time there, see the spatial distribution of residence times in fig. 2(c).

To assess the impact of sensing and motility noise, we computed memory kernels $w(\tau )$ for different values of the binding rate constant λ and the rotational diffusion coefficient D_rot, see fig. 3(a) (as well as fig. S3 in the SM for the corresponding test scores). All kernels display a characteristic shape comprising a lobe of positive values followed by a lobe with negative values, with total integral close to zero, ideally suited to compute a smoothed time derivative [22]. This characteristic shape of the kernel is markedly different from the optimal kernel $w(\tau ) \approx 1-2\tau /T_{\textrm{max}}$ with approximately linear shape that one would obtain in the limit of zero motility noise $D_{\textrm{rot}} = 0$.

Zoom In Zoom Out Reset image size

We define an effective measurement time T_eff by the simple condition $w(T_{\textrm{eff}}) = 0.05\, w_{\textrm{min}}$, where w_min is the minimum of the kernel. Figure 3(b) shows that T_eff follows an empirical power law with exponent α:

$$\begin{equation} T_{\textrm{eff}} \sim D_{\textrm{rot}}^{-\alpha},\quad \alpha \approx 0.2, \ldots,0.3. \end{equation} \tag{ 3 }$$

The exponent α depends weakly on λ and approaches $\alpha \approx 1/3$ for $\lambda \rightarrow \infty$. We can simulate this limit case of zero sensing noise, corresponding to the limit $\lambda \rightarrow \infty$, directly, by letting $s(t)=\lambda \,c(\textbf{r}(t))$. A different definition of T_eff, e.g., $w(T_{\textrm{eff}}) = 0.1\, w_{\textrm{min}}$ gave similar results. Likewise, a larger initial radial distance R₀ gave virtually the same power law, see fig. S4 in the SM.

The empirical power law reported in eq. (3) formalizes the trade-off choice between precision (requiring long T_eff to average out sensing noise), and accuracy (requiring short T_eff to reflect the current value of the rapidly changing swimming direction).

We quantified chemotactic performance as a function of both sensing and motility noise, see fig. 4. Remarkably, a small amount of motility noise can increase the chemotactic index, presumably because it facilitates course corrections for agents moving almost perpendicular to the gradient. This parallels the observation in fig. 2(a) that for ideal agents typical FPTs are shorter for D_rot = 0.1 s⁻¹ compared to the case $D_{\textrm{rot}} = 0$.

Zoom In Zoom Out Reset image size

We used supervised machine learning, applying support vector machines to a suitable training dataset comprising short trajectories, to determine optimal memory kernels for run-and-tumble chemotaxis in the presence of sensing and motility noise. Chemotactic agents employing these learned kernels for binary decision making in continuous time can find hidden targets releasing chemoattractant molecules efficiently. The chemotactic performance of these machine-learning agents is close to that of an ideal agent without sensing noise. The corresponding positive chemotactic index implies a class imbalance between swimming directions pointing up-gradient or down-gradient, which breaks a symmetry between false-positive and false-negative decisions. By computing memory kernels for sub-ensembles with homogeneous levels of sensing noise, we can dissect the dependence of an effective measurement time T_eff of gradient sensing by temporal comparison on sensing and motility noise, which revealed an empirical power law $T_{\textrm{eff}} \sim D_{\textrm{rot}}^{-\alpha }$ with α in a tight range $\sim 0.2,\ldots ,0.3$. In contrast, in a limit of high sensing noise, Strong et al. previously predicted the different power law $T_{\textrm{eff}} \sim D_{\textrm{rot}}^{-1}$ [21]. Our empirical power law for gradient sensing by temporal comparison is also different from the analytical result $T_{\textrm{eff}} \sim D_{\textrm{rot}}^{-1/2}$ previously obtained for tracking a time-varying concentration [18], or estimating a time-varying concentration gradient by spatial comparison [19]. Lastly, we showed that chemotactic performance becomes maximal at a small, but non-zero level of motility noise.

Previous experiments estimated the rotational diffusion coefficient of E. coli bacteria during runs as $D_{\textrm{rot}} \approx 0.1\,\textrm{s}^{-1}$ [26], which matches the predicted value for a Brownian particle of same size [21]. Although, this value is close to the optimal rotational diffusion coefficient that maximizes the time-averaged chemotactic index in our simulations, it seems unlikely that bacteria would optimize D_rot by adjusting their size.

Using a diffusion coefficient of chemoattractant D = 10³ μm²/s [27], we estimate a maximal binding rate constant λ = 10⁷ μM⁻¹ s for the idealized case that every molecule diffusing to the surface of spherical bacterial cell of radius a = 1 μm becomes absorbed. Previous experiments showed that chemotactic performance is maximal at a concentration gradient c ∼ 0.1 μM/μm [25], which corresponds to the gradient strength right at the target in our simulations. Our choice of radial concentration profile probes concentrations in a range of $c(R_{\textrm{target}}) = 10\,\mathrm{\mu }{\textrm {M}}$ to $c(R_{\textrm{max}}) = 0.1\,\textrm {nM}$, spanning from a gradient strength close to the source that gave maximal chemotactic performance in previous experiments [26], to a concentration just below the presumed sensitivity threshold for bacterial chemotaxis [28].

Experimentally measured chemotactic response kernels to the chemoattractant aspartate have a shape very similar to the optimal memory kernels predicted here [29]; a similar result was previously obtained using a minmax-strategy for an ensemble of different concentration gradients [22].

While we assumed a total randomization of swimming direction after instantaneous tumbling events, random re-orientations are only partial for E. coli and tumbling takes about 0.15 s. Biological cells performing chemotaxis exhibit sensory adaption, i.e., respond to relative changes of concentration [30]. In our case, the symmetry of the training problem implies symmetric memory kernels with total integral approximately zero; thus, decisions are insensitive to any constant offset of the concentration signal, see also [22]. In future work, one could envision more sophisticated navigation strategies, with agents that estimate local signal-to-noise ratios and adapt their decision kernels appropriately [31], with different costs for false-positive and false-negative decision errors.

Our work contributes to a recent interest to employ machine learning to theoretically understand optimal navigation strategies in the presence of noise using various learning algorithms [31–35]. Insights from these studies are already been put to use to improve navigation capabilities of artificial microswimmers [32].

We thank Julian Rode and all members of the Biological Algorithms group, as well as Marc Timme for stimulating discussions. ROR was supported by the Ministry of Science and Art of the Federal State of Saxony, Germany through grant 100400118 to Marc Timme and BMF, financed with tax funds on the basis of the budget adopted by the Saxon State Parliament (Forschungsprojektförderung Titelgruppe 70 des Sächsischen Staatsministerium für Wissenschaft und Kunst). BMF is supported by the DFG through FR3429/4-1 (Heisenberg grant), and under Germany's Excellence Strategy - EXC-2068 - 390729961. ROR and BMF acknowledge support through the Center for Advancing Electronics Dresden (cfaed).

Data availability statement: The data that support the findings of this study are openly available at the following URL/DOI: https://doi.org/10.5281/zendo.7540545