help button home button Biophys. J.
HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS

Originally published as Biophys J. BioFAST on December 22, 2006.
doi:10.1529/biophysj.106.084897
This Article
Right arrow Abstract Freely available
Right arrow Full Text (PDF)
Right arrow All Versions of this Article:
biophysj.106.084897v1
92/6/1900    most recent
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Right arrow reprints & permissions
Citing Articles
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Guirao, B.
Right arrow Articles by Joanny, J.-F.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by Guirao, B.
Right arrow Articles by Joanny, J.-F.
Biophysical Journal 92:1900-1917 (2007)
© 2007 The Biophysical Society

Spontaneous Creation of Macroscopic Flow and Metachronal Waves in an Array of Cilia

Boris Guirao and Jean-François Joanny

Laboratoire Physico-chimie Curie, Institut Curie, Paris, France

Correspondence: Address reprint requests to Boris Guirao, Laboratoire Physico-chimie Curie (UMR 168), Institut Curie, 26 rue d'Ulm, 75248, Paris cedex 05, France. Tel.: 33-1-42-34-64-71; E-mail: boris.guirao{at}m4x.org.


    ABSTRACT
 TOP
 ABSTRACT
 INTRODUCTION
 SPONTANEOUS ALIGNMENT OF AN...
 AXONEMAL BEATING
 LEFT-RIGHT BEATING SYMMETRY...
 DISCUSSION AND CONCLUDING...
 APPENDIX I: AVERAGE FORCE...
 APPENDIX II: FOURIER MODE...
 ACKNOWLEDGEMENTS
 REFERENCES
 
Cells carrying cilia on their surface show many striking features: alignment of cilia in an array, two-phase asymmetric beating for each cilium, and existence of metachronal coordination with a constant phase difference between two adjacent cilia. We give simple theoretical arguments based on hydrodynamic coupling and an internal mechanism of the cilium derived from the behavior of a collection of molecular motors to account qualitatively for these cooperative features. Hydrodynamic interactions can lead to the alignment of an array of cilia. We study the effect of a transverse external flow and obtain a two-phase asymmetrical beating, faster along the flow and slower against the flow, proceeding around an average curved position. We show that an aligned array of cilia is able to spontaneously break the left-right symmetry and to create a global average flow. Metachronal coordination arises as a consequence of the internal mechanism of the cilia and their hydrodynamic couplings, with a wavelength comparable to that found in experiments. It allows the cilia to start beating at a lower adenosine-triphosphate threshold and at a higher frequency than for a single cilium. It also leads to a rather stationary flow, which might be its major advantage.


    INTRODUCTION
 TOP
 ABSTRACT
 INTRODUCTION
 SPONTANEOUS ALIGNMENT OF AN...
 AXONEMAL BEATING
 LEFT-RIGHT BEATING SYMMETRY...
 DISCUSSION AND CONCLUDING...
 APPENDIX I: AVERAGE FORCE...
 APPENDIX II: FOURIER MODE...
 ACKNOWLEDGEMENTS
 REFERENCES
 
Many cells and bacteria have cilia or flagella on their surfaces. Examples are sperm cells that have one flagellum used for propulsion, the green alga Chlamydomonas that uses two flagella, and the much studied protozoan Paramecium, which is covered by ~4000 cilia that produce a very efficient motion with a velocity of order 1 mm/s in water, corresponding to 10x the Paramecium size per second. Humans have ciliated cells in several organs: in the brain (ependymal cells that create cerebrospinal fluid flow (1Go)), the retina (photoreceptor connective cilia), the respiratory tract (epithelial cells), the ear (hair bundles), the Fallopian tube, or the kidney (2Go).

Cilia have two major roles: i), detection (sensory cilia or flagella), for example, in the retina, the ear, and the kidney; and ii), propulsion or creation of fluid flow (motile cilia or flagella) as for Paramecium or in the respiratory tract where the fluid flow is used to move away the mucus.

The common structure of most cilia and flagella is an axoneme wrapped by the plasma membrane. The (9 + 2) axoneme is made of nine microtubule doublets arranged on a circle around a central pair of microtubules (3Go). The cilium or flagellum is attached to the cell membrane by a basal body made of nine microtubule triplets that has a structure very similar to that of a centriole. The basal body is attached to the cell membrane by anchoring fibers (4Go). Typically the radius of an axoneme is 0.1 µm. The main structural difference between cilia and flagella is their length (10 µm for cilia, up to 100x longer for flagella (5Go)). Dynein molecular motors are attached to the nine microtubule doublets. Upon consumption of adenosine-triphosphate (ATP), dynein motion generates forces that induce a sliding between adjacent microtubules. Because the whole structure is attached at its base, this sliding motion induces the bending of the cilium or flagellum and its beating (3Go,6Go).

We here focus on ciliated cells creating fluid flow. These are cells with cilia on their surface, beating in one preferred direction in a coordinated way. One central feature of cilia beating is the existence of two phases with a broken symmetry. Each beating can be decomposed into an effective stroke that propels the fluid and a recovery stroke where the cilium is coming back against the flow. In the example of Paramecium in water, the effective stroke lasts typically 9 ms whereas the recovery stroke lasts 26 ms. The typical beating frequency in water is 30 Hz (7Go). The beating of Paramecium cilia is three-dimensional but for some species like Opalina or Chlamydomonas, it remains essentially planar (8Go). In this work, we discuss the role of an external velocity field in this left-right symmetry breaking between the effective stroke and the recovery stroke for planar beating.

One of the most striking features of an assembly of beating cilia is that they all beat in the same direction: the surrounding fluid can only be propelled efficiently if all the beatings have the same orientation. In mature ciliated cells, the beating direction is defined by the anchoring of the basal foot on the basal body. Only newly formed or developing cilia are randomly oriented (9Go). When they start beating, they tend to spontaneously align to finally beat in the same direction. One of the questions addressed in this article is the nature of the parameters that control this orientation.

Another important feature of ciliated cells is the existence of waves propagating all along the surface. These are called metachronal waves and might be due to the coordination of adjacent cilia, for example, via hydrodynamic interactions. Experimentally metachronal waves are observed to propagate in all possible directions: in the direction of the effective stroke (symplectic metachronal waves), in the opposite direction (antiplectic), or even in a perpendicular (laeoplectic or dexioplectic) or oblique direction. The origin of these waves and the mechanisms controlling their formation are not well understood. We show in this article that metachronism can arise naturally from the hydrodynamic couplings between cilia. Using a two-state model for the dynein motion as an internal mechanism of the cilia, metachronism appears to be a local minimum in the oscillation threshold of the motors (10Go,11Go). It decreases the amount of ATP necessary for beating, but it also increases the beating frequency and induces a rather stationary flow that accounts for the regular swim of ciliated organisms like Paramecium. This last effect could be the major interest of metachronal coordination.

The local [Ca2+] concentration has a strong influence on the beating pattern of cilia or flagella. For example detergent-treated Paramecium are able to swim forward at low [Ca2+] concentration (<10–6M) and backward at high [Ca2+] concentration (>10–6M) because of ciliary reversal: the directions of effective and recovery strokes are switched (12Go,13Go). In any case, the wild-type Paramecium can have a very efficient backward motion monitored by calcium tanks in its body. We only discuss here qualitative aspects of the role of calcium.

In this article, we address the question of the spontaneous alignment of an array of beating cilia and the possibility of a spontaneous symmetry breaking in the beating that leads to the appearance of a macroscopic fluid flow. The internal mechanism of the cilia is described by the model of references (10Go,11Go) that is based on a two-state model to describe the cooperative effects between dynein motors and only considers the relative sliding of two microtubules in the axoneme. The coordination between the cilia is due to hydrodynamic interactions that are discussed in details in a coarse-grained description where the effect of the cilia on the flow is replaced by an effective force.

The outline of the article is as follows. In the next section, we give a simple model for the alignment of beating cilia. In the "Axonemal beating" section, we discuss the beating of one cilium following the model of Jülicher and Camalet (10Go,11Go). Finally, in "Left-right beating symmetry breaking and metachronal coordination", we discuss the spontaneous breaking of the left-right symmetry of the beating due to the flow created by the cilia themselves and the existence of metachronal waves. The interest of our approach is to study the role of hydrodynamic effects on some characteristic features of ciliated cells; proteins, Ca2+ waves, or more generally chemical signals may not be the only answers to all these questions.


    SPONTANEOUS ALIGNMENT OF AN ARRAY OF CILIA
 TOP
 ABSTRACT
 INTRODUCTION
 SPONTANEOUS ALIGNMENT OF AN...
 AXONEMAL BEATING
 LEFT-RIGHT BEATING SYMMETRY...
 DISCUSSION AND CONCLUDING...
 APPENDIX I: AVERAGE FORCE...
 APPENDIX II: FOURIER MODE...
 ACKNOWLEDGEMENTS
 REFERENCES
 
Experimental results
In an assembly of cilia covering the surface of a mature cell, cilia are beating in a preferred direction, and only newly formed or developing cilia are randomly oriented (9Go). We first discuss the experiments showing how this preferred orientation is chosen.

As mentioned before, the ciliary axoneme grows from a basal body analogous to a centriole. Two basal body appendages, the basal foot and the striated rootlet, located in the plane of the effective stroke, confer an asymmetrical organization to the basal body. The basal foot is laterally associated with two consecutive triplets and points in the direction of the effective stroke (14Go,15Go). The striated rootlet, associated with the proximal end of the basal body, sinks into the cytoplasm in the opposite direction (14Go). These two appendages define therefore an orientation of a cilium before beating motion.

During ciliogenesis, newly formed basal bodies migrate toward the cell membrane where they anchor with no apparent order. Anchoring induces axoneme assembling, and cilia grow in random orientations. While cilia are growing, they do not beat immediately. A reorientation by rotation of the basal bodies in a common direction occurs at the final stage of ciliogenesis, when mature cilia beat (16Go). The preferred direction of the assembly is then well defined. In primary ciliary dyskinesia (PCD) also known as immotile cilia syndrome, axonemes are incomplete, and the ciliary activity is abnormal or absent: the fluid is poorly or not propelled. At the cell level, the basal bodies are randomly oriented (17Go).

These experimental facts suggest that the beating and orientation of cilia are closely related. Our working hypothesis is that the alignment of an assembly of beating cilia is mostly due to hydrodynamic coupling between cilia. The global flow created by the other cilia tends to orient a given cilium and above a certain beating amplitude, all cilia orient on average in the same direction. We now give a simple modeling of this cooperative alignment.

Alignment transition
We assume in the following that the beating is planar. It is the case for Opalina, for example, but not exactly for Paramecium where the recovery stroke is not in the plane of the effective stroke. Near the top of the ciliary layer, observations show that the velocity is time independent and uniform (18Go). Consequently, we average the beating over one time period and replace each cilium of length L (and its effective and recovery stroke) by a single force (stokeslet) Formula, parallel to the surface, created in the fluid of viscosity {eta} at height h < L above the membrane, as sketched in Fig. 1.


Figure 1
View larger version (10K):
[in this window]
[in a new window]

 
FIGURE 1  (a) Beating pattern of a cilium. The fluid is efficiently propelled during the effective stroke. During the recovery stroke, the cilium comes back close to the surface, minimizing the viscous effects. (b) Effective force in the fluid Formula applied at a height h above the cell membrane, to mimic the cilium beating.

 
We choose the x axis in the direction of the effective stroke and the z axis perpendicular to the cell surface that we approximate by a plane. The stokeslet along the x axis (Formula) is located at point Formula. The velocity created at point Formula by this stokeslet with a no-slip boundary condition on the plane z = 0 is given by:

Formula
where the response tensor G is given in Blake (19Go) and reads:

Formula 1(1)
where Formula 1, Formula 1 and {alpha} = x, y.

To simplify the hydrodynamic problem, we assume in all the following that the cilia are far away from each other. In this asymptotic limit, it is consistent to describe the effect of the cilium in the fluid by a stokeslet. We introduce the two-dimensional vector along the cell surface Formula 1 and consider the limit Formula 1. We are interested in the velocity in the vicinity of the ciliary layer, typically 0 < z < 1.5L. At lowest order in z/r, the velocity reads:

Formula 2(2)
where we have defined a dimensionless height Formula 2. In the following, we use mostly the velocity Formula 2. At this order, the flow field is a radial flow centered on the cilium-stokeslet. Note that because of the no-slip boundary condition on the surface, the force appears in the velocity field (Eq. 2) in the combination fh homogeneous to a momentum.

We consider now a regular array of cilia; this is a reasonable assumption for Paramecium, which shows a beautiful and very regular array of cilia on its surface (20Go). We assume that this array is in an infinite plane. Cilium i is defined by its position in the xy plane (the cell surface) by a vector Formula 2 and by the angle of its plane of beating with the x axis, {phi}i as displayed on Fig. 2. The total velocity at the cilium i at height z, Formula 2 is the sum of all the velocities Formula 2 created by the other cilia j != i:

Formula 2


Figure 2
View larger version (10K):
[in this window]
[in a new window]

 
FIGURE 2  Hexagonal lattice of cilia with a distance d between two neighboring cilia. Cilium j exerts on the fluid a force f in the direction {phi}j; Formula 2, where Formula 2 is the unit vector from cilium j to cilium i.

 
We single out the Formula 2 dependence of Formula 2 writing Formula 2 with Formula 2 and

Formula 3(3)
where Formula 3 is a unit vector from cilium j to cilium i and Formula 3, as shown on Fig. 2. We now use a mean field approximation, replacing the velocity Formula 3 by its average over the directions Formula 3 given by:

Formula 4(4)
where P({phi}) is the probability for a cilium to make an angle {phi} with the x axis. The velocity at cilium i is Formula 4. The mean field approximation assumes that the fluctuations around a given angle determining the direction of the flow are small. We choose the x axis in the direction of the flow without any loss of generality, so that the probability P({phi}) is peaked around {phi} = 0.

To determine the probability P({phi}), we write a stationary Fokker-Planck equation {partial}{phi}J = 0. The probability current J = P{partial}t{phi}Dr{partial}{phi}P is the sum of two terms: a convection term and a diffusion term where Dr is a rotational diffusion coefficient. The beating plane can fluctuate due to thermal fluctuations. Because of the flow, if the beating plane of one cilium is at an angle {phi} with the flow direction, the cilium is subject to a torque Formula 4 along the z axis that tends to align it in the direction of the flow, where Formula 4 is the velocity of the global flow and {alpha} a viscosity coefficient involving the geometry of the cilium. A rotating cilium is also subject to a viscous torque Mviscous opposing the rotation Formula 4 where {zeta} is the rotational friction constant. The total torque on the cilium vanishes (Formula 4) and the probability distribution satisfies the Fokker-Planck equation:

Formula 4

We define the effective temperature as Dr{zeta} = kBT. We introduce the dimensionless velocity Formula 4 and we impose the normalization condition Formula 4. We obtain:

Formula 5(5)
where Formula 5 is the modified Bessel function defined in Gradshteyn and Ryzhik (21Go). The velocity Formula 5 can then be self-consistently determined by calculating Formula 5, and summing over all the lattice sites. We obtain:

Formula 6(6)
where Formula 6 is a constant depending on the nature of the lattice, d is the lattice constant (the distance between cilia), and Formula 6 is a modified Bessel function (21Go). For a hexagonal lattice, the explicit calculation yields:

Formula 7(7)
(for a square lattice Formula 7). The self-consistent Eq. 6 for the flow velocity can be discussed by expanding the integrals I0 and I1 in the vicinity of Formula 7: there are two solutions Formula 7, and a solution at a finite velocity that exists only within a certain range of parameters:

Formula 8(8)

This solution exists only if:

Formula 9(9)

Within the mean field approximation Eq. 9 defines a dynamical phase transition between a nonmoving fluid with randomly oriented cilia and a moving fluid with a global flow Formula 9 given by Eq. 8 where all cilia are spontaneously aligned in the same orientation. This dynamical phase transition is second order (with a continuous velocity at the transition) and it is associated to a spontaneous breaking of the initial O(2Go) symmetry.

The influence of some of the parameters can be directly analyzed on Eq. 9. A decrease of the distance d between two cilia favors the alignment, increasing the hydrodynamic coupling. A decrease of the temperature T also favors alignment as the random thermal motion opposes it. An increase of the effective hydrodynamic force of one cilium f is associated to an increase in the interactions between cilia and leads to a better alignment. The same effect occurs for {alpha} and the cilium length L. Finally, increasing h helps to create a global flow, since the velocity on the membrane vanishes and the higher the force is exerted, the more efficient.

An analysis only involving the geometrical and kinetic parameters of the cilium requires the estimation of the parameters f, {alpha}, and h. Even if the exact motion of the cilium is necessary for a quantitative study, it is worth going further qualitatively in the analysis. The height h is of the order of the cilium size h ~ L. The calculation of {alpha} is given in Appendix I for a general beating (Eq. 48). It turns out that {alpha} is linked to the difference of the areas covered during the effective and recovery strokes. Here we approximate Formula 9, where {xi}{perp} is the perpendicular friction constant per unit length of the cilium (10Go,22Go) and Formula 9 the amplitude of the tip movement.

To give a simple estimation of the effective hydrodynamic force, we consider that the friction coefficient takes the perpendicular value {xi}{perp} during the effective stroke, and the parallel value {xi}|| during the recovery stroke. Introducing the beating frequency {omega}, we estimate Formula 9 (see Eq. 48). A more precise calculation of these two quantities as a function of the beating patterns is given in Appendix I. The numerical factor is such that Formula 9. Consequently, we obtain:

Formula 9
The difference between the two local drag coefficients {xi}{perp} and {xi}|| is the key to an efficient beating. Increasing the amplitude Formula 9 or the frequency {omega} of the beating favors the alignment of cilia, as could be expected. Increasing the viscosity of the medium could also promote the transition at first sight because it seems to increase the coupling between cilia (since {xi}{perp}, {xi}{perp} {propto} {eta}). Nevertheless, this last argument requires more attention because the beating pattern (force exerted by the cilium, amplitude) and its frequency depend on the medium viscosity. Indeed, in Machemer (23Go) measures a decrease of the frequency when {eta} increases, as we get in "Axonemal beating". Moreover, we can also expect a decrease of Formula 9 with {eta}. Thus, increasing viscosity may lower the coupling between cilia (as suggested in Gheber et al. (24Go)) and disadvantage their alignment.

Equation 9 is satisfied even for rather large d. Indeed, taking Formula 9, {xi}{perp} ~ 10{eta}water (25Go), {omega} ~ 102rad s–1, kBT ~ 10–21J we get:

Formula 9

Thus, Eq. 9 remains true even for d ~ 10L: hydrodynamic effects should be sufficient to induce alignment in several cases since for many ciliated cells, d < L. In the following, we take a closer look at the internal beating mechanism of one cilium and then at the beating of an array of cilia to obtain a more precise and quantitative description.


    AXONEMAL BEATING
 TOP
 ABSTRACT
 INTRODUCTION
 SPONTANEOUS ALIGNMENT OF AN...
 AXONEMAL BEATING
 LEFT-RIGHT BEATING SYMMETRY...
 DISCUSSION AND CONCLUDING...
 APPENDIX I: AVERAGE FORCE...
 APPENDIX II: FOURIER MODE...
 ACKNOWLEDGEMENTS
 REFERENCES
 
In this section we discuss the beating mechanism of a single cilium. We follow closely the work of Camalet and Jülicher (11Go), which mimics the cilium by two microtubule filaments sliding along one another under the action of the dynein motors and uses a two-state model to describe the collective motion of the dyneins. We use as boundary conditions for the motion those introduced recently (I. Riedel-Kruse, A. Hilfinger, J. Howard, and F. Jülicher, unpublished data; (26Go,27Go)) that seem to have good experimental support (28Go). In the next section, we use the same model to discuss the coordination between cilia.

Equation of motion
Each microtubule doublet within the axoneme can be described effectively as an elastic rod. Deformations of this rod lead to local sliding displacements of neighboring microtubules. Here, we only consider planar deformations in the (xz) plane. In this case the geometrical coupling between bending and sliding can be captured by considering two parallel elastic filaments (corresponding to two microtubule doublets) with a constant separation a along the whole length of the rod (see Fig. 3). At one end, which corresponds to the basal end of an axoneme, the two filaments are elastically attached and are allowed to slide with respect to each other, but not to tilt (26Go,29Go). The basal connection is characterized by an elasticity k and a frictional drag {gamma}. The configurations of the axoneme are described by the shape of the filament pair given by the position of one filament Formula 9 at arclength s. The shape of the other filament is then given by Formula 9, where Formula 9 is the filament normal. In the following, we describe the filament conformation by the angle {psi} between the local tangent vector and the z axis or by the deformation h in the transverse direction defined in Fig. 12.


Figure 3
View larger version (11K):
[in this window]
[in a new window]

 
FIGURE 3  Two filaments (full curves) Formula 9 and Formula 9 at constant separation a are connected to the membrane at the bottom end where s = 0. Internal forces f(s) are exerted in opposite directions, tangential to the filaments. The displacement {Delta} at the tip s = L is indicated.

 

Figure 12
View larger version (8K):
[in this window]
[in a new window]

 
FIGURE 12  Sketch of a beating cilium in a plane at an angle {phi} with the direction of the external flow V.

 
The energetics of the filament pair is due to the bending elasticity. In addition to filament bending, we also take into account internal stresses due to the active elements (dyneins). We characterize them by the force per unit length f(s) acting at position s in opposite directions on the two filaments. This force density corresponds to a shear stress within the cilium that tends to slide the two filaments with respect to each other. The local curvature is C = {partial}s{psi} (see Fig. 12). The sliding displacement, {Delta}(s, t) is related to the sliding displacement at the base {Delta}0(t) by {Delta}(s, t) – {Delta}0(t) = a{psi}(s, t), because we impose the boundary condition {psi}(0, t) = 0.

A configuration of a filament pair of length L is associated to the free energy functional:

Formula 9

Here, {kappa} denotes the total bending rigidity of the filaments. The inextensibility of the filaments is taken into account by the Lagrange multiplier {Lambda}(s), which enforces the constraint Formula 9. The first term of this expression is the elastic energy due to the basal sliding occurring with a connection of elasticity k. The tangent component of the integrated forces acting on the filament between s and L is denoted by {tau}(s). Assuming that there is no external force applied at the end s = L of the cilium:

Formula 9

We assume for simplicity that the hydrodynamic effects of the surrounding fluid can be described by two local friction coefficients {xi}{perp} and {xi}|| for normal and tangential motion. The total friction force per unit length exerted by the cilium on the fluid is then Formula 9. The force balance at arclength s can then be written as:

Formula 10(10)
which leads to:

Formula 11(11)
where the derivatives with respect to arclength have been denoted by a dot.

The beating of the filaments is very sensitive to the boundary conditions imposed at its ends. As the force density Formula 11 is equilibrated by the density of friction force exerted by the fluid on the system (see Eq. 10), the boundary contributions coming from the free energy variation {delta}G are equilibrated by external forces Formula 11 and torques Formula 11 applied at the ends. At the free end of the cilium, both the external force and the external torque vanish and:

Formula 12(12)

At the base, s = 0, the boundary conditions are:

Formula 13(13)

The external torque and force are chosen in such a way that the base is fixed (Formula 13) and that the cilium remains perpendicular to the surface (Formula 13 or {psi}(0) = 0). The final boundary condition is associated to the basal sliding

Formula 14(14)

The determination of the cilium motion requires a model for the shear force created by the dyneins that we calculate using a two-state model.

Two-state model for the cilium
Following Camalet and Jülicher (11Go), we now introduce the two-state model of coupled molecular motors (30Go,31Go) to describe the internal mechanism of the cilium. This model allows the calculation of the shear force f due to the dyneins.

Each motor has two different chemical states, a strongly bound state, 1, and a weakly bound state, 2. The interactions between a motor and a filament in both states are characterized by potential energy landscapes W1(x) and W2(x), where x denotes the position of a motor along the filament. The potentials have the filament symmetry: they are periodic with period l, Wi(x) = Wi(x + l) and are, in general, spatially asymmetric, Wi(x) != Wi(–x).

In the presence of ATP, the motors undergo transitions between states with transition rates {omega}1 and {omega}2. Introducing the relative position {xi} of a motor with respect to the potential period, (x = {xi} + nl with 0 ≤ {xi} < l and n an integer), we define the probability Pi({xi}, t) for a motor to be in state i at position {xi} at time t. The relevant Fokker-Planck equations are:

Formula 14

Formula 14
where {nu} = {partial}t{Delta} = a{partial}t{psi}(s) is the sliding velocity between the two filaments.

The simplest choice of the two potentials Wi, is a saw-tooth potential (with barrier height Formula 14) representing a strongly bound state for W1, and a flat potential W2 representing a weakly bound state. Here, for simplicity, we use the symmetric potentials:

Formula 14

Formula 14

Although this choice is somehow arbitrary, we checked that the final results only depend qualitatively on the actual shape of the potentials and of the rates defined below.

When a number of motors act together to propel a filament, the direction of motion is a collective property: the filament might move in either direction (31Go). The absence of asymmetry in the potentials implies that an individual motor is not able to move directionally. It is not the case for an assembly of motors: even with a symmetric potential, provided that detachment can only take place at a localized position near the bottom of a potential well, oscillations can occur.

We define the distance from equilibrium {Omega}:

Formula 14

{Omega} is related to the chemical potential difference between ATP and its hydrolysis products, {Delta}µ = µATPµADP µP. At equilibrium, {Delta}µ = 0 and {Omega} = 0. We assume for simplicity that the binding rates {omega}2 and the detachment rate {omega}1 are given by:

Formula 14

Formula 14

Note that, with this choice the sum {omega}1 + {omega}2 = {nu}(1 + {Omega}) does not depend on {xi}, and that if {Omega} = 0, {omega}1 = 0 and no directional movement is possible. Here {nu} is a constant transition rate. If we assume that the motors are uniformly distributed along the filaments with a density {rho}, the probabilities P1 and P2 satisfy the relationship P1 + P2 = {rho}. The Fokker-Planck equations reduce then to a single equation for P = P1:

Formula 15(15)

This model leads to an expression for the shear force per unit of length f(s, t) created by the dyneins and driving the cilium beating. Using the results of Prost et al. (30Go), and the fact that W2 is a constant:

Formula 16(16)
where K is an elastic stiffness per unit length mimicking the influence of the nexins that are proteins acting as springs in the axonemal structure, and {lambda} is an internal friction coefficient per unit length modeling the friction encountered by the motors. Equations 11 and 15 allow in principle a complete calculation of the beating motion.

In the following, we assume that the beating occurs with a "small" amplitude, which means that both Formula 16 and Formula 16. A quick look at the beating pattern of a cilium of Paramecium shows that the beating occurs with a large amplitude. Nevertheless, this approximation allows us to extract interesting information on the parameters controlling the beating. Moreover, the work of Hilfinger shows that larger amplitude beating patterns are very similar to small amplitude patterns (26Go). We must however keep in mind that our approach is valid only if the system stays close to an oscillation bifurcation, which is consistent with the fact that we consider only small movements.

We use the deformation h, rather than {psi} or {Delta} to describe cilium motion and work at second order in |h| so that Formula 16. In the absence of any external flow, the equation of motion 11 projected on Formula 16 imposes that Formula 16 (10Go,11Go). The projection of the equation of motion on Formula 16 then yields:

Formula 17(17)

The nonlinear terms are not important here but they will turn out important in the following section. Indeed, experiments show that the beating is clearly asymmetrical (ES and RS), and we must expand at least to the next order if we want to capture this phenomenon. With this variable and at this order, the boundary conditions read:

Formula 18(18)

The explicit solution of the equation of motion 17 is obtained by Fourier expansion in time Formula 18. The explicit derivation of the equation satisfied by the Fourier components is given in Appendix II. At linear order the effect of the motors is characterized by a susceptibility:

Formula 19(19)

Using dimensionless variables, Formula 19, Formula 19, Formula 19, and Formula 19, the equation of motion for the Fourier component n reads:

Formula 20(20)

Note that Formula 20 where Sp is the "sperm number" introduced in Lowe (32Go). The boundary conditions at the base Formula 20 are:

Formula 21(21)

At the free end of the cilium, Formula 21:

Formula 22(22)
with Formula 22 where we have introduced the dimensionless parameters Formula 22 and Formula 22.

Beating pattern
In the absence of any external flow the beating is symmetric and the Fourier components n = 0 and n = 2 of h vanish for symmetry reasons. Close to the oscillation bifurcation threshold, cilium beating is dominated by the first Fourier component n = 1. The solution of the linear equation of motion 20 is written as the sum of four exponentials:

Formula 23(23)
with

Formula 23
where the two inverse decay lengths are given by:

Formula 24(24)

The boundary conditions are explicitly discussed in Appendix II. The condition for existence of nonvanishing solutions is given by Eq. 57 of Appendix II. This is a complex equation that gives therefore two conditions that determine both the critical value of the distance from equilibrium where the oscillations start {Omega}c, and the reduced oscillation frequency Formula 24. The critical value {Omega}c is a Hopf bifurcation threshold: there are no oscillations if {Omega} ≤ {Omega}c and cilium beating is only possible if {Omega} ≥ {Omega}c. At the bifurcation threshold, the amplitude of the oscillations vanishes. It is not possible to calculate the amplitude of the oscillations above the bifurcation threshold with the linear theory presented here. This requires a complete determination of the third order terms in the equation of motion that goes far beyond the scope of this work. This very complex problem is addressed in the work of Hilfinger and Jülicher (26Go).

Determination of the Hopf bifurcation threshold and the beating frequency at this threshold does not seem to be possible analytically. We therefore rely on a numerical solution, choosing reasonable values of the various parameters. We study a cilium of length L = 12 µm, which is the length of a Paramecium cilium. It moves in a fluid of viscosity {eta} ~ 4{eta}water = 410–3 Pa s, which is higher than the water viscosity to take into account the proteins above the cell body. We estimate {kappa} = 410–22Nm2 corresponding to 20 microtubules and to what is measured in Ishijima and Hiramoto (33Go). The other parameters are l = 10 nm, a = 20 nm, U = 10 kT, Formula 24, {rho} = 5108m–1, and {lambda} = 1 Pa s very similar to those of Camalet (10Go). To match the typical frequency observed in Paramecium, and to obtain a realistic pattern of the beating, we take {xi}{perp} = 35{eta} = 1410–2 Pa s and we choose k = 6.5 Nm–1, {gamma} = 7.3{eta}water, and {nu} = 600 s–1. The value of {xi}{perp} is rather high, but it must include the hydrodynamic interactions of the cilium with the cell surface (34Go), which are not taken into account if we use the classical friction per unit length of a cylinder (35Go).

The numerical resolution of Eq. 57 then yields Formula 24. This corresponds to a critical beating frequency Formula 24, which is the typical value for Paramecium (7Go). The calculated beating pattern of the cilium is shown on Fig. 4. It corresponds to a superposition of waves propagating from the base of the cilium to the free tip as observed experimentally.


Figure 4
View larger version (25K):
[in this window]
[in a new window]

 
FIGURE 4  Approximate cilium deformation Formula 24 at different times steps (corresponding to different gray levels) during a beating period with Formula 24. The beating is symmetrical with respect to the vertical axis. Deformations are propagating from base to tip.

 
A detailed study shows that the equation giving the bifurcation threshold and the beating frequency has several solutions Formula 24 with Formula 24. A first guess would be that the axoneme starts beating at the lowest threshold Formula 24. However, it is known experimentally that during the beating the deformation waves propagate from the base to the tip and not from the tip to the base (36Go). The first two oscillating modes correspond to waves propagating in the opposite direction. To be consistent with the observations we do not consider them here. A better choice of the transition rates {omega}1 and {omega}2 would perhaps allow to justify this choice. The direction of propagation of the wave is extremely sensitive to the boundary conditions. We have allowed here basal sliding as suggested by some experiments and we have imposed that the cilium is clamped at its base with an angle {psi} = 0. This also seems consistent with some experiments analyzed in Hilfinger (26Go). The other extreme limit of a completely free cilium (a vanishing external torque at the base) leads to a wave propagating form the base to the tip for the first mode. We have tried to use an intermediate boundary condition where the torque at the base is an elastic torque and varying the related stiffness, however, we were not able to obtain a beating pattern looking like the experimental one. We therefore proceed, considering only the third beating mode Formula 24.

On the other hand, it is not surprising that the lowest beating threshold corresponds to a tip to base propagation because the tip is free to move whereas the base is not. The same phenomenon has been observed in Dreyfus et al. (25Go): the authors construct an artificial flagellum with magnetic particles. One end is attached to a red blood cell and the other end is free to move. The beating is driven by an applied transverse oscillating magnetic field. As a result, a bending wave propagates from the free end, which is more mobile. The easiest way to move a filament being to start with its free end, our result may not be an artifact of the model. The cell must then impose an ATP concentration such that Formula 24 so that the third mode is selected.

The beating frequency fc varies with the viscosity of the medium. Experimentally, when methyl-cellulose is added in water, the viscosity increases significantly. We predict here a decrease of fc with increasing external viscosity as observed in the experiments of Machemer (23Go) and in numerical simulations (37Go) (see Table 1). We observe an approximate linear decrease of the beating frequency when plotted against log({eta}/{eta}w) as in the simulations performed in Gueron and Levit-Gurevich (37Go) (we find similar values for the frequency).


View this table:
[in this window]
[in a new window]

 
TABLE 1  Viscosity and beating frequency

 
The effect of [Ca2+] on the beating can also be studied qualitatively. As mentioned in the introduction, [Ca2+] has a strong influence on the beating pattern. Calcium concentration variations are at the basis of the shock responses of many organisms, changing the ciliary-type beating into a flagellar-type beating in Chlamydomonas, or switching the directions of effective stroke and recovery stroke in Paramecium, or in reversing the direction of the "wave" propagation on the flagellum and thus reversing the direction of the movement in Chritidia (38Go).

This last example can be explained qualitatively within our approach. In Chritidia, both directions are possible for the deformation wave propagation. Calcium may affect the chosen mode of beating, allowing the system to choose Formula 24 and its tip-to-base pattern instead of Formula 24.

Calcium is also likely to change the attachment/detachment rate (and thus change the parameter {nu}) or the boundary conditions at the base of the cilium (and thus change k). In Chlamydomonas, calcium has a contractile effect on the striated fibers connecting the basal bodies of the two flagella (39Go). These changes induce a change in the beating pattern, and may result in a switch from base-to-tip to tip-to-base wavelike propagation.


    LEFT-RIGHT BEATING SYMMETRY BREAKING AND METACHRONAL COORDINATION
 TOP
 ABSTRACT
 INTRODUCTION
 SPONTANEOUS ALIGNMENT OF AN...
 AXONEMAL BEATING
 LEFT-RIGHT BEATING SYMMETRY...
 DISCUSSION AND CONCLUDING...
 APPENDIX I: AVERAGE FORCE...
 APPENDIX II: FOURIER MODE...
 ACKNOWLEDGEMENTS
 REFERENCES
 
In the presence of a transverse external flow, the beating can no longer be symmetrical as sketched on Fig. 5. The cilium tends to beat faster and quite straight in the direction of the flow, whereas it comes back slower and more curved against the flow. This looks like a two-phase beating with an effective and a recovery stroke.


Figure 5
View larger version (12K):
[in this window]
[in a new window]

 
FIGURE 5  Effect of an external flow Formula 24 on the beating of a single cilium. (a) Symmetrical beating. (b) Broken symmetry due to the external flow.

 
If the beating is asymmetrical, the cilium exerts a force in the fluid that can itself produce a flow. In a certain range of parameters, one can therefore expect that a continuous flow is spontaneously generated by hydrodynamic interactions between cilia: an assembly of cilia, beating symmetrically, is able to break spontaneously this left-right symmetry of the beating to create a global flow. This idea of a spontaneous breaking of the left-right symmetry has already been suggested (40Go) with a more abstract system (called rowers) having two internal energy states.

In this section, we first study the effect of an external velocity imposed by the experimentalist on the beating symmetry of a single cilium. We then consider an array of aligned cilia and determine the conditions under which this assembly of cilia breaks its left-right symmetry and generates a global flow. Metachronal coordination between cilia naturally emerges from hydrodynamic coupling as a local minimum of the oscillation threshold {Omega}c.

External breaking of the beating symmetry: cilium submitted to an external flow
We impose an external flow Formula 24 along the x axis for simplicity, and satisfying the no-slip boundary condition on the cell membrane as in "Spontaneous alignment of an array of cilia: a simple model", with the same notations (Formula 24). It is found experimentally that the velocity above the cilia sublayer is time independent (18Go), justifying our choice. This flow is in this first part externally fixed and we consider the limit of vanishingly small flows. The force per unit length exerted by the cilium on the fluid Formula 24 depends on the external velocity Formula 24.

Formula 25(25)

The equation of motion 11 reads then:

Formula 26(26)

The boundary conditions are the same as in the absence of the neighboring cilia and are given by Eqs. 1214. Following the same procedure as for a cilium in the absence of flow, we find the equation of motion for the deformation of the cilium h:

Formula 27(27)

The introduction of the external flow breaks the h ->h symmetry (or left-right symmetry) introducing in Eq. 17 terms of zeroth and second order in h in the equation of motion. The boundary conditions do not depend on the external flow. Because the height z where we consider the flow depends on the arclength s, we write (expanding Eq. 47):

Formula 28(28)

Note that we cannot simply write Formula 28 and get rid of the second term because it is of order |h|2. As above, we expand the deformation of the cilium h in Fourier components in time. Using the same notations as before, the equation of motion of the Fourier components can be written as:

Formula 29(29)
for n = 0, 1, 2 and where we have introduced the new dimensionless parameters:

Formula 29

In the limit of small external velocities, we have neglected terms of order Formula 29, and we keep this approximation throughout our calculations.

The equation for the first mode is identical to Eq. 55, with the same boundary conditions. At this order in U, the fundamental mode is not affected by the external flow. Consequently, the oscillation threshold and the beating frequency are the same as in the absence of flow and the Fourier component Formula 29 is given by Eq. 23.

The zeroth Fourier component Formula 29 gives the average deformation of the cilium. It is a solution of:

Formula 30(30)
with the same boundary condition as before. Nevertheless, Formula 30 does not vanish at first order in velocity because of the broken symmetry due to the external flow that is reflected in the right-hand side of Eq. 30. We write the solution as the sum of two contributions: Formula 30, corresponds to the curvature of the cilium under the flow V(z) at equilibrium, i.e., in the absence of any beating (Formula 30), and Formula 30, corresponds to the corrections to this equilibrium deformation due to the beating when there is enough ATP in the medium: Formula 30. If as above, we ignore the elasticity of the nexins (Formula 30):

Formula 31(31)
where Formula 31 is the amplitude of the first Fourier mode of the oscillation defined in Eq. 23 and Formula 31 is a function calculated numerically only depending on Formula 31 (see Eq. 23). In the limit U = 0, Formula 31 as expected. The average deformation of the cilium is plotted on Fig. 6 a, which shows the bent shape under the action of the external flow.


Figure 6
View larger version (29K):
[in this window]
[in a new window]

 
FIGURE 6  (a) Average position of a cilium that is curved in the direction of the flow, Formula 31. (b) Second Fourier component of the deformation Formula 31 at different times during a beating period. The scale is dilated: Formula 31 with the parameters Formula 31, Formula 31, and Formula 31.

 
The second Fourier component gives the asymmetry of the beating. It is obtained from the equation of motion:

Formula 32(32)

We write the solution of Eq. 32 as:

Formula 32
where Formula 32 is calculated numerically from Formula 32 as well. Here also, in the limit U = 0, Formula 32. The plot Formula 32 at different times on Fig. 6 b, leads to a complicated pattern.

The total deformation of the cilium Formula 32. is plotted at different times equally spaced on Fig. 7. To stress the fact that the beating is easier and faster in the direction of the flow, and more difficult and slower against the flow, we have chosen rather large values of the parameters, Formula 32, Formula 32, and Formula 32, and we plot Formula 32 for Formula 32 on Fig. 7.


Figure 7
View larger version (22K):
[in this window]
[in a new window]

 
FIGURE 7  Beating pattern at the base of the cilium (Formula 32) with the parameters Formula 32, Formula 32, and Formula 32: the cilium beats faster in the direction of the flow and slower in the opposite direction around a curved average position.

 
The external flow thus breaks the left-right symmetry in two ways. First the average position of the cilium is not the vertical axis but a cilium curved in the direction of the flow. Second, the beating itself is no longer left-right symmetric: the cilium goes faster in the direction of the flow and comes back slower against the flow. The beating pattern looks like a two-phases beating with an effective stroke and a recovery stroke. This effect increases with Formula 32 and U and is consequently even more significant in real ciliary beatings. The external flow may therefore be an important factor in the asymmetry of the beating.

Another important result is that, because the beating propagates a base-to-tip deformation, the curved cilium exerts a finite average force in the fluid in the direction of the flow. Thus, if an external flow breaks the left-right beating symmetry, the cilia create a force in its direction and can amplify this flow. This is the basis of the left-right spontaneous symmetry breaking that we discuss in the next section.

The external flow is not always the only source of symmetry breaking. If it were so, Paramecium would always go in the same direction once it started moving. This is not the case because this organism is able to go backward when it bumps into an obstacle, thanks to the release of calcium that reverses the beating. Nevertheless, some ciliary assembly never reverse their beating and create a flow in a unique direction. This is the case for the respiratory, ependymal and Fallopian cilia for instance. In those examples, the flow is transverse and always in the same direction. Our approach may be relevant for this category of cilia for which an external flow might take a significant part in the breaking of the beating symmetry.

Spontaneous breaking of the beating symmetry: array of aligned cilia
We now consider a regular array of cilia on a cell body, beating all in the same direction. Starting from a symmetrical beating, we show that the left-right symmetry is spontaneously broken within a certain range of the parameters controlling the beating due to the hydrodynamic couplings between cilia.

Equations of motion
For a cilium located in the xy plane at position Formula 32, we call Formula 32 the velocity created by the other cilia at the point Formula 32 of arclength s. The equation of motion of the cilium is similar to that obtained previously with an external flow field and we write up to the third order in h as:

Formula 33(33)
where the projection of the local external velocity on the cilium normal is:

Formula 34(34)

The boundary conditions for the motion are the same as in the previous section. The velocity Formula 34 created at arclength si of the cilium i by a cilium j is given by:

Formula 34
where G is the second order hydrodynamic tensor given by Eq. 1 and Formula 34 is the force per unit of length created by the beating of the cilium j. The total velocity at the arclength si of the cilium i is thus given by:

Formula 34
As in "Spontaneous alignment of an array of cilia: a simple model", we consider the limit Formula 34 and we only keep terms of the second order in s/r, r being the distance between two cilia in the array so that:

Formula 34

This means that only the velocity along the x axis created by the component fx of Formula 34 plays a role and that we can ignore the other component Vz of the velocity appearing in Eq. 34. Using the notations of Fig. 2, and noting that Formula 34, we obtain:

Formula 35(35)

As in the previous sections, we expand the velocity, the force, and the cilium deformation in Fourier modes in time. For simplicity, we only consider here the first two Fourier components and do not look at h2 that characterizes the asymmetry of the beating. The Fourier components of the velocity are related to those of the force by:

Formula 35

The Fourier components of the force f0 = < fx > and f1 are calculated using the expression of fbeat and its average over one time period given by Eq. 48 in the small movements approximation:

Formula 36(36)
where Formula 36 is the imaginary part of a complex number. We assume that all cilia are identical, and that they all beat with the same pattern. The only difference in the beating patte