Filopodial-Tension Model of Convergent-Extension of Tissues

In convergent-extension (CE), a planar-polarized epithelial tissue elongates (extends) in-plane in one direction while shortening (converging) in the perpendicular in-plane direction, with the cells both elongating and intercalating along the converging axis. CE occurs during the development of most multicellular organisms. Current CE models assume cell or tissue asymmetry, but neglect the preferential filopodial activity along the convergent axis observed in many tissues. We propose a cell-based CE model based on asymmetric filopodial tension forces between cells and investigate how cell-level filopodial interactions drive tissue-level CE. The final tissue geometry depends on the balance between external rounding forces and cell-intercalation traction. Filopodial-tension CE is robust to relatively high levels of planar cell polarity misalignment and to the presence of non-active cells. Addition of a simple mechanical feedback between cells fully rescues and even improves CE of tissues with high levels of polarity misalignments. Our model extends easily to three dimensions, with either one converging and two extending axes, or two converging and one extending axes, producing distinct tissue morphologies, as observed in vivo.

The observed asymmetry of filopodial protrusion led us to propose a filopodial-tension mechanism for CE based on anisotropic filopodial pulling forces between cells. We explicitly model the number of cell-cell connections, their range, angular distribution, strength, and frequency of formation and breakage. We define an appropriate set of metrics to quantify both the effects of model parameters and planar-polarization defects (such as misalignments and the passive cells) on the dynamics of tissue-level CE. Since our filopodial-tension model extends naturally to three dimensional tissues, we discuss the two types of 3D CE and their corresponding tissue morphologies.

Anisotropic Filopodial-Tension Model of Convergent Extension
Experiments show that long filopodia continuously form and retract during CE in epithelial sheets and that these filopodia preferentially form in-plane along angles near the axis of tissue contraction. Each model cells therefore extends and retracts filopodia (which we represent using the model concept of a link) distributed within a range of angles around the directions perpendicular to the cell's planar-polarity axis. To simulate the observed binding of filopodial tips to other cells and the roughly length-independent pulling forces which retracting filopodia generate, in our model, an extending link binds to the cell it contacts, then generates a constant (length independent) tension force between the cells it connects [20,25,29]. We then test whether this tension force is sufficient to explain observed local cell intercalation and global tissue CE.
In the filopodial-tension model (Fig 2, S1 Movie) cells form and eliminate links representing filopodia with a defined set of neighboring cells (terms in boldface identify model objects). Each cell carries a polarization vector (perpendicular to its planar-polarity axis) (Fig 2, red arrow) that defines its preferred direction of filopodial protrusion (Fig 2, blue horizontal line). We simplify the model by having the links connect the centers-of-mass of cells rather than connecting the actin cortex of one cell to the actin cortex of the contacted cell, as do real filopodia. Because filopodia typically form in a pair of growth cones roughly along the convergence axis and with a typical maximal length, we allow a cell to form links within a range of angles ±ϑ max around this axis on either side of the cell with a maximum length of approximately r max . Specifically, a cell can form a link only with those cells whose centers-of-mass lie within a distance r max from its center of mass and within an angle ±ϑ max of its polarization axis (Fig 2, blue horizontal line). A cell can have at most n max links to other cells at any time (including links formed and received) and only one link is allowed between any pair of cells. The actual number of links a cell forms may be less than n max . Each link between a pair of cells exerts a tension force of magnitude λ force along the line connecting the cells' centers-of-mass. To model the finite lifetimes of filopodia, we define a relaxation time, t interval , after which we remove the links of all cells and create new ones. In a simulation in which the links form and then persist indefinitely, the cells only move a few microns (lattice sites) from their original locations and the tissue does not converge or extend.
We implement the filopodial tension model using the Cellular Potts model (CPM, also known as the Glazier-Graner Hogeweg model, GGH), where each cell is represented as a collection of lattice sites with the same cell index. An effective-energy cost function, H, specifies the In active convergent-extension, the cells in the tissue generate deforming forces due to anisotropic adhesion or pulling forces between cells (red arrows), while in passive convergent-extension, the surrounding environment deforms the tissue (blue arrows). Cell intercalation occurs in types of CE, but the axis of cell elongation is typically perpendicular to the axis of elongation in active CE and parallel in passive CE.  Filopodial Model of CE cell's properties (see supplemental material). The tension force along a link between a pair of cells is independent of its length and acts along the vector between their centers-of-mass. In the GGH/CPM formalism, the tension has the form: where the sum is over all pairs of linked cells, λ force is the strength of the pulling force between cells σ and σ', l σ,σ' is the current distance between the cells, and the term H 0 aggregates all the other GGH/CPM cost function terms. The GGH/CPM simulations evolve stochastically from random lattice-site updates subjected to the effective-energy cost function, H. The time unit is the Monte Carlo Step (MCS), defined as the rate of lattice-site updates (see supplemental material for more details on the GGH/CPM formalism).
The filopodial-tension model has five intensive parameters (λ force , t interval , r max , n max , ϑ max ) and one extensive parameter (N, the number of cells), making a complete sensitivity analysis computationally costly. We therefore fixed all parameters to reference values that are within the ranges observed in vivo and produced biological plausible convergent-extension (Table 1), then studied the effects of varying each intensive parameter one-at-a-time. The biological parameters proposed by the model can be directly measured experimentally, but since the concept of a filopodial-based CE is new and applies more readily to CE of deep tissues, which are not as easily visualized as epithelial sheets, appropriate experimentally-derived values are harder to find. The most studied cases are chicken limb-bud mesenchymal intercalation [30] (t interval = 2.2 hours; r max = 3 cell diameters; n max = 11; ϑ max = 45°), Xenopus gastrulation and notochord formation [31][32][33] (t interval = 2.0-2.7 min; r max = 1.5 cell diameters; n max = 8-9; ϑ max = 60°), and Xenopus Keller explants [23,34] (t interval = 0.5-1.0 hour; r max = 1.5 cell diameters; n max = 8-9; ϑ max = 30°).

Metrics
All simulations start with a mass of identical cells uniformly distributed inside a rough circle. Each cell has the same planar-polarization vector (V). To quantify the degree of tissue deformation we calculate the inverse aspect ratio between the length of the minor (L -) and major (L + ) axes of the tissue (Fig 3B, green line). Initially the aspect ratio is close to 1 and decreases in time to a final value κ (Fig 3B, dashed red line) that depends on the filopodial tension parameters (λ force , t interval , r max , n max , ϑ max ), the number of cells in the tissue (N) and the surface tension of the tissue γ (defined below).
The final inverse aspect ratio quantifies the maximum elongation of the tissue, but does not convey how fast the tissue elongates. To quantify the elongation rate, we define the elongation time (τ) the time an initially isotropic tissue takes for its major axis (L + ) to double the length of its minor axis (L -), which is equivalent to the time when the inverse aspect ratio (L -/L + ) first decreases to 0.5 ( Fig 3B, dashed blue lines). We consider CE to fail if L -/L + never reaches 0.5.
Since both the filopodial-tension model and the GGH/CPM are stochastic, we average the value of the elongation time (τ) over 10 simulation replicas. Because the tissue inverse elongation ratio converges to the same value independent of the simulation seed or initial conditions (S1 Fig), unless specified otherwise, we calculate the final inverse aspect ratio κ for a single simulation replica, with the standard deviation indicating the fluctuations in κ around its final value for that replica.

Surface Tension vs. Filopodia Tension
Successful CE depends on the ability of intercalating cells to generate forces stronger than the internal and external forces that oppose tissue deformation. Here, the opposing forces come from the superficial tension (γ) between the cells and the external medium, defined as [35]: where J c,c is the contact energy between cells and J c,M is the contact energy between cells and medium (see supplemental material). When the filopodial tension is weak compared to the surface tension (λ force < 2γ), cells do not intercalate and CE fails. For larger filopodial tensions, the elongation time (τ) decreases as a power of λ force (τ / λ force -1.25±0.03 ) ( Fig 4A, red line). The final inverse aspect ratio (κ) decreases monotonically with increasing λ force (Fig 4B). Increasing γ shifts the κ vs. λ force curve to the right and decreasing γ shifts the κ vs. λ force curve to the left (Fig 4B, inset). Normalizing the filopodial tension by the surface tension (λ force /γ) collapses the κ vs. λ force curves (Fig 4B), showing the linear relationship between λ force and γ. The surface tension (γ), however, has little effect on the elongation time (τ), which depends on λ force , but is relatively insensitive to surface tension ( Fig 4A). The κ vs. λ force /γ curve is sigmoidal on a log-log scale (Fig 4B), because the shape of the tissue changes little for weak filopodial tensions and because the total number of   cells limits κ for strong filopodial tensions (see Fig S3A). At the inflection point of κ vs. λ force /γ, the tensions of the links (λ force ) balances the external surface-tension forces that oppose tissue elongation (λ force /γ~6). Near this inflection point κ varies as an approximate power law of λ force (κ / λ force -1.51±0.08 ).

Parameter Sensitivity
Next we studied how the remaining filopodial tension parameters affect CE, specifically, the mean lifetime of the filopodia, modeled as the time interval between link formation and breakage (t interval ); the maximum length of the filopodia, modeled as the maximum distance of interaction between the cells' centers-of-mass (r max ); the maximum number of filopodial interactions per cell (n max ); and the maximum angle between the filopodial direction and the cells' convergence axis (ϑ max ). Fig 5A shows that, for the reference parameter values ( Table 1) the lifetime of filopodia, t interval , has no effect on τ or κ for t interval ≲ 200 MCS. For the reference parameter values, 200 MCS corresponds to the typical time the cells require to rearrange their positions in response to a given set of neighbors interactions. Increasing filopodial lifetimes above 200 MCS slows cell intercalation (increasing the elongation time) and increases the tissue's final inverse aspect ratio (corresponding to less deformation). The maximum range (r max ) of filopodia interaction has different effects on the final inverse aspect ratio (κ) and elongation time (τ). For r max < 2 cell diameters, κ decreases as a power law in r max (κ / r max -3.5±0.2 ), then saturates for r max ! 2 cell diameters, while the elongation time (τ) decreases monotonically with increasing r max ( Fig 5B). The same effect is seen with respect to the maximum number of links (n max ): κ decreases as a power law in n max for n max < 4 (κ / n max -1.5±0.03 ) and saturates for n max > 4 (this saturation makes sense since the cell typically has 4 neighbors within the range of its filopodia for r max = 2 and ϑ max = 45°) while τ decreases monotonically with increasing n max ( Fig 5C). Thus r max and n max have affect the rate of cell intercalation more than the final inverse aspect ratio, while the tissue's surface tension affects only the final inverse aspect ratio and not the rate of cell intercalation (Fig 4). Both κ and τ are concave with respect to the maximum angle of filopodial protrusion (ϑ max ), since for small ϑ max the number of cell center-of-mass within the cones defined by ϑ max and r max is very small, while for ϑ max = 90°the forces on the cell are symmetric since it extends filopodia uniformly in all directions. In both limits CE fails (Fig 5D). Since the net intercalation force is the difference between the tension forces parallel and perpendicular to the Filopodial range (r max ): κ decreases with increasing r max for r max < 2 cell diameters and is constant for r max > 2 cell diameters. τ decreases monotonically with increasing r max . CE fails for r max < 1.5 cell diameters. (C) Number of filopodial interactions (n max ): κ and τ decrease monotonically with increasing n max , however κ decreases more slowly for n max > 4. (D) Angular range of filopodia (W max ): for n max = 3 (reference value), τ (blue squares) and κ (red dots) decrease monotonically with increasing W max for small W and increase monotonically with increasing W max for large W, with minima at W max = 30°and W max = 40°, respectively. CE fails for W max > 70°. For n max = 7, τ vs. W max (blue line) and κ vs. W max (red line) are also concave curves with minima at W max = 40°and W max = 50°, respectively. CE fails for W max > 80°. convergence axis (roughly R y max 0 ðcosðyÞ À sinðyÞÞdy ), we might expect the force to be greatest (and thus κ and τ to be smallest) when ϑ max = 45°and for their values to increase symmetrically away from ϑ max = 45°. The curves, however, have different minima and are not symmetric: the smallest final inverse aspect ratio (κ) is around ϑ max = 40° (Fig 5D, red dots) and the smallest elongation time (τ) is around ϑ max = 30° (Fig 5D, blue squares).
This asymmetry is caused by the limited number of neighbors with which a cell can form a link. Both the maximum number of links per cell (n max ) and the number of cells within the link interaction range (r max ) can limit the actual number of links a cell forms. If the maximum number of links per cell is lower than the number of cell neighbors within a cone of range r max (e.g. n max = 3) and angle ϑ < ϑ max , increasing ϑ max leads to more links with cells at larger ϑ and thus reduces the net tension force applied along the direction of the convergence axis. In effect, large ϑ max causes the cell to waste its limited number of filopodia. For large n max , links form to all cells within the cone of range r max and small θ regardless of the value of ϑ max .Thus, for large n max (e.g. n max = 7), the κ and τ curves are roughly symmetrical around their minima at ϑ max~4 5°(blue and red lines in Fig 5D).

Contact-Mediated Pulling
The filopodial tension model assumes that cells can extend filopodia, contact and pull other cells that lie within a given distance, even if they do not touch each other before filopodial extension. An example would be the formation of adhesion junctions between cells which coupled to a contractile stress fiber in both cells. To model these cases, we defined a contact-mediated cell tension model, which is identical to the filopodial tension model except that the maximum link length r max in the filopodial tension model is replaced with the condition that cells must be in touch before pulling on each other (Fig 6A).
The qualitative results for the contact-mediated cell tension model do not differ much from the filopodial tension model. The κ x λ force curve is sigmoidal on a log scale, τ decreases with a power law (κ / λ force -1.18±0.06 ) and CE fails for λ force < 20 (Fig 6B). The dependence of κ on the number of filopodial interactions (n max ) is still a power law (κ / n max -1.5±0.03 ) and saturates when n max = 4. The elongation time (τ), however, does not keep decreasing as it does for the filopodial tension model, but also saturates around n max = 4 links (Fig 6C), as few cells have more than 4 neighbors with centers near the convergence plane. The (κ, τ) x ϑ max curves have minima at ϑ max = 40°and ϑ max = 35°, respectively, but are less skewed than in the filopodial tension model (compare Figs 6D and 5D). CE fails for ϑ max < 10°and ϑ max > 70°.

Polarization Misalignment
Convergent-extension requires cells to have consistent planar polarity throughout an extensive region of tissue. This correlated orientation might result from a long-range bias from a morphogen gradient, cellular or intercellular differences in protein expression [36], or from a boundary-relay mechanism [37,38]. In our previous simulations we assumed that all cells had perfectly aligned polarization vectors (Fig 2, red arrows), i.e., they all pointed in the same direction with the same magnitude, and they maintained their internal orientation throughout the simulation. To study the effect of polarization misalignment on CE we added a zero-mean Gaussian distributed displacement angle to the cells' polarization vectors and varied the standard deviation of the distribution (σ) while keeping the mean direction (here, the vertical axis) constant. Since the final elongation ratio is sensitive to the distribution of polarization vectors, the values of κ were averaged over 5 simulations.
The filopodial tension model tolerates small polarization misalignments, with a tissue with a displacement angle of σ = 10°reaching the same final inverse aspect ratio as in the perfectly aligned case with little decrease in elongation rate (an 11% increase in τ). The tissue remained aligned with the mean direction of cell polarization (the vertical axis) for small misalignments (σ < 40°, Fig 7B), but bent at around σ = 50° (Fig 7C). For polarization misalignments with σ > 60°, CE fails and the tissue breaks its symmetry, acquiring more complex shapes such as the caltrop (see Fig 7D). Both metrics are exponential functions of the variance σ 2 (Fig 7A).

Mechanical Feedback Rescue of CE
So far we have assumed that the polarization vector of the cells remains constant throughout the entire process. In reality, however, cells are constantly communicating with their neighbors either through signaling or through mechanical interactions. During CE cells establish and maintain their polarity through the planar cell polarity (PCP) pathway, however there is Cells only pulls neighbors (here, 3) that share a common surface area (shown in red) and that lie inside a maximum angle with respect to the convergence plane (here the horizontal axis). (B) Dependence of τ and κ with λ force is qualitatively the same as before (Fig 4). (C) Dependence with n max is reversed, with the speed of intercalation (τ -1 ) saturating after n max = 3 and κ still decreasing. (D) The (κ τ,) x W max curves are more symmetric, but the tissue still elongates more and faster at lower angles. Filopodial Model of CE growing evidence that mechanical feedback may also play a role in the maintenance of global tissue polarity during development [39][40][41][42].
The presence of a mechanical feedback mechanism may rescue CE in tissues with high polarization alignment defects. In order to investigate this we developed a simple model of mechanical feedback and applied it to the misalignment polarization cases that have been described in the last section (Fig 7). The feedback model assumes that the pulling forces on a cell due to filopodial interactions affect its polarization vector. We implement a simple phenomenological version of such an interaction by calculating the line of tension from the sum of all filopodial interactions of the cell with its neighbors (Fig 8A). From this line of tension we extract an orthogonal vector T which is averaged with the previous cell polarization vector V in the following way: where V t+Δt is the polarization vector of a cell at time t+Δt, V t Ã is the normalized polarization of the same cell at time t, T t Ã is the normalized tension vector of the cell at time t, and w is a feedback weighting factor ranging from 0 (no feedback) to 1 (no memory) (Fig 8A). This iterative processes repeated at discrete time intervals set equal to filopodia lifetime (Δt = t interval ). For simplicity, we do not distinguish between the pulling forces generated by the cell from the pulling forces that their neighbors exert on it. The normalized tension vector of a cell is calculated from the vector that maximizes the sum of all projections of the normalized lines of force from all the cells' neighbors: where the sum is over all the cell's neighbors that pull on it, ϕ i is the angle of the line of force between the cell and the pulling neighbor i, and ϕ T is the angle that defines the tension vector T Ã = (cos ϕ T , sin ϕ T ).
For tissues with a high starting level of polarization misalignment (σ ! 40°), addition of this mechanical feedback mechanism usually leads to a lower final elongation ratio κ, as long as the feedback factor w is below 0.1. For these cases, tissue elongation times (τ) decrease with higher feedback levels (Fig 8C, green and red lines), while k usually decreases with lower feedback levels. For tissues with a low starting level of polarization misalignment (σ 20°), addition of this mechanism leads to lower final elongation ratios only for small levels of feedback (w 0.001) (Fig 8B, blue and black lines), while the time of tissue elongation (τ) remains relatively unchanged with respect to the case with no feedback (Fig 8C, blue and black lines). In all simulations where the addition of mechanical feedback rescues CE, the cells established a global polarization axis emergently. We chose to implement the feedback update iteratively, which leads to fast destabilization of the tissue for high levels of feedback, as is expected for a case with no memory. Although a continuous model would be relatively more robust, we expect the same destabilization effect when the weighting factor w approaches 1.

CE in the Presence of Non-active Cells
Next we varied the number of intercalating cells in the tissue to check if there is a minimum number of active cells needed to drive CE and how this change tissue dynamics. We defined two types of cells without filopodia: passive cells, which lack filopodia but can be pulled by the filopodia of other cells; and non-responsive, or refractory cells, which cannot be pulled by the filopodia of other cells. The former would correspond to cells whose surface adhesion molecules were compatible with those of the cells extending filopodia and the latter to cells with incompatible adhesion molecules. The parameters for cells which produced filopodia were the same as in Table 1. Since the final elongation ratio is sensitive to the distribution of active/nonactive cells, the values of κ were averaged over 5 simulations.
For tissues with a mixture of active and passive cells, both κ and τ decrease monotonically with the percentage of active cells in the tissue (Fig 9, red dots). However, even a fraction of Filopodial Model of CE active cells can drive CE. For 40% or more active cells (S2 Movie), the tissue deforms almost as much as a tissue composed entirely of active cells (Fig 9B, red dots), though the elongation time increases with the percentage of passive cells up to twice that for a tissue of all active cells (Fig 9A, red dots). For higher fractions of passive cells the final inverse aspect ratio increases significantly with the fraction of passive cells (Fig 9B). E.g., for 90% passive and 10% active cells (Fig 9C and S3 Movie), the tissue's final inverse aspect ratio never drops below 0.3 ( Fig  9B) and the elongation time τ is more than ten times that for a tissue of all active cells (Fig 9A). In all simulations, the active cells migrate towards the midline of the elongating tissue, leaving the passive cells at the lateral margins (Fig 9D).
The presence of relatively high CE despite the presence of only 10% of active cells can be explained by the relative tissue area covered by the filopodia. Every active cell can pull on neighbors that lie up to a distance r max from its center of mass and within an angle ϑ max on each side of its convergence plane (see Fig 2), thus covering an area of 2ϑ max r max 2 . For the reference parameters (ϑ max = π/2 and r max = 2 cd, see Table 1) and a population of 10% of active cells (N/10, where N is the number of cells in the tissue) this amounts to~1.2N(cd) 2 , which more than covers the whole area of the tissue (N(cd) 2 ). Refractory cells have a stronger effect on CE than passive cells. CE fails when the percentage of refractory cells is above 60% (Fig 9B, blue squares), while with passive cells it only fails for percentages higher than 90% (Fig 9B, red squares). For higher fractions of active cells, the two populations sort out, with the active cells extending normally and the refractory cells displaced to both sides of the elongating tissue (Fig 9F). Surface tension between the cells and the surrounding medium causes the refractory cells to form droplet-like clusters which bend the extending active-cell tissue into a wavy bar (Fig 9F and S4 Movie).

3D Versions
The 2D filopodial tension model is a reasonable description of cells within epithelial sheets, where cell movement is confined to a plane. However, in many situations cell intercalation occurs in 3D. That is the case in radial intercalation during epiboly of the developing Xenopus Laevis embryo, where cells in a multilayered epithelium intercalate and converge perpendicular to the plane of the sheet [43]. The filopodial tension model can be easily extended to three dimensions, but due to the extra degree of freedom, it breaks in two versions, depending on which axis is rotated: In equatorial or extensional intercalation, obtained by rotating the 2D model around the polarization vector (the red arrow in Fig 2), the cells pull on all neighbors that lie in a convergence plane (Fig 10A). At the tissue level, equatorial intercalation results in the convergence of the tissue along the two directions perpendicular to the polarization vector and its extension along the polarization vector ( Fig 10A' and 10A'').
In bipolar or convergent intercalation, obtained by rotating the 2D model around the convergence line (the blue line in Fig 2), the cells pull on all neighbors that lie along a convergence axis (Fig 10B). At the tissue level, the bipolar intercalation results in the convergence of the tissue along the axis of convergence and its expansion in the other two directions (Fig 10B' and  10B").
Beginning with a spherical tissue with all the cells polarized in the same vertical direction, the 3D equatorial model produces a tissue resembling a prolate spheroid (cigar shaped, Fig  10A", S5 Movie), while the bipolar model produces a tissue resembling an oblate spheroid (lentil shaped, Fig 10B").
The bipolar model has more biological correspondence than the equatorial model: cells with unipolar or bipolar protrusive activity are much more common during development than cells with equatorial protrusive activity, and the resulting tissue shape from the 3D bipolar model corresponds to the thinning and expansion associated with radial intercalation.
For both versions of the 3D model, the dependence of the parameters κ and τ with λ force , r max , n max and t interval are qualitatively the same as in the 2D model. The results only differ qualitatively with respect to ϑ max . For the same values of r max and n max , the 3D convergence model is slightly less skewed than the 2D version, with the best value for κ around ϑ max = 45°a nd the best value for τ around ϑ max = 35° (Fig 11B). The 3D extension model, however, presents a more drastic change in the (κ and τ) vs. ϑ max curve when compared to the 2D. While the  Fig 5D), which in turn has a range of optimal values slightly lower than in the 3D convergence model (B).
The reasons for the asymmetry is that the final shape of the tissue in the 3D extension model-a two cell diameter tube orthogonal to the convergence plane-is more sensitive to perturbations than the lentil shape tissue obtained by the 3D bipolar model. Two cells pulling each other along a convergence axis leads them to be aligned in a plane perpendicular to the direction of the pulling. This plane fully coincides with the 3D bipolar model extension plane (Fig 10B"), but only partially with the extension plane of the 3D extension model (Fig 10A"). In the simulations results shown in Fig 11A, the value of ϑ max = 30°represents the optimal maximum angle value where the pulling forces in the 3D extension model are still able to align the tissue without destabilizing it.

Discussion
Here we developed a new 2D model of for active cell intercalation. The model drastically differs from previous existing models by its explicit use of pulling forces between cells rather than anisotropy on adhesion energies or surface tensions on the cell surface. A recent experimental work by Pfister et al. [28] supports our force-driven model. The use of forces makes the model easily adaptable to three dimensions, but the extra degree of freedom gives two ways in which this can be achieved: either by rotating the model around the polarization vector or around the convergence line (see Fig 10A and 10B).
The model was implemented in CompuCell3D using the Cellular Potts formalism. The core of the model, however, is independent of the mathematical formalism and can be easily implemented on other types of agent-based formalisms such as cell-center or vertex models, as long as they provide a volume exclusion mechanism. Although some quantitative results might differ, we expect the same qualitative results. The advantages of implementing and simulating our model using the Cellular Potts formalism include the ease to manipulate and study the effects of the surface tension on the tissue dynamics (Fig 4) and the addition of the model as an integrated part of the CompuCell3D simulation package [44] that allows for it to be immediately reusable by others.
Model validation can be done at two different levels quantitative and semi-quantitative: when both microscopic and macroscopic measurements are available for a specific tissue, using model parameters that agree with those measured for the cells in the tissue should result in tissue-level model output (here, the rate of convergence and final aspect ratio) agreeing quantitatively with that of the experimental tissue. This agreement should persist under different experimental conditions. This form of validation shows that the hypothesized mechanisms included in the model are sufficient to reproduce the experiment quantitatively. Note that a model can only show the sufficiency of modeled mechanisms, not their necessity, since a different set of mechanistic hypotheses might yield the same results. In our case, since we are not modelling a specific tissue, model validation can only be semi-quantitative. We show that changes in the properties of the modeled cells (such as average number of filopodia, length and angular distributions) and tissue properties (such ratio of active and non-active cells) predict relative changes in the rate of CE and final tissue aspect ratios in real tissues that agree with those in experiments when the corresponding parameter and conditions are similarly modified. In this case, agreement demonstrates the plausibility of the hypothesized mechanisms, but detailed quantitative validation requires additional experimental measurements. We hope that our semi-quantitative analysis of the filopodial tension model of CE will inspire the additional experimental measurements that a more detailed quantitative validation requires.
Our model predicts that external forces, such as surface tension and pressure, can only affect the final degree of tissue elongation (κ) (Fig 4), whereas the internal parameters that regulate cell-intercalation can affect both tissue dynamics (as measured by τ) and the final tissue shape (Figs 4 and 5). This can be easily tested experimentally by either changing the properties of the external environment (the surrounding cells/matrix) of the intercalating tissue or culturing it ex vivo.
Of the five cell intercalation parameters, the time interval between link formation/breakage (t interval ) had negligible effects as long as it is below the typical time that the cells take to rearrange positions and/or shapes in response to a given set of external forces. This might be different if a refractory time interval between pulls is added to the model. We expect that in the presence of such refractory time, an increased frequency in link formation/breakage would slow down the speed of intercalation and reduce the final elongation ratio.
The model also predicts that the maximum range of cell interaction (r max ) and the maximum number of links per cell (n max ) had no effect on the final tissue elongation after r max = 2 and n max = 3, but the time of elongation kept decreasing for higher values of r max and n max ( Fig  5B and 5C). It was not possible to increase those parameters indefinitely to check if the speed of intercalation would also saturate because the simulated cells start to fragment past r max ! 6 or n max ! 7. The current implementation of the model, however, does not allow for more than one active link between cell pairs, which would likely decrease elongation time.
All cell intercalation models assume some type of increased cell activity along the convergence axis, which is often translated, as is the case here, into the assumption that the cells are bipolar (one exception being our 3D equatorial/extensional filopodial tension model). This however is not necessarily true and we expect the model to also work in cases where the simulated cells are either monopolar in opposite directions of the same convergence axis or randomly alternate being monopolar in each direction of the convergence axis.
Our model also suggests that CE can be successfully achieved even in the presence of relatively high degrees of polarization defects. We predict that tissues containing polarization misalignments of up to ±10°will be practically indistinguishable to the optimally aligned case. Even severe misalignments (about ±50°) would still lead to some CE, although to a much lesser degree of final elongation ratio and with longer tissue elongation times (Fig 8). Experimental disruptions of the PCP pathway that alter the global alignment of cells in a dose-dependent manner would provide a way to test some of these predictions. Addition of a simple mechanical feedback mechanism by which the cells readjust their polarization in response to the pulling forces from the neighbors does not have major effects on the speed of tissue elongation ( Fig  8C), but can fully rescue and even improve on the final elongation ratio of tissues with low or even severe polarization misalignments (Fig 8B). We choose to implement a minimal phenomenological feedback mechanism to explore the general response and self-organization of the tissue, but the formulation of the model allows the replacement of this generic mechanism with more detailed feedback models that reflect a specific tissue.
In cases where some cells fail to polarize, the severity of the effects on CE will depend on the type of interaction between the polarized (intercalating) cells and the unpolarized (nonintercalating cells). If the polarized (or active) cells can still pull on the unpolarized cells, then CE still happens even in a situation where the vast majority of cells (95%) are not active (Fig  9), although at a great reduction in both speed and final elongation ratio. If, on the other hand, the unpolarized cells are non-responsive and cannot be pulled, then the reductions in speed and final elongation ratio are much more sensitive to the presence of unpolarized cells (Fig 9) and CE completely fails when the population of active cells falls below 25% (Fig 9C).
Another prediction of the model is the separation between the intercalating cells and the non-responsive/refractory cells (Fig 9D). Such defects could be induced experimentally by randomly distributed knock-out of intercalating cells, e.g. using electroporation of tissues with a dominant negative or RNAi, and would provide a way to further test model predictions.
Finally, the model reduces to the more common implementations when the maximum range of interaction is replaced by the common contact area condition. In this case, instead of contracting (or increasing the tension of) the cell's surfaces that are aligned with the polarization vector or the global direction of convergence, we pull the neighbors that are not aligned with it (Fig 6D). In both cases active CE is achieved by the same principle, promoted cell-cell activity along one axis and inhibition along the other.