## ABSTRACT

The body of vertebrate embryos forms by posterior elongation from a terminal growth zone called the tail bud. The tail bud is a source of highly motile cells that eventually constitute the presomitic mesoderm (PSM), a tissue that plays an important role in elongation movements. PSM cells establish an anterior-posterior cell motility gradient that parallels a gradient associated with the degradation of a specific cellular signal (FGF) known to be implicated in cell motility. Here, we combine the electroporation of fluorescent reporters in the PSM with time-lapse imaging in the chicken embryo to quantify cell diffusive movements along the motility gradient. We show that a simple microscopic model for random cell motility induced by FGF activity along with geometric confinement leads to rectified tissue elongation consistent with our observations. A continuum analog of the microscopic model leads to a macroscopic mechano-chemical model for tissue extension that couples FGF activity-induced cell motility and tissue rheology, and is consistent with the experimentally observed speed and extent of elongation. Together, our experimental observations and theoretical models explain how the continuous addition of cells at the tail bud combined with lateral confinement can be converted into oriented movement and drive body elongation.

## INTRODUCTION

Most vertebrate species exhibit an elongated body axis. This characteristic pattern is established during embryogenesis as the tissues progressively form in an anterior-to-posterior direction. Microsurgical ablation of the posterior presomitic mesoderm (PSM), which contains the precursors of skeletal muscles and the axial skeleton, severely reduces posterior elongation movements, indicating that this tissue plays a major role in the control of posterior elongation of the embryonic axis. Analysis of cell motility in the chicken embryo PSM (Bénazéraf et al., 2010) shows that there is an anterior-posterior gradient in the activity of cells. However, locally, the motility inside the PSM of chicken and zebrafish embryos is manifested by random, undirected, Brownian-like cellular motion (Bénazéraf et al., 2010; Lawton et al., 2013). These random, diffusive movements contrast with the oriented cell intercalation movements controlling elongation of the anterior parts of the embryo. Within the PSM, this motility gradient is downstream of a chemical gradient of the secreted fibroblast growth factor FGF, as shown schematically in Fig. 1A. FGF is known to play an important role in cell motility (Delfini et al., 2005); indeed, increasing FGF concentration in the PSM causes cell motility to increase without any orientation preference (Bénazéraf et al., 2010), reducing the speed of elongation. Earlier qualitative models focused on the role of the motility gradient as a driving force for elongation (Bénazéraf et al., 2010), but how this motility gradient is related to elongation is still unclear.

The observation that at a local level cell motility is random and undirected needs to be reconciled with the emergence of global body elongation. As new cells enter the PSM, they are exposed to a high concentration of FGF and become highly motile, but do not move in an oriented manner. Owing to the confinement of motility driven by the presence of relatively immobile and stiff lateral tissues (Bénazéraf et al., 2017), the cells generate an effective mechanical pressure. Thus, a potential mechanism for the rectification that leads to the observed directional cell velocity and body elongation could be this mechanical pressure. Because the expression of FGF is highest at the posterior PSM and decreases away from it, there is a natural reduction in motility anteriorly, consistent with observations. Eventually, the effects of anterior adhesion dominate over motility and this causes cellular condensation into somites.

To quantify these processes, here we use experimental observations to measure the effective diffusivity of cells as well as their advection speed as a function of their location relative to the last formed somite. Our observations suggest a minimal microscopic cellular description of a zone of proliferating cells with high motility, which we use to develop both a quantitative cellular model and an equivalent macroscopic continuum theory, thus providing a theoretical and computational framework for body elongation. These complementary approaches yield simple expressions for the speed and scale of body elongation that are consistent with our experimental measurements, with potential implications for our understanding of outgrowth-driven morphogenesis in other settings in which FGF-driven gradients in motility and contractility drive the extension of the limb bud (Gros et al., 2010) and the gut (Nerurkar et al., 2019).

## RESULTS

### Experimental observations

Our experiments were carried out with chicken embryos at Hamburger-Hamilton (HH) stages 10-11 (Hamburger and Hamilton, 1992), corresponding to the period when the elongation of the embryo is most substantial (Denans et al., 2015). A time series of an elongating PSM is shown in Fig. 1B (see also Movie 1). To measure the elongation rate, we registered the movement of the embryo with respect to the last-formed somite at the beginning of the experiment (*t*=0) as depicted by the black dashed line in Fig. 1B, and tracked the advancement of the Hensen's node, *L*_{e}(*t*) as a function of time. In Fig. 1C, we see that *L*_{e}(*t*) increases linearly with time, with a mean elongation rate *V*=(2.8±0.3)×10^{−2} μm/s (averaged from five embryos).

To evaluate the role of cell motility on overall body elongation, we then examined the movement of cells by electroporating the PSM cells with fluorescent reporters that label cell nuclei specifically (Hatakeyama and Shimamura, 2008), as shown in Fig. 1D (see also Movie 1). As a first approximation, the motion of the cells can be considered two-dimensional in the anterior-posterior and medio-lateral directions because the relative dorso-ventral depth of the PSM is small. In the reference frame of the last-formed somite at *t*=0, during a fixed acquisition time of 4 h, posterior cell trajectories showed a larger net displacement than the anterior ones, consistent with prior experiments (Bénazéraf et al., 2010). To quantify the variations of cell motility along the body axis, we divided the PSM into three regions (anterior, middle and posterior), as shown schematically in Fig. 1E, and obtained the mean square displacement (MSD) of the cells given by < Δ*r*^{2}(*t*)>, where defines the distance that the cell travels in a lag time *t*, as shown in the inset of Fig. 1F (Qian et al., 1991; Wirtz, 2009).

*D*is the effective cell diffusivity and

*v*is the local population drift velocity (Qian et al., 1991). This model for cell dynamics is in accordance with previous findings for chicken (Bénazéraf et al., 2010) and zebrafish (Lawton et al., 2013) PSM elongation. Fig. 1F shows the mean MSD curves of the three regions (anterior, middle, posterior) defined in Fig. 1E and the fit to Eqn 1 (black dashed curve). From the fits, we obtain

*D*

_{Post}=(3.5±0.7)×10

^{-2}μm

^{2}/s,

*D*

_{Mid}=(2.1±0.4)×10

^{-2}μm

^{2}/s,

*D*

_{Ant}=(1.4±0.3)×10

^{-2}μm

^{2}/s,

*v*

_{Post}=(2.0±0.2)×10

^{-2}μm/s,

*v*

_{Mid}=(1.1±0.1)×10

^{-2}μm/s and

*v*

_{Ant}=(0.6±0.1)×10

^{-2}μm/s (mean±s.e.m.). These estimates confirm the presence of a motility gradient of cells along the AP axis. [Note that one might argue that the ‘drift’

*v*could be the result of a persistent random walk rather than a drift. However, this kind of dynamics will not show a gradual change in the slope of 〈Δ

*r*

^{2}(

*t*)〉 as is observed here.]

Our observations quantify how the cells move in space and time along the PSM. When new cells are added close to the tail bud (TB), they are highly motile but move in random directions; as they move further away from the zone of proliferation in the TB, on average they gradually slow down and stop moving.

Next, we provide microscopic and macroscopic chemo-mechanical mathematical and computational models based on the observed gradient in motility and the confinement due to the more rigid lateral tissues, and show that they explain the observed rectified cell motility and the resulting posterior body elongation.

### Microscopic cellular model

*r*^{i}(

*t*), assuming that inertial effects are negligible, is:where γ is the viscous friction coefficient experienced by each cell, is a potential that prevents cell-cell overlap, and

*ζ*_{i}(

*t*) is random force. The viscous friction is a result of the interaction of cells with their environment and with each other. The repulsive interaction between cells of diameter

*a*guarantees that two cells cannot occupy the same position; a simple form that suffices for this is given by:andNote that in our simulations the cell-cell excluded volume interactions are conservative and thus collisions are elastic but overdamped. The random force

*ζ*_{i}(

*t*) is assumed to have zero-mean and normally distributed with Gaussian statistics such that:andwhere is the single-cell activity/motility and

*ζ*

_{i,α}are the

*x*or

*y*components of

*ζ*_{i}.

*τ*is the slowest time scale associated with kinetics of degradation of FGF, which initiates once the cell enters the TB. We will assume that the initial FGF level and

*τ*are both constant, which is approximately true in the experimental setup for the observed times of 4 h. At the posterior-end corresponding to the TB, we assume that cells move as the body elongates owing to a constant rate at which they are added in the space that is not occupied by other cells (Fig. 2A,B). The movement of the cells from the TB to the PSM region is limited by the available space in the PSM. In our model, if there is no free space, no cells are added, so the rate of adding cells is limited by the motion of the TB, the diffusion and advection of cells away from the TB and the maximal cell packing density

*ρ*

_{0}[the ratio of the area covered by the cells (disks) to the total PSM area]. We do not include cell division in our model, because it has been shown that interfering with cell division in the PSM does not affect the tailbud elongation considerably during the time scale of our analysis (Bénazéraf et al., 2010).

The TB boundary is modeled as a wall of immobile cells attached to their neighbors by elastic springs with spring constants, *k*_{chain}=2*k*, and allowed to move posteriorly as a result of the mechanical pressure exerted by the motile cells anterior to it. Our simulations show that after an initial transient state, the motion of the wall reaches a steady state whereby cells added at a constant rate cause the wall to move at a constant velocity (Movie 2). The quasi-one-dimensional dynamics manifest in the simulations arises because of the impermeable boundary conditions in the direction perpendicular to the direction of elongation.

The fraction of motile cells is much larger near the moving wall, and as one moves anteriorly, this density falls off quickly because of the decrease in FGF activity. Changing FGF activity by changing *τ*, and thus varying the probability with which each individual cells stops moving (and therefore the total fraction of moving cells), changes the velocity and motility profiles (see Fig. 2C,D).

*x*=

*a*,

*y*∈[0,

*Y*

_{max}]) only when there is a gap of order of the size of one cell

*a*(see Fig. 2B). Thus, our model is a cellular analog of the Brownian ratchet introduced in the context of molecular polymerization (Peskin et al., 1993). In the PSM region where the internal tissue resistance to cell motion dominates any external resistance, as in our simulations, the rate of elongation is limited by the waiting time for a gap to open that allows for the addition of a cell of size

*a*, i.e. . Here, is the activity-driven diffusivity of a single cell (consistent with Eqn 1). In this limit of diffusion-limited elongation (

*T*

_{a}≪

*τ*), the speed of elongation scales as

*a*/

*T*

_{a}, i.e.:Similarly, the length scale over which the fraction of motile cells falls off exponentially is:We note that the other limit, i.e.

*τ*≪

*T*

_{a}corresponding to adding cells that are not active, will lead to a jammed state and is not relevant here. In Fig. 2C,D, we show the decay in the speed of the cells and the effective diffusivity as a function of the distance from the TB for different values of

*τ*. In the inset, we see that the rescaled elongation velocity and cell diffusivity using the relations shown in Eqns 9 and 10 as a function of location from the wall is consistent with our simple scaling arguments. However, as can be noticed from the ordinate in the inset of Fig. 2C, the wall velocity in our simulations, , is a consequence of the collective effect of elongation driven by the addition of multiple cells. To test that this is indeed the case, we varied both cell diffusivity and maximal cell packing

*ρ*

_{0}and confirmed that the overall elongation rate (see Fig. S3) is consistent with experimental results (Fig. S4).

### Macroscopic continuum theory

*ρ*(

*x*,

*t*) and the velocity field

*v*(

*x*,

*t*) of motile cells as a function of location

*x*in a fixed lab frame. A hydrodynamic description of the diffusion, advection and degradation of motile cells can then be written in terms of the equations for mass and momentum balance (in the viscously dominated limit) as:

The first equation describes the variations in the density of actively motile cells that also gradually decay at a rate *τ*, whereas the second equation characterizes how the stresses associated with the active pressure generated by the random motility of the cells causes the tissue to respond and push the TB. Here, *η* is the viscosity of the PSM and *ξ* is the viscous friction associated with motion of the elongating PSM relative to the surrounding tissues (endoderm, ectoderm, neural tube and lateral plate). We note that these friction coefficients define a characteristic viscous screening length *l*_{ξ}∼(*η*/*ξ*)^{1/2} that we assume is of the order of a few cell sizes, i.e. *l*_{ξ}∼*a*.

The equations need a closure relation for the pressure; a minimal relation assumes that the active pressure is proportional to the density of motile cells *p*∼*αρ*, consistent with both our microscopic model, with , and also consistent with earlier simulations of active Brownian particles (Mallory et al., 2014). Although other relations of the form *p*∼*αρ*^{q} are plausible, we choose to stay with a minimal model here.

*x*∈[

*s*(

*t*), ∞), with the boundary of the TB being

*s*(

*t*). Motile cells enter the domain at a rate proportional to the difference between their local density and the maximum cell density

*ρ*

_{0}, so that the cell flux is , where is the injection rate of cells from the progenitor zone into the TB. This flux must be balanced by diffusion and advection of cells from the boundary, i.e.which is reminiscent of a generalized Stefan-like condition in moving boundary problems in solidification. It is also necessary to satisfy force balance at the moving boundary so that:where

*F*is the resisting pressure exerted by the tissue ahead of the TB, and . In the cellular simulation, this force is a result of the dynamic friction between the wall cells and the substrate and thus depends on the velocity; here, we have assumed a more general form. Far from the TB, we assume that due to degradation of FGF, the density and velocity of motile cells vanishes so that:Together, Eqns 11 and 12, along with the above boundary conditions (Eqns 13-15), determine the spatiotemporal evolution of the density and velocity fields in the elongating embryo as well as the speed of elongation of the embryo itself.

To understand the dependence of the solution of Eqns 11 and 12 on the problem parameters, we can rewrite the equations in a form that depends on five dimensionless variables:

Here Π_{1} is the ratio of the cell size *a* and the elongation length scale , Π_{2} is the ratio of the maximum active stress *αρ*_{0} and the external viscous friction *ξVL*∼*ξV*^{2}*τ*, and Π_{3} is the ratio of the internal viscous stress *ηV*/*L* and the active stress *αρ*_{0}. We note that the parameter , where *l*_{ξ} is the viscous screening length defined ealier (after Eqn 12). Finally, there are two dimensionless parameters associated with the boundary conditions Π_{4}, which is the ratio of the external stress at the PSM boundary *F* and the active stress *αρ*_{0}, and Π_{5}, which is the scaled ratio of the internal cell diffusivity and the rate of addition of cells at the PSM boundary .

*V*, all the lengths by elongation length scale

*L*=

*Vτ*, and the density by the maximal close packing density

*ρ*

_{0}, the dimensionless form of Eqns 11-15 can be written in terms of the dimensionless variables , and as:andalong with the boundary conditions:andFor parameter values, we assume that , the typical degradation time scale

*τ*≈1−2×10

^{4}s, the viscosity of the PSM

*η*∼10

^{4}−10

^{5}Pa · s and the friction coefficient

*ξ*≈10

^{12}−10

^{13}Pa · s/m

^{2}measured using the micropipette aspiration technique (Guevorkian et al., 2010) (see section 1 in the supplementary information). This yields Π

_{1}≈0.005−0.08 and Π

_{2}Π

_{3}=0.03−10.

Using the following parameter values Π_{1}=0.01, Π_{2}=2, Π_{3}=1, Π_{4}≈0.005, Π_{5}=0.003, we solve the initial boundary value problem (Eqns 17-21) using a finite difference method. In Fig. 3B,C, we see that our calculated profiles for the velocity *v*(*x*) and diffusivity compare well with the experimentally observed profiles, and are suggestive of a simple exponential law (see also Fig. S5). In Fig. 3B, we show that the system evolves into a steady state following a short transient (insets), with an active zone near the TB that moves at constant speed. Both the cellular motility (Fig. 3B) and diffusivity (Fig. 3C) decay towards the anterior where the FGF degrades. Comparing our numerical results for the moving average of cell velocity 〈*v*〉 and the rescaled diffusivity with the experimentally determined values as a function of position relative to Hensen's node (see Materials and Methods) we see that the profiles compare well. Furthermore, on converting the results to dimensional values, we can estimate the active stress *αρ*_{0}≈5−50 Pa, qualitatively similar to measurements on amniote embryos (Zhou et al., 2009; Serwane et al., 2017). From the simulations, we can also estimate the typical length associated with motility decay *L*∼1200 μm and the length scale *l*_{ξ}∼1800 μm associated with friction, which are of the same order as the experimental estimates, *L*=*Vτ*∼260−520 μm and , 3-15 times the cell size.

A sensitivity analysis of our continuum model (see section 3 in the supplementary information) shows that the parameters Π_{1}, the ratio of the cell size to the elongation length scale that arises from activity, and Π_{3}, the ratio of viscous stresses to the active stress, affect the velocity and motility profile significantly. The other parameters, Π_{2}, Π_{4} and Π_{5}, have a negligible effect on our results over a range relevant for our problem (for the range of parameter values used, see Table S3; for the associated results, see Figs S6-S8). The dependence of the length scales corresponding to the typical motility and velocity profiles as a function of the dimensionless parameters Π_{1} and Π_{3} is shown in Fig. S10.

## DISCUSSION

In vertebrate embryos, posterior structures are formed sequentially by a combination of cell proliferation and cell motility that together leads to body elongation. Although it has long been observed that the elongation process involves the posterior displacement of the TB, the physical mechanisms underlying this have not been clearly elucidated. Here, we have quantified this process experimentally and theoretically and shown that it occurs as a result of two effects: the addition of motile cells at a boundary (the TB) and confinement both laterally and anteriorly. These two effects lead to the rectification of random cell diffusivity, which generates the forces underlying elongation and the emergence of characteristic velocity and activity/diffusivity length scales. Simple cellular models and a one-dimensional continuum framework capture the essence of the behavior in terms of a pair of fundamental dimensionless parameters that characterize the ratio of the active driving stress to the internal and external viscous stresses, and are consistent with a narrow boundary layer of activity near the TB, in agreement with observations. Natural extensions of the models include adding the effects of dimensionality, such as the lateral growth in the absence of FGF decay (Bénazéraf et al., 2010), two-dimensional compressive stresses that generate vortical flows near the TB (Xiong et al., 2020), and the gradients in viscosity from the anterior to the posterior region (Mongera et al., 2018). Such multidimensional models could also be modified to describe other embryonic outgrowths processes that exhibit graded diffusive behavior of cells downstream of FGF signaling, such as the budding of vertebrate limbs (Gros et al., 2010) and the onset of gut elongation (Nerurkar et al., 2019).

## MATERIALS AND METHODS

### Chicken embryo preparation and electroporation

Fertilized chicken eggs were obtained from a commercial provider (Les Couvoirs de l'Est, Willgottheim, France) and incubated at 37°C in a humidified incubator. After 24 h, HH stage 4-5 embryos were mounted on filter paper and transferred ventral side up to 35 mm agar/albumen petri dishes for injection (Chapman et al., 2001). Electroporation of the PSM was performed using H2B-mCherry or H2B-Venus nuclear markers as described previously (Bénazéraf et al., 2010). Electroporated embryos were returned to incubator and left to grow to HH stage 10-11 before imaging.

### Time-lapse imaging and track analysis

The imaging procedure used here was similar to a previously described procedure (Czirók et al., 2002; Rupp et al., 2003). Briefly, the embryos were transferred to custom-made six-well observation chambers containing agar/albumin gel, and positioned ventral side up. The time-lapse imaging was performed at 37°C using a motorized upright microscope (Leica DMR, Leica Microsystems) with a 10× objective (N.A. 0.3) and a CCD digital camera (QImaging Retiga 1300i) at a rate of 10 frames/h. At each time point, brightfield and fluorescent images of the embryo were taken at three various fields to cover the total length of the PSM. Post-acquisition processing was performed on images as described previously (Czirók et al., 2002; Rupp et al., 2003) to obtain a 2D time series. Cell tracking and trajectory analysis was performed on fluorescent images using custom-made MATLAB (MathWorks) routines. Trajectories were analyzed on 4-h long movies. To account for cell position modification due to tissue drift, we divided our movies to 1 h segments. The location of each cell was defined with respect to the node at the beginning of each 1 h trajectory. For each cell trajectory, the MSD was calculated and adjusted with Eqn 1 to obtain *D* and *v*. Further, the moving average filter, Smooth (MATLAB), was used for the plots of *D* and *v* as a function of position to obtain Fig. 3.

### Discrete cell simulations

The discrete cell simulations were performed by discretizing Eqn 2 using the Euler–Maruyama method and integrating in time. In the initial state, all the particles are arranged in a square lattice (particles are slightly perturbed from the lattice) in a region that has dimensions of 80*a*×160*a*. At a given time, the particles can be at *x*∈[0,*X*_{max}], *y*∈[0,*Y*_{max}], where *X*_{max} is the position of the rigid boundary and *Y*_{max}=160*a* is the lateral width. We use reduced, dimensionless units in which all energies are given in terms of a typical energy *w*: *U**=*U*/*w*, all lengths in terms of the typical cell size *a*, *r**=*r*/*a* and masses in terms of a mass unit *M*, *m**=*m*/*M*. Using these fundamental units, we can also reduce the effective temperature, i.e. the active random motility , time , fraction *ρ**=*Na*^{3}/*V* (*V* is the volume) and any other physical quantity of interest, where and *ρ** are all dimensionless.

## Footnotes

**Author contributions**

Conceptualization: O.P., L.M.; Methodology: I.R., K.G., L.M.; Software: I.R., A.G.; Validation: A.G., L.M.; Formal analysis: I.R., K.G., A.G., L.M.; Investigation: I.R., K.G., L.M.; Resources: O.P., L.M.; Data curation: K.G.; Writing - original draft: I.R., K.G., L.M.; Writing - review & editing: A.G., O.P., L.M.; Visualization: I.R., A.G.; Supervision: O.P., L.M.; Project administration: L.M.; Funding acquisition: O.P., L.M.

**Funding**

This work was partially supported by the French Agence Nationale de la Recherche (ANR-14-CE32-0009-01 to K.G.), the Human Frontier Science Program (RGP0051/2012 to O.P.), and grants and fellowships from the Schlumberger Foundation and the MacArthur Foundation, and the National Institutes of Health (1R01HD097068 to O.P. and L.M.). Deposited in PMC for release after 12 months.

## References

*Development*

*Nature*

*Nat. Phys.*

*Dev. Dyn.*

*J. Microsc.*

*Proc. Natl. Acad. Sci. U.S.A.*

*Elife*

*Curr. Biol.*

*Phys. Rev. Lett.*

*Dev. Dyn.*

*Dev. Growth Differ.*

*Development*

*Physical Review E*

*Nature*

*Nature*

*Biophys. J.*

*Biophys. J.*

*BioTechniques*

*Nat. Methods*

*Annu. Rev. Biophys.*

*Dev. Cell*

*Development*

**Competing interests**

The authors declare no competing or financial interests.