Cell mechanics, environmental geometry, and cell polarity control cell–cell collision outcomes

Abstract

Interactions between crawling cells, which are essential for many biological processes, can be quantified by measuring cell–cell collisions. Conventionally, experiments of cell–cell collisions are conducted on two-dimensional flat substrates, where colliding cells repolarize and move away upon contact with one another in “contact inhibition of locomotion” (CIL). Inspired by recent experiments that show cells on suspended nanofibers have qualitatively different CIL behaviors than those on flat substrates, we develop a phase field model of cell motility and two-cell collisions in fiber geometries. Our model includes cell–cell and cell–fiber adhesion, and a simple positive feedback mechanism of cell polarity. We focus on cell collisions on two parallel fibers, finding that larger cell deformability (lower membrane tension), larger positive feedback of polarization, and larger fiber spacing promote more occurrences of cells walking past one another. We can capture this behavior using a simple linear stability analysis on the cell–cell interface upon collision.

---

1 Introduction

Collective cell behaviors play crucial roles in many physiological processes such as embryonic development and morphogenesis and pathologies such as cancer metastasis. These collective behaviors are facilitated and regulated by intercellular interactions.¹ The simplest system to investigate cell–cell interactions is a pair of cells, where cell behavior after cell–cell collision can be tracked, allowing experimenters to determine the biological and physical factors regulating cell–cell interactions. One example of a stereotypical outcome of two-cell collision is contact inhibition of locomotion (CIL), during which a cell reverses its direction of motion upon contact with another cell so that both cells migrate away from each other.² CIL is a common behavior in cell–cell collisions for many cell types including fibroblasts and neural crest cells when they are studied on two-dimensional flat surfaces.^2,3 CIL is important for cells' directional collective migration^4,5 and pattern formation in embryogenesis.^6,7 During CIL, the cell polarizes itself away from contact, e.g. by redistributing the activity of proteins including Rho-family GTPases such as Rac1 and RhoA.^8,9 In order to increase the frequencies of cell collisions to study the rules of CIL in a more controllable way, experimentalists have used one-dimensional microchannels¹⁰ and 1D micropatterned fibronectin stripes for cells to travel along, constraining their movement and collision in straight lines.^11–13 The outcomes in these collisions are dependent on various cellular properties including cell–cell adhesion, and can include reversing away from a head-on collision (e.g. classical CIL), cells managing to squeeze past one another to continue in their original direction (“walk-past”), or sticking together and/or forming a train of pairs of cells (“training”).^11–14 Cell–cell collisions also depend strongly on the geometry of the cell's environment – as do many other cell behaviors.¹⁵ Cells on suspended nanofibers,¹⁶ which mimic single fibers from a biologically realistic extracellular matrix (ECM)¹⁷ have drastically different collision behaviors than cells on featureless 2D substrates.¹⁸ These nanofibers more closely resemble native ECM^19,20 than macroscale fibers (e.g. those with radii similar to cell sizes), which have also been studied.²¹ On nanofibers, collision outcomes between cells crawling along a single fiber are predominantly walk-past, while collisions between cells attached to two parallel fibers predominantly result in training (or walk-past if they have recently divided).¹⁸ Large rates of cells walking past other cells have also been observed for cancerous cells on micropatterns.^13,22

What determines collision outcomes between cells in various environments? What cellular properties and environmental factors control cell motion and interactions, giving rise to a certain cell–cell collision result? To explore these questions, we model cell behaviors during collisions in suspended fiber geometries using numerical simulations, studying the effects of mechanical properties of cells and the fiber geometry on collision outcomes. In particular, our focus is on the “walk-past” behavior where cells manage to push past one another, and capturing the role of the large cell shape changes that are required in this process. The complicated deformation involved in such processes cannot be captured by our original self-propelled particle model,¹⁸ which describes cells as spherical beads. To model a deforming cell with a complex shape, we use a phase field approach to simulate cell motion and two-cell collisions. Phase field models have been broadly used to study cell migration at the single-cell^23–26 or collective level,^27–33 including studying cell collisions on flat surfaces³⁴ or micropatterned substrates^26,35 as well as hydrodynamically-interacting cells.³⁶ The phase field approach has also been extended to 3D,^27,37,38 including simulations of cells on macroscale fibers with significant curvature effects.³⁸ Our model is most closely related to those of ref. 26 and 31. We describe the cell's polarity – the biochemical difference between its front and back – with a simple field driven by positive feedback, akin to ref. 39 and 40. This simple model is distinct from previous modeling of single cells^41–43 and cell–cell collisions using more detailed Rho GTPase polarity models.^35,44 We also include cell–cell and cell–fiber adhesion (primarily facilitated by E-cadherin junctions and focal adhesions, respectively¹) which can have significant effects on collective cell migration.⁴⁵ We find we can control whether cells walk past each other or form trains with two cells moving in the same direction by tuning cell deformability, the strength of positive feedback driving cell polarity, and the distance between fibers. We capture these behaviors qualitatively with a simple linear stability analysis focused on whether cells can break symmetry at the collision to walk past one another. Our simulation results and linear stability theory, though based on a minimal modeling approach, help to reveal factors playing key roles in controlling cell–cell collisions.

2 Theoretical model and simulation method

We use a two-dimensional phase field model to simulate pairs of cells moving along a single fiber or two parallel fibers (Fig. 1 middle and right columns, compared with experimental images from ref. 18 in the left column). A cell's shape is represented by a dimensionless scalar field ϕ(r,t) where ϕ transitions continuously between one inside the cell and zero outside the cell, implicitly defining the cell boundary by ϕ = 0.5 (black contours overlaid on the phase fields of Fig. 1). We will generally work in simulation units, measuring lengths relative to the numerical lattice on which we solve our partial differential equations (with lattice size l₀), along with a time scale τ₀, and an energy scale ε₀. Roughly matching cell sizes and velocities as well as fiber diameters, we estimate that these scales correspond to a physical scale of l₀ of 0.25 μm and τ₀ of 0.5 min. In this section, we first introduce free energies of cells in a fiber geometry, defining cell–substrate and cell–cell interactions, and then specify the mechanisms of cell motility and polarity that regulate cell dynamics, before briefly describing the simulation method of cell collisions. Default values of simulation parameters can be found in Table 1 and details of the numerical method can be found in Appendix B.

Fig. 1 Examples of cell fluorescence images (left column, images of fibroblasts reprinted from ref. 18) and cell phase fields (ϕ, middle and right columns) on suspended fibers (white dashed lines in phase fields). Cells expand along fibers and become either spindle-like ((a), on a single fiber) or parallel-cuboidal ((b), on two parallel fibers). Fluorescence images of both spindle- and parallel-shaped cells (on fibers of various diameters as labeled in each image) in the left column have been previously published in ref. 18, with actin stained in red and nuclei in blue and scale bars 50 μm (see details in ref. 18). In the phase fields of the middle and right columns, the cell boundary (ϕ = 0.5) is marked by a thin black contour. Middle column shows cells with no active forces. Right column shows the steady-state shapes of single cells with the full dynamic simulation. In the right column, thick magenta line segments indicate the locations where motility forces are largest, where the absolute value |[Doublestruck P]ϕ²(1 − ϕ)²χ| > 0.025. The polarity fields [Doublestruck P] associated with these cells are shown in Fig. 6.

Table 1 Summary of default parameter values used in simulations

Parameter	Value	Location
Length unit (lattice size) l0	0.25 μm
Time unit τ0	0.5 min
Cell boundary width d		Eqn (1)
Cell initial radius R0	40l0
Cell area constraint strength λ	1000ε0
Cell–fiber adhesion strength E0	2000ε0
Cell adhesion reference area Aref	1l0^2
Cell–cell area exclusion strength g	1ε0/l0^2	Eqn (4)
Cell–cell adhesion strength σ	5ε0/l0^2
Transport coefficient 1/γ	50l0^2/(ε0τ0)	Eqn (6)
Friction coefficient η	1ε0τ0/l0^4	Eqn (7)
Cell velocity noise strength ṽ0	0.1l0/τ0
Motility force strength A	200ε0/l0^3	Eqn (8)
Polarization decay time constant τ	100τ0	Eqn (9)
Simulation time step Δt	0.002τ0	Eqn (B.1)

---

2.1 Free energy functionals

We first define cell free energy functionals that control cell shapes in our 2D phase field. We model a cell as having a characteristic tension opposing extension of its membrane (a line tension in our 2D model), a constraint on area, and an adhesive energy promoting its spreading along the nanofibers. As a result, the free energy functional for a phase field ϕ_i(r,t) of a single cell i is:

The first term in eqn (1) is a double-well potential similar to a free energy functional usually studied in Cahn–Hilliard models.^26,30,31 The parameter α, which is related to cell membrane tension, has units ε₀/l₀². The double-well free energy enforces the phase field to have ϕ_i = 1 inside cell i and ϕ_i = 0 outside the cell, and penalizes the formation of a longer cell boundary. The characteristic length d sets the cell interface width. The line tension of the cell boundary is Γ (units of energy/length, ε₀/l₀), given by³⁰

The line tension Γ is directly proportional to α, so in this work, we vary α to control line tension. The second term of eqn (1) is a soft constraint on the cell area , which penalizes deviations from the area of a circle with radius R₀. The cell radius is R₀ = 40l₀ = 10 μm, roughly reasonable for many eukaryotic cells. If only line tension and cell area constraint were taken into account, a free, static cell would evolve to a circle of radius R₀.

The third term in eqn (1) models the adhesive energy between the cell body and the fiber substrate with an adhesion strength E₀.^46,47 This negative energy promotes the expansion of the area of adhesive contact between cell i and the fiber A_ad(ϕ_i). This adhesion energy saturates to its minimum value −E₀ as the adhesive area exceeds the scale A_ref = 1l₀² reflecting the finite concentration of adhesive molecules (e.g., proteins like integrin that mediates the formation of focal adhesions^1,48). This saturating form is necessary, because a simpler definition such as −E₀A_ad(ϕ_i) without A_ref can make the cell phase field tend to spread around the entire fiber substrates, leading to cell breaking apart or numerical instability. We note that different definitions of cell–substrate adhesion, such as treating adhesion as a frictional force⁴⁹ or considering the adhesion dynamics with a deformable substrate,^50–52 have been used in other models. We define the cell–fiber adhesive area in eqn (1) as:

where we have introduced a separate phase field χ(r) that indicates the locations of fibers, which are present wherever χ = 1 and absent elsewhere (χ = 0). χ is nonzero over a width of 500 nm (two lattice spacings l₀), which matches the diameter of nanofibers used in the experiments of ref. 18. The expression in eqn (3) measures the overlapping area of cell–fiber contact where both ϕ_i and χ are nonzero. We use the somewhat nonobvious definition of eqn (3) to avoid numerical issues. This form of A_ad has a functional derivative δA_ad/δϕ_i = 6ϕ_i(1 − ϕ_i)χ – we can see that changes in the cell boundary, where ϕ_i(1 − ϕ_i) is nonzero, change the adhesive area, as we would expect. Our initial idea was using a simpler form ; we found this choice to lead to issues where ϕ_i has a low nonzero value along the fibers even far from the cell boundary.

The combined effect of our line tension and cell–fiber adhesion is that single cells that are initially circular and have no active driving forces (middle panel of Fig. 1) spread along the fibers (indicated by dashed white lines), but come to a steady state determined by the competition of line tension, which tends to promote a circular shape, and cell–substrate adhesion, which tends to elongate the cell along the fibers.

Cell–cell interactions are described by an additional free energy functional for cells i and j:^26,28

which includes the effects of both cell–cell repulsion (the exclusion of cell area overlap where both ϕ_i and ϕ_j are one, with an excluding strength g), and the adhesion with strength σ at cell–cell interface where cell boundaries (ϕ(1 − ϕ) ≠ 0) touch. The total free energy functional is:where i, j = 1 or 2.

2.2 Dynamics of the phase field model

Here we introduce the motility mechanisms that drive cell motion. The phase field ϕ_i evolves by:³¹which includes a velocity field v_i(r,t) arising from the cell's self-propulsion, which advects the cell boundary, as well as the second term on the right hand side, which describes the relaxation of cell shape to the local minimum of F.

The velocity field v_i for each cell is created by an active motility force density f_i(r,t) (force per unit area, with units ε₀/l₀³) arising from actin polymerization and myosin contraction on cell i. We assume that the relationship between these arises from the overdamped limit of viscous fluid flow,³¹ leading to:

The first term on the right hand side comes from a force balance between the friction force density −ηv_i (with a friction coefficient η) and the combination of the passive force density (δF/δϕ_i)∇ϕ_i (due to cell shape relaxation) and active motility force density f_i. Because cell–cell collisions are stochastic, we have added a velocity noise term ṽ_i, which is sampled and implemented every t_sample onto the cell velocity field. This term is uniform in space, modeling fluctuations in the whole cell's velocity. We assume the added noise is Gaussian, ṽ_i = ṽ₀(Xx̂ + Yŷ) drawing independent random variables X ∼ [scr N, script letter N](0, 1) and Y ∼ [scr N, script letter N](0, 1) at each sampling time point t_sample. This random velocity is held constant for a period of t_sample until the next sampling time point. We choose t_sample = 0.002τ₀ which happens to be equal to the simulation time step Δt used in our algorithm, but note that it should not vary if Δt is changed; see Appendix B. With our typical range of parameters (those in Table 1 and used in the plots of Fig. 3 and 4), we obtain average cell centroid speeds around 0.2 to 2 μm min⁻¹, depending on the specific parameters chosen. This range of cell speeds is roughly consistent with that seen in experiments.¹⁸ We note the mean cell speeds are much larger than the velocity noise strength ṽ₀ = 0.1l₀/τ₀ = 0.05 μm min⁻¹ – our noise here is a small perturbation, meant to break symmetry rather than realistically represent the degree of variability in cell speeds observed in ref. 18.

2.3 Cell polarization field

Where the cell exerts protrusive vs. contractile force will depend on the cell's polarity – what part of the cell is the “front.” We characterize this polarity using the dimensionless scalar field [Doublestruck P]_i(r,t) of cell i. This single field [Doublestruck P]_i is a summary of multiple biochemical processes and molecular details that characterize a cell's polarity, including the distributions of Rho GTPases like Rac1 and RhoA, as well as actin and myosin, centrosome orientation, etc.,¹ phenomenologically modeling their force-generating effects in the cell. Our convention is that protrusion occurs where [Doublestruck P]_i > 0 (cell front, where Rac1 facilitates actin polymerization) and contraction occurs where [Doublestruck P]_i < 0 (cell rear, with RhoA activity and myosin-driven cytoskeleton contraction).^8,9,35 We assume that cell–substrate contact is necessary for force generation,^53,54 so force density f_i is localized at the overlaps between cell boundary and fiber substrates (also the typical locations of focal adhesions⁵⁵). To localize force to this region, we need both χ = 1 and ϕ_i²(1 − ϕ_i)² to be nonzero, suggesting:

f_i = A[Doublestruck P]_iϕ_i²(1 − ϕ_i)²χn_i (8)

where A sets the strength of protrusion and n_i = −∇ϕ_i/|∇ϕ_i| is the outward normal unit vector to cell i.

We want to capture persistent cell migration. We use a very simple positive feedback mechanism, akin to those proposed by ref. 39 and 40 in which successful protrusions promote further protrusion (increasing polarity). Within our phase field framework, we thus promote increases in polarity [Doublestruck P]_i where the cell has already started to protrude – the regions where ∂ϕ_i/∂t > 0, and decrease [Doublestruck P]_i where contractions occur ∂ϕ_i/∂t < 0. The simplest implementation of this idea is:

where β is a dimensionless parameter that sets the positive feedback strength, and τ sets the timescale over which polarity decays. We adopt, for simplicity, the convention that polarity evolves in the reference frame of the substrate – appropriate if polarity is tightly associated with the basal membrane or with focal adhesions – so that there is no velocity advecting [Doublestruck P]_i. An initially static cell will have ∂ϕ_i/∂t = 0 everywhere, which tends to make [Doublestruck P]_i disappear. To break this symmetry, we set an initially unbalanced distribution of polarity [Doublestruck P]_i inside the cell to start dynamic simulations of moving cells (Section 2.4 and Appendix B). A large enough random perturbation in polarity also should break the initial symmetry, though we do not explore this in depth. While the polarization field [Doublestruck P]_i is created only interior to the cell, where ϕ_i > 0, it may persist transiently after the cell leaves, decaying with time constant τ; see Fig. 6 in the Appendix. In these images, for the sake of transparency, we show the polarity field in the whole area, while only the polarity near the cell boundary actually leads to any effect in the model. We note that this transient nonzero [Doublestruck P]_i behind a moving cell is a nonphysical side effect of the numerical model we choose which does not include a term to explicitly degrade polarity once the cell leaves the substrate-localized field behind. When evolving the polarity dynamics of eqn (9), we set [Doublestruck P]_i = 0 far from the cell (where ϕ < 0.001); see Appendix B. This is similar in spirit to the approach of ref. 56 for solving reaction-diffusion equations in an arbitrary geometry.

2.4 Simulation method

We solve the phase field and polarity equations eqn (6) and (9) using a finite difference method (Appendix B). To start a simulation trajectory, we initialize the phase field system by setting a circular cell of radius R₀ on the fibers and allowing the cell to spread over the fibers without active forces, obtaining a static “relaxed” cell shape. These are shown in the middle panel of Fig. 1. We then simulate the free motion of a single isolated cell by driving the obtained static cell to move along the fibers while allowing the polarity field to evolve normally by eqn (9) starting with an asymmetric polarization field. The cell in motion reaches a steady state with a stable shape after long enough time (e.g. in the right panel of Fig. 1). We use such dynamic steady states as initial cell conditions when simulating two-cell collisions, ensuring that the cells involved in a collision are generated from identical parameters and already stabilized on fibers (having traveled for sufficiently long time and not affected by any particular initial state we choose for the simulation) prior to the collision. We simulate two-cell collisions under various conditions, whose outcomes depend on cellular properties and characteristics of the fiber environment. More simulation details can be found in Appendix B.

3 Results

If we simulate single cells on a single fiber, we obtain spindle-shaped cells (Fig. 1(a)), consistent with the experiments of ref. 18. Cells simulated on two parallel fibers take on a near-rectangular shape (e.g. in the middle and right panels of Fig. 1(b) where fibers are separated by spacing L = 90l₀ or 22.5 μm, which is comparable to fiber spacings in ref. 18) – we call these cells “parallel” or “parallel-cuboidal” as in the experimental work of ref. 18. See comparisons between experimental images and simulations in Fig. 1.

What happens if we collide two cells together? We find that if we collide cells on a single fiber (spindle cells), they almost always form a train (nearly 100% probability at our default parameters) with one cell moving in its original direction of motion and forcing back the other cell which repolarizes and turns around (Fig. 2(a), Movie 1). Experimentally, we earlier found that spindle cells most commonly walk past one another¹⁸ – but we would not expect to see this in our two-dimensional model, because for cells attached to single fibers, the walk-past involves cells moving past each other in the third dimension (which could be potentially treated using 3D extensions of phase field modeling like those in ref. 38 and 57). Instead, we see training, which was the second-most-common experimental outcome for spindle collisions in ref. 18. Because of this expected discrepancy, we do not focus on spindle-spindle collisions in this paper. Parallel-cuboidal cells are more planar in 3D real space and are better captured by the 2D phase field, so they are our primary interest here.

Fig. 2 Possible collision outcomes of cells on fibers. Directions of motion are marked by white arrows. (a) In our simulations, the collision of two spindle-shaped cells on a single fiber almost always results in cells sticking together and moving in the same direction (“training”). (b) The collision of two parallel-cuboidal cells on two fibers can also result in training, or possibly reversing from each other akin to classical CIL, which is most likely to occur in simulations without cell–cell adhesion. In addition, cells can walk past each other while keeping their original directions of motion. During a walk-past, both cells shrink to a single fiber and become teardrop-like. In other collision outcomes (training and reversing), one or both of the cells can also shrink to single-fiber shape (not shown in this figure). Cell–cell collision movies for common outcomes and a few unusual cases can be seen in Movies 1–8.

If we simulate collisions of parallel cells attached to two fibers (also separated by L = 22.5 μm in Fig. 2(b)), we find, depending on the parameters, multiple distinct outcomes, including training, walk-past, and rarely, cell reversals. This is consistent with the behavior of fibroblasts on parallel nanofibers, which we earlier found to either form trains with near-100% probability (cells that have not recently divided) or most commonly walk past one another (cells that have recently divided, which are faster, and cells in low calcium media).¹⁸ In a training event, cell pairs spontaneously break symmetry and end up moving together either leftward or rightward (see Movie 2). In a walk-past, in order to keep their original opposite directions of motion and avoid each other, the two cells initially attached to both fibers will continue protruding on different single fibers along their respective directions upon contact (as in the last row of Fig. 2). Each cell gradually shrinks to attach only onto the fiber where it continues protruding and detaches from the other fiber, and becomes teardrop-like in shape on single fiber (as spindle cells in motion, see Movie 3). The cells are in steady motion on single fibers after collision and do not resume the initial two-fiber parallel shapes. In the least likely outcome of reversing, both cells repolarize and turn around as in a classical CIL (see Movie 4). In most cases, cells in training or reversing keep their parallel shapes throughout the process, but in a few cases one cell (see Movies 5 and 6) or both cells (see Movies 7 and 8) could shrink to single fiber (similar to cells experiencing walk-past).

What factors control the outcomes of cell–cell collisions? We study the effects of varying three parameters on collision outcomes between parallel cells on two fibers: the line tension parameter α (which is proportional to line tension Γ, characterizing the softness or deformability of a cell), the positive feedback strength β of cell polarity, and the spacing L between parallel fibers. We sweep the parameters (α, β) over the ranges of 0.1ε₀/l₀² ≤ α ≤ 0.2ε₀/l₀² with an interval of Δα = 0.02ε₀/l₀² and 10 ≤ β ≤ 50 with an interval of Δβ = 4, and simulate three fiber spacings L = 90l₀, 80l₀ and 70l₀ (22.5 μm, 20 μm and 17.5 μm, respectively, similar to experimental values from ref. 18). For each of these parameter sets, we perform 96 independent simulation trajectories and obtain the statistics of the number of each possible collision consequence, identifying the dominant outcome and analyzing both qualitative and quantitative effects of the parameters. These results are shown in Fig. 3(a) (left panel) and (b) (right panel), where the phase diagrams in (a1) and (b1) illustrate the dominant collision outcomes. We call an outcome dominant if it occurs at least 53 times out of 96 simulations (or 55% of the time), while we also show the coexistence of significant presence of both walk-past and training (labeled by superposition of symbols in phase diagrams) when the difference of frequencies between the two outcomes is less than 10 (roughly 10%).

Fig. 3 Phase diagrams showing dominant cell collision outcomes on two parallel fibers separated by a spacing L = 90l₀ (22.5 μm) in (a1) and L = 70l₀ (17.5 μm) in (b1). On top of each phase diagram is a schematic showing a pair of cells traveling toward each other before collision along parallel fibers of corresponding spacing. (a2) and (a3) show the frequencies of training and walk-past results, respectively, calculated from 96 independent simulation trajectories at each point in phase diagram (a1), while (b2) and (b3) show the same frequencies for each point in (b1). The blue dashed line in all the figures is a global fit to the linear stability analysis predicting the transition between training and walk-past. Note that when L decreases from (a) to (b), more points are dominated by training and only a few are dominated by walk-past, and the linear stability prediction of walk-past threshold (dashed blue line) also moves upward accordingly.

3.1 Effect of line tension, positive feedback strength, and fiber spacing on cell–cell collision outcomes

We show the results of our simulations on the widest fiber spacing L = 90l₀ in Fig. 3(a) (left panel). Most outcomes of these collisions are either walk-past events (stars) or training events (crosses) in Fig. 3(a1). There is a transition between primarily training outcomes (at high α – i.e. high line tension – and low positive feedback strength β) and primarily walk-past (high β and low α). More quantitatively, the probability of training and walk-past at each parameter set is shown in Fig. 3(a2) and (a3). This transition from training to walk-past at high β and low α makes intuitive sense. The force driving the cells to walk past one another depends on the positive feedback β – more polarized cells are more persistent and more likely to walk past. However, to walk past, cells must deform, which is resisted by the cell's line tension Γ ∼ α. We will later formalize this intuition in a linear stability theory in Section 3.3, which predicts the shape of the dashed line in Fig. 3. We also find that at sufficiently small α (≤ 0.14ε₀/l₀² in the range of the phase diagram) and large enough β, parallel cells are themselves not stable. A cell initially attached to two fibers will become unstable and unable to keep hanging onto both fibers under the strong driving force, and shrink to a single fiber even without interacting with another cell; these cases are indicated by blue dots in Fig. 3(a1). An example of this behavior is shown in Movie 9. With even smaller α or larger β, cell phase fields break apart and fail to simulate cell behaviors normally (as labeled by cyan diamonds). In these cases, cells experience large-scale deformations, such as being stretched along a fiber to a large extent while narrowing in the perpendicular direction so that the cell interior thickness is no longer larger than the finite boundary width d used in our phase field model. This sort of behavior would be expected in the limit of zero line tension, when the cell would prefer to completely wet the substrate.

How do cell–cell collisions change when we change the fiber spacing L? We show outcomes of cell–cell collision simulations at narrower fiber spacing (L = 70l₀) in Fig. 3(b) (right panel). Comparing with the left panel, we see that training is far more common at smaller fiber spacing – the transition between training and walk-past moves toward smaller α and larger β. This shift in the transition line is predicted by our linear stability analysis in Section 3.3, given by the blue dashed line in Fig. 3(b1). Given our earlier intuition, we think of decreasing the fiber spacing L as effectively making cells harder to deform, with similar effects as increasing the line tension (or α). Unsurprisingly, we find that collision outcomes with spacing L = 80l₀ are intermediate between L = 70l₀ and L = 90l₀ (Fig. 7 in the appendix.)

3.2 Collision outcomes of two cells on parallel fibers with zero cell–cell adhesion

Experimentally, fibroblasts in low-calcium media (which have low cell–cell adhesion due to the calcium-dependence of cadherin interactions) were found to almost always walk past one another in parallel collisions.¹⁸ However, our physical intuition for our model suggests a different outcome – cell–cell adhesion will decrease the energy per unit length of the cell–cell contact, making it easier for one cell to walk past another. We test this idea by simulating collisions between cells with σ = 0, removing cell–cell adhesion entirely (Fig. 4). We find that outcomes are mostly similar to the results with nonzero cell–cell adhesion (compared with Fig. 3 left panel), but can be qualitatively different in some regions of the phase diagram. Most dramatically, we see that cell reversal (classical CIL) emerges as the dominant outcome of collision (labeled by plus signs) with large line tension parameter α and small β. In this parameter range, cells are least likely to deform and training would be the dominant result with nonzero σ. This suggests that cell–cell adhesion helps keep cells stuck together until they come to a consensus direction and form a train. We also suspect that, similarly, if we disrupted cell–cell adhesion in experiments of fibroblast collisions in the spindle case, we could convert relatively common training events into reversals.

Fig. 4 (a) Phase diagram showing dominant cell collision outcomes without cell–cell adhesion, on two parallel fibers separated by a spacing L = 90l₀ (22.5 μm). The absence of cell boundary adhesion leads to the dominance of both cells reversing from each other upon collision at large α and small β. The linear stability prediction of walk-past threshold (dashed line) moves slightly to the left compared to the one at nonzero cell–cell adhesion in Fig. 3 left panel. (b), (c) and (d) show the frequencies of training, walk-past, and reversing results, respectively, calculated from 96 independent simulation trajectories at each point in the phase diagram.

What happens to the transition boundary between training and walk-past events? At smaller α and/or larger β across the central part of the diagram, the collision outcomes are dominated by training, similar to Fig. 3(a1), but with slightly more training and fewer walk-past results in the upper left half of the diagram – the absence of cell–cell adhesion shifts the transition line up and to the left on the graph. This small shift of the transition between walk-past and training is consistent with the relatively small effect of cell–cell adhesion on the interfacial energy, and can be incorporated into our linear stability analysis (see Section 3.3 and Appendix C). The shift of the blue dashed line can thus be predicted again using the linear stability analysis. Subplots in Fig. 4(b), (c) and (d) display the distributions of probabilities of training, walk-past, and reversing respectively, also calculated from 96 independent simulation results at each point.

We found that decreasing cell–cell adhesion increases training and decreases walk-past – in conflict with experimental observations of fibroblasts in a low-calcium environment which reduces cell–cell adhesion and leads to walk-past almost 100% of the time after collisions,¹⁸ but consistent with our physical intuition. We discuss potential explanations for this discrepancy later in the Discussion (Section 4).

3.3 Linear stability analysis

In an attempt to quantify the locations of transition between walk-past and training in the phase diagrams in Fig. 3 and 4, we develop a simple linear stability theory for when walk-past is possible. Our core idea is that a necessary – though not necessarily sufficient – condition for walk-past to occur is that upon an initially symmetric cell–cell contact, the symmetric state should be unstable, leading to the cell–cell interface tilting (illustrated in Fig. 5). The growth of this tilted interface is penalized by factors like the line tension, which should suppress the development of the instability at cell–cell contact.

Fig. 5 Illustration of two cells on parallel fibers approaching each other and forming a cell–cell interface upon contact. Before the collision, the shape of both cells 1 (left) and 2 (right) represented by phase fields in (a) can be approximated by a rectangle of height L (the fiber spacing) and length l (along fibers) in (b). Upon collision, with the formation of interface between cell boundaries (thick green dashed line in (b)), the length of one side of a rectangular cell on a fiber extends to l + u and its other side on the other fiber shortens to a length of l − u, contributing to the increase of perimeter of each cell.

We first consider the energetic cost for the deformation of a cell attached onto parallel fibers separated by spacing L, assuming this energy is dominated by the line tension Γ. Initially, the cell shape resembles a rectangle of height L and length l (the left images of Fig. 5). The perimeter of this rectangular cell is 2(L + l). Upon collision with another cell, we see that the boundary of the cell–cell contact tilts and elongates. We approximate this by assuming a near-linear contact shape, with the length of one horizontal boundary segment along one fiber slightly elongating to l + u and the horizontal segment along the other fiber shortening to l − u. The length of the cell–cell contact thus extends from its initial length L to . The change of perimeter compared to an undeformed cell is . At the onset of collision we suppose u ≪ L, and when given the line tension Γ in eqn (2), we have the extra energetic cost H of forming the stretched interface approximated as (for each cell):

Our full phase field model has an overdamped dynamics in the absence of motility. If there were no active forces, the boundary would simply relax back to its equilibrium shape. Because of this, we expect an overdamped relaxation of the deformation u, i.e. with a > 0 a constant prefactor. Competing against this relaxation is the polarity model, which is driven by positive feedback – successful protrusion leads to more protrusion. This suggests that there should be a term that drives positive feedback in u, which promotes the growth of the cell–cell interface. We assume this term can be approximated by, at small u, some positive function f(β)u controlled by positive feedback strength β. This is a simple physical intuition in the first order approximation, in which ∂u/∂t = f(β)u will promote more ∂u/∂t when u > 0. Combining the two competing effects leads to

When f(β) < 4aΓ/L, any initial deformation of the interface u shrinks with time – u = 0 is the stable solution, making a walk-past unlikely to happen as an outcome of the collision. When f(β) > 4aΓ/L, the straight interface deformation u = 0 is linearly unstable, and we see that cells will break symmetry at the cell–cell contact. This is a necessary ingredient for one cell to go to the top fiber and one to the bottom, as occurs in a walk-past event. If the function f(β) is linear, we would then predict that walk-past is only possible once β > β_crit, with β_crit a linear function of Γ/L, or equivalently α/L.

Our linear stability theory predicts that walk-past occurs if β is above a critical value, and that this critical value increases with Γ (which is proportional to α) and decreases with fiber spacing L – this matches our physical intuition, and is qualitatively consistent with our results in Fig. 3. Can we make this more quantitative? Because we have not linked a or the linear function f(β) to our simulation parameters, we must set these from our simulation results. We find that the value β_crit

provides a good global description of the transition between walk-past and training in Fig. 3 (and also Fig. 7) – this is the blue dashed line. We write this equation in terms of the nondimensionalized parameters and . The natural scale for α arises from thinking about the total boundary energy for a cell of perimeter 2πR₀, which is , with (eqn (2)). Our threshold in eqn (12) quantitatively captures the transition between dominances of training and walk-past outcomes, and explains why small α, large β and large L could more predominantly cause a walk-past.

Within our linear stability theory, what is the role of cell–cell adhesion? Until now, we have assumed that the energy of the cell's boundary is solely controlled by the line tension Γ, without addressing the quantitative effects of the weak cell–cell adhesion in our model (σ). However, our cell–cell adhesion term in eqn (4) will lower the energy of the cell–cell boundary, effectively changing the cell–cell interfacial tension to Γ − cσ, with c some constant (see, e.g. ref. 28). We compute this value of c (cσ turns out to be small relative to Γ) and show how this alters our predicted β_crit in Appendix C. The effect of the absence of cell–cell adhesion is relatively minor, shifting the transition line leftward to smaller α by an amount of 0.0134ε₀/l₀² in the phase diagrams from Fig. 3(a1) to Fig. 4(a) as a result of changing σ = 5ε₀/l₀² (in eqn (4)) to σ = 0. We note that this calculation of c is fairly rough, so we see some disagreements with the predicted line in, e.g.Fig. 8 with σ = 0 and smaller fiber spacing L where some points dominated by walk-past are below the shifted line. However, our model correctly captures the direction and overall scale of the shift in the phase diagram due to cell–cell adhesion.

There are some subtleties in our linear stability analysis that we should mention. First, the symmetry-breaking we are describing here must take place in competition with other processes we have not described – e.g., a repolarization away from the cell–cell contact which could lead to a reversal. Our theory will thus never capture the transition to reversals in Fig. 4. Secondly, we note that the same basic competition between positive feedback driving asymmetry and line tension suppressing it can also be present in a single cell, and not just a cell–cell collision, though there will be less positive feedback, since only one cell will be pushing at the interface. We think this explains the cases where cells are unstable on both fibers (blue dots in Fig. 3), which will also be more prominent at low α and high β. We also note that our specific fit parameters, in eqn (12), are tuned to the transition between walk-past and training being predominant, not the first appearance of walk-past. We believe this point is easier to identify in a finite, noisy data set.

4 Discussion

Our simulations and linear stability analysis show that both mechanical factors like the cell tension and biochemical factors like the positive feedback strength, along with the fiber spacing, can control the outcome of cell–cell collisions in the parallel geometry, switching cells between walk-past and training. Experimentally, NIH 3T3 cells in the parallel geometry almost always train, unless they have recently divided, in which case they most commonly walk past one another.¹⁸ Cells that have recently divided have higher speeds¹⁸ – and in our model, cells with higher positive feedback values β also have higher cell speeds (see Fig. 9 for a plot of cell speed changing with β). Our model is in reasonable agreement with this central finding of ref. 18. Our approach also has significant advantages in simulating collisions in parallel geometry over our simple self-propelled particle model. Though the active particle model in ref. 18 can generate walk-past events even on parallel fibers, these events are somewhat unnatural, requiring either an extreme deformation of the cells or cells being pushed off the fibers. In particular, we could not find parameters of the particle model that allowed walk-past nearly 100% of the time, as seen in the low-calcium experiments of ref. 18. In contrast, our phase field model that simulates cell shape changes is able to generate a robust transition between training 100% of the time to predominantly walk-past as observed in experiments (e.g. by simply adjusting β). Our model in this paper allows more realistic study of cell shape mechanics and fiber geometry than ref. 18. We note that although the 2D phase field model cannot be readily mapped to the 3D particle model for direct comparison, efforts have been made to relate phase field approaches to simplified minimal models under certain assumptions (such as sharp-interface limit and quasi-circular approximation for cell shapes),^42,50,52 and the applicability of these simplifying methods to cell collisions on fibers can be explored in future endeavors.

Our simulation results provide some possible directions for experimental explorations. The idea of controlling cell speeds and collision outcomes through positive feedback of cell polarity may be more strongly experimentally tested by looking at, e.g. cells with up-regulation of Rac1 activity as a reasonable guess for the relevant positive feedback. The activation of Rac1 can promote stronger cell polarization,⁵⁸ leading to larger forces and higher cell speeds, which are modeled in this work by increasing β, indicating Rac1 regulation could be a controlling factor of the positive feedback mechanism. This regulation can be experimentally realized by controlling molecular pathways for Rac1 activation^58,59 or through photo-activation.⁶⁰ Our model could be further challenged by modification of cell tensions^59,61,62 to change the deformation energy of the boundary. Cell membrane tension may be directly modified through lipid composition or osmotic pressure,⁶³ or via Rac inhibition,⁵⁹ or comparisons could be made between, e.g. ectoderm, mesoderm, endoderm, which have systematically different tensions⁶¹ – though each of these approaches naturally will have potential confounding factors. We also predict a relatively large effect of the fiber spacing – a small decrease of inter-fiber separation L (about 22%) can drastically alter the distribution of dominant outcomes in phase diagrams from Fig. 3(a1) to (b1). This is possible to test, though not straightforward, as depending on the stiffness of the fibers, pairs of parallel fibers may be pulled together.^16,18 To reduce fiber deflection and keep spacing uniform, fibers with higher stiffness (and also larger diameter in the scale of μm rather than nm) could be used to study the effects of fiber spacing, though we note changing fiber diameter also can change cell protrusions and shapes (cells form shorter protrusions and become less elongated along fibers with larger diameters¹⁸). Through manipulating cell tensions, positive feedback, and fiber spacing in experiments, our simulation results of how cell collision outcomes depend on these factors (including the linear stability theory) can be tested in a more realistic and controllable way.

One of the most striking common findings in both Singh et al.¹⁸ and our numerical studies is the near absence of the classical CIL – the mutual reversal, either on one fiber or two fibers. Within our model, reversals are largely present only when the cell–cell adhesion is zero and the cell's positive feedback strength β is low and line tension is high – suppressing walk-past. However, this is also the largest contrast between our results and the experimental results in Singh et al.¹⁸ – we predict relatively minimal effect of cell–cell adhesion on the transition between training and walk-past (in fact, the absence of adhesion promotes slightly more training). By contrast, parallel NIH 3T3s exhibit training with normal cell–cell adhesion but would switch to walk-past in calcium-free media with low adhesion. This inconsistency suggests either (1) our model is oversimplified, (2) disruption of E-cadherin also disrupts cell–cell signaling that is necessary for contact interactions,⁸ and this also promotes walk-past, or (3) the low-calcium media perturbation used in ref. 18 also changed another key variable, e.g. changing the cell's tension or positive feedback. Calcium naturally plays a large role in many factors of cell migration,^64,65 and the low-calcium media can therefore affect more than just cell–cell adhesions. To more precisely manipulate cell–cell adhesion in experiments, genetic modification and other methods can be used to create mutants lacking adhesion⁶⁶ or cells with quantitatively different adhesion levels⁵⁹ to test whether the adhesion promotes training or walk-past. To model the biochemical signaling between cells that can affect collision outcomes, more complicated polarity models that explicitly simulate cell biochemistry can be used (like those modeling Rho GTPases in ref. 28 and 35) with direct interactions between cell polarities. This can also include a contact repolarization mechanism that leads to CIL on a flat surface and controls cell–cell collision,^26,28,35 which could potentially promote walk-past on double fibers when cells can more easily break symmetry and push each other onto different fibers upon contact, and also help to reveal the role of cell–cell adhesion in the process.⁶⁷

Our model provides a key new idea – that walk-past may be controlled by linear stability of the symmetric cell–cell interaction, and thus managed by cell mechanics and environment geometry. These results could be brought into closer connection with experimental work in a few different directions, either via expanding our current “bottom-up” model, or via a data-driven “top-down” approach.^22,68 It would be interesting to see if our linear stability hypothesis also holds in more realistic bottom-up models of cell polarity.^25,28,35 However, extension of these models to the fiber geometry may not be straightforward. Our model reproduces spindle and parallel cell shapes similar to experimental observations, but real cell bodies are much more stretched and elongated along fibers.^16,55 More detailed models of protrusion formation and stabilization, e.g. following ref. 69 and 70, might allow for improvement in this direction. Another potential path to closer contact with experiment is to connect to data-driven models for cell collision learned from experimental measurements.^22,68 The measurements of ref. 22 show that noncancerous cells are characterized by more cell–cell repulsion and cell–substrate friction, while metastatic cells feature more cell–cell attraction and anti-friction, leading to more sliding events (walk-past). Connecting these quantitative details to the underlying simple linear stability idea could be a fruitful area for future work.

Author contributions

YL, ASN, and BAC conceptualized the research. YL developed software and ran simulations. YL, ASN, and BAC analyzed data.

Conflicts of interest

There are no conflicts to declare.

Data availability

Supplementary information is available. See DOI: https://doi.org/10.1039/d5sm00359h. The supplementary information includes 9 movies, whose captions are provided in Appendix D. Data for this article, including sample codes of numerical simulation and statistical results, are available at https://doi.org/10.5281/zenodo.16923125.

Appendices

A Supplementary figures

Additional figures mentioned in the main text, including illustrations of polarity fields associated with the cells in Fig. 1 (Fig. 6), simulation outcomes of cell collisions at intermediate fiber spacing L = 80l₀ (Fig. 7) and at L = 70l₀ and L = 80l₀ without cell–cell adhesion (Fig. 8), and the dependence of model cell speed on positive feedback strength β (Fig. 9), are shown in this section.

Fig. 6 Polarization fields of cells in motion (both on a single fiber and on double parallel fibers) shown in the right column of Fig. 1.

Fig. 7 (a) Phase diagram showing dominant cell collision outcomes on two parallel fibers separated by a spacing L = 80l₀ (20 μm). (b) and (c) show the frequencies of training and walk-past results calculated from 96 independent simulation trajectories at each point in the phase diagram.

Fig. 8 Phase diagrams showing dominant cell collision outcomes without cell–cell adhesion on two parallel fibers separated by spacings L = 70l₀ (17.5 μm) in (a1) and L = 80l₀ (20 μm) in (b1). (a2)–(a4) show the frequencies of training, walk-past, and reversing results calculated from 96 independent simulation trajectories at each point in (a1), while (b2)–(b4) show the same frequencies for each point in (b1).

Fig. 9 Example of cell centroid speed dependent on the positive feedback strength β of cell polarity when α = 0.2ε₀/l₀² and fiber spacing L = 90l₀. Each data point is an average speed from 20 simulations and the standard error is smaller than the marker size.

B Simulation details

In our numerical modeling, the 2D simulation space used for the cell phase field ϕ_i and its associated polarization field [Doublestruck P]_i is discretized with lattice size l₀ = 0.25 μm, and the simulation time step is Δt = 0.002τ₀ = 0.001 min (see Table 1). We use a periodic boundary condition for the simulation space, which has a dimension of 200l₀ × 200l₀ for a single cell (Fig. 1) or 400l₀ × 200l₀ for two cells (Fig. 2). At the beginning of each simulation trajectory, ϕ_i is initially in a circular shape of radius R₀ = 40l₀, which is then placed on top of the fiber substrate field χ(r) (consisting of either a single fiber or double parallel fibers whose width is 2l₀) and relaxed without active forces. We use a forward in time, central in space (FTCS) scheme to calculate the spatiotemporal evolution of phase field ϕ_i:where the time derivative ∂ϕ_i/∂t is calculated at r = (x, y) (position on the numerical lattice) and time t (using eqn (6)) from the instant evaluation of:for a single cell i without interacting with other cells, andall at time t. The calculation of v_i(r,t), which is needed for calculating ∂ϕ_i/∂t, is done by using eqn (7), which also uses the calculation of eqn (B.2) and (B.4) and implements the random noise ṽ_i every t_sample = 0.002τ₀ (with the noise generated at each sampling time point and kept being applied for a period of t_sample each time). Here t_sample is chosen to be equal to Δt, so that the velocity noise is sampled and added at each simulation time step and updated at the next time step in our algorithm. We note that changing the time step Δt should not vary the noise application timescale t_sample, and would require changing the frequency and duration of the application of the velocity noise. For instance, if Δt is changed to Δt_new = Δt/n then the noise should be sampled every t_sample = nΔt_new, after which point the same noise also persists for nΔt_new until being updated at the next sampling time point. This gives a consistent perturbation onto cell velocity. We mention this only for completeness – we have presented no results where the time step is changed.

When we initialize our simulations, we first relax the cell's shape by evolving it with no active motility force, i.e.f_i = 0. This cell shape relaxation lasts for a duration of 40τ₀, resulting in a static spindle shape (on one fiber) or a parallel shape (on two fibers) as in the middle column of Fig. 1.

After the initial static relaxation period, an asymmetric polarization field is applied to the cell, kick-starting the positive feedback and providing a driving force for cell motion. For example, for a cell to start moving rightward, we set [Doublestruck P]_i = 1 in the right half of the initially static, symmetric cell, and [Doublestruck P]_i = 0 in the left half of the cell (and also outside the cell). The polarization field [Doublestruck P]_i evolves with time through:

in which Θ(ϕ_i − ϕ_min) is a Heaviside step function equal to one wherever ϕ_i > ϕ_min = 0.001 and zero where ϕ_i < ϕ_min = 0.001. This sets [Doublestruck P]_i = 0 far from the cell. The time derivative ∂[Doublestruck P]_i/∂t is evaluated at r and t by using eqn (9), which requires the calculation of ∂ϕ_i/∂t using eqn (6) (which involves eqn (B.2)–(B.4) as we have illustrated above). Here the nonzero cell polarity gives rise to a motility force field f_i(r,t) calculated by eqn (8) using the current ϕ_i(r,t) and [Doublestruck P]_i(r,t). f_i(r,t) is further used in eqn (7) to obtain v_i(r,t), which is substituted into eqn (6) for the evaluation of ∂ϕ_i/∂t which is then used in both eqn (9) and (B.1) (and then (B.5)) to obtain ϕ_i(r, t + Δt) and [Doublestruck P]_i(r, t + Δt).

As the dynamic simulation proceeds, the freely migrating cell takes an asymmetric shape (as in the right column of Fig. 1) and the polarization field takes the value of [Doublestruck P]_i > 0 near the cell front and [Doublestruck P]_i < 0 near the rear (as in Fig. 6). Both cell shape and polarization field become stabilized in steady cell motion before the end of this free motion period of 400τ₀ (unless the phase field breaks down or parallel cell becomes unstable and shrinks to single fiber as in the upper left corner of each phase diagram in Fig. 3 and 4, for which collision is not simulated). In the last simulation period afterwards (cell collision), we place two identical cells moving in opposite directions (already in steady motion, with phase and polarization fields generated in identical conditions from the previous dynamic stabilization period) together in the same simulation space (as in the left column of Fig. 2) and make them continue moving and experience a head-to-head collision. The presence of cell–cell interaction requires using the functional derivative of the full free energy δF/δϕ_i (not just δF_i/δϕ_i for single cell in eqn (B.2)):

for j ≠ i. We perform the simulation for another 400τ₀, and observe the collision outcomes by keeping track of both cells using their centers of mass, which we estimate by using the method introduced in ref. 72 in a periodic boundary condition. We classify trajectories either as training or walk-past, and all other trajectories as reversals, which we have verified is reasonable by visualizing trajectories. See code for details.

C Effect of cell–cell adhesion in the linear stability theory

To estimate the quantitative effect of the adhesion between cell boundaries in touch, specifically how much it affects the location of the transition line in the phase diagrams between dominantly walk-past and training as collision outcomes predicted from the linear stability analysis, we start by considering the simplified cell shapes made of straight-line boundary segments in the upper panel of Fig. 10 (similar to those in Fig. 5(b)) of two cells coming into contact with one another at horizontal position x = 0. When crossing a cell boundary in our simulations, the 2D phase field of a model cell transitions between ϕ = 1 and ϕ = 0 in a smooth way (across a boundary width ) that can be described by a hyperbolic tanh function. In the example of Fig. 10, which shows the 1D phase field profiles of cells 1 and 2 along the x direction (as cross sections of their 2D phase fields) in the lower panel, functions and are used to represent the right boundary of cell 1 and the left boundary of cell 2, respectively, which are in contact at x = 0 (where θ₁(0) = θ₂(0) = 0.5). We suppose the interface formed between the cells has a length of l_int, and using the second half of eqn (4), we calculate the contribution of cell–cell adhesion to the total free energy (lowering of free energy) for each cell to be:

Fig. 10 Example of cell–cell contact in 2D (upper panel) and 1D (lower panel).

This means that – from the standpoint of our linear stability theory – that the cell–cell adhesion will act as a negative line tension. Though we derived this for two cells with a straight cell–cell boundary, we expect this expression can be used to estimate the adhesion energy between phase field cells with an interface of length l_int in any shape as long as the interface thickness d times the curvature of the boundary is sufficiently small.⁷³ For the same interface, the energetic contribution from line tension Γ is according to eqn (2).

To quantify the amount of shift of the location of linear stability prediction (dashed threshold) in the phase diagrams, we combine both line tension and adhesion effects and find that by having ΔE_adh = ΔΓ_adhl_int (or ), the cell–cell adhesion effectively changes the line tension to Γ + ΔΓ_adh = Γ − cσ where c = 0.008l₀, and contributes an amount of Δα_adh = −0.0134ε₀/l₀² to the parameter α given other parameters in Table 1. When the adhesion is turned off by changing σ = 5ε₀/l₀² to σ = 0, the threshold moves leftward horizontally to smaller α by 0.0134ε₀/l₀².

D Movie captions

Movie 1: A typical collision between spindle-shaped cells on a single fiber, which almost always results in training. Here we use parameters α = 0.1ε₀/l₀² and β = 10.

Movie 2: A collision between parallel-cuboidal cells on two parallel fibers, which results in a training. Here we use parameters α = 0.1ε₀/l₀² and β = 10, and fiber spacing L = 22.5 μm.

Movie 3: An example of walk-past between parallel cells upon collision on two fibers. Here we use parameters α = 0.1ε₀/l₀² and β = 14, and fiber spacing L = 22.5 μm (here walk-past is a possible but not dominant outcome). Both cells shrink to different single fibers (with each cell randomly shrinking to one of the two fibers) and become teardrop-shaped during the walk-past.

Movie 4: An example of two-cell reversing (like a classical CIL) upon collision on two parallel fibers. Here we use parameters α = 0.1ε₀/l₀² and β = 10, and fiber spacing L = 22.5 μm (here reversal is a possible but not dominant outcome). Both cells keep their parallel shapes after collision.

Movie 5: An example of collision between parallel cells on two fibers, which results in a training during which one cell shrinks to a single fiber. Here we use parameters α = 0.12ε₀/l₀² and β = 14, and fiber spacing L = 22.5 μm.

Movie 6: An example of collision between parallel cells on two fibers, which results in a reversal during which one cell shrinks to a single fiber. Here we use parameters α = 0.14ε₀/l₀² and β = 14, and fiber spacing L = 22.5 μm (here reversal is a possible but not dominant outcome).

Movie 7: An example of collision between parallel cells on two fibers, which results in a training during which both cells shrink to single fibers. Here we use parameters α = 0.1ε₀/l₀² and β = 14, and fiber spacing L = 22.5 μm.

Movie 8: An example of collision between parallel cells on two fibers, which results in a reversal during which both cells shrink to single fibers. Here we use parameters α = 0.1ε₀/l₀² and β = 14 as well as zero cell–cell adhesion, and fiber spacing L = 22.5 μm (here reversal is a possible but not dominant outcome).

Movie 9: A single cell initially attached to two parallel fibers becomes unstable and shrinks to a single fiber on its own under strong driving force. Here we use parameters α = 0.1ε₀/l₀² and β = 30, and fiber spacing L = 22.5 μm.

Acknowledgements

The authors acknowledge support from the National Science Foundation award MCB 2119948 (YL and BAC) and MCB 2119949 (ASN). This work was carried out at the Advanced Research Computing at Hopkins (ARCH) core facility (rockfish.jhu.edu), which is supported by the National Science Foundation (NSF) grant number OAC 1920103. We thank Pedrom Zadeh for useful conversations, including feedback on the polarity model we use, originally developed for ref. 71, and thank Emiliano Perez Ipiña and Wei Wang for a careful reading of the draft.

References

1. B. Ladoux and R.-M. Mège, Mechanobiology of collective cell behaviours, Nat. Rev. Mol. Cell Biol., 2017, 18, 743–757.
2. B. Stramer and R. Mayor, Mechanisms and in vivo functions of contact inhibition of locomotion, Nat. Rev. Mol. Cell Biol., 2017, 18, 43–55.
3. C. Carmona-Fontaine, H. K. Matthews, S. Kuriyama, M. Moreno, G. A. Dunn, M. Parsons, C. D. Stern and R. Mayor, Contact inhibition of locomotion in vivo controls neural crest directional migration, Nature, 2008, 456, 957–961.
4. E. Theveneau, L. Marchant, S. Kuriyama, M. Gull, B. Moepps, M. Parsons and R. Mayor, Collective chemotaxis requires contact-dependent cell polarity, Dev. Cell, 2010, 19, 39–53.
5. B. A. Camley, J. Zimmermann, H. Levine and W.-J. Rappel, Emergent collective chemotaxis without single-cell gradient sensing, Phys. Rev. Lett., 2016, 116, 098101.
6. J. R. Davis, C.-Y. Huang, J. Zanet, S. Harrison, E. Rosten, S. Cox, D. Y. Soong, G. A. Dunn and B. M. Stramer, Emergence of embryonic pattern through contact inhibition of locomotion, Development, 2012, 139, 4555–4560.
7. J. R. Davis, A. Luchici, F. Mosis, J. Thackery, J. A. Salazar, Y. Mao, G. A. Dunn, T. Betz, M. Miodownik and B. M. Stramer, Inter-cellular forces orchestrate contact inhibition of locomotion, Cell, 2015, 161, 361–373.
8. A. Roycroft and R. Mayor, Molecular basis of contact inhibition of locomotion, Cell. Mol. Life Sci., 2016, 73, 1119–1130.
9. R. Mayor and C. Carmona-Fontaine, Keeping in touch with contact inhibition of locomotion, Trends Cell Biol., 2010, 20, 319–328.
10. B. Lin, T. Yin, Y. I. Wu, T. Inoue and A. Levchenko, Interplay between chemotaxis and contact inhibition of locomotion determines exploratory cell migration, Nat. Commun., 2015, 6, 6619.
11. R. A. Desai, S. B. Gopal, S. Chen and C. S. Chen, Contact inhibition of locomotion probabilities drive solitary versus collective cell migration, J. R.Soc., Interface, 2013, 10, 20130717.
12. E. Scarpa, A. Roycroft, E. Theveneau, E. Terriac, M. Piel and R. Mayor, A novel method to study contact inhibition of locomotion using micropatterned substrates, Biol. Open, 2013, 2, 901–906.
13. D. F. Milano, N. A. Ngai, S. K. Muthuswamy and A. R. Asthagiri, Regulators of metastasis modulate the migratory response to cell contact under spatial confinement, Biophys. J., 2016, 110, 1886–1895.
14. W. Shi, S. Gupta, C. Copos and A. Mogilner, Collective mechanics of small migrating cell groups, Semin. Cell Dev. Biol., 2025, 166, 1–12.
15. G. Charras and E. Sahai, Physical influences of the extracellular environment on cell migration, Nat. Rev. Mol. Cell Biol., 2014, 15, 813–824.
16. K. Sheets, S. Wunsch, C. Ng and A. S. Nain, Shape-dependent cell migration and focal adhesion organization on suspended and aligned nanofiber scaffolds, Acta Biomater., 2013, 9, 7169–7177.
17. S. Meehan and A. S. Nain, Role of suspended fiber structural stiffness and curvature on single-cell migration, nucleus shape, and focal-adhesion-cluster length, Biophys. J., 2014, 107, 2604–2611.
18. J. Singh, A. Pagulayan, B. A. Camley and A. S. Nain, Rules of contact inhibition of locomotion for cells on suspended nanofibers, Proc. Natl. Acad. Sci. U. S. A., 2021, 118, e2011815118.
19. A. Jana, A. Tran, A. Gill, A. Kiepas, R. K. Kapania, K. Konstantopoulos and A. S. Nain, Sculpting rupture-free nuclear shapes in fibrous environments, Adv. Sci., 2022, 9, 2203011.
20. A. Jana, A. Sarkar, H. Zhang, A. Agashe, J. Wang, R. Paul, N. S. Gov, J. G. DeLuca and A. S. Nain, Mitotic outcomes and errors in fibrous environments, Proc. Natl. Acad. Sci. U. S. A., 2023, 120, e2120536120.
21. N. D. Bade, R. D. Kamien, R. K. Assoian and K. J. Stebe, Curvature and Rho activation differentially control the alignment of cells and stress fibers, Sci. Adv., 2017, 3, e1700150.
22. D. B. Brückner, N. Arlt, A. Fink, P. Ronceray, J. O. Rädler and C. P. Broedersz, Learning the dynamics of cell–cell interactions in confined cell migration, Proc. Natl. Acad. Sci. U. S. A., 2021, 118, e2016602118.
23. D. Shao, W.-J. Rappel and H. Levine, Computational model for cell morphodynamics, Phys. Rev. Lett., 2010, 105, 108104.
24. D. Shao, H. Levine and W.-J. Rappel, Coupling actin flow, adhesion, and morphology in a computational cell motility model, Proc. Natl. Acad. Sci. U. S. A., 2012, 109, 6851–6856.
25. E. Perez Ipiña, J. d'Alessandro, B. Ladoux and B. A. Camley, Deposited footprints let cells switch between confined, oscillatory, and exploratory migration, Proc. Natl. Acad. Sci. U. S. A., 2024, 121, e2318248121.
26. P. Zadeh and B. A. Camley, Picking winners in cell–cell collisions: wetting, speed, and contact, Phys. Rev. E, 2022, 106, 054413.
27. M. Nonomura, Study on multicellular systems using a phase field model, PLoS One, 2012, 7, e33501.
28. B. A. Camley, Y. Zhang, Y. Zhao, B. Li, E. Ben-Jacob, H. Levine and W.-J. Rappel, Polarity mechanisms such as contact inhibition of locomotion regulate persistent rotational motion of mammalian cells on micropatterns, Proc. Natl. Acad. Sci. U. S. A., 2014, 111, 14770–14775.
29. W. Wang, R. A. Law, E. Perez Ipiña, K. Konstantopoulos and B. A. Camley, Confinement, jamming, and adhesion in cancer cells dissociating from a collectively invading strand, PRX Life, 2025, 3, 013012.
30. B. Loewe, M. Chiang, D. Marenduzzo and M. C. Marchetti, Solid-liquid transition of deformable and overlapping active particles, Phys. Rev. Lett., 2020, 125, 038003.
31. G. Zhang, R. Mueller, A. Doostmohammadi and J. M. Yeomans, Active inter-cellular forces in collective cell motility, J. R.Soc., Interface, 2020, 17, 20200312.
32. R. Mueller, J. M. Yeomans and A. Doostmohammadi, Emergence of active nematic behavior in monolayers of isotropic cells, Phys. Rev. Lett., 2019, 122, 048004.
33. G. Zhang and J. M. Yeomans, Active forces in confluent cell monolayers, Phys. Rev. Lett., 2023, 130, 038202.
34. J. Löber, F. Ziebert and I. S. Aranson, Collisions of deformable cells lead to collective migration, Sci. Rep., 2015, 5, 9172.
35. D. A. Kulawiak, B. A. Camley and W.-J. Rappel, Modeling contact inhibition of locomotion of colliding cells migrating on micropatterned substrates, PLoS Comput. Biol., 2016, 12, e1005239.
36. W. Marth and A. Voigt, Collective migration under hydrodynamic interactions: a computational approach, Interface Focus, 2016, 6, 20160037.
37. Y. Cao, E. Ghabache, Y. Miao, C. Niman, H. Hakozaki, S. L. Reck-Peterson, P. N. Devreotes and W.-J. Rappel, A minimal computational model for three-dimensional cell migration, J. R. Soc., Interface, 2019, 16, 20190619.
38. B. Winkler, I. S. Aranson and F. Ziebert, Confinement and substrate topography control cell migration in a 3D computational model, Commun. Phys., 2019, 2, 82.
39. F. Thüroff, A. Goychuk, M. Reiter and E. Frey, Bridging the gap between single-cell migration and collective dynamics, eLife, 2019, 8, e46842.
40. I. Niculescu, J. Textor and R. J. de Boer, Crawling and gliding: a computational model for shape-driven cell migration, PLoS Comput. Biol., 2015, 11, e1004280.
41. W. Marth and A. Voigt, Signaling networks and cell motility: a computational approach using a phase field description, J. Math. Biol., 2014, 69, 91–112.
42. B. A. Camley, Y. Zhao, B. Li, H. Levine and W.-J. Rappel, Crawling and turning in a minimal reaction-diffusion cell motility model: Coupling cell shape and biochemistry, Phys. Rev. E, 2017, 95, 012401.
43. J. Algorta, J. P. Town, O. D. Weiner and L. Edelstein-Keshet, Investigating local negative feedback of Rac activity by mathematical models and cell motility simulations, bioRxiv, 2025, preprint DOI:10.1101/2025.05.05.651928.
44. B. Merchant, L. Edelstein-Keshet and J. J. Feng, A Rho-GTPase based model explains spontaneous collective migration of neural crest cell clusters, Dev. Biol., 2018, 444, S262–S273.
45. Y. Zhang, E. E. Bastounis and C. Copos, Emergence of multiple collective motility modes in a physical model of cell chains, npj Syst. Biol. Appl., 2025, 11, 52.
46. I. B. Bischofs, F. Klein, D. Lehnert, M. Bastmeyer and U. S. Schwarz, Filamentous network mechanics and active contractility determine cell and tissue shape, Biophys. J., 2008, 95, 3488–3496.
47. P. J. Albert and U. S. Schwarz, Dynamics of cell shape and forces on micropatterned substrates predicted by a cellular Potts model, Biophys. J., 2014, 106, 2340–2352.
48. B. J. Dubin-Thaler, G. Giannone, H.-G. Döbereiner and M. P. Sheetz, Nanometer analysis of cell spreading on matrix-coated surfaces reveals two distinct cell states and STEPs, Biophys. J., 2004, 86, 1794–1806.
49. E. L. Barnhart, K.-C. Lee, K. Keren, A. Mogilner and J. A. Theriot, An adhesion-dependent switch between mechanisms that determine motile cell shape, PLoS Biol., 2011, 9, e1001059.
50. F. Ziebert and I. S. Aranson, Effects of adhesion dynamics and substrate compliance on the shape and motility of crawling cells, PLoS One, 2013, 8, e64511.
51. J. Löber, F. Ziebert and I. S. Aranson, Modeling crawling cell movement on soft engineered substrates, Soft Matter, 2014, 10, 1365–1373.
52. M. S. Mizuhara, L. Berlyand and I. S. Aranson, Minimal model of directed cell motility on patterned substrates, Phys. Rev. E, 2017, 96, 052408.
53. C. E. Chan and D. J. Odde, Traction dynamics of filopodia on compliant substrates, Science, 2008, 322, 1687–1691.
54. N. Xia, C. K. Thodeti, T. P. Hunt, Q. Xu, M. Ho, G. M. Whitesides, R. Westervelt and D. E. Ingber, Directional control of cell motility through focal adhesion positioning and spatial control of Rac activation, FASEB J., 2008, 22, 1649–1659.
55. A. Jana, I. Nookaew, J. Singh, B. Behkam, A. T. Franco and A. S. Nain, Crosshatch nanofiber networks of tunable interfiber spacing induce plasticity in cell migration and cytoskeletal response, FASEB J., 2019, 33, 10618–10632.
56. J. Kockelkoren, H. Levine and W.-J. Rappel, Computational approach for modeling intra- and extracellular dynamics, Phys. Rev. E:Stat., Nonlinear, Soft Matter Phys., 2003, 68, 037702.
57. A. Moure and H. Gomez, Phase-field model of cellular migration: Three-dimensional simulations in fibrous networks, Comput. Methods Appl. Mech. Eng., 2017, 320, 162–197.
58. N. Ma, E. Xu, Q. Luo and G. Song, Rac1: a regulator of cell migration and a potential target for cancer therapy, Molecules, 2023, 28, 2976.
59. S. Moazzeni, K. Kyker-Snowman, R. I. Cohen, H. Wang, R. Li, D. I. Shreiber, J. D. Zahn, Z. Shi and H. Lin, N-Cadherin based adhesion and Rac1 activity regulate tension polarization in the actin cortex, Sci. Rep., 2025, 15, 4296.
60. L. Rossetti, S. Grosser, J. F. Abenza, L. Valon, P. Roca-Cusachs, R. Alert and X. Trepat, Optogenetic generation of leader cells reveals a force-velocity relation for collective cell migration, Nat. Phys., 2024, 20, 1659–1669.
61. J.-L. Maître, H. Berthoumieux, S. F. G. Krens, G. Salbreux, F. Jülicher, E. Paluch and C.-P. Heisenberg, Adhesion functions in cell sorting by mechanically coupling the cortices of adhering cells, Science, 2012, 338, 253–256.
62. S. Yamada and W. J. Nelson, Localized zones of Rho and Rac activities drive initiation and expansion of epithelial cell–cell adhesion, J. Cell Biol., 2007, 178, 517–527.
63. A. Diz-Muñoz, D. A. Fletcher and O. D. Weiner, Use the force: membrane tension as an organizer of cell shape and motility, Trends Cell Biol., 2013, 23, 47–53.
64. S. SenGupta, C. A. Parent and J. E. Bear, The principles of directed cell migration, Nat. Rev. Mol. Cell Biol., 2021, 22, 529–547.
65. W.-C. Hung, J. R. Yang, C. L. Yankaskas, B. S. Wong, P.-H. Wu, C. Pardo-Pastor, S. A. Serra, M.-J. Chiang, Z. Gu and D. Wirtz, et al., Confinement sensing and signal optimization via Piezo1/PKA and myosin II pathways, Cell Rep., 2016, 15, 1430–1441.
66. E. Bastounis, B. Álvarez-González, J. C. del Álamo, J. C. Lasheras and R. A. Firtel, Cooperative cell motility during tandem locomotion of amoeboid cells, Mol. Biol. Cell, 2016, 27, 1262–1271.
67. B. A. Camley and W.-J. Rappel, Physical models of collective cell motility: from cell to tissue, J. Phys. D: Appl. Phys., 2017, 50, 113002.
68. D. B. Brückner and C. P. Broedersz, Learning dynamical models of single and collective cell migration: a review, Rep. Prog. Phys., 2024, 87, 056601.
69. D. Caballero, R. Voituriez and D. Riveline, Protrusion fluctuations direct cell motion, Biophys. J., 2014, 107, 34–42.
70. S. Lo Vecchio, R. Thiagarajan, D. Caballero, V. Vigon, L. Navoret, R. Voituriez and D. Riveline, Collective dynamics of focal adhesions regulate direction of cell motion, Cell Syst., 2020, 10, 535–542.
71. P. Zadeh and B. A. Camley, Inferring nonlinear dynamics of cell migration, PRX Life, 2024, 2, 043020.
72. L. Bai and D. Breen, Calculating center of mass in an unbounded 2D environment, J. Graphics Tools, 2008, 13, 53–60.
73. K. R. Elder, M. Grant, N. Provatas and J. M. Kosterlitz, Sharp interface limits of phase-field models, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2001, 64, 021604.

---