Contributed by Herbert Levine, November 21, 2012 (sent for review September 9, 2012)
Recent experiments have shown that spreading epithelial sheets exhibit a long-range coordination of motility forces that leads to a buildup of tension in the tissue, which may enhance cell division and the speed of wound healing. Furthermore, the edges of these epithelial sheets commonly show finger-like protrusions whereas the bulk often displays spontaneous swirls of motile cells. To explain these experimental observations, we propose a simple flocking-type mechanism, in which cells tend to align their motility forces with their velocity. Implementing this idea in a mechanical tissue simulation, the proposed model gives rise to efficient spreading and can explain the experimentally observed long-range alignment of motility forces in highly disordered patterns, as well as the buildup of tensile stress throughout the tissue. Our model also qualitatively reproduces the dependence of swirl size and swirl velocity on cell density reported in experiments and exhibits an undulation instability at the edge of the spreading tissue commonly observed in vivo. Finally, we study the dependence of colony spreading speed on important physical and biological parameters and derive simple scaling relations that show that coordination of motility forces leads to an improvement of the wound healing process for realistic tissue parameters.
One of the many remarkable properties of multicellular tissues is their ability to regenerate, even from severe damage, to a state very similar to their original form. This wound healing process is not only crucial for regaining basic tissue functionality, but also critical for restoring protection from infection, for example, by bacteria that can invade the organism via breaches in surface tissues. Depending on the severity of damage to the tissue, wound healing can involve several stages, including inflammation, angiogenesis, the regeneration of extracellular matrix and basement membranes, and reepithelialization (1⇓–3). During the latter process, the surrounding epithelium covers the wound by cell division and migration. Epithelial tissues are confluent arrangements of tightly adhesive cells in single or multiple layers that present the foremost barrier of the body against invasion.
In recent years, the epithelialization phase of wound healing has been studied in the laboratory, using a variety of in vitro models. Typically, epithelial cells are grown on a substrate to form a cohesive, monolayered sheet and a wound is created by scratching, by laser ablation, or by removing agarose blocks (4). It was discovered that the leading edge of the epithelial tissue often does not move uniformly when spreading over the substrate but exhibits long finger-like protrusions that move faster than the surrounding epithelial cells (4⇓–6). Furthermore, these fingers give rise to large-scale flow patterns within the tissue (7). Large-scale flows were also observed away from the tissue edge by Angelini et al. (8, 9). They observed that at low cell densities before forming mature epithelial sheets, cells in the bulk of the tissue exhibit spontaneous, large-scale, swirl-like flow patterns. These swirls have a complex dependence on cell density: Their typical size increases with increasing cell density, whereas their velocity decreases.
Further insight into the spontaneous motility of epithelial cells has been gained from traction force measurements of cohesive colonies (10, 11). In a seminal work, Trepat et al. (11) studied the spreading of large monolayered colonies of epithelial cells several millimeters in diameter and constructed 2D force maps of the tissue by averaging the forces exerted by cells at different distances from the edge. Surprisingly, they discovered that even cells in the middle of the epithelial sheet, many cell diameters away from the boundary, exert active motility forces. This finding is in sharp contrast to the standard picture of wound healing in which only the cells adjacent to the edge of the tissue exert motility forces due to their loss of contact inhibition and pull the rest of the sheet behind them. By integrating the forces from the tissue edge over the colony, Trepat et al. determined the tension in the sheet generated by the motility forces. Over 80% of this tension originated from cells more than 50 μm from the boundary of the sheet. Therefore, bulk motility forces are relevant from a biological perspective, as they constitute the major driving force for the spreading process of the epithelial sheet.
These experimental studies raise several important questions: How do cells in the center of the tissue know the direction of the edge and what is the mechanism for the orientation of motility forces in general? Can the mechanical properties of spreading epithelial sheets such as the high tensile stress, in combination with cohesive expansion and cell division, be understood from a simple model? How does the speed of spreading depend on fundamental properties of the system, including the magnitude of motility forces and their level of stochasticity, the expansion pressure of the cells within the sheet, and the friction with the underlying substrate? In what regime do spontaneous swirls occur and what leads to the formation of finger-like protrusions at the tissue edge? And to what extent does alignment of motility forces improve the process of wound healing?
A number of studies have addressed these questions. Gov (12) suggested a coupling between motility forces of neighboring cells as an explanation for the coordination of motility forces. In this picture, each cell has a planar cell polarity, which has a tendency to align with the polarity of its neighbors. Lee and Wolgemuth studied the same mechanism in a continuum model (13, 14). In a one-dimensional model, Puliafito et al. (15) studied mechanics and cell proliferation in spreading colonies. In their model, cells are described by springs with a preferred length with friction and active motility forces with the substrate. Cells grow when they are stretched beyond their preferred length and division is implemented as a function of cell size. In their model, cell motility is described by a Gaussian random noise with an imposed outward bias close to the tissue edges. More recently, Serra-Picamal et al. (16) investigated temporal stress patterns in spreading tissues, using a similar model and including sequential fronts of cytoskeletal reinforcement and fluidization. Bulk motility effects such as the spontaneous swirls observed by Angelini et al. (8, 9) were investigated using a particle-based model by Henkes et al. (17), based on earlier studies of collective motion (18, 19). Pressure-based spreading dynamics of monolayered tissues were studied in detail using particle-based models by Drasdo et al. and Galle et al. (20, 21).
In this work, we show that a very simple class of models that does not invoke cell signaling, planar cell polarity, or direct interactions between the polarization of neighboring cells can explain a large number of observations related to wound healing and the spreading of epithelial colonies.
We propose that the central property underlying the coordination of motility forces is the tendency of a single cell to align its motility force with its velocity. From a biological perspective this type of interaction can easily be imagined. Motility forces are generated by lamellipodia, thin protrusive sheets filled with actin filaments, that can be hidden beneath neighboring cells (22) and a lamellipodium exerting a force in a direction other than the cell’s velocity will be quickly realigned with the cell’s direction. In fact, it is well known that the orientation of lamellipodia of motile cells is highly sensitive to physical influences like gradients in substrate elasticity, substrate friction, or the location of the nucleus (23⇓–25). Therefore, we propose that the motility forces of cells in the tissue tend to align with the flow. Spontaneous flow patterns arise if a group of cells begins to move in a given direction, dragging other cells along, which in turn align their motility forces with the flow.
This mechanism can lead to an enhanced spreading speed of colonies and give rise to tension within them: Outwardly oriented motility forces contributing to the expansion of the colony are favored, as the tissue readily expands under tension and cell density is constantly replenished via cell division. On the other hand, flows resulting from motility forces pointing inward are resisted by increases in cell density. In this picture, cell division gives rise to a ratchet-like mechanism where expansion is followed by cell division that, in turn, prevents a reversal of flow direction. This leads to a long-range bias in the orientation of motility forces, which thereby contribute to the expansion of the colony.
Note that the mechanism we suggest is part of a general class of models used to study collective motility like flocking of birds, herding of wildebeest, or microorganism vortices (26⇓–28). Thus far, however, the effectiveness of this simple mechanism as the driving force of wound healing and spreading of cohesive epithelial tissues has not been studied.
To demonstrate the effectiveness of such a mechanism, we use a mechanical simulation of an epithelial sheet. The specifics of our simulation are described in ref. 29. In previous studies, this type of simulation was used to study the properties of 3D tissue aggregates like their growth, competition, rheology, surface tension, cell sorting, and the diffusion of cells within the tissue (29⇓–31). Here, to study monolayered epithelial tissues, we use a 2D version of this simulation. In our simulation model, individual cells are represented by two particles that interact via a repulsive force 
, giving the cell the tendency to expand in size. Here, 
is the unit vector pointing toward the other particle constituting the cell and r is the distance between the two particles, whereas B is the expansion parameter. If the cell size represented by the distance between these two particles surpasses a threshold 
, the cell divides at a rate 
. The details of this division process are also described in previous work (29); in brief, two new spheres are placed a small distance 
away from the spheres constituting the old cell. In addition to friction with the substrate and hydrodynamic dissipation between cells, intracellular dissipation between the particles constituting a cell is implemented as a resistance to abrupt changes in cell volume with a friction coefficient 
. As the rate of cell death is usually very low in the experiments we aim to explain, it is set to zero in our simulations.
Galilean invariance is broken by the presence of a friction force proportional to the cell velocity 
, which represents the interaction of cells with the substrate. (The main results of ref. 11 are independent of substrate elasticity as shown in their work. Therefore, we do not model the substrate explicitly in our simulation.) Volume exclusion and adhesion of cells are described via a radial force of the form 
, which is repulsive at short and attractive at intermediate distances and acts between particles constituting different cells. Here, once again r is the distance between the interacting particles and 
is the unit vector pointing from the particle on which the force is acting toward the interacting neighboring particle, whereas 
, 
are coefficients. The interacting particle feels a force of the same magnitude in the opposite direction. Note that the interaction radius of a particle is normalized to 1 and the force vanishes outside a range 
. Dissipative particle dynamics (DPD), which locally conserve momentum (32), are used to model friction and fluctuations within the tissue. In this framework, friction forces oppose the relative motion of particles within a range 
of each other with a friction coefficient 
. Details about DPD are readily available in the literature (32) and our implementation is described in detail in ref. 29.
In addition to friction forces with the substrate, cells in our simulation exert active motility forces against the substrate to propel them in a given direction. We use a simple stochastic model to describe the orientation of the motility force and the influence of the cell velocity on this orientation (Fig. 1): We assume that each cell switches between a motile and a nonmotile state. The cell switches from the nonmotile to the motile state at a constant rate 
. In the motile state, the cell exerts a motility force of a fixed magnitude m, in a randomly chosen direction. The rate of switching from the motile to the nonmotile state, on the other hand, depends on the orientation of the motility force relative to the cell velocity. In general, this dependence can be a complex function of the velocity 
and the motility force 
. As a minimalistic model, we choose just two rates 
and 
, depending on whether 
and 
are parallel or antiparallel, where typically 
. In other words, if the component of the motility force along the cell velocity is positive, the reorientation rate is given by 
, which is smaller than the reorientation rate when the component is negative, 
. Note that in the simulation, both particles constituting a cell have the same motility force and the reorientation rate is determined from their averaged velocity. (In our simulation, the effective cell velocity relaxes toward the instantaneous cell velocity with a relaxation time τ, averaging out fast velocity fluctuations on timescales much shorter than τ that can occur in DPD dynamics.) Schematically, the equation of motion of a particle in our simulation can be summarized as
where 
is the intracellular friction force between particles constituting the same cell (coefficient 
) and 
represents the friction forces between neighboring particles constituting different cells (coefficient 
). Furthermore, 
is a momentum-conserving noise force between particles.
Illustration of the motility–velocity alignment mechanism used in the simulations. Each cell can either be in a nonmotile state, in which it exerts no motility force, or in a motile state, in which it exerts a motility force of a fixed magnitude m. The rate at which cells go from the nonmotile into a motile state 
is constant and the direction of the motility force is chosen randomly when the cell transitions into the motile state. However, the stability of the motile state depends on the alignment of the motility force 
with the cell velocity 
. If the motility force has a positive component along the cell velocity 
, the transition rate into the nonmotile state is given by 
. Otherwise, the rate is given by 
. Hence, for 
there is a positive coupling between cell velocity and the orientation of the motility force, whereas for 
the motility force is oriented randomly.
To facilitate a systematic analysis of the parameters, we define a standard parameter set and explore parameter space by varying this set. Parameter values are presented relative to this standard set and are denoted by an asterisk; e.g., 
means 
. Further details and the parameters used can be found in Supporting Information and Table S1. All results are presented in simulation units.
We begin our analysis by studying the spreading dynamics and the mechanical properties of the tissue in situations that closely correspond to the experimental configuration in refs. 11 and 15. Namely, we simulate the growth of colonies in an unbounded geometry on a substrate, initialized either from a single cell or from a group of cells located in a small region of the substrate. Fig. 2 and Movie S1 show snapshots of a colony grown from an initial 50 cells, spreading over the substrate. The colony spreads cohesively and assumes the shape of a disk. Fig. 3A shows the cell density and the velocity field for this simulation. Simulations initialized from a single cell and from 500 cells are presented in Movies S2 and S3, respectively, and show similar qualitative behavior.
Snapshots of a growing colony at different time points, representative of the different regimes of spreading dynamics, with the simulation time and the number of cells indicated. The large tension in the tissue can lead to holes, which form and subsequently close again (Movie S1).
(A) Cell density (color coded) and cell velocity (vectors) in the spreading colony at 
. (B) Traction force map in the colony at 
. The vector field represents the local traction forces exerted on the substrate, whereas the color code represents the y component of this force. The traction forces are computed by calculating the vector sum of the friction forces and motility forces. Motility forces align with the velocity field and dominate over substrate friction forces at the time point presented. This leads to tension in the colony. Movie S6 shows the evolution of traction forces in the spreading colony with time.
To probe the mechanical state inside the tissue, we use a similar approach to that used in refs. 11 and 33. Using the traction forces exerted on the substrate, we reconstruct the stress field in the tissue. In Fig. 3B, the traction force map of the colony is shown at a time point corresponding to Fig. 3A. The traction forces felt by the substrate are computed as the vector sum of the motility forces and background friction forces, as both the motility forces and the background friction imposed in our simulation are balanced by the substrate in reality. From Fig. 3B, it is apparent that although the traction force field is highly disordered, velocity–motility coupling leads to global alignment of motility forces with the cell velocity. As a result, the y component of the traction force has predominantly negative values in the upper half of the colony and predominantly positive values in the lower half of the colony. This, in turn, gives rise to tension in the center of the colony.
In general, it is not possible to reconstruct all three components of the stress tensor in two dimensions from the traction force field without making detailed assumptions on the rheology of the tissue. However, the 1D stress map can be determined by averaging the traction forces in the colony in one direction. Fig. 4 shows the temporal evolution of this 1D pressure field in the tissue, as the colony expands. The stress determined in this manner is independent of tissue rheology and can be compared directly with the experimental results obtained by Trepat et al. (11). The corresponding profiles for different initializations are presented in Figs. S1 and S2. Note that by averaging over a finite segment through the tissue in the x direction as in ref. 11, rather than over the entire colony, the tension profile presented in Fig. 4 is enhanced and exhibits large regions of tension even for very large colony sizes, as shown in Fig. S3. This is because parts of the colony close to the edge that are under pressure do not contribute to the stress in the center. However, the integral of the traction forces over the entire y direction is guaranteed to vanish only if the whole colony in the x direction is taken into account.
One-dimensional pressure field in the simulation presented in Fig. 2 as a function of time. Positive values (red) represent pressure and negative values (blue) correspond to tension. The temporal evolution of the stress field is characterized by a phase during which the tissue is under tension due to radial motility forces followed by a final phase where the colony is under pressure.
The result presented in Fig. 4 is representative of the generic picture for the evolution of the stress field within the model tissue. It can be categorized as follows: After a short-lived transient phase where the colony is small and where motility forces are oriented in random directions, the velocity field aligns the motility forces radially and gives rise to tension throughout the colony. This phase is similar to the experimentally observed stress profile reported in ref. 11. As the colony continues to grow, cells in regions close to the edge collectively move at the maximum velocity possible from the motility force. This results in a small velocity gradient close to the edge of the colony and in turn leads to a buildup of pressure in peripheral regions of the tissue, as cells continue to divide. In the subsequent phase, as cell division dominates over expansion, it eventually leads to a pressure-limited, high cell density state in the bulk of the tissue.
The collective motility of cells and the tension it generates can lead to an inversion of the density profile, where cell density becomes lowest in the center of the tissue as recently reported experimentally (34). Very large tension can even lead to the formation of holes in the colony. Whether this is a physiological parameter regime can be told only by experiments. Holes appearing in some published material (15) suggest that this is possible. Movies S4 and S5 present simulations with an increased division rate 
. In these simulations, holes do not arise in the tissue and the tension profile in the colony is somewhat less pronounced (Fig. S4).
Note that the mechanism of velocity–motility coupling we suggest here does not always give rise to tension within spreading colonies. Rather, we find tension in a parameter regime, where growth pressure is small compared with motility forces. Nevertheless, even for parameters where the entire colony is under pressure, coordinated motility can enhance spreading speed if the cell division rate is low and the friction with the substrate is high.
To obtain a better understanding of colony expansion, we analyze the spreading process after perturbing important model parameters. In particular, we are interested in the role of the global coordination of motility forces. In Fig. 5, the square root of the area covered by the colony is plotted as a function of time for different parameters. Whereas small colonies expand exponentially in time, limited by cell division, for large colonies a constant speed of expansion arises from friction with the substrate, which balances the combination of motility and expansion forces in the tissue. For large colonies, cell division is no longer limiting for colony growth, but the maximum spreading speed is the main indicator of spreading efficiency.
Square root of the area covered by the colony as a function of time. Large colonies radially expand at a constant velocity, which leads to a quadratic dependence of the covered area on time. The parameters varied in the subplots with the rest of the parameters identical to the standard tissue, from the fastest to the slowest spreading speed, are given by (A) the degree of coupling between cell velocity and motility given by the ratios 
, 10.0, 1.0, and 0.1, with 
fixed; (B) substrate frictions 
, 1.0, 4.0, and 10.0; and (C) rates of division for cells surpassing the size threshold, 
, 1.0, 0.2, and 0.1.
Fig. 5A shows the effect of coupling between cell velocity and motility forces, as given by the ratio 
. Indeed, we find that a positive coupling given by 
significantly enhances the asymptotic spreading velocity. This can be understood as follows: In the absence of velocity–motility coupling, cellular motility forces point in random directions and contribute to the pressure in the tissue. Cells in the bulk can divide only by pushing layers of cells at the edge over the substrate. The balance of bulk pressure and the substrate friction force determines the spreading velocity. The presence of velocity–motility coupling alleviates this problem. Instead of being passively pushed over the substrate by expansion pressure and hindering growth of the colony, cells close to the tissue edge collectively move toward the edge and thereby create the space required for cell division and reducing pressure in the sheet. On the other hand for 
, velocity–motility coupling counteracts the expansion of the colony (Fig. 5A).
The effect of substrate friction on spreading is illustrated in Fig. 5B. It is apparent that friction with the substrate plays an important role in the speed of the wound healing process. For fixed motility forces, substrate friction determines the maximum spreading velocity. Finally, in Fig. 5C, the division rate 
for cells surpassing the size threshold is varied. Increasing this rate leads to faster spreading of small colonies, where expansion is limited by the rate of cell division. On the other hand, the asymptotic spreading velocity for large colonies in the motility-dominated regime remains unchanged.
As shown in our simulation, velocity–motility coupling can enhance the spreading velocity of monolayered tissues and give rise to tension within much of the expanding colonies. In this section, we address two questions related to these dynamics: First, in which parameter regime does coordination of motility forces provide an advantage for spreading tissues and is this the case for real monolayered epithelial sheets? Second, our simulation shows a transition from large-scale tension to pressure in the bulk of spreading colonies as the colony expands. At what size does this transition take place and at what colony radius do we expect this transition for biological tissues?
For the first of these questions, we study a simple one-dimensional model and compare the spreading speed of a colony spreading purely due to tissue expansion pressure with that of a colony spreading via optimally aligned motility forces. By comparing the respective spreading velocities for large colonies, we obtain a condition that gives the parameter regime in which motility-based spreading is advantageous compared with pressure-based expansion. We first consider purely pressure-based spreading: As a simplification, we assume incompressibility and a dependence of cell division on pressure that can be described by a first-order expansion (30, 35). Hence, the continuity equation is given by 
, where κ and 
are expansion parameters. Furthermore, we neglect effects due to the internal viscosity of the tissue, assuming that the dominant mode of dissipation is friction with the substrate, in which case the force balance condition is given by 
, where ξ is a friction coefficient. With the edges of the tissue located at 
, the boundary conditions impose a vanishing stress at the tissue edges 
. For symmetry reasons the velocity vanishes at the origin 
. The spreading velocity of the colony is given by the cell velocity at the tissue edge 
. For large colonies 
, the spreading velocity converges to 
.
On the other hand, in the case of colony expansion based on an optimal alignment of motility forces in the absence of expansion pressure, the maximum spreading speed is simply given by the maximum velocity of an individual cell. Note that consistent with this assumption, experimentally, the spreading speed of epithelial sheets is comparable to the velocity of single motile cells, 
μm/h (7).
For spreading based on motility alignment to be advantageous, 
must be large compared with the maximum pressure-based spreading velocity 
. For real tissues, many parameters can be inferred only indirectly. To estimate the friction coefficient ξ, we use the force balance condition 
. From the work of Trepat et al. (11), we know the typical magnitude of traction forces on the substrate, 
Pa. As a length scale, we use the typical size of such stress patches, 
μm. Together with the typical cell velocity 
, this yields a friction coefficient 
Pa h/μm2. Furthermore, we assume a typical cell division rate 
of one cell division per day based on ref. 11. From these estimates, we obtain the condition 
Pa on the expansion pressure of the tissue for motility-based spreading to be advantageous. There are no precise experimental values for this expansion pressure available. The analysis of Trepat et al. (11), however, did not detect a pressure contribution in the stress profile of the colony, suggesting that 
is significantly smaller than the condition derived above. Hence, unless the expansion pressure in epithelial sheets is on the order of several kilopascals, a mechanism for the global alignment of motility forces provides a distinct advantage for spreading and wound healing.
We now turn to the second question regarding the transition from bulk tension to pressure. One of the conclusions of our work is that the spreading regime observed in ref. 11, in which the colony is under tension, is a transient effect. Both for very small and for very large colonies, the tissue is under pressure. The size, below which the colony is under pressure, depends on initial conditions and on the details of the coupling mechanism. On the other hand, the colony size above which the expanding tissue is under pressure arises when exponential cell division overtakes the linear growth in colony diameter. Hence, the motility-based spreading velocity 
μm/h must be compared with the spreading speed of an exponentially growing colony of diameter d with a uniform cell division rate 
. Once again, assuming a division rate of about one cell division per day, 
, we expect the transition from global tension to pressure to take place around a colony diameter given by 
mm. We conclude from this analysis that even colonies several millimeters in size can be under tension, consistent with the observations of Trepat et al. (11). Applying this argument to the simulation presented in Fig. 2, we measure the effective cell division rate for small colonies, 
. The asymptotic spreading velocity can be determined from 
or can alternatively be measured from Fig. 6, 
. Together this yields a diameter 
(in simulation units) for the onset of pressure within the colony, consistent with the results of Fig. 4.
Evolution of the edge of spreading tissues at different time points after removal of a hard wall in the y direction that was in place during initialization. Both tissues have identical parameters except that the adhesion force is set to zero in the tissue on the right-hand side, 
. Whereas the adhesive tissue shows pronounced finger-like protrusions, the tissue without adhesion spreads much more uniformly, consistent with our analysis (Movies S7 and S8). The parameters deviating from the standard tissue are given by 
, 
, 
, 
, and 
.
The coupling between cell velocity and motility forces in our simulations can give rise to finger-like protrusions of the tissue edge. To mimic wound healing experiments, the simulation was initialized by growing the tissue in a long narrow compartment with a periodic boundary condition in one direction and a hard wall constraining it in the other direction. After setting the cell velocities to zero and converting the cells into the nonmotile state, the simulation was started by replacing the hard walls with periodic boundary conditions and increasing the compartment size, effectively creating a wound that can be closed by the tissue from both sides. An example of a simulation is shown in Fig. 6 for the model with (Fig. 6A) and without (Fig. 6B) adhesion. For the adhesive case, the straight boundary of the adhesive tissue is clearly unstable, leading to pronounced finger-like protrusions. A linear stability analysis of a continuum version of our model in the absence of cell division shows that a planar tissue interface is unstable in the presence of velocity–motility coupling and demonstrates that the absence of adhesion (which appears as a decrease of surface tension) tends to stabilizes the interface.
Experiments often display large-scale velocity correlations in the bulk of spreading epithelial tissues (8, 9). Our model is able to capture these correlations as shown in Fig. 7A where we present snapshots of the velocity field computed using our model for three different densities. These snapshots clearly show large-scale collective motion in the form of swirls. To quantify this motion, we compute the velocity correlation function. Before calculating this correlation function, we subtract the mean velocity vector across the field of view, which is always orders of magnitude smaller than the mean squared velocity, and calculate the autocorrelation function of the difference 
. Then the velocity correlation function at a distance r is given by
where brackets denote an average over all directions, making the correlation function depend solely on the scalar distance r. In Fig. 7B we plot 
as a function of r for four different densities, including the ones used in Fig. 7A. On short length scales, the velocity correlation function decays exponentially with a single characteristic length scale that increases for increasing cell density (Fig. 7B, Inset). The length scales extracted from the fits are 1.2, 1.3, 1.45, and 1.8, for, respectively, densities 1, 1.25, 1.5, and 2 of the reference tissue density. These should be compared with a characteristic cell length of 0.7 for the reference tissue. Thus, our length scales are comparable to the experiments in Angelini et al. (8), which reported decay lengths of ∼3 cell lengths.
(A) Snapshots of the velocity field for low, intermediate, and high confluent densities. The typical size of the swirl patterns increases with cell density. Movie S9 shows this effect. (B) Velocity correlation function for different densities, at a ratio 
. We set the density at confluence to 1. The correlation length of the velocity correlation function increases with cell density. Left Inset shows a zoom on the region of the minima. Right Inset shows the decay of the velocity correlation function on a logarithmic scale. (C) Average cell speed as a function of density. This velocity decreases significantly with increasing density. (D) Phase diagram of aligned vs. swirling motion for different degrees of velocity–motility coupling 
and for different cell densities. The line indicates the critical cell density for the phase transition.
As in the experimental results of ref. 8, 
does not decay monotonically but shows a small negative minimum (Fig. 7B, Inset). This negative value indicates the presence of swirls due to the antiparallel velocity vectors on opposite sides of the swirl pattern and its location determines the length scale of the swirl (8). This length scale increases for increasing density, consistent with the experiments. Moreover, the average speed of cells (plotted in Fig. 7C) decreases consistently for increasing densities, matching the experimental results of ref. 9.
Above a critical density, a phase transition occurs and unconstrained cells start to move together at a uniform velocity, in agreement with the results obtained by Szabo et al. (19), or alternatively organize in a solid-like phase as in ref. 17. To gain a better understanding of this phase transition, we have investigated its dependence on cell density and the degree of coupling between velocity and motility given by 
. A tissue is considered “uniformly moving” if the order parameter for the alignment of the speed is at least 0.5. As shown in our phase diagram (Fig. 7D), decreasing the degree of coupling 
or decreasing cell density favors disorganized swirling motion. For higher degrees of coupling 
, the critical density for the phase transition to uniform motion decreases.
In this study, we show that a broad range of experimental effects observed for monolayered epithelial sheets can be understood with a simple model of single-cell behavior. We study this model using detailed mechanical simulations. A key element in our model is the alignment of the motility apparatus of the cell and its velocity. Similar alignment is known to occur in other contexts, including durotaxis or haptotaxis (23⇓–25).
The model is able to qualitatively reproduce the mechanical properties of spreading colonies observed in ref. 11. In particular, our simulations display the experimentally observed buildup of tension throughout spreading colonies (Fig. 4). We show that this tension is caused by a long-range alignment of motility forces. In our simulation, this occurs in a dynamical manner due to the behavior of individual cells and we do not need to impose special properties of cells close to the edge or particular boundary conditions as in other studies (13⇓–15). Our model predicts that the observed tension is a transient phenomenon and that, eventually, pressure from cell growth and from randomly oriented motility forces dominates over tension generated by aligned motility forces. Note, however, that in real tissues this pressure may be very small compared with the tensile stresses and occur only for very large colonies, as cells in the center of the tissue undergo contact inhibition and lose their spontaneous motility.
In addition, the simulations reveal the emergence of finger-like protrusions at the leading edge of an advancing tissue sheet (Fig. 3), similar to the ones observed in recent experiments. Finally, our 2D, particle-based simulations are able to mimic the large-scale collective motion observed in experiments (Fig. 7). This is accomplished without invoking complex cell–cell interactions for the global coordination of motility forces, as suggested in other approaches (12⇓–14). In particular, our model reproduces the characteristic increase in swirl size and decrease in swirl velocity with increasing cell density observed experimentally.
There are important biological questions related to the wound healing process of monolayered epithelial tissues that are beyond the scope of the model introduced here. Whereas motile cells in the bulk of the tissue generate most of the stress that drives the spreading process, cells in mature epithelial tissues do not display motility forces. The signaling cascade that leads to this switch in phenotype for bulk cells is not known and may be more complex than the current biological picture of contact inhibition, based purely on cell–cell contact, for which several molecular models have been proposed (36, 37). Mechanical effects like the buildup of tension in the tissue may play an important role for both division and motility of cells in the tissue (15).
Possible biological mechanisms for the influence of cell velocity on the alignment of motility forces include the mechanical interaction with the substrate, which can destabilize lamellipodia that are not aligned with the cell velocity, or an influence of the position of the nucleus relative to the lamellipodium. Although there are some experimental indications that motility forces of cells tend to align with cell velocity, there are currently no experiments establishing this effect directly. A test of this idea would be to exert an external force on a single cell or cells in a tissue to study its effect on the orientation and the magnitude of motility forces. One possibility to implement such an external force experimentally would be the utilization of a fluid flow to exert a viscous drag on cells (38).
A more indirect verification of our hypothesis could be achieved via statistical analysis of spatial and temporal correlations between velocities and traction forces in spreading epithelial sheets. Ideally, agarose blocks would be removed to stimulate the spreading of an initially resting, confluent, monolayered epithelium. Such an experiment could buttress our suggested mechanism by showing that the velocity field precedes the alignment of motility forces or falsify it by showing the reverse. Some promising results in this direction were obtained in a very recent work by Serra-Picamal et al. (16), which analyzes the evolution of both velocity and traction force fields in time in spreading epithelial colonies and appears to show that traction forces temporally trail the velocity field in the tissue in a wave-like pattern.
We acknowledge useful discussions with Eshel Ben-Jacob and Assaf Zaritsky. This work was supported by National Science Foundation (NSF) Grant DMS-1068869 and by the NSF Center for Theoretical Biological Physics (Grant NSF PHY-0822283).
This contribution is part of the special series of Inaugural Articles by members of the National Academy of Sciences elected in 2011.
Author contributions: M.B., J.E., E.H., W.-J.R., and H.L. performed research; M.B., J.E., and E.H. analyzed data; and M.B., W.-J.R., and H.L. wrote the paper.
The authors declare no conflict of interest.
This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1219937110/-/DCSupplemental.
