A stochastic model for bacteria-driven micro-swimmers
Christian Esparza Lopez, Albane Thery, Eric Lauga

TL;DR
This paper develops a stochastic fluid dynamic model for bacteria-driven micro-swimmers, analyzing their motion, diffusion, and chemotaxis, with implications for designing biomedical micro-robots.
Contribution
It introduces a novel analytical and computational model for bacteria-driven micro-swimmers, including effects of run-and-tumble motion and chemotaxis, validated by simulations.
Findings
Mean squared displacement scales with swimming speed and number of bacteria at short times.
At long times, MSD depends on particle size, not bacteria count.
Micro-swimmers can perform chemotaxis with drift velocity proportional to chemical gradient.
Abstract
Experiments have recently shown the feasibility of utilising bacteria as micro-scale robotic devices, with special attention paid to the development of bacteria-driven micro-swimmers taking advantage of built-in actuation and sensing mechanisms of cells. Here we propose a stochastic fluid dynamic model to describe analytically and computationally the dynamics of microscopic particles driven by the motion of surface-attached bacteria undergoing run-and-tumble motion. We compute analytical expressions for the rotational diffusion coefficient, the swimming speed and the effective diffusion coefficient. At short times, the mean squared displacement (MSD) is proportional to the square of the swimming speed, which is independent of the particle size (for fixed density of attached bacteria) and scales linearly with the number of attached bacteria; in contrast, at long times the MSD scales…
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Taxonomy
TopicsMicro and Nano Robotics · Molecular Communication and Nanonetworks · Microfluidic and Bio-sensing Technologies
\makeFNbottom
[TABLE]
††footnotetext: Department of Applied Mathematics and Theoretical Physics, University of Cambridge, Cambridge CB3 0WA, United Kingdom.††footnotetext: ‡ [email protected]
1 Introduction
Miniaturisation of actuators and efficient power sources are two of the biggest technical challenges in the design and fabrication of microscopic robots 1, 2. As is often the case, Nature can offer insight into overcoming these challenges. Flagellated bacteria, such as the well-studied E. coli, are known to be efficient swimmers with intricate sensing capabilities 3 and they have inspired scientists to integrate living cells and synthetic components into bio-hybrid devices 4, 5, 6, 7, 8, 2, 9, 10, 11, 12. Bacteria-driven micro-swimmers have received special attention due to their potential applications in medicine such as targeted drug delivery 13, 14, 15, 16, 8, 17, 18, 19. Commonly, a micro-swimmer consists of a synthetic substrate such as a micro-plate 7, 20 or a micro-bead 8, 16, 21, 17 with just a few or many bacteria attached to their surfaces (see experimental example in Fig. 1a). The attachment of the cells to the surface is often random but can be controlled by patterning the substrate using chemical or physical techniques 5, 22. The dynamics of a passive bead driven by surface-attached bacteria is a result of the precise behaviour of each cell. When free swimming, a peritrichous bacterium such as E. coli moves in a series of nearly straight paths (‘runs’), interrupted by quick changes in orientation (‘tumbles’). This so-called run-and-tumble motion consists therefore of two stages: (i) the running phase, in which all the helical flagellar filaments rotate counterclockwise (when measured from behind the cell) and form a bundle which aligns with the cell body and whose rotation in the fluid propels the cell forward; (ii) the tumble phase, in which at least one flagellar filament rotates in the clockwise direction and leaves the bundle resulting in a reorientation of the whole cell 3, 23. During the running phase, the bacterium experiences a constant propulsive force and torque along the same direction. Both, force and torque, balance with the corresponding hydrodynamic translational and rotational drags, so that the cell is overall force- and torque-free. In contrast, in the tumbling phase the bacterium experiences no propulsive force but a torque which reorients the cell body 15.
The running and tumbling events for E. coli are known to be Poisson distributed 25 with each run lasting for approximately one second while tumbles last just a tenth of a second 23. These times can be altered by chemical concentration gradients in the environment surrounding the cell. Bacteria such as E. coli can measure temporal differences in concentration of certain chemicals during their locomotion, for example aspartate, and reduce their tumbling rate if moving up the gradient 26, 27. As a result, the trajectory of a swimming bacterium is an isotropic random walk in homogeneous environments and a biased random walk when in the presence of a chemical gradient (if the cell is sensitive to it). The bias in the direction of motion is known as chemotaxis and when displayed by bacteria attached to passive particles it can potentially provide chemotactic abilities to the particles themselves 13, 28. Indeed, there are many other taxis techniques which can be used to control the trajectories of synthetic micro-swimmers, such as phototaxis 29 and magnetotaxis 30. However, for medical applications chemotaxis appears to be the most natural, as not only it does not require the use of external force fields but can also be linked to the ubiquitous presence of chemical gradients in the human body.
Helical trajectories of particles with surface-attached bacteria have been observed in experiments with small beads on which only one or two bacteria can be attached 15, 31. This suggests that the cells, despite being attached to a surface, still follow a run-and-tumble dynamics, and each bacterium applies a constant force and constant torque to the bead. Following this observation, a number of mathematical models have been proposed to describe the motion of bacteria-driven micro-swimmers 7, 14, 15. The swimmers are modelled as spherical particles actuated by random forces and torques taken as the sums of the individual propulsive forces and torques applied by each bacterium. Given the small velocity scales, the motion occurs at low Reynolds number and is governed by Stokes laws which relate linearly the applied forces and torques to the linear and angular velocities of the bead. Although these models have been validated numerically against experimental results 13, 14, 15, 16, 8, very few analytical expressions have been derived. For example, it is known that for large number of surface-attached bacteria , the swimming speed increases as 16, 14 but the dependence on the size of the micro-swimmer is unclear.
A detailed mathematical analysis of the motility properties of the micro-swimmers is essential for their optimal design and fabrication. In this study, we thus develop a stochastic model for bio-hybrid spherical micro-swimmers with E. coli bacteria attached to their surface. We derive analytical expressions for the rotational diffusion coefficient and the mean squared displacement (MSD), based on the following assumptions: (i) low-Reynolds number flow, (ii) a large number of uniformly attached bacteria, (iii) negligible thermal noise compared to the random activity of each cell and (iv) run-and-tumble dynamics for each bacterium. We first consider the situation where the chemical environment is homogenous and next investigate the effects of external chemical gradients, assuming a linear response of each bacterium. In particular, we show rigorously that the individual capabilities of E. coli cells to move up concentration gradients are inherited by the micro-swimmers on which they are attached. Our analytical results are validated against numerical simulations (using a standard Brownian dynamics framework) and past experimental results 13, 14, 15, 16, 8.
2 Mathematical model and numerical simulations
2.1 Micro-swimmers in homogeneous environments
Following previous studies 15, 14, we model the micro-swimmers as passive spherical particles of radius with a number of uniformly distributed bacteria attached to their surface. Reported densities of attachment in experiments range from one bacterium per up to one per 16. The bacteria are assumed to be fixed in position and orientation with respect to the surface of the particle. This condition is satisfied in practice by using strong chemical binding such as streptavidin-biotin interactions 15. Previous numerical investigations suggest that the direction of the flagellar bundle is unaffected by the fluid flow in the vicinity of the swimmer 14, therefore the flagellar bundle is assumed to align with the orientation of the cell body with no change as the swimmer moves.
Each cell is assumed to perform its own, independent run-and-tumble dynamics. The reaction torques in the running and tumbling states have magnitudes and respectively 15. These torques can be neglected if we focus on sufficiently large particles such that , where is the average propulsive force exerted by each bacterium 16, 15 and the angle defined in Fig. 2 denotes the orientation of the bacterium flagella bundle with respect to the radial direction. In other words, for large particles, the reorientation of the swimmer is dominated by the moment-arm torque induced by the propulsive forces of the cells. Finally, thermal noise can also be neglected as it induces typical forces three orders of magnitude smaller than the propulsive forces from the bacteria. We model therefore each bacterium as a two-state machine which exerts a force of magnitude when running and no force when tumbling. The transition between the running and the tumbling states is modelled as a continuous time Markov chain with transition rates and from the tumbling to the running state and vice-versa, respectively (see Fig. 3). The values for E. coli in homogeneous environments (no chemical gradients) are s*-1* and s*-1*. Furthermore, the swimming motion occurs at low Reynolds number 25. The linear and angular velocities of the micro-swimmer, denoted by and , are obtained by force and moment balance using Stokes law
[TABLE]
where is the dynamic viscosity of the fluid, and are the total force and torque applied by the surface-attached cells, while and are unit vectors which determine the position and orientation of the -th bacterium (see Figs. 1b and 2). Specifically, each bacterium is assumed to push on the fluid along the direction of its bundle of flagellar filaments and thus to exert a force on the particle along . In Eqs. (1) and (2), the stochastic variable determines the state of the -th bacterium with when running and when tumbling. The position and orientation of the micro-swimmer evolve in time according to
[TABLE]
where the dot denotes a derivative with respect to time and is a body-fixed unit vector. We integrate Eq. (3) numerically using an Euler method and Eq. (4) using the mid point method presented in Refs. 32, 33 with a time step s. We consider different configurations, sampling the geometrical angles , , and , which describe the location and orientation of the cells on the bead surface (see Fig. 2), from uniform distributions in the intervals , , and respectively.
The maximum and minimum deviation angles from the radial direction are and following Ref. 14.
2.2 Run-and-tumble dynamics
If we denote by and the probabilities of finding a particular bacterium in the running and the tumbling states respectively, the master equation for the run-and-tumble dynamics is given by
[TABLE]
where and are the transitions rates defined earlier. The system in Eq. (5) can be easily integrated by considering the normalisation condition , with the result
[TABLE]
where is a constant determined by the initial state of the bacterium. In particular, if the bacterium is running at time we have , hence
[TABLE]
where we denote by , the probability for a bacterium to be running or tumbling at time given that it was running at time . Furthermore, in the steady state, the results of Eq. (6) reduce to
[TABLE]
In our numerical simulations, we assume that we have waited long enough so that the system has reached steady state. We therefore determine the state of the variables by sampling a pseudo random number from a uniform distribution and comparing it with for the first step and , for the following steps. Here is obtained from Eq. (6) by taking . In the steady state, the mean of the variables is
[TABLE]
while the time autocorrelation is given by
[TABLE]
If the bacteria behave independently from each other, the probability of the -th bacterium being running at time , given that the -th bacterium is also running at time , is for and otherwise. Hence
[TABLE]
We will use Eqs. (9) and (11) to evaluate the MSD bellow.
2.3 Chemotaxis of a single bacterium
Bacteria such as E. coli are able to navigate through chemical gradients by using an inhibition of the mechanism that inverts the polarity of their bacterial motor when the cell swims in a favourable direction 3. Polarity inversion is responsible for tumbling events, and therefore runs are extended when swimming up (resp. down) the chemoattractant (resp. repellent) gradient 34. To investigate whether the individual chemotaxis of single cells contributes collectively to chemotaxis of the micro-swimmer, we use the simplest mathematical model that describes chemotaxis of a single bacterium. Specifically, we model each organism as equipped with an internal sensor that responds to temporal variations in the chemical concentration of its environment by modifying the tumbling rate, , i.e. the expected rate at which running cells stop running and transition to a tumble 27. For weak chemical concentrations sensed by the bacterium, the variation in is captured by a linear response as 35
[TABLE]
where is the mean tumbling rate when is constant. The function is a memory kernel which describes the response and has been measured for E. coli by Segall et al. 26. For a constant concentration gradient, the variation in the concentration experienced by the bacterium will depend on its position as follows
[TABLE]
where is a constant background concentration and and are, respectively, the position and velocity of the bacterium at time measured from the beginning of the last running event. The memory kernel satisfies the condition of adaptability 27
[TABLE]
which ensures when is kept constant. Inserting Eqs. (13) into Eq. (12) we obtain the tumbling rate as
[TABLE]
For computational convenience, we consider in our numerical simulations the simplest possible adaptive kernel. Specifically, we take only two impulses localised at and with , and therefore write , where is the magnitude of the response. This models a cell which compares concentrations at two different, specific times in the past, and changes its tumbling rate accordingly. Substitution of this kernel into Eq. (15) yields
[TABLE]
Note that, since any function can be expressed as a superposition of impulses, it is possible to treat more general chemotactic responses (such as those proposed by Clark et al. 36 or Celani et al. 37) in a similar manner by approximating by a sum of delta functions, , where ’s are delay times and where the ’s determine the intensity of the response subject to .
2.4 Chemotaxis of bacteria-driven micro-swimmers
Chemotaxis is included numerically in our model by modifying the tumbling rate of each bacterium according to Eq. (16). We make the simplifying assumption that the particle is permeable to the chemical so that the chemical gradient is not perturbed by its presence. In general, we should evaluate in Eq. (16) as the position of each cell. However, if we consider shallow gradients such that , then it is appropriate to neglect variations in concentration along the surface of the particle and substitute the position of the centre of the micro-swimmer for in Eq. (16).
We run the simulation from the previous sections for an ensemble of beads, now evaluating the tumbling rate at each step using Eq. (16) with taken to be the position of the centre of the bead. We further take the position of the two impulses in the kernel to coincide with the maximum and minimum of the response kernel measured by Segall et al. 26, i.e. and . The intensity of the response is taken as 37. As we require the second term in brackets in Eq. (16) to remain small and , where is the swimming speed of the micro-swimmers 16, we choose , which corresponds to concentration gradients , so that and the tumbling rate is reduced by approximately one tenth. The chemical gradient is set up along the direction.
3 Theoretical and numerical results
We start this section with the description of the three-dimensional trajectories of the micro-swimmers as obtained numerically and we observe that for long time-scales these become three-dimensional random walks. We then introduce a coarse-grained model to describe the diffusive behaviour of the micro-swimmers. We define the rotational diffusion coefficient, , the effective speed, , the effective diffusion coefficient, and derive analytical expressions for each one of them. Next, we show that in the presence of a non uniform chemical concentration field, to which the bacteria respond chemotactically, the micro-swimmers perform a biased random walk and we quantify their response in terms of a drift speed , for which we derive an analytical expression. Throughout this section, we validate the analytical results against our numerical simulations.
3.1 Three-dimensional trajectories of bacteria-driven micro-swimmers
We present first the results in isotropic environments, i.e. constant concentration field and thus constant tumbling rate. Typical trajectories for the centre of the particle, resulting from numerical integration of equations (3) and (4), are shown in Fig. 4. At short time scales, the applied force and torque are nearly constant and the resulting trajectories are noisy helices. For a given distribution of attached bacteria, the average force and torque acting on the particles in the steady state are given by
[TABLE]
where denotes the average over the probability distribution Eq. (6). Defining the vectors and , we can express the average radius and pitch of the trajectories as follows
[TABLE]
In contrast, at longer times the trajectories become three-dimensional random walks and their diffusive behaviour determines the motility properties of the micro-swimmers. It is well known in the theory of Brownian motion that for a random walk governed by rotational diffusion, the MSD is given by 38
[TABLE]
where the angle brackets denote ensemble average, is the effective speed at which the bead moves and is the orientation correlation time, which is the time scale of decay for the orientation correlation function, . In other words, is the time it takes the micro-swimmer to forget its initial orientation. For short time scales, that is , Eq. (21) reduces to , which represents ballistic motion. On the other hand, for time scales such that the MSD is linear in , a dependence typical in diffusion processes. The constant of proportionality is the effective diffusion coefficient which is defined in three dimensions by
[TABLE]
Since is a macroscopic property of the micro-swimmers, it is likely that a simplified description of the trajectories ignoring the microscopic details would still lead to the same result. Furthermore, ignoring the fine structure of the driving mechanism will render the chemotaxis analysis more tractable (see Sec. 3.3). Inspired by this we propose a coarse grained model as described in the next section. We will show that the micro-swimmers MSD is described by Eq. (21) and we will derive analytical expressions for , and .
3.2 Coarse-grained modelling
We start by observing that the trajectories of the micro-swimmers such as those illustrated in Fig. 4 and Fig. 5a consist of almost undisturbed helical paths interrupted by sudden changes in direction. We thus construct a coarse-grained model replacing every helical trajectory by a straight path along an average axis denoted by for the th path (see Fig. 5b). Each straight path is travelled at a constant speed for an average time . The resulting trajectory is therefore a chain of equal-length links. The polar and azimuthal angles between consecutive paths and are denoted by and respectively. We further take as uniformly distributed, due to homogeneity of the space, but we will assume that there is some persistence in the direction of motion, as can be observed in Fig. 4.
In polymer physics, a mathematically-identical setup is used in the freely-rotating chain model of a polymer and it is thus a textbook result that the velocity correlation function satisfies 39, 40, 41
[TABLE]
where is the average cosine of the angle between consecutive paths. Since a change in direction on a trajectory is the result of several consecutive tumbling events of the bacteria, which are Poisson distributed, the probability of observing turning events in an interval of time must be Poisson distributed as well. Considering this and the fact that the turning rate is we obtain
[TABLE]
The MSD follows from a double integration in time of Eq. (3.2), leading to the result
[TABLE]
Comparing with Eq. (21) we observe that the correlation time and the average flight time are related by
[TABLE]
The average flight time can also be obtained from the definition of the Kuhn length as 40, 41
[TABLE]
where is the maximum length of the trajectory (see Appendix A). Substitution of Eq. (3.2) and Eq. (26) into Eq. (27) yields
[TABLE]
and therefore
[TABLE]
From rotational diffusion 42, 43, the correlation time is related to the rotational diffusion coefficient , in spatial dimensions, by . Therefore for our micro-swimmer . The value of depends on the distribution of the torque that acts on the micro-swimmer, for example, reorientation due to thermal noise yields from the Stokes-Einstein relation 43. In the next section we will calculate the value of explicitly, considering the variance of the torque given by Eq. (2).
3.2.1 ** Rotational diffusion coefficient**
We now calculate the rotational diffusion coefficient for the particle by considering the variance of its angular displacement. The orientation of the micro-swimmer at time can be obtained from its initial orientation by applying the rotation matrix , where is a matrix such that for any vector we have , with the angle vector
[TABLE]
Given an orthonormal basis , the components of on satisfy , which in the Gaussian limit (thermal noise 42) reduces to , where denotes the variance. For large we may assume that is spherically uniform, thus we define the rotational diffusion coefficient as
[TABLE]
Substituting Eq. (2) in Eq. (30) leads to
[TABLE]
where Eq. (32) defines the vectors . The variance of the angle vector for a fixed configuration of attached bacteria is then given by
[TABLE]
For large we may replace the sum by times the average over configurations
[TABLE]
where the subscript denotes that the average is taken over the distributions of and , which is equivalent to the average over the angles , , and . As the positions and orientations of the bacteria are identically and independently distributed, the correlation function on the right hand side of Eq. (34) does not depend on , hence
[TABLE]
Remembering that defines the direction of motion in our coarse-grained model, we have
[TABLE]
from Eq. (3.2). On the other hand, since is the angle between and and the bacterial distributions are independent and identical, then . Therefore
[TABLE]
Note that, in doing this approximation, we are ignoring the periodic component of the correlation function of the vectors . As it can be seen in Fig. 7 this decays on a time scale shorter than and, since we are interested on the diffusive behaviour of the micro-swimmer, we can ignore it***The torque fluctuates around its mean on a time-scale in the order of which as seen below is much smaller than .. A substitution of Eq. (37) into Eq. (32) yields
[TABLE]
We expect the rotational diffusion time scale to be larger than the running and tumbling times, i.e. to be in the limit , which will be verified later. In this limit, Eq. (3.2.1) simplifies to
[TABLE]
We note that vanishes for , and indeed if all the bacteria are oriented in the radial direction they produce no moment-arm torque from their propulsive force and there is no change in the orientation of the particle. To treat this case one must include the effect of the reaction torques, the calculations are however more complicated since the direction of varies in time with respect to the direction of . On the other hand, the effects of thermal diffusion (which had been neglected) can be easily included by adding to the MSD and by replacing by where and are the thermal translational and rotational diffusion coefficients respectively.
In the particular case of a micro-swimmer immersed in water we have and with the parameters , , and , we obtain
[TABLE]
Therefore Eq. (39) is valid as long as we have . For a uniform distribution of angle of attachment in the range , , hence the condition is accomplished provided that the particle is sufficiently large, . We now proceed to calculate the micro-swimmers effective speed.
3.2.2 ** Effective speed**
We define the effective speed as the ensemble average projection of the velocity along the angular velocity, that is
[TABLE]
where denotes the average over configurations. For a large number of bacteria , we may assume that the vector is uniformly distributed on the sphere and therefore . Substitution in Eq. (41) yields
[TABLE]
where we have assumed that the bacteria behave independently of each other, that is . Notice that is the fraction of time that the bacteria exert a force on the sphere, hence the effective speed is simply the average swimming speed of each bacterium multiplied by the average number of bacteria pushing along the axis of the helical path . Notably, for a fixed density of attached bacteria, is independent of the size of the micro-swimmer as , this is in agreement with the numerical results of Arabagi et al. 14. On the other hand, is proportional to the square root of the number of bacteria, which is in agreement with the experimental results of Behkam et al. 16.
3.2.3 ** Long-time diffusion of bacteria-driven micro-swimmers**
We now use the results above to obtain an analytical expression for the effective diffusion coefficient, , defined by Eq. (22). Substitution of Eqs. (39) and (42) into Eq. (22) with leads to
[TABLE]
Our result makes two predictions with important consequences for experiments: (i) for a setup with a uniform density of surface-attached bacteria, we predict that the long-time diffusion coefficient increases with the square of the size of the particle; (ii) the diffusion constant is, perhaps surprisingly, independent of the value of the propulsive force exerted by each bacterium. Recall however that diffusion is obtained in the “long-time" limit, which is defined as and therefore a variation in the value of changes this limit. Note that when for all , we obtain a maximum value for the rotational diffusion and therefore a minimum value of the linear diffusion. In contrast, when for all swimmers we obtain a singular expression for the diffusion coefficient, due to the fact that there is no reorientation as discussed in Sec. 3.2.1. This limit is treated in Appendix A.
We next validate our analytical results against numerical simulations. In Fig. 6 we show a log-log plot of the MSD as a function of time for a micro-swimmer of radius m, with density of attached bacteria , a propulsive force pN and transition rates s and s. Time is measured in units of . At very short times , the force and torque remain constant and therefore the MSD is simply the instantaneous mean squared velocity multiplied by the time squared, that is (dashed magenta line). The solid blue line is Eq. (21) and, for times-scales shorter than , is ballistic with a MSD approximately . On the other hand, for time-scales larger than , the MSD approaches the asymptotic limit represented by the dashed dotted green line, with given by Eq. (43).
On Fig. 7 we illustrate the exponential decay of the velocity correlation function (solid red line) for the same micro-swimmer. The dashed line is the theoretical prediction from Eq. (3.2) with and given by Eqs. (39) and (42) respectively. The shaded region is the (numerical) standard deviation. The coarse-grained model captures accurately the late time behaviour of the computational results. Finally, Fig. 8 shows the dependence of with the bead radius , for a micro-swimmer with the same set of parameters , , and . The red circles and error bars represent the numerical mean and the standard deviation of , which where obtained from a linear fit to the MSD in the regime . The dashed blue line represents our theoretical prediction Eq. (43), which agrees very well with the numerical results. The inset in Fig. 8 confirms that is independent of the propulsive force of the individual bacteria, .
3.3 Chemotaxis of bacteria-driven micro-swimmers
In the previous sections, we analysed the motion of the micro-swimmers in chemically-homogeneous environments. Here we consider the case where instead the surrounding environment has a (weak) gradient of a solute to which the cells respond chemotactically. The question then emerges, whether the chemotaxis of individual cells translates to chemotaxis of the particle to which they are attached? The answer is not obvious a priori since bacteria in diametrically opposite positions on the particle can modify its motion in a symmetric fashion and as a result the chemotactic response can be vanishing. Experiments show however that the micro-swimmers perform indeed chemotaxis 13, 15. Here we rationalise this result using the coarse grained model developed in the last section and we give an analytical expression for the chemotactic drift, which is validated against numerical simulations.
3.3.1 ** Chemotactic drift speed**
As illustrated in the inset of Fig. 9, when the tumbling rate varies according to Eq. (12), the trajectories are biased in the direction of increasing concentration (positive ) and the drift increases with the magnitude of the chemical gradient. This result can be rationalised observing the dependence of the average displacement on the magnitude of the chemical gradient. Using the coarse-grained model developed in Sec. 3.2 and following the classical de Gennes analysis for a single cell 35, 44, 45, 46, we will compute the chemotactic drift speed . The drift speed is obtained by averaging the mean distance travelled in a particular run, over the possible directions that the run can take and dividing the result by the average flight time. As the tumbling rate does not change when the bead moves in the plane perpendicular to , the average displacement is non vanishing only in the direction of the gradient (i.e. ). Hence the drift velocity is
[TABLE]
where is the probability of tumbling between and and the subscript denotes that averages are taken over the directions of the velocity. Let us denote by the probability of not stopping within seconds. By the same argument we used to compute in Sec. 3.2 we know that is exponentially decaying and decays with rate , i.e.
[TABLE]
and . The probability of stopping after seconds is therefore given by
[TABLE]
Substitution of Eqs. (45) and (46) into Eq. (44) and using integration by parts gives the result
[TABLE]
where is the velocity of the particle in the direction of the gradient. Remembering that from Eq. (29), and using Eq. (12) we have
[TABLE]
this means that for shallow gradients, the tumbling rate for the micro-swimmer decreases in the same fashion as the tumbling rate for an individual bacterium, as long as . We use Eq. (48) to expand up to linear order in the gradient as follows
[TABLE]
where we have ignored the constant which gives terms that do not contribute to the drift, when we consider adaptive kernels with . Inserting Eq. (49) in Eq. (47) and expressing as a superposition of impulses we obtain the drift velocity as
[TABLE]
Here we have assumed that the distribution of is still uniform at first order in , hence and the velocity correlation along is given by
[TABLE]
where we have used the fact that is spherically symmetric. Substitution of Eq. (51) into Eq. (50) yields the drift as
[TABLE]
It is known that vanishes for 26. Furthermore, decreases quadratically with increasing , therefore for sufficiently large radius of the particle we have for all relevant values of (for example, when and , ) and we may then expand the exponential in the last line of Eq. (3.3.1) to first order in the rotational diffusion coefficient. This leads to the general prediction for drift velocity
[TABLE]
where we have used Eqs. (29) which leads to . Our prediction for is seen to be independent of the radius of the particle and is linear in the concentration gradient.
For the particular choice of used in our computations, we obtain
[TABLE]
These results are validated against our simulations in Figs. 9 and 10. In Fig. 9 we first plot as a function of for fixed , , , s and the same set of parameters , and as in the isotropic case. Correspondingly, we plot in Fig. 10 the drift speed as a function of for fixed radius m and the same set of parameters. Once again, the red circles and error bars represent the numerical values of the mean and standard deviation of . These values were obtained from a linear fit to the average displacement in the regime . The dashed blue lines represent the prediction of Eq. (54) and we can observe that the agreement with the numerical results is excellent.
4 Discussion
In this paper, we developed a coarse-grained model that allows to accurately predict the diffusive behaviour of bacteria-driven micro-swimmers, both in homogeneous environments and in the presence of (weak) chemical gradients. The analytical expressions for the effective diffusion coefficient (Eq. 43) and the chemotactic drift speed (Eq. 54) were validated against numerical simulations using the stochastic model presented in Section 2.1 with Eq. (12) for the chemotactic change of the tumbling rate of the bacteria. We found that the effective swimmer velocity, , is independent of the swimmer size, but increases with the square root of the number of surface-attached bacteria , when is large. On the other hand, the rotational diffusion coefficient, , increases linearly with and is inversely proportional to the particle radius , in contrast with the thermal rotational diffusion coefficient which decreases as . For shallow chemical gradients , we found that the micro-swimmers respond chemotactically if we assume that each bacterium on its surface does, and that the resulting chemotactic drift speed is proportional to and independent of , the latter being a consequence of the bacteria being constrained by the bead and each bacterium experiencing the same change in concentration.
The expression we obtained for the effective velocity in this paper can be directly compared with experimental results. Behkam et al. reported for micro-swimmers of radius and densities of attachment between one bacterium per and with 16. In comparison, our model predicts , in good agreement. The same study also reported that the swimming speed increases with the square root of the number of bacteria , for , as predicted theoretically. The same scaling result was reported by Arabagi et al. 14 and Zhuang et al. 15. Furthermore, both of these studies reported swimming velocities independent of the particle radius in homogeneous environments as well as in the presence of a chemical gradient, again confirming our theoretical results. Unfortunately there are no past experimental measurements of the effective diffusion coefficient. In their numerical simulations however, Arabagi et al. 14 reported the value for a bead of radius and , in good agreement with Eq. (43) which predicts . Regarding the chemotactic drift of particles with surface-attached bacteria, Zhuang et al. 13, 15 found that is proportional to the magnitude of the chemical gradient and that it is independent of the particle size, both of which are predicted by Eq. (54).
The analytical results in this paper will be useful for the practical design of bacteria-driven micro-swimmers. The model can be extended to include other steering strategies such as phototaxis, pH-taxis and external electromagnetic fields 29, 20. The assumption of a chemically permeable particle could also be relaxed, and the perturbation to the concentration field by the presence of the particle (and the cells) should be solved for. Another possible extension is the study of deformable random walkers, inspired by experiments with bacteria-driven micro-swimmers in which E. coli are attached to red blood cells 47. One of the ultimate applications of bacteria-driven micro-swimmers is targeted drug delivery, which requires a suspension of these particles to move collectively. A mathematical description of the collective behaviour of multiple micro-swimmers is therefore needed and this requires an explicit consideration of their hydrodynamic interactions. The interactions depend critically on the distribution of attached bacteria on the surface of each particle, which can be controlled using micro-manipulation 48 and nanoprinting methods 5, 49. A theory including hydrodynamic interactions, taking clues from classical work on active fluids 50, 51 will therefore allow to move bacteria-driven micro-swimmers closer to applications.
Conflicts of interest
There are no conflicts to declare.
Acknowledgements
This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement 682754 to EL).
Appendix A Singular case
When all bacteria point in the radial direction, there is no applied torque and every vector remains constant in time. In this case, the velocity correlation function is easily computed as follows
[TABLE]
where we have used the fact that the orientations and positions of bacteria are uncorrelated. Integrating twice, we obtain the MSD
[TABLE]
This is essentially ballistic and for we have with given by Eq. (42).
This result might appear counter intuitive, as one would expect that a large number of uniformly distributed radial forces acting on a sphere should account for no net thrust. However, this is only true if the time autocorrelation of the forces decays to zero. As the autocorrelation of the epsilon variables decays to a finite value, from Eq. (11), so does the velocity correlation function, hence the ballistic motion. Tilting the forces away from the radial direction allows for reorientation of the sphere, reducing the probability of finding a given force pointing in the same direction at two different times accounting for a vanishing autocorrelation at large time-scales, hence the long-time diffusive behaviour for the case . Notice that is the maximum average length that a trajectory can have (with not necessarily null). Indeed, the arc length of a helix is the same as that of a straight line traversed at the same speed, therefore can be obtained by considering the square root of Eq. (41) in the limit when is parallel to and multiplying the result by , which is exactly .
Notes and references
- Sitti 2009
M. Sitti, Nature, 2009, 458, 1121.
- Martel 2012
S. Martel, Biomed. Microdevices, 2012, 14, 1033–1045.
- Berg 2004
H. C. Berg, E. coli in Motion, Springer New York, 2004.
- Hosseinidoust et al. 2016
Z. Hosseinidoust, B. Mostaghaci, O. Yasa, B.-W. Park, A. V. Singh and M. Sitti, Adv. Drug Deliver. Rev., 2016, 106, 27 – 44.
- Carlsen and Sitti 2014
R. W. Carlsen and M. Sitti, Small, 2014, 10, 3831–3851.
- Ceylan et al. 2017
H. Ceylan, J. Giltinan, K. Kozielski and M. Sitti, Lab Chip, 2017, 17, 1705–1724.
- Julius et al. 2009
A. Julius, M. Sakar, E. Steager, U. Cheang, M. Kim, V. Kumar and G. Pappas, 2009 IEEE International Conference on Robotics and Automation, 2009, pp. 1004–1009.
- Darnton et al. 2003
N. Darnton, L. Turner, K. Breuer and H. C. Berg, Biophys. J., 2003, 86, 1863–1870.
- Katuri et al. 2017
J. Katuri, X. Ma, M. M. Stanton and S. Sánchez, Accounts Chem. Res., 2017, 50, 2–11.
- Schwarz et al. 2017
L. Schwarz, M. Medina-Sánchez and O. G. Schmidt, Appl. Phys. Rev., 2017, 4, 031301.
- Wang and Pumera 2018
H. Wang and M. Pumera, Adv. Funct. Mater., 2018, 28, 1705421.
- Bastos-Arrieta et al. 2018
J. Bastos-Arrieta, A. Revilla-Guarinos, W. E. Uspal and J. Simmchen, Front. Robot. AI, 2018, 5, 97.
- Zhuang and Sitti 2016
J. Zhuang and M. Sitti, Sci. Rep-UK, 2016, 6, 32135.
- Arabagi et al. 2011
V. Arabagi, B. Behkam, E. Cheung and M. Sitti, J. Appl. Phys., 2011, 109, 114702.
- Zhuang et al. 2017
J. Zhuang, B.-W. Park and M. Sitti, Adv. Sci., 2017, 4, 1700109.
- Behkam and Sitti 2008
B. Behkam and M. Sitti, Appl. Phys. Lett., 2008, 93, 223901.
- Schauer et al. 2018
O. Schauer, B. Mostaghaci, R. Colin, D. Hürtgen, D. Kraus, M. Sitti and V. Sourjik, Sci. Rep., 2018, 8, 9801.
- Cho et al. 2012
S. Cho, S. J. Park, S. Y. Ko, J.-O. Park and S. Park, Biomed. Microdevices, 2012, 14, 1019–1025.
- Nguyen et al. 2016
V. D. Nguyen, J.-W. Han, Y. J. Choi, S. Cho, S. Zheng, S. Y. Ko, J.-O. Park and S. Park, Sensor. Actuat. B-Chem., 2016, 224, 217 – 224.
- Steager et al. 2011
E. B. Steager, M. S. Sakar, D. H. Kim, V. Kumar, G. J. Pappas and M. J. Kim, J. Micromech. Microeng., 2011, 21, 035001.
- Behkam and Sitti 2007
B. Behkam and M. Sitti, Appl. Phys. Lett., 2007, 90, 023902.
- Stanton et al. 2015
M. M. Stanton, J. Simmchen, X. Ma, A. Miguel-López and S. Sánchez, Adv. Mater. Interfaces, 2015, 3, 1500505.
- Lauga 2016
E. Lauga, Annu. Rev. Fluid Mech., 2016, 48, 105–130.
- Kim et al. 2012
D. Kim, A. Liu, E. Diller and M. Sitti, Biomed. Microdevices, 2012, 14, 1009–1017.
- Block et al. 1982
S. M. Block, J. E. Segall and H. C. Berg, Cell, 1982, 31, 215–226.
- Segall et al. 1986
J. E. Segall, S. M. Block and H. C. Berg, P. Natl. Acad. Sci. USA, 1986, 83, 8987–8991.
- Schnitzer 1993
M. J. Schnitzer, Phys. Rev. E, 1993, 48, 2553–2568.
- Park et al. 2013
S. J. Park, S.-H. Park, S. Cho, D.-M. Kim, Y. Lee, S. Y. Ko, Y. Hong, H. E. Choy, J.-J. Min, J.-O. Park and S. Park, Sci. Rep. UK, 2013, 3, 3394.
- Steager et al. 2007
E. Steager, C.-B. Kim, J. Patel, S. Bith, C. Naik, L. Reber and M. J. Kim, Appl. Phys. Lett., 2007, 90, 263901.
- Carlsen et al. 2014
R. W. Carlsen, M. R. Edwards, J. Zhuang, C. Pacoret and M. Sitti, Lab Chip, 2014, 14, 3850–3859.
- Edwards et al. 2013
M. R. Edwards, R. W. Carlsen and M. Sitti, Appl. Phys. Lett., 2013, 102, 143701.
- Simo et al. 1995
J. C. Simo, N. Tarnow and M. Doblare, Int. J. Numer. Meth. Eng., 1995, 38, 1431–1473.
- Zupan and Saje 2011
E. Zupan and M. Saje, Adv. Eng. Softw., 2011, 42, 723 – 733.
- Sourjik and Wingreen 2012
V. Sourjik and N. S. Wingreen, Curr. Opin. Cell Biol., 2012, 24, 262 – 268.
- de Gennes 2004
P.-G. de Gennes, Eur. Biophys. J., 2004, 33, 691–693.
- Clark and Grant 2005
D. A. Clark and L. C. Grant, P. Natl. Acad. Sci. USA, 2005, 102, 9150–9155.
- Celani and Vergassola 2010
A. Celani and M. Vergassola, P. Natl. Acad. Sci. USA, 2010, 107, 1391–1396.
- Chandrasekhar 1943
S. Chandrasekhar, Rev. Mod. Phys., 1943, 15, 1–89.
- Lovely and Dahlquist 1975
P. S. Lovely and F. Dahlquist, J. Theor. Biol., 1975, 50, 477 – 496.
- Rubinstein and Colby 2003
M. Rubinstein and R. H. Colby, Polymer Physics (Chemistry), Oxford University Press, 2003.
- Yamakawa 1997
H. Yamakawa, Helical Wormlike Chains in Polymer Solutions, Springer Berlin Heidelberg, 1997.
- Saragosti et al. 2012
J. Saragosti, P. Silberzan and A. Buguin, PLOS ONE, 2012, 7, 1–6.
- Doi and See 1996
M. Doi and H. See, Introduction to Polymer Physics, Oxford University Press, USA, 1996.
- Locsei 2007
J. T. Locsei, J. Math. Biol., 2007, 55, 41–60.
- Taktikos et al. 2013
J. Taktikos, H. Stark and V. Zaburdaev, PLoS ONE, 2013, 8, e81936.
- Pankratova et al. 2018
E. V. Pankratova, A. I. Kalyakulina, M. I. Krivonosov, S. V. Denisov, K. M. Taute and V. Y. Zaburdaev, PLoS ONE, 2018, 13, e0190434.
- Alapan et al. 2018
Y. Alapan, O. Yasa, O. Schauer, J. Giltinan, A. F. Tabak, V. Sourjik and M. Sitti, Sci. Robotics, 2018, 3, eaar4423.
- Barroso et al. 2015
Á. Barroso, S. Landwerth, M. Woerdemann, C. Alpmann, T. Buscher, M. Becker, A. Studer and C. Denz, Biomed. Microdevices, 2015, 17, 26.
- Rozhok et al. 2005
S. Rozhok, C.-F. Shen, P.-L. Littler, Z. Fan, C. Liu, C. Mirkin and R. Holz, Small, 2005, 1, 445–451.
- Marchetti et al. 2013
M. C. Marchetti, J.-F. Joanny, S. Ramaswamy, T. B. Liverpool, J. Prost, M. Rao and R. A. Simha, Rev. Mod. Phys., 2013, 85, 1143.
- Jülicher et al. 2018
F. Jülicher, S. W. Grill and G. Salbreux, Rep. Prog. Phys., 2018, 81, 076601.
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1Sitti 2009 M. Sitti, Nature , 2009, 458 , 1121.
- 2Martel 2012 S. Martel, Biomed. Microdevices , 2012, 14 , 1033–1045.
- 3Berg 2004 H. C. Berg, E. coli in Motion , Springer New York, 2004.
- 4Hosseinidoust et al. 2016 Z. Hosseinidoust, B. Mostaghaci, O. Yasa, B.-W. Park, A. V. Singh and M. Sitti, Adv. Drug Deliver. Rev. , 2016, 106 , 27 – 44.
- 5Carlsen and Sitti 2014 R. W. Carlsen and M. Sitti, Small , 2014, 10 , 3831–3851.
- 6Ceylan et al. 2017 H. Ceylan, J. Giltinan, K. Kozielski and M. Sitti, Lab Chip , 2017, 17 , 1705–1724.
- 7Julius et al. 2009 A. Julius, M. Sakar, E. Steager, U. Cheang, M. Kim, V. Kumar and G. Pappas, 2009 IEEE International Conference on Robotics and Automation, 2009, pp. 1004–1009.
- 8Darnton et al. 2003 N. Darnton, L. Turner, K. Breuer and H. C. Berg, Biophys. J. , 2003, 86 , 1863–1870.
