Motility of Colonial Choanoflagellates and the Statistics of Aggregate Random Walkers

We illuminate the nature of the three-dimensional random walks of microorganisms composed of individual organisms adhered together. Such aggregate random walkers are typified by choanoflagellates, eukaryotes that are the closest living relatives of animals. In the colony-forming species Salpingoeca rosetta we show that the beating of each flagellum is stochastic and uncorrelated with others, and the vectorial sum of the flagellar propulsion manifests as stochastic helical swimming. A quantitative theory for these results is presented and species variability discussed.

Active microparticles, self-propelled by stored energy or that available from the environment, typically exhibit directed motility combined with rotational diffusion, leading to random walks that at large times are statistically similar to their equilibrium counterparts.For artificial swimmers such as Janus particles [1], powered by inhomogeneous surface chemical reactions, the source of randomness is the same thermal fluctuations that translate Brownian particles, but here rotate them [2].In biology, several paradigms for stochastic locomotion exist.For single-celled organisms, stochastic beating leads to noisy swimming paths [3], and active processes such as flagellar bundling/unbundling by bacteria [4] and synchronization/desynchronization in algae [5] enhances this stochasticity.'Obligate' eukaryotic polyflagellates such as the ciliate Paramecium [6] and the alga Volvox carteri [7], exhibit large-scale flagellar coordination, and increased regularity of motion.
Here we study motility in an important example of a 'facultative' colonial organism, the choanoflagellate Salpingoeca rosetta (Fig. 1), which exhibits uni-and multicellular forms with variable cell number.Single cells of S. rosetta, like other microorganisms, are random walkers (see SI).We report three main experimental results: (i) individual flagella of the constituent cells beat stochastically, (ii) flagella on a given colony display negligible cross-correlation, and (iii) the swimming trajectories of colonies are stochastic helices.These results suggest a hitherto unrecognized class of microorganisms, here called aggregate random walkers (ARWs): those built by stitching together individual random walkers [8].We construct a minimal model to explain this motility.
Choanoflagellates are the closest unicellular relatives of animals [9].They filter feed by using their flagellum to drive fluid through an eponymous funnel-shaped collar [10].This beating also confers motility.Fig. 1 shows a rosette colony of S. rosetta which is held together by an extracellular matrix, filopodia, and intercellular bridges [11].Colonies form by cell division, not aggregation [12].The evolutionary advantage of the colonial form is not fully understood, but it is triggered by certain bacteria [11,13], and theory suggests that chain-like colonies have enhanced nutrient uptake [14].
S. rosetta (obtained from Dr. Barry Leadbeater, University of Birmingham, UK) were cultured in artificial seawater (36.5 g/L Marin Salts (Tropic Marin, Germany)).To provide a food source for prey bacteria, organic enrichment (4 g/L Proteose Peptone (Sigma-Aldrich, USA), 0.8 g/L Yeast Extract (Fluka Biochemika)) was added to the cultures at 15 µl/ml.Cultures were grown at 22 • C and split every 4 − 7 days.To study the flagella beat dynamics, colonies were stuck to poly-L-lysine (0.01%, Sigma) treated microscope slides and flagella beats imaged at 500 fps (Fastcam SA3, Photron, USA) in bright field.Image template matching was employed to track the motion of the only slightly moving colonies, and in the local frame of the organism, the bases of the flagella were tracked by techniques similar to that of the active contour model [15], yielding as a readout of beating the angle θ(t) defined in Fig. 1b.
Figs. 2a,b show θ(t) from two flagella on the same colony, and it is clear that they have distinct frequencies.In general, the beating frequencies f , found by Fourier transforming θ(t), show a surprisingly high variability  (Fig. 2e).The normalized autocorrelation C θ (∆t) = θ(t)θ(t + ∆t) t / θ(t) 2 t for a single flagellum is plotted in Fig. 2c.Similar to the function discussed [3,16] in the context of flagellar beating in Chlamydomonas, the data are consistent with C θ = exp(−|t|/τ ) cos(2πf t), the envelope of which is shown in the figure.The decay time τ also shows a very high degree of variability (Fig. 2f), but all are < 1 s suggesting high stochasticity.Within colonies, the cross correlation between flagella C θ,θ (∆t) = θ(t)θ (t+∆t) t / θ(t) 2 t θ (t) 2 t (Fig. 2d) is negligible (only a very slight signal can be made out, which we attribute to the overall wiggling of the colonysee Supplemental Video 1).All cross-correlation signals were found to be less than 0.05 (Fig. 2g).The lack of correlation between beating flagella in colonies makes S. rosetta an ideal model organism for ARWs.
In studying the swimming trajectories of S. rosetta, ensemble averages taken over many colonies will eliminate features related to colony-specific morphology (cell location and flagella orientation).To overcome this lack of 'ergodicity', we obtained long tracks of 36 individual colonies.In-house software logged and synchronized the position of the xy-stage (MS-2000, ASI, USA) to a camera (Imaging Source, Germany) filming in bright field at 15 fps.This enabled tracking of colonies moving in three dimensions at distances much longer than the field of view.To track the particles a Gaussian-mixture model [17] was applied to estimate the moving background and subsequently the tracks were manually controlled.Fig. 3a shows three examples, all ∼ 20 minutes in length.On close inspection (inset of Fig. 3a) we observe that the trajectories are noisy helices.The mean squared displacement ∆r 2 xy = [x xy (s + t) − x xy (s)] 2 s (Fig. 3b), shows an early time ballistic ∼ t 2 behavior (for t < 1 s) and late time diffusive ∼ t form (inset) similar to that of Janus particles [2].However, comparing these curves to those of conventional active Brownian particles (see SI), we observe a different intermediate time behavior.These bumps (Fig. 3b) appear precisely because some of the constituent cells may beat off center and induce internal (effective) torques producing stochastic helical trajectories.To highlight the underlying regularity of this helical swimming, we calculate the velocity autocorrelation 3c) which oscillates at the frequency of the induced rotation and decays on a timescale of several oscillation periods.
Active random walks have been the attention of much research [18,19], but only recently have rotational torques been incorporated.External torques appear on e.g.magnetotactic bacteria in the presence magnetic fields [20] and gyrotactic organisms such as certain algae in gravitational fields [21] and can be treated analytically [22].However, the present internal torques can be treated analytically only in 2D [23] and numerical [24] or approximative [25] methods are needed in 3D.Below we develop an approximate 3D theory with the goal of simple analytical functions that can be used to extract physical quantities and interpret the data.
The diffusion of a random walker can be described by the Langevin equation dx , where D is the translational diffusion constant and W (t) a standard vector Wiener process with dW i (t) dW j (t ) = δ ij δ(t − t ).The case v(t) = 0 is a passive particle and leads to the projected mean squared displacement ∆r 2  xy ≡ ∆x 2 + ∆y 2 = 4 Dt.Building an ARW from passive particles leads to no new behavior, but motile particles also have a stochastic velocity term.In the simplest case in two dimensions, the speed v is constant and v(t) = v (cos θ(t), sin θ(t)) evolves stochastically through θ(t).The choice dθ = √ 2D r dW r (t) leads to the conventional result ∆r 2 xy = (2v 2 /D 2 r ) D r t + e −Drt − 1 + 4 Dt, which behaves ballistically, ∆r 2 ∼ t 2 , at early times, but diffusively, ∆r 2 ∼ t, at longer times (see SI) with an enhanced diffusion constant D ∞ = D + v 2 /2D r , describing well the motion of Janus particles [2].This is a 2D result.Contrary to passive random walkers, active random walkers' effective diffusion constant can vary with dimension.
The corresponding 3D result is D ∞ = D + v 2 /3D r given a similar definition of 3D rotational diffusion.Typically, D v 2 /D r and passive diffusion can be ignored.
The Reynolds number for S. rosetta is Re ∼ 10 −4 .At such low Reynolds numbers, inertia is negligible and the fluid dynamics becomes governed by the linear Stokes equation.Accordingly, self-propelled choanoflagellates are both force-and torque-free.We assume that S. rosetta are sphere-like such that couplings between translations and rotations can be ignored.Heuristically, the velocity of a colony v(t) is approximately a linear sum of the velocities that the constituents would have had swimming independently, v(t) ≈ η v i (t), the factor η accounting for the change in drag with the radius a of the colony, as η ∼ a −1 .If some of the walkers comprising the colony, placed at positions {r i }, beat off center, an angular velocity ω(t) ≈ η r r i × v i will also be induced, where η r ∼ a −3 .Since {v i } and {r i } are given in the local coordinate system of the particle, they must be rotated along with the particle.For a two dimensional ARW, this motion is described by , and D r is an effective rotational diffusion constant which can be calculated if the individual stochastic processes are prescribed.With v(t) constant, constant ω(t) yields circles in the absence of noise.In three dimensions, such motion leads to helices, making (2D-projected) three dimensional ARWs behave very differently from 2D ones and necessitating a full 3D theory.
With the goal of a minimal model, we take the swimming speed to be constant and let the direction of the velocity evolve with the rotation of the organism according to dv(t) = dΩ(t) × v(t), where dΩ(t) = ω(t) dt + √ 2D r Λ(t) dW r (t).Here, Λ(t) is a random unit vector orthogonal to v(t), uncorrelated in time.To limit further the number of model parameters we assume the magnitude of ω(t) is constant, while its direction obeys dω(t) = √ 2D r Λ(t) × ω(t) dW r (t).This last update makes the system analytically quite intractable and thus we shall seek an approximate solution.As motivation, consider the case D r = 0 with specified initial conditions v(0) = v 0 , ω(0) = ω 0 .This system can be solved exactly to yield where and R is some orthogonal matrix.Eq. ( 1) describes a helix of radius v ω /ω 0 and mean speed v p (averaged over 2π/ω 0 ).The form of (1) inspires an approximative solution in the presence of noise in which the deterministic helix parameters define a continuous-time random walk with helix-like steps, the matrix R becoming a stochastic matrix process.As an effective description we assume , where the matrix factors are rotations around the x, y, z axes.α, β, γ are taken independent and identically distributed with dα = √ 2D r dW α (t).This approximation makes the system much more manageable.While the approach breaks x−y symmetry, simulations show it to be an overall good approximation for the statistics of interest In the stationary limit we find (see SI): The xy-projected mean squared displacement becomes where the constant a 0 enforces ∆r 2 xy (t = 0) = 0.As t → ∞ we obtain ∆r 2 xy = 4D ∞ t, where These results have been verified by simulations using the Euler-Maruyama method.It has previously been shown that reciprocal swimming enhances diffusion [19], and the last terms of ( 13), which are major contributions to the diffusion constant, embody this phenomenon.Equations ( 10) and ( 11) describe the approximate functions corresponding to the data of Figs.3c and 3b, respectively.The diffusion constant D ∞ can be extracted from the linear late-time behavior of ∆r 2 xy (dashed gray in inset of Fig. 3b), and can be used in (13) to fix one of the model parameters in terms of the others.The remaining three are fitted simultaneously to the curves of Fig. 3b and 3c.The experimental data are well-described by the model as shown by the dashed lines in the figures.The relative magnitudes of the extracted velocities, v p and v ω , reveal how much energy the organisms spent on effective (v p ) and circular (v ω ) swimming, for example the blue curve in Fig. 3 has v p = 10.6 µm/s and v ω = 14.2 µm/s.While not producing the precise morphologies of the colonies, the fitted velocities combined with the extracted frequency ω 0 , do constrain the possible configurations.Using the fitted velocities and a colony radius a ∼ 5 µm, we find an effective translational force of ∼ 1 pN, and using ω 0 , an effective torque ∼ 4 pN• µm: the small residual forces that propel and rotate a colony are on the order of that of a single cell.
Just as flagella beating in S. rosetta varies between cells, morphology varies between colonies as a result of the cell division process [12].This stochasticity enables two colonies of similar size to swim very differently.To quantify this, we used in-house software to track ∼ 750 colonies of varying size swimming in quasi-2D between two cover slips, and when a colony was in focus the area of an ellipse fitted to its outline served as an estimator of size (see Supplemental Video 2).This method, while introducing uncertainty in area, enables high throughput.To obtain model parameters, long tracks are needed.The parameters for 36 such tracks are given in the SI, and the speed v 2 p + v 2 ω of those tracks is shown in Fig. 4 as green circles.There is a slow increase in speed with colony size.This trend can be explained by simple ad-hoc models such as random orientation of cells in a sphere-like structure: drag scales linearly with radius a but maximum propulsive force (the case where all propulsive forces point in one direction) scales like a 2 .However, there is an intriguing lack of very slow swimmers which would be predicted by such a model.Indeed giving cells an orientation more parallel with its location would only yield slower swimming speeds.More importantly, Fig. 4 shows just how different colonies of similar size are: the stochastic processes underlying colony formation have high variances.From fits of the long tracks this stochasticity seems to apply to all model parameters (SI).This is contrary to e.g.bacterial clumps where rotation rate clearly decreases with size [8].Contrary to the phototactic response of Chlamydomonas and Volvox in which the time-scale of rotation is matched to inner chemistry [26], or the chemotactic response of sperm cells in which curvature and torsion of swimming paths are directly manipulated by the single beating flagellum [25], due to this stochastic morphology of S. rosetta, knowledge of the overall colony morphology and motion (e.g.ω 0 ) is arguably not available at the single-cell level, rendering 'deterministic' chemotactic strategies difficult.Thus one of the most important issues is the possibility of chemotaxis in aggregate random walks through suitable modulation of the independent constituents [27].
A fundamental operation in the theory of stochastic processes is their summation to yield a single effective process.The corresponding operation for random walkers, 'stitching' them together, yields ARWs.As we have shown, there is a crucial complexity for random walkers: the underlying flagellar beating can also yield rotations, so the 'summation rules' differ.Our results suggest that for simple random walkers the ARWs can be described approximately through four numbers: v ω , v p , ω 0 , and D r .The question of the correct 'summation rules' for general random walkers (e.g.anisotropic, hydrodynamically translation-rotation coupled, or non-identical con-stituents) remains open.Likewise, the transition, via e.g.self-assembly or flagella growth, from high to low stochasticity in ARWs with non-independent constituents is intriguing.The present exemplar, S. rosetta, is a very good approximation to what one might call an ideal biological ARW: independent constituents and a roughly spherical shape.Its mode of swimming raises many interesting questions about the evolution of multicellularity and on the nature and origin of noise, both internal and environmental.

Single cells
Slow-swimmer S. rosetta unicells, similar in morphology to the individual cells that comprise a colony, clearly exhibit random walk behaviour.Fig. 5 shows the average square displacement of 32 S. rosetta slow-swimmers each filmed for ∼1.5 minutes.The behaviour is well-described by the equation for conventional active random walkers, The active rotational diffusion constant for both single cells and colonies (main text) are on the order of 0.1 s −1 .With a beat frequency of ∼ 40 Hz this corresponds to a ∼ 3 • rotational deviation per beat.The thermal rotational diffusion constant D thermal r = k B T /8πµa 3 ranges from 0.012 to 0.0014 s −1 for radii 2.5 to 5.0 µm, at least an order magnitude below the active one.

Comparison of fit parameters
Fig. 6 shows scatter plots of fit parameters of the model to 36 different S. rosetta colonies, indicating the high variances of all parameters.We note, however, that the determination of some parameters are difficult in certain regions.For instance, ω 0 and v ω are hard to determine when either one becomes small, and accordingly we have forced them to zero in these cases and plotted them in blue.Naturally, these cases will have a higher v p as is clear in the two plots in the left part of Fig. 6.
Applying the same area estimator as in Fig. 4 of the main text, the parameters can also be plotted as a function of size.Just as with swimming speed, Fig. 7 shows that the model parameters have very high variances and no clear dependency on size.For a subset of the short tracks we were able to fit the model well enough to estimate ω 0 and these are shown as blue crosses.However, the short track colonies for which good estimated could be obtained are biased towards high ω 0 (and v ω ).Nonetheless, there is no clear tendency for larger colonies to rotate slower as is the case for e.g.bacterial clumps [1].For an interesting example of a big fast-spinning colony see end of Supplemental Video 2 in which a colony has formed as a dumb-bell shape.

FIG. 2 .
FIG. 2. (color online) Flagellar beating dynamics.(a,b) Time series of the base angle θ(t) on two flagella within a single colony.(c) Autocorrelation function of θ for one flagellum, with fit of the envelope to an exponential decay (dashed red).(d) Cross correlation of θ between two flagella on the same colony.(e) Peak frequencies of n = 23 tracked flagella.(f) Decay time of autocorrelation in single flagella.(g) Magnitude of cross correlations between flagella in same colonies.
FIG. 4.(colour online) Speed versus size for ∼ 750 colonies.Colony size is estimated by median xy-projected area.Colonies are blue and green circles, green speed being v 2 p + v 2 ω from model fits to the long tracks.Green crosses are single-celled fast swimmers[11].Colored background indicates running mean and standard deviation.