Soft Matter
View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

PAPER

Cite this: Soft Matter, 2024,
20, 2739

View Journal | View Issue

Structure formation induced by non-reciprocalcell–cell interactions in a multicellular system†
Biplab Bhattacherjee,

Masayuki Hayakawa

and Tatsuo Shibata

Collective cellular behavior plays a crucial role in various biological processes, ranging fromdevelopmental morphogenesis to pathological processes such as cancer metastasis. Our previous
research has revealed that a mutant cell of Dictyostelium discoideum exhibits collective cell migration,
including chain migration and traveling band formation, driven by a unique tail-following behavior at
contact sites, which we term ‘‘contact following locomotion’’ (CFL). Here, we uncover an imbalance offorces between the front and rear cells within cell chains, leading to an additional propulsion force inthe rear cells. Drawing inspiration from this observation, we introduce a theoretical model that
Received 25th December 2023,
Accepted 10th February 2024

incorporates non-reciprocal cell–cell interactions. Our findings highlight that the non-reciprocal

DOI: 10.1039/d3sm01752d

the observed collective cell migrations. Furthermore, we present a comprehensive phase diagram,showing distinct phases at both low and intermediate cell densities. This phase diagram elucidates arsc.li/soft-matter-journal

specific regime that corresponds to the experimental system.

interaction, in conjunction with self-alignment interactions, significantly contributes to the emergence of

1 Introduction
The coordinated movement of biological cells as a group is afundamental process in several biological phenomena, including morphogenesis,1 wound healing,2 and tumor progression.3
This cooperative migration involves intricate intercellular interactions mediated by secretory chemicals4 and cell–cell contacts.
The release of chemical factors can induce chemotactic responses,
which promote directional movement of the population.5,6 Moreover, populations of cells that are bound to each other by tightintercellular adhesions, typically found in epithelial cells, achieve
coordinated migration by acquiring motility.7–12 Recent studies
have shown that transient intercellular contacts, such as contact
inhibition of locomotion (CIL)13 and contact following locomotion (CFL)14 or contact following (CF), are also crucial for inducing
cooperativity in cell populations.15,16 CIL is a phenomenon where
migrating cells avoid moving on top of other cells and migrateaway from each other upon contact.17 It has been shown that CIL
plays a key role in the formation of supracellular polarity, which
determines the orientation of movement as a group.18 On the
Laboratory for Physical Biology, RIKEN Center for Biosystems Dynamics Research,
2-2-3 Minatojima minamimachi, Chuo-ku, Kobe 650-0047, Japan.
E-mail: biplab.bhattacherjee@riken.jp, masayuki.hayakawa@riken.jp,
tatsuo.shibata@riken.jp† Electronic supplementary information (ESI) available: Movie 1: chain formation
at low density r = 0.01 with (a, k) = (0.04, 0.014). Movie 2: traveling bandformation at high density r = 0.2 with (a, k) = (0.04, 0.014). Movie 3: dynamics
of cell-pairs of KI cells at low density. See DOI: https://doi.org/10.1039/
d3sm01752d

This journal is © The Royal Society of Chemistry 2024

other hand, CFL is a tail-following behavior observed in pairs of
migrating cells, typically found in the slime mold Dictyostelium
cells.14,19–21 CFL is responsible for the formation of a stream ofcells during its development, and more recently, has been foundto be involved in a collective traveling band formation of mutant
Dictyostelium cells (KI cells) that lack chemotactic activity.14,22
Various theoretical models have been proposed to understandthe role of these transient intercellular contacts in the emergenceof collective order of migrating cells.15,16,18,23–26
The collective cell migration with various types of intercellular contacts can be modeled by the framework of selfpropelled elements. These elements move in the direction of
their own polarity by consuming internal energy, and interactions
between them induce collective eﬀects. Since these systems areout of equilibrium, the reciprocity between the elements is often
broken.27–31 Reciprocity in the context of interacting elements
illustrates the property of pair interactions between two elements
A and B, where the action on A due to the presence of B
accompanies an equal and opposite reaction on B. Nonreciprocal interactions between mesoscopic parts of a system
can often arise due to the complex environments surrounding
the elements, as well as diffusiophoretic forces between them.32
Furthermore, living cells possess their own internal degrees of
freedom, and their states can differ through their interactions
with one another. As a result, even a pair of cells can have nonreciprocal effective interaction forces. In a recent study, the front–
rear interaction symmetry was found to be broken during thecollective migration of three different types of cells: humanumbilical vein endothelial cells (HUVECs), kidney epithelial

Soft Matter, 2024, 20, 2739–2749 | 2739

View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

Paper

Fig. 1 Typical snapshot of chain formation (left) and band formation(right) observed in KI cells.14 The white scale bars in the images mark thelength scale of 50 mm (left) and 100 mm (right), respectively.

sheets (MDCK cells) and metastatic breast cancer cells.33 Nonreciprocal interactions have also been proposed for the pedestrian
dynamics.34 The statistical mechanics of the systems with nonreciprocal interactions has also been discussed.35
A mutant cell of Dictyostelium discoideum that lacks chemotactic activity to cAMP under starvation conditions, called the
KI cell, exhibits organized coordinated motions, such as traveling band formation (Fig. 1 (right)) and chain formation (Fig. 1(left)) at intermediate and low cell densities, respectively.14,22
These collective behaviors are attributable to the contact following locomotion (CFL), without relying on extracellularchemical guidance cues. In fact, a knockout mutant that fails
to express the cell–cell adhesion molecule TgrB1 exhibited
reduced CFL activity and did not form traveling bands.14 In
this study, we develop a theoretical model that captures thesecollective behaviors and try to understand the mechanisms
behind the formation of spatially asymmetric structures andspecifically the role of non-reciprocal forces in the formation ofchains and bands.
In this study, we first show an analysis of experimental data,in which we find non-reciprocal interactions for the front andthe rear cells moving as a pair. We then introduce a theoreticalmodel that includes the non-reciprocal interaction observed
experimentally. We find that the model exhibits a variety ofphenomena. At low cell density, transient chain structures with
their cell–cell interaction lifetime 4300 s are observed. Such an
asymmetric structure formation is induced by the non-reciprocalinteraction. When the adhesive forces are dominant, the chain
structures become stable over time and forms stable elongatedclusters. In the intermediate cell density, high density band
structures are formed mediated by both the adhesion force of
CFL and the self-alignment interactions. Finally, the parameter
range corresponding to the experimental system was identified as
the range of adhesion and alignment strength in which bothchains and traveling bands are formed at low and intermediate
densities, respectively.

2 Experimental observations of CFL
induced traveling band formation
In our previous experimental study,14 a mutant strain ofDictyostelium discoideum called the KI cell, which lacks chemotactic activity, was observed to exhibit collective cell migration,
forming structures such as traveling bands and chains. The

2740 | Soft Matter, 2024, 20, 2739–2749

Soft Mattertypical size se of the cells is of the order of 10 to 20 mm, and theaverage speed ve of cells is around 0.3 to 0.4 mm s1. Isolated
single KI cells exhibit a persistent random walk with a persistence time tp of about 320 s, obtained from the temporalcorrelation function of the velocity.14
Within a particular cell density regime of about 105 cells per
cm2 (intermediate density), KI cells form band-like structures
as shown in Fig. 1 (right) (see also Fig. 1(c) and (d) in ref. 14).
These structures are periodically arranged in space with aninterval of around 1 mm, equivalent to 50se to 100se. The band
structures propagate in a sea of low-density cells, perpendicular
to the direction in which they span with a speed of about
0.5 mm s1, which is faster than the cell speed of about 0.3 to
0.4 mm s1, indicating that there is turnover of cells in thetraveling band. The cells enter the band from its front and leave
it from its tail. Thus, what is propagating is the density of cells
rather than a particular group of cells. The average polar order
within the band is greater than 0.6, indicating that the band isin a polar-ordered phase.14 In contrast, the low-density region
outside the band lacks an overall direction, with an averagepolar order of less than 0.3.
The formation of such patterns should be induced by theinteractions between cells, which are short ranged. Such cell–
cell interactions can be identified in the low density situation.
In such a situation, cells are not completely isolated withoutcollective direction, rather they are associated with each otherforming a chain structure, in which cells follow the cell in frontof them as shown in Fig. 1 (left) (see also Fig. 4(a) in ref. 14).
This dynamical motion is called ‘‘contact following locomotion(CFL)’’ or ‘‘contact following (CF)’’. This chain structure is notstable, but the lifetime tl-exp of a contact between a pair of cells
is typically of the order of tl-exp = 300 s (Fig. 4(c) in ref. 14).
Therefore, we call it transient chain formation. When a cell thathas separated from other cells makes contact with another cell,
it often reforms a chain structure. Thus, the chain structure
repeatedly forms and collapses.
What kind of force between the cells in a pair maintains the
CFL? Upon observing a pair of cells, we noticed interesting
features in their motion. The rear cell exhibited faster and more
fluctuating motion compared to the front cell. Fig. 2 (top) showsa histogram of the diﬀerence between the rear cell speed vr andthe front cell speed vf. A statistical test indicates that the averagevalue of the speed diﬀerence is significantly positive, revealingthat the speed of rear cells is faster on average than the front cell.
Comparing both sides of the distribution (Fig. 2 (bottom)), we
found that the frequency of vr 4 vf was greater than in the
opposite case. This suggests that the force acting on the rear cellis eﬀectively greater than that acting on the front cell, tomaintain the CFL. We speculate that an eﬀective force is induced
upon contact that propels the rear cells in the direction of thefront cells, while the front cells are less aﬀected by the rear cells.
This force imbalance could maintain the CFL and the chain-likestructure. The specific type of motion observed during the
dynamical evolution of a cell pair depicting CFL, resulting in
non-isotropic speed distribution while maintaining the pair, isdiscussed in Appendix A.

This journal is © The Royal Society of Chemistry 2024

View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

Soft Matter

Fig. 2 (top) Histogram of the diﬀerence of rear cell speed vr and front cell
speed vf, vr  vf, in a pair of cells. Pairs with relative angle between velocity
vectors less than p/2 were considered. The total number of samples in thehistogram was 1124. The average speed of front and rear cells was 0.343(mm s1) and 0.357, respectively. The average diﬀerence in the speeds was
vr  vf ¼ 0:0140. Using a one-tailed t-test, it was found that with statisticallysignificant P-value = 0.00535, the mean value of the speed diﬀerence vr 
vf is positive, suggesting that the speed of the rear cell was statisticallyfaster than that of the front cell. bottom comparison of both sides (positive
and negative ranges) of the histogram illustrating that the frequency ofvr 4 vf was greater than in the opposite case.

3 Model of collective motion inducedby contact following locomotion
Here, we present a theoretical model aimed at elucidating therole of CFL in driving collective motion of KI cells. The cells
possess polarity that is induced by spontaneous symmetry
breaking,36–38 which causes them to move in the direction of
their polarity. This motion can be modulated by intracellular
perturbations resulting from stochastic processes within thecell. One such process includes signal reactions that lead tospontaneous polarity formation. Consequently, the direction ofcell polarity can fluctuate, resulting in the motion of cells being
random over a long time scale.37 Our previous analysis has
indicated that the correlation time of cell motion is approximately 320 seconds.14 Thus, the motion of individual cells canbe described using Langevin equations.
We consider two kinds of interactions between cells: selective contact interaction, CFL, and volume exclusion interaction.
When a pair of cells come close and move almost in the same
direction, the rear cell head is attached to the front cell tail, thereby
inducing selective contact. Conversely, when cells come close but
move in opposite directions, only mechanical collisions areobserved without CFL. The selective contact mechanism is supported by trans-membrane proteins, TgrB1 and TgrC1, expressed
during the development of Dictyostelium discoideum cells, which

This journal is © The Royal Society of Chemistry 2024

Paper
can bind with each other and mediate heterophilic interactions forinter-cellular contact.20,39,40 A previous study reported that TgrC1
operates at the rear-end of the front cell, while TgrB1 operates at
the front-end of the rear cell.20 A mutant lacking tgrb1 indeed did
not exhibit CFL activity or tail-following behavior,14 indicating therole of these proteins in CFL. Consequently, the CFL term depends
on the polarity direction of the cells in contact.
Cells possess polarity that arises from the formation ofasymmetric patterns of molecules within the cell. Once thepolarity direction is established, cells move in the direction of
their polarity. However, due to intercellular interactions, such as
CFL, the velocity direction of cells may deviate from their polarity
direction. To compensate for this deviation, cells may reorganize
their internal polarity through processes that may be induced by
the cell–cell adhesions that lead to CFL or mechanical forces
resulting from interactions between neighboring cells. These
processes give rise to a self-alignment between the directions ofcell polarity and cell velocity, promoting the coordinated motionof cells.
Based on the above considerations, we can now derive theequation of motion for an individual cell i in the system. The
system comprises N elements in a two-dimensional space withdimensions A = Lx  Ly. Cell i is located at position ri(t) = (xi(t),yi(t)) with a polarity direction yi(t). Then, the evolution equa:tions for the velocity ri = vi and polarity yi are given by:XX
ni þ avi ¼ v0 ^Caij þ m
Fij
(1)
j2Rc

j2R0


 pﬃﬃﬃﬃﬃﬃﬃ
y_ i ¼ k sin yvi  yi þ 2Dzi

(2)where cell i moves in the direction of polarity, n̂i = (sin yi, cos yi),with a constant speed v0. The second and third terms on the righthand side of eqn (1) are the CFL and volume exclusion interactions, respectively, and the first term on the right hand side ofeqn (2) is the self-alignment term. These terms are eﬀective with
the interaction between cells. With zi being a zero mean Gaussian
white noise with unit variance, isolated cells exhibit a persistent
random walk, with the correlation time given by 1/D.
For the CFL term, we construct a force Caij that is exerted oncell i by cell j when it is in contact with cell j, given by
Caij ¼^rij 
ni  ^rij ;1þ^ni  ^
nj 1 þ ^4(3)

with r̂ij = (rj  ri)/|rj  ri|. This term describes the forces induced
by selective contact, which were observed in the experiment,
where the broken reciprocity arises. This force is eﬀective
within the distance Rc = 1.9. Here, we choose the unit of length
such that the diameter of cell s is unity, i.e., s = 1. The strengthof CFL is specified by a. The first parenthesis in the expressionof Caij indicates the eﬀect of relative orientation on the force.
When the two cells in a pair are oriented in the same direction(n̂in̂j = 1), it is maximum. In contrast, when they are in the
opposite direction (n̂in̂j = 1), the force Caij vanishes. The
second parenthesis is the tail-following term. When the polarityof cell i, n̂i, and the unit vector from cells i to j, r̂ij, have the sameorientation, i.e., cell i is oriented to the direction of cell j, theSoft Matter, 2024, 20, 2739–2749 | 2741

View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

Papermaximum following force is induced in cell i in the direction ofcell j. Such a situation is expected when cell i is located in the
tail of cell j and oriented to cell j. In contrast, when cell i is thefront cell and is oriented towards the opposite direction from
cells i to j, the force Caij vanishes.
For the volume exclusion interaction, we consider a soft
repulsive core with the assumption that the cells mostly have acircular shape with a diameter of s. The corresponding force
acting on any cell i in the vicinity of another cell j is given by
Fij = riU(rij), obtained from the WCA potential,41
8 n


 o
< 4e s=rij 12  s=rij 6 ; if rij  R0Uðrij Þ ¼
(4):0;
otherwisewith rij = |rij| being the distance between the centers of the cells.
In the third term on the right hand side of eqn (1), theparameter m is the strength of this repulsion force. Throughout
this paper, we use e = 1 and s = 1, which set the units of energy
and length, respectively, and m = 1. The cut-oﬀ R0 = 21/6 sets theupper bound of the repulsion force.
The first term on the right hand side of eqn (2) gives the selfalignment42,43 eﬀect between the polarity direction yi and thevelocity direction yvi . The velocity direction is defined by vi/|vi| =(sin yvi , cos yvi ). The coeﬃcient k determines the strength of theinternal force that reorganizes the cell polarity to self-align withthe velocity direction, and sin(yvi  yi) corresponds to the firstterm in the Fourier series, satisfying the symmetry with respect
to yi - yi + 2p, yvi - yvi + 2p and (yi, yvi ) - (yi, yvi ).
The density r = (N/A)(p/4)s2 plays a crucial role in theexperiment. In this paper, we consider two density values: r =0.01 and r = 0.2, corresponding to low and intermediate
densities in the experiment, respectively. The diﬀusivity D =0.003125 s1 is derived from the observed persistence timescale of a single cell in the experiment, s2/D = 1/D = 320 s. Basedon the values of ve and se, KI cells take around 30 to 60 seconds
to travel across the cell. Thus, the speed v0 is estimated to be
between 1/60 = 0.016 s1 and 1/30 = 0.033 s1 to becomecomparable to the experimental situation. An average speedof v0 = 0.025 s1 is taken for the simulation. The dimensionlessPéclet number Pe = v0s/D = 8 sets the activity of the cells in thissystem. We use the Euler integration scheme with an integration time step of dt = 103ts to perform the numerical simulations. We take N = 1024 for the low-density analysis, N = 512 toshow configurations in Fig. 3(F) (for better visualization), andN = 4096 for the intermediate-density analysis. The parameters
we will explore in this paper are the strengths of the adhesioninteraction a and the self-alignment interaction k. The dimension (Lx  Ly) of the system is determined by r and N. For the
low and intermediate densities, (Lx  Ly)(r,N) is given by (304.74 
263.91)(0.01,1024) and (136.284  118.025)(0.2,4096), respectively.

4 Results
Before comparing the simulation behavior of the present modelwith the experimental observation, we explore the phenomena
observed in the model. For the two diﬀerent densities (r = 0.01

2742 | Soft Matter, 2024, 20, 2739–2749

Soft Matterand r = 0.2), by changing a and k, this model exhibits a variety
of dynamic pattern formation. To characterize individual
phases, we use the polar order parameter and the mean
eccentricity of cell clusters. The polar order parameter hPi isdefined as

+
*
N

1 X

~
hPi ¼
vi  ;(5)
Nv0  i¼1 
where hi indicates an average taken over multiple steady state
configurations.
A cluster of cells is defined as a set of cells which areinterconnected via the interactions. Here, any two cells which
are interacting with each other, i.e., within a distance less thanRc, are considered to be connected. Clusters characterized by amonolayer arrangement of cells, or at most a bi-layered structure, which are elongated along their direction of motion, arereferred to as chains. To identify spatial symmetry breaking inthe structures such as chains and clusters, which are mostly
elongated along their direction of motion, we introduce themean eccentricity at any phase point (a, k) in the phase space asvﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃu" 2 #ubieða; kÞ ¼ t 1 ; where bi and ai are the semi-minor
ai iand semi-major axis of the cluster i, calculated from its Gyration tensor. The averaging hii is taken over all such clusters i, inaround 5000 well separated steady state configurations. An
elongated structure can be identified by its e. For symmetricstructures, e remains close to zero. As the structure elongates, eshifts toward unity.
4.1

Chain formation at low density

We first study the case when the density is low (r = 0.01). In thesmall (a, k) regime, the system resides in the isotropic phase(Fig. 3(A) and (F)-(i)). In this phase, the polar order parameterhPi is close to zero. e is also observed to be close to zero. The
typical average lifetime tl of a pair of cells (hereafter called pairlifetime) in this phase is quite small (few tens of seconds). Thus,
the cells move in random directions and the cell distribution is
spatially homogeneous without a specific spatial pattern. With
increasing (a, k), tl keeps on increasing. As it reaches at a = a1(k),marked by the red line with pentagons in Fig. 3(A), tl becomeslarger than tl-exp = 300 s and a transient chain formation phaseemerges as shown in Fig. 3(F)-(ii) (see also Movie S1, ESI†). The
formed chain structures move along their direction of elongation (see Appendix B). Within the range of chain formation, we
can observe the growth of e beyond 0.4. The polar orderparameter of the entire system remains close to zero. This phaseis considered to explain the formation of transient chain-likestructures observed at low density in the experiment.
With the further increase of (a, k) at a = a2(k) marked by theblue line with diamond in Fig. 3(A), the elongated structures
stabilize over time with e 4 0.6 (elongated cluster phase),
depicting the formation of infinitely long lived clusters(tl - N), elongated mostly along their direction of motion(Fig. 3(F)-(iii)). The whole system of cells then forms a very few

This journal is © The Royal Society of Chemistry 2024

View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

Soft Matter

Paper

Fig. 3 Low density (r = 0.01): (A) phase diagram in the a–k plane with the mean eccentricity e(a, k) color coded. The red pentagons depicts a1(k) for
which tl(a, k) B tl-exp, identifying the onset of the transient chain formation phase. The green filled triangles denote the order–disorder transition points inthe (a, k) phase plane and consequently identify the onset of the stable elongated cluster phase. The solid lines connecting the points depict thecorresponding phase boundaries. (B) The
order parameter P(a, k) is shown for different values of k as a function of a. A disorder–order transition is
ﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
ppolarobserved and the maxima in DPða; kÞ ¼ hP2 i  hPi2 (errorbar) identify the transition points ac(k). (C) The data collapse of the hP(a, k)i curves are shownwith the transformation: a - [a  ac(k)]k2.0. The transition point ac(k) seems to decay as a power law, ac(k) B k1.2 (inset). (D) The order–disordertransition is associated with the pair-lifetime tl(a, k) becoming infinitely large. The variation of the tl(a, k) vs. a is shown for different values of k. (E) The datacollapse of tl(a, k) vs. a curves for different k values is shown for the transformation: a - ak0.9 and tl(a, k) - tl(a, k)k0.2. The growing exponential
Bexp(4200ak0.9) fits well with the data collapse. (F) Typical configurations in the three different phases are shown (isotropic phase in (i), transient chain
formation phase in (ii) and stable elongated cluster phase in (iii)) for k = 0.3. Color indicates the direction of motion.

polar clusters with overall polarisation hPi becoming close to 1.
Thus, a disorder–order transition happens as we go from thetransient chain formation phase to the elongated cluster phase.
The variation of the polar order hPi across the transition isshown in Fig. 3(B). The critical transition points ac(k) (green
line with filled triangles in Fig. 3(A)) are identified by the
maxima in the fluctuation of the order parameter as DPða; kÞ ¼
pﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
hPða; kÞ2 i  hPða; kÞi2 (error bars in Fig. 3(B)). The hP(a, k)icurves show data collapse by the transformations, a - [a 
ac(k)]k2.0 (Fig. 3(C)), with ac(k) following a power law decay of
the form ac(k) B k1.2 (Fig. 3(C) inset). Data collapse of theorder parameter growth curve indicates that the growth of orderparameter follows the same law close to the transition within
the range of a and k we studied and for a given N.44
Note that the disorder–order transition line indicated by ac(k)and the boundary line a2(k) marking the onset of long-lived pairs

This journal is © The Royal Society of Chemistry 2024

overlap closely. This suggests that these two events are notindependent, rather the formation of the ordered phase stabilizes
pairs, making their pair-lifetime infinitely large. The typical
growth of the tl(a, k) as a function of a are shown in Fig. 3(D)for different values of k. It shows an initial exponential growthfollowed by a saturation due to the finiteness of the simulation
time. A data collapse with the transformation tl(a, k) - tl(a, k)k0.2
and a - ak0.9 is observed (Fig. 3(E)). The collapsed data fits an
exponentially growing function of the form Bexp(Aak0.9). This
indicates an exponentially fast divergence to the long living pairformation phase near the order–disorder transition.
4.2

Traveling band formation at intermediate density

We next study the case when the density is intermediate (r =0.2). We observe four diﬀerent phases, i.e., isotropic phase,microscopic cluster phase, macro polar cluster phase and band

Soft Matter, 2024, 20, 2739–2749 | 2743

View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

Paper
formation phase. The phase diagram with the observed phases
is shown in Fig. 4A. Each of these phases is identified by theirunique spatial distribution of the cells and is quantified bytheir density distributions shown in Fig. 4C. To observe the
distribution of the local density, the whole system is divided
into 20  20 square grids of dimensions lx  ly = 6.8  5.9. The
local densities of any grid i is estimated as the number of cells(nil) per unit area, i.e. ril = nil/(lx  ly). The distribution Dp(rl) isthen calculated over all such spatial grids and over a minimumof 50 000 well separated steady state configurations.
When both a and k are suﬃciently small, the cellular
dynamics is dominated by diﬀusion and hence the systemstays in the isotropic phase (black triangles in Fig. 4A), depicting a homogeneous spatial distribution of cells (Fig. 4B(a)). The
velocity vectors of the cells are observed to be pointing inSoft Matterrandom directions, resulting hPi to be close to zero. The distribution of the local cell density Dp(r) is found to fit a Gaussianwith its peak at the mean density (r = 0.2) of the system(Fig. 4C(a)).
With the increase of a, keeping k to be small, we observe a
microscopic cluster phase (green circles in Fig. 4A), where thesystem consist of microscopic polar clusters with the lengthscales of these clusters much smaller than that of the system. A
typical configuration is shown in Fig. 4B(b), with the polar
direction marked with color-code. The value of the polar orderparameter hPi in this phase is close to 1/2. The associateddensity distribution is shown in Fig. 4C(b). The density distribution Dp(rl) can be fitted to a left truncated Gaussian
distribution. As shown in Fig. 4B(b), the spatial distributionof cells is non-homogeneous with vacant spaces (zero density).

Fig. 4 Intermediate density (r = 0.20): (A) phase diagram in the a–k plane, at r = 0.2, depicts the isotropic phase (black triangle), microscopic cluster
phase (green circle), macro polar cluster phase (blue half-filled circle) and the band formation phase (red square). Typical configurations of the different
phases (at phase points marked by open black circles) are shown in (B). The density distribution of the phases in (B) is depicted in (C). The corresponding
phase points (a, k) are mentioned on top of each plot. The red line in (A) marks the onset of longer pair lifetime tl 4 300 s at low density. The red squares
above the red line in (A) identify the experimental regime. Four points are marked as blue diamonds in (A), within the experimental regime at which thecorresponding low (r = 0.1) and intermediate density (r = 0.2) configurations are shown in Fig. 5. (C) Density distribution Dp(rl) of (a) isotropic,(b) microscopic cluster, (c) macro polar cluster and (d) band formation phases identified by (a) Gaussian distribution with a peak close to the mean density(r = 0.2), (b) left truncated Gaussian distribution, (c) bimodal distribution with one peak at a very high density and another peak at zero identifies, and(d) bimodal distribution with both the peaks at non-zero densities, respectively. Color in (B) indicates the direction of motion.

2744 | Soft Matter, 2024, 20, 2739–2749

This journal is © The Royal Society of Chemistry 2024

View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

Soft Matter
When both a and k are large, the system shows a compact
macro polar cluster phase (blue half filled circle in Fig. 4A),
with the length scale corresponding to the size of the cluster
equivalent to the order of the system size. A typical configuration is shown in Fig. 4B(c). The density distribution Dp(rl) of
this phase (Fig. 4C(c)) is a bimodal distribution, with one peakat high density corresponding to the large high density polar
cluster and another peak at zero, corresponding to the larger
vacant spaces. This phase is a completely ordered phase with
hPi observed to be close to 1.0. In this phase, the cells move
mostly as a connected cluster. However, with evolution of time,
the shape of the cluster and the direction of propagation keep
on changing, seemingly independently.
For suﬃciently large self-alignment interaction k, with comparatively smaller adhesion strength a or even no adhesion,
band formation arises. Here, the term ‘‘band’’ corresponds to amoving high density patch propagating perpendicular to its
length within a low density sea of cells. A regime of bandformation phase is depicted in Fig. 4A (red square). A typicalband structure (Fig. 4B(d)) moves perpendicular to its spanning
direction with a sharp front and a comparatively fluctuating tail(Movie S2, ESI†). The density distribution Dp(rl) of a bandformation phase is shown in Fig. 4C(d). The distribution Dp(rl)has a peak at low density corresponding to the low density ofsea of cells, and a high density shallow peak corresponding to
the band. Compared to the macro polar cluster phase, in thisphase, the low density sea of cells is dominating compared to
the large vacant spaces available in the macro polar cluster
phase. hPi is observed to be of the order of 1/2.
Note that, the phase points not identifiable by their statistical characteristics within our observation time are marked asasterisks.

Paper
4.3 Parameter region that explains both transient chain andtraveling band formations
In the experiment, the band formation and the transient chain
formation has been observed by just changing the cellulardensity.14 Such a parameter region can be identified by observing
the overlap between the transient chain formation phase at
r = 0.01 and the band formation phase at r = 0.2. In Fig. 4A,
the red line indicates the onset of transient chain formation
phase in the low density (r = 0.01). The red squares above this redline identify the experimental phase that depicts the chain
formation at r = 0.01 and the band formation at r = 0.2.45 Within
this phase, typical spatial configurations of cells at four diﬀerent
parameter sets (marked by blue diamonds in Fig. 4A) are shownat both r = 0.2 (Fig. 5(i)) and at r = 0.01 (Fig. 5(ii)). In thissubsection, we show that the characteristics of the traveling band
observed experimentally can be reproduced in our simulations.
The spatial profile of the local density along the direction ofband propagation rl(x) indicates that the high density band has a
sharp front and a decaying tail (Fig. 6(a)). We first measured theband propagation speed vb as the rate of displacement of the sharphigh density front. The mean speed of the band is higher than thatof the cellular speed (v0 = 0.025) (Fig. 6(b)), consistent with theexperimental observation. Actually, in the simulation, the cells
enter the band from its front and leave it from its tail. The band
speed vb depends on the parameter value. As shown in Fig. 6(b),the band speed vb increases as the CFL strength a increases.
From Fig. 6(a), the bandwidth w is estimated to be aroundw = 29.4, which corresponds to about 300 mm to 600 mm (forse = 10 mm and 20 mm, respectively), which is consistent withthe experimental observation.14 For the cell size (se = 12 mm)
used for the simulations, the average bandwidth corresponds
to around 350 mm. From Fig. 6(b), vb is between 0.03 andFig. 5 Typical configurations at (i) intermediate density (r = 0.2) depicting band formation and at (ii) low density (r = 0.01) depicting transient chain
formation are shown at four phase points (marked as blue diamonds in Fig. 4A) (a) (0.04, 0.014), (b) (0.01, 0.04), (c) (0.006, 0.08) and (d) (0.00081, 0.8)
identified by the set of parameters (a, k) within the phase corresponding to the experimental system. Color indicates the direction of motion depicted in
the color-wheel. Note that the amplitude of the velocity vectors in (ii) are multiplied by 1.5 for better visibility.

This journal is © The Royal Society of Chemistry 2024

Soft Matter, 2024, 20, 2739–2749 | 2745

View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

Paper

Soft Matter

Fig. 6 (a) A typical spatial profile of the local density rl(x) in space along the direction of band propagation shown at (a, k) = (0.0006, 0.06). To calculate
rl(x), we divided the system into 30 strips of dimension (lbin  Ly) = (4.54s  118.03s) along the x-axis. The cyan dashed line marks the hrli = 0.3 line andthe solid blue arrow identifies the bandwidth, which is of about B29s. Note that the local density is above 0.3 within the band. (b) The variation of theband propagation speed vb(a, k) against a is shown (blue circle) for a fixed value of k = 0.06. The cellular speed v0 = 0.025 is shown as the red dashed line,
for comparison. The variation fits a linear equation of the form vb(a) B 1.057a + 0.028 (blue line). (c) The average local polarisation, hPl(rl)i, of cells as afunction of local density is shown for a fixed k = 0.06. Some a values are chosen within the band formation phase at a = 0.0001 (cyan pentagon), 0.0002(blue diamond), 0.0006 (green triangle), 0.001 (pink circle) and 0.002 (red square). The values a = 0.004 (brown cross) and a = 0.01 (black asterisk)correspond to the regime just below and above the experimental regime, respectively.

0.04 (s1), which corresponds to about 0.36 mm s1 to 0.48 mm s1.
Thus, the typical time interval up to which cells stay within bands
is given by w/(vb  v0) and thus ranges within 2000 to 5000 s. In
the experiment, the bandwidth is observed to be in the range from
200 to 700 mm and the band speed is observed to be approximately
0.5 mm s1. Thus the time spent by cells within the band is within
the range of 1000 to 3500 s. Hence the observed average durationof the cell within the band in the simulations is similar to theexperimentally observed result.14
We next study the polar order inside and outside the bands.
To identify the band area, we introduce a local density threshold of rl = 0.3, above which areas are considered as the band
region. As there are density fluctuations at the tail of the band,the threshold is taken to be 0.3, which is greater than thesystem density 0.2. In Fig. 6(c) we show the average of the local
order parameter Pl(rl) as a function of the local density for k =0.06. We find that the band area with rl greater than rl = 0.3
exhibits the local polar order parameter value larger than 0.63.
In contrast, the area outside of the band, whose density is close
to rl = 0.1 as found in Fig. 6(a), exhibits the local polar orderparameter which is close to 0.3. Thus, the traveling band is thepolar ordered phase, which phase-separates with the disordered background. These results in the simulation is consistent
with the experimental observation.14
Why is the local polar order parameter in the band area notclose to 1, but distributed around 0.63? As observed in ref. 14, it was
found that the migration direction was widely distributed in therange of about 60 degrees with the standard deviation to be about
30 degrees. Such heterogeneity in the direction of cell motion within
the band can also be seen in the present model. In Fig. 7(a) and (b),we show two typical configuration depicting band formation. As wezoom in a certain part of the bands, shown in Fig. 7(c) and (d), we
observe that there are some cells moving in directions which arealmost 90 degrees separated from the other cells within the band.
The standard deviation in the distribution of cell orientation within
the band is about 16 degrees. Because of this heterogeneity, thelocal order parameter in the ordered traveling band is much lessthan 1. However, the origin of such large heterogeneity in theordered traveling band remains to be elusive.

5 Discussion

Fig. 7 Formation of bands at the phase point (a, k) = (a) (0.01, 0.04) and(b) (0.04, 0.014). Both the bands move along the negative x-axis. The areas
marked as a square box in (a) and (b) are magnified in (c) and (d),respectively. Color indicates the direction of motion depicted in thecolor-wheel.

2746 | Soft Matter, 2024, 20, 2739–2749

In this article, we presented a theoretical model to reproduce thecollective behavior of KI cells, a mutant of Dictyostelium cell thatlacks chemotaxis activities.14,22 The KI cells are known to exhibit
distinctive patterns of aggregation, forming chain and traveling band
structures in low and intermediate cellular densities, respectively(Fig. 1).14 Our model successfully reproduced these phenomena andshed light on the underlying mechanisms driving these patternformations.

This journal is © The Royal Society of Chemistry 2024

View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

Soft Matter

Paper

We first showed that when cells are moving in a chain, there is
an imbalance of force between the front and rear cells (Fig. 2).
This observation prompted us to introduce a ‘‘Contact Following
Locomotion’’ (CFL) term into our model. The CFL term induced
non-reciprocal interaction between cells, creating an additional
propulsion force in the rear cells, thereby facilitating the formation of chain structures during migration. Moreover, themodel demonstrated a wide range of behaviors, including bothchain and traveling band formations (Fig. 3 and 4). We note thatin the previous model,14 non-reciprocal intercellular interactions
were considered to describe contact following in the evolution
equation for the cell polarity, corresponding to a signallingmechanism induced by contact following. Since the speed wasfound to depend on the cell position relative to the other cells in apair, in the present model, the non-reciprocity induced by thecell–cell contact is implemented as a mechanical force in theevolution equation for the velocity.
Our results suggest that the non-reciprocal eﬀect in the CFL
term primarily contributes to the formation of asymmetricalstructures, such as chains. In contrast, the self-alignment eﬀect
might be suﬃcient for the emergence of traveling bands (Fig. 8).
However, it is essential to recognize that the parameters governing the CFL term and the self-alignment term may not be entirely
independent in real cells. In fact, self-alignment of cell polarity
direction may occur through the cell–cell adhesion responsible
for the CFL.20 Therefore, it is plausible that in actual cells, theparameters of the CFL term and the self-alignment term areinterdependent and cannot be modified independently.
In the tgrb1 mutant, which lacks the expression of theadhesion protein TgrB1, the formation of chains is severely
impaired, and traveling bands do not form.14 This suggests thatthe cells were in the state of an isotropic phase due to the lossof both eﬀects of CFL and self-alignment terms.
While we have identified parameter regions that correspondto the experimental observations, the parameter space allowing

both chain and traveling band formations is still wide. Specifically, in regions with low values of the parameter k, the selfalignment eﬀect is weak, yet the CFL term induced band formation, emphasizing the importance of the CFL term in promoting this behavior. Future research should aim to pinpoint the
precise parameter values that best explain experimental results.
Our model demonstrates that the tail-following behavior
induced by the CFL term and the self-alignment eﬀect work inconcert to not only facilitate the formation of chain and travelingband structures but also to produce a variety of other collective
behaviors. These two eﬀects are also present in a variety ofphenomena, including flocking in birds46,47 or schooling infish,48 although their physical origins may diﬀer. Therefore, thecollective behaviors induced by the combination of these eﬀects
may be universal in a wide range of active matter systems.

Fig. 8 Typical configuration showing band formation at r = 0.2, when theadhesive force a is set to zero. The k is fixed at 0.05. The color wheel
corresponds to the polarity direction of individual cells.

Fig. 9 Three types of motions of a cellular pair, i.e. front-rear dynamics,
intermittent sliding dynamics and side-by-side motion are schematicallyshown.

This journal is © The Royal Society of Chemistry 2024

Author contributions
The project was conceptualised by T. S. The experiments were
designed by M. H. and T. S., and implemented by M. H. The
theory was designed by B. B. and T. S. The simulation codes
were developed by B. B. Data analysis was performed by B. B.
and T. S. The draft was written by B. B. and T. S. and criticallyexamined by M. H.

Conflicts of interest
There are no conflicts to declare.

Appendices
Appendix A:

Speed asymmetry of front and rear cells

According to our observation of the movement of a pair of KI cells,
the arrangement of the two cells can be broadly divided into two
diﬀerent types, (i) front-rear and (ii) side-by-side. In the later case,to maintain the pair, the cells need to move with a similar speed.
Thus the side-by-side arrangement (Fig. 9 (bottom)) has no contribution in the speed asymmetry. In the front-rear arrangement,

Soft Matter, 2024, 20, 2739–2749 | 2747

View Article Online

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

Paper

Soft Matter

however, the speed symmetry is generally broken. During CFL, thefront–rear arrangement is observed to be prevalent.14 During theevolution of a cell-pair in a front–rear arrangement, the head of therear cell is generally adhered to the tail of the front cell. However,intermittently, the head of the rear cell moves along the lateralsurface of the front cell. Though the cells are still in contact, thetail–head contact may or may not be disrupted. This we term as
the intermittent sliding (Fig. 9 (middle)).
Now we try to understand the mechanism which gives rise tounequal average speeds of cells moving together in a front–reararrangement. If we follow the trajectory of a typical pair (shown in
Fig. 10) we observe that the rear cell trajectory has larger fluctuations compared to the front cell trajectory. Consider the velocitiesof the front and the rear cell to be vf and vr respectively. We candecompose the velocities in parallel and perpendicular components as (vf8, vf>) and (vr8, vr>) w.r.t. the tangent of the trajectory,t̂. Since the motion of the front-cell defines the trajectory direction, the direction of its motion is always along t̂, implyingvf> = 0. Thus for the front cell in a pair, at every instance,

and using eqn (8) we re-write eqn (9) as(D
E ¼ hvf i; if vr? ¼ 0
2
2 1=2
hvr i ¼ vfk þ vr?
4 hvf i; if vr? 4 0(10)Thus, in the front–rear arrangement, the intermittent sliding
mechanism allows the rear cell to have a larger average speed
compared to the front cell.
In Fig. 10, we show the typical trajectories of the front cell(red) and the rear (blue) cell during their evolution as a pair.
The starting point of the trajectory is marked by the asterisk.
Note that, at the initial part of their dynamics as a pair, they
were in side-by-side arrangement, and the blue and red lines
almost follow each other. After certain time, the blue trajectory
falls behind the red, as the cellular pair keeps on moving as afront–rear pair. After the transformation to the front–rearsituation, the rear cell (blue) shows larger fluctuations in its
trajectory. The same dynamics can be observed in Movie S3,
ESI.†

vf (t) = vf8(t),(6)

Appendix B:

hvfi = hvf8i,(7)

To characterize the direction of the elongated structures and/or
chains, we first obtain the eigen direction yei
i of cluster i, whichcorresponds to the largest eigenvalue of the gyration tensor of thecluster. A typical cluster is shown in Fig. 11 (left). The average velocity
2P
3
sin yvj
j2i
6
7
direction of such a cluster is given by yvi ¼ tan1 4P
5. Thuscos yvj

where the averaging is taken over time along the whole trajectory. Now for a pair of cells moving together, the average speedof the cells along t̂ has to be the same. Thus,
hvf8i = hvr8i.

(8)chain elongation

j2i

However for the rear cell, the perpendicular component is nonzero, which we observe from the intermittent sliding dynamics.
Thus for the rear cell, we can write the average speed as,
hvri = h[vr82 + vr>2]1/2i(9)Fig. 10 The trajectories of the front cell (red) and the rear cell (blue)
during the evolution of a cellular pair, with their starting points marked withasterisk. The trajectory corresponding to the rear cell depicts larger
fluctuations. Note that the dots correspond to time points separated by
15 seconds, which is the frame rate of the experimental observation.

2748 | Soft Matter, 2024, 20, 2739–2749

v
we define Dy = h|cos dyi |ii, where dyi = (yeii  yi ). The averaging hii
is taken over all such clusters i in the same configuration andmany (B5000) such well separated configurations in the steadystate. In Fig. 11 (right) we show the phase diagram in the low
density (r = 0.01) with Dy color coded. When the cluster ismoving along the direction of its elongation, di is close to zero(or p), correspondingly Dy is close to unity. On the other hand if
the cluster moves perpendicular to its elongation di is close top/2 making Dy to be close to zero. Note that, within the chain
formation region and beyond in the phase diagram, Dy becomes
greater than 0.7, implying that the angle di becomes less thanp/4, identifying structures elongated along their direction of
motion.

Fig. 11 (left) A typical cluster is shown. The eigen direction ~yei (blue arrow)corresponding to the largest eigenvalue of the gyration tensor and thevelocity direction ~yv (red arrow) of the cluster are depicted. (right) The a–kphase diagram in the low density r = 0.01, with the cos Dy color coded also
shown.

This journal is © The Royal Society of Chemistry 2024

View Article Online

Soft Matter

Acknowledgements

Open Access Article. Published on 04 March 2024. Downloaded on 10/20/2024 5:36:54 AM.
This article is licensed under a Creative Commons Attribution 3.0 Unported Licence.

We acknowledge T. Ishibashi, R. Nishizawa, G. Ogita, M.
Tarama, S. Tarama and V. Semwal for simulating discussions
and important insights regarding the project. This was was
supported by the core funding at RIKEN Center for Biosystems
Dynamics Research, Kakenhi grants 22H05170 (TS) and
19K12770 (MH), and JST CREST Grant JPMJCR1852 (TS).

Notes and references
1 M. Dicko, P. Saramito, G. B. Blanchard, C. M. Lye, B. Sanson
and J. Étienne, PLoS Comput. Biol., 2017, 13, 1–30.
2 A. Grada, M. Otero-Vinas, F. Prieto-Castrillo, Z. Obagi and
V. Falanga, J. Invest. Dermatol., 2017, 137, e11–e16.
3 P. Friedl, J. Locker, E. Sahai and J. E. Segall, Nat. Cell Biol.,2012, 14, 777–783.
4 D. Dormann and C. J. Weijer, Curr. Opin. Genet. Dev., 2003,
13, 358–364.
5 A. Haeger, K. Wolf, M. M. Zegers and P. Friedl, Trends Cell
Biol., 2015, 25, 556–566.
6 M. Chuai and C. J. Weijer, Curr. Opin. Genet. Dev., 2009, 19,343–349.
7 E. Theveneau and R. Mayor, Cell. Mol. Life Sci., 2013, 70,
3481–3492.
8 P. Friedl and R. Mayor, Cold Spring Harbor Perspect. Biol.,2017, 9, a029199.
9 E. Scarpa and R. Mayor, J. Cell Biol., 2016, 212, 143–155.
10 K. Sato, T. Hiraiwa, E. Maekawa, A. Isomura, T. Shibata andE. Kuranaga, Nat. Commun., 2015, 6, 10074.
11 K. Sato, T. Hiraiwa and T. Shibata, Phys. Rev. Lett., 2015,
115, 188102.
12 T. Yamamoto, T. Hiraiwa and T. Shibata, Phys. Rev. Res.,
2020, 2, 043326.
13 B. Stramer and R. Mayor, Nat. Rev. Mol. Cell Biol., 2017, 18,
43–55.
14 M. Hayakawa, T. Hiraiwa, Y. Wada, H. Kuwayama andT. Shibata, eLife, 2020, 9, e53609.
15 T. Hiraiwa, Phys. Rev. Lett., 2020, 125, 268104.
16 T. Hiraiwa, Phys. Rev. E, 2019, 99, 012614.
17 C. Carmona-Fontaine, H. K. Matthews, S. Kuriyama, M. Moreno,
G. A. Dunn, M. Parsons, C. D. Stern and R. Mayor, Nature, 2008,
456, 957–961.
18 E. Theveneau, L. Marchant, S. Kuriyama, M. Gull, B. Moepps,
M. Parsons and R. Mayor, Dev. Cell, 2010, 19, 39–53.
19 C. J. Weijer, J. Cell Sci., 2009, 122, 3215–3223.
20 T. Fujimori, A. Nakajima, N. Shimada and S. Sawai, Proc.
Natl. Acad. Sci. U. S. A., 2019, 116, 4291–4296.
21 D. Li and Y.-l Wang, Proc. Natl. Acad. Sci. U. S. A., 2018, 115,
10678–10683.
22 H. Kuwayama and S. Ishida, Sci. Rep., 2013, 3, 2272.

This journal is © The Royal Society of Chemistry 2024Paper
23 B. A. Camley, Y. Zhang, Y. Zhao, B. Li, E. Ben-Jacob,
H. Levine and W.-J. Rappel, Proc. Natl. Acad. Sci. U. S. A.,
2014, 111, 14770–14775.
24 B. A. Camley, J. Zimmermann, H. Levine and W.-J. Rappel,
Phys. Rev. Lett., 2016, 116, 098101.
25 D. A. Kulawiak, B. A. Camley and W.-J. Rappel, PLoS Comput.
Biol., 2016, 12, e1005239.
26 B. Smeets, R. Alert, J. Pešek, I. Pagonabarraga, H. Ramonand R. Vincent, Proc. Natl. Acad. Sci. U. S. A., 2016, 113,
14621–14626.
27 J. Dzubiella, H. Löwen and C. N. Likos, Phys. Rev. Lett., 2003,
91, 248301.
28 E. A. Lisin, O. F. Petrov, E. A. Sametov, O. S. Vaulina,
K. B. Statsenko, M. M. Vasiliev, J. Carmona-Reyes andT. W. Hyde, Sci. Rep., 2020, 10, 13653.
29 B. Sabass and U. Seifert, Phys. Rev. Lett., 2010, 105, 218103.
30 I. Sriram and E. M. Furst, Soft Matter, 2012, 8, 3335–3341.
31 F. Jiménez-Ángeles, K. J. Harmon, T. D. Nguyen, P. Fenterand M. Olvera de la Cruz, Phys. Rev. Res., 2020, 2, 043244.
32 R. Soto and R. Golestanian, Phys. Rev. Lett., 2014, 112, 068301.
33 J. LaChance, K. Suh, J. Clausen and D. J. Cohen, PLoS
Comput. Biol., 2022, 18, 1–28.
34 D. Helbing and P. Molnár, Phys. Rev. E: Stat. Phys., Plasmas,
Fluids, Relat. Interdiscip. Top., 1995, 51, 4282–4286.
35 A. V. Ivlev, J. Bartnick, M. Heinen, C.-R. Du, V. Nosenko and
H. Löwen, Phys. Rev. X, 2015, 5, 011035.
36 M. Nishikawa, M. Hörning, M. Ueda and T. Shibata, Biophys. J., 2014, 106, 723–734.
37 T. Hiraiwa, A. Nagamatsu, N. Akuzawa, M. Nishikawa andT. Shibata, Phys. Biol., 2014, 11, 056002.
38 Y. Arai, T. Shibata, S. Matsuoka, M. J. Sato, T. Yanagida andM. Ueda, Proc. Natl. Acad. Sci. U. S. A., 2010, 107, 12399–12404.
39 S. Hirose, R. Benabentos, H.-I. Ho, A. Kuspa and G. Shaulsky,
Science, 2011, 333, 467–470.
40 S. Hirose, B. Santhanam, M. Katoh-Kurosawa, G. Shaulsky
and A. Kuspa, Development, 2015, 142, 3561–3570.
41 J. D. Weeks, D. Chandler and H. C. Andersen, J. Chem. Phys.,1971, 54, 5237–5247.
42 C. A. Weber, T. Hanke, J. Deseigne, S. Léonard, O. Dauchot,
E. Frey and H. Chaté, Phys. Rev. Lett., 2013, 110, 208001.
43 O. Dauchot and V. Démery, Phys. Rev. Lett., 2019, 122,068002.
44 The analysis presented here does not specify the nature of
the transition and a detailed analysis is underway.
45 The animations of chain and band formation are provided
for two such phase points in the ESI†.
46 M. Nagy, Z. Ákos, D. Biro and T. Vicsek, Nature, 2010, 464,
890–893.
47 M. Yomosa, T. Mizuguchi and Y. Hayakawa, PLoS One, 2013,8, 1–6.
48 D. Bumann and J. Krause, Behaviour, 1993, 125, 189–198.

Soft Matter, 2024, 20, 2739–2749 | 2749