Is Cell segregation just like oil and water?

Florian Franke1,2*, Sebastian Aland3,2, Hans-Joachim Böhme1,2, Anja Voss-
Böhme1,2, Steffen Lange1,2

1 DataMedAssist, HTW Dresden, 01069 Dresden, Germany

2 Faculty of Informatics/Mathematics, HTW Dresden - University of Applied
Sciences, 01069 Dresden

3 Faculty of Mathematics and Computer Science, TU Freiberg, 09599 Freiberg

*<EMAIL_ADDRESS>
## Abstract

Understanding the segregation of cells is crucial to answer questions about
tissue formation in embryos or tumor progression. Steinberg proposed that
separation of cells can be compared to the separation of two liquids. Such a
separation is well described by the Cahn-Hilliard (CH) equations and the
segregation indices exhibit an algebraic decay with exponent $1/3$ with
respect to time. However, experimental cell segregation data also reveals
other scaling exponents and even slow logarithmic scaling laws. These
discrepancies are commonly attributed to the effects of collective motion or
velocity-dependent interactions. By using a cellular automaton (CA) model
which implements a dynamic variant of the differential adhesion hypothesis, we
demonstrate that it is possible to reproduce biological cell segregation
experiments from Méhes et al. with just adhesive forces. The segregation in
the CA model follows a logarithmic scaling, which is in contrast to the
proposed algebraic scaling with exponent $1/3$. However, within the less than
two orders of magnitudes in time which are observable in the experiments, a
logarithmic scaling may appear as a pseudo-algebraic scaling. In particular,
we demonstrate that the CA model can exhibit a range of exponents $\leq 1/3$
for such a pseudo-algebraic scaling. Additionally, we reproduce the
experimental data with segregation of the CH model, although the CH model
matches in an intermittent regime, where it exhibits an algebraic scaling with
an exponent $1/4$, which is smaller than the asymptotic exponent of $1/3$.
This corroborates the ambiguity of the scaling behavior, when segregation
processes are only observed on short time spans.

## Introduction

Pattern formation of cells and cell segregation are complex and crucial
processes, in particular in the context of embryogenesis. When different types
of cells are intermixed, they start to segregate into homogeneous domains [1,
2, 3, 4, 5]. This behavior has been shown for many different cell types in
several species, for instance hydra [6, 7], zebra fish [8] and chicken [9,
10]. Why and how cells rearrange themself in a certain way is still not fully
understood, and various theories and hypotheses have been formulated to
explain the process of cell segregation [11, 12, 13, 14, 15, 16, 8, 17, 18,
19, 20, 21, 5].

One of the most well-known theories in the context of cell segregation is the
differential adhesion hypothesis of Steinberg [15, 22], which focuses on the
impact of adhesion on cell segregation. He proposed that the sorting behavior
of cells results from differences in the adhesion strengths between different
cell types, which implies that sorting is driven by the minimization of the
surface energy. Additionally, he suggested that a mixed cell population will
always minimize its total adhesive free energy and conjectured that cells
segregate like demixable fluids, e.g., water and oil.

The separation of fluids is theoretically well studied. The kinetics of this
separation can be modeled with the Cahn-Hilliard Navier-Stokes equations [23,
24, 25]. The level of segregation is typically quantified by segregation
indices, the interface length between clusters of different type or by the
average cluster diameter. For a narrow cluster size distribution, the average
cluster diameter scales inverse-proportional to the interface length and
segregation indices, see SI text. An increase of the level of segregation
corresponds to a decrease of the former two measures and, accordingly, an
increase of the latter one, the average cluster diameter. For the Cahn-
Hilliard Navier-Stokes model, it is well known that during segregation the
interface length exhibits an algebraic decay over several orders of magnitude
in time. The exponent of this algebraic scaling depends on the flows in the
model, which is influenced, among others, by the length scale of the system
[26], ranging from $1/3$ [27, 26, 28] for the diffusive regime, into which the
scenario of segregating biological cells falls, to $2/3$ for the laminar or
turbulent regime [26]. Note that, on the temporal scale, these exponents are
only reached asymptotically and can be preceded by exponents down to $1/6$ in
an intermittent regime [28]. In either case, the average cluster diameter is
inverse-proportional to the interface length, that is both the cluster
diameter and the interface length scale algebraically with exponents that are
equal in absolute value but have opposite signs.

In contrast to fluid segregation, not only one but a variety of agent-based
models have been used to simulate the segregation of biological cells [29, 11,
30], since there is a variety of cell-based mechanisms, such as active cells
or cell interaction mechanisms beside adhesion, which have potential influence
on the segregation and need to be studied. While algebraic scalings of the
segregation indices over time can be observed in most of these models, the
corresponding exponents vary over a wide range of $1/40-1/3$, depending on
which segregation mechanisms are incorporated and which models are used [17,
18, 19, 31]. One of the earliest attempts of simulating cell segregation is
the Cellular-Potts-Model (CPM) of Glazier and Graner [13, 14], in which
segregation results from differential adhesion. While the observed segregation
indices display a logarithmic decay, successive studies concluded that the
segregation indices actually follow a logarithmic decay only initially and
settle to an algebraic one for longer times [18, 17, 19]. Nakajima and
Ishihara [18] used the CPM to study the effects of even and uneven cell type
ratios on the segregation process. They found the exponent of the algebraic
scaling to decrease for increasingly asymmetric mixtures of cells, with
exponents ranging from $1/3$ for a $50/50$ ratio down to $1/4$ for a $90/10$
ratio. In any case, they observed the average cluster diameter to be inverse-
proportional to the segregation indices. Belmonte et al. [17] modeled
segregation by a self-propelled particle model with velocity alignment to
study the influence of collective motion. They also observed algebraic scaling
with an exponent of maximal $0.18$ concluding that even weak collective motion
accelerates cell segregation. Beatrici and Brunnet [19] studied a specific
particle system incorporating velocity differences between cell types, the
boids model, and concluded that velocity differences are sufficient to
generate algebraic segregation even without collective motion. Depending on
the chosen velocities and cell ratios between fast and slow cell types, they
observe both logarithmic and algebraic scaling, the latter with exponents
ranging from $0.18$ to $0.22$. The latter finding is supported by a study of
Strandkvist et al. [31] who found an algebraic scaling with exponents ranging
from $0.025$ to $0.17$ with a particle system incorporating velocity
differences between cell types.

Concerning data, several experiments have been conducted on cell segregation.
Rieu and Sawada [6], Schötz et al. [32] and Beysens et al. [33] conducted
experiments with hydra cells and zebra fish cells. They noticed similarities
of cell behavior to fluids by comparing characteristics of cell segregation
with those expected for viscous fluids according to hydrodynamic laws. For
instance, they compare the ratio of viscosity to surface tension and the time
course of relaxation to the equilibrium and the characteristics of the reached
equilibria. Cochet-Escartin et al. [34] studied hydra cells both in
experiments and in CPM simulation to determine whether differences in tissue
surface tension are sufficient for segregation. They found algebraic scaling
with exponent $0.74$ for the experiments and $0.5$ for the simulations.
However, they only measured cell segregation in the experiments over half an
order of magnitude in time. In contrast, Méhes et al. [20] studied the
influence of collective motion in experiments with fish and human cells and
measured algebraic cell segregation indices with exponent of $0.31$ for more
than two orders of magnitude in time. They further measured the average
cluster diameter, with an algebraic increase with exponent of $0.74$. This
means that the cluster diameter was not inverse-proportional to the
segregation indices, indicating that the cluster size distribution is not
narrow. They suspected that this behavior was a result of collective motion,
which they concluded to be a segregation promoting effect.

In summary, in the context of cell segregation, an algebraic scaling with an
exponent that differs from $1/3$, the value expected for fluid segregation,
has been attributed to additional intercellular interaction besides
differential adhesion [16, 11, 35]. Such mechanisms include collective motion
[20, 36] or velocity-dependent interaction of the cells [17, 19, 31]. However,
for both, experiments and numerical simulations, the algebraic decay of the
segregation indices is usually only observed during the last two orders of
magnitude of time [19, 20, 18, 31] or on an even shorter time interval [34].

We use a cellular automaton according to Voss-Böhme and Deutsch [21], which
solely incorporates adhesive forces between cells, to reproduce cell
segregation experiments from Méhes et al. [20]. We find a good match between
experiments and simulations over the whole time span of the experiments. This
is surprising, since our model always generates logarithmic scaling of the
segregation indices over time. The match between the model and the proposed
algebraic scaling with exponent $1/3$ in the experiments is possible since the
experimental observation is limited to less than two orders of magnitude in
time. To make this point more pronounced we will use the term pseudo-algebraic
scaling for such behavior in the following. Depending on the model parameters
and the considered time interval, we observe this pseudo-algebraic scaling
with a range of exponents $\leq 1/3$. In the light of such possible
misinterpretations, experimental segregation may actually scale
logarithmically. We propose that, while additional effects like collective
motion might be segregation promoting, the main factor that rices cell
segregation may still be adhesive forces. Moreover, we also find an match
between the experimental data and the Cahn-Hilliard model. For this
comparison, we develop a mapping between the length scales of the cellular
automaton and the Cahn-Hilliard model, such that only a single parameter of
the CHM, the mobility constant which sets the time scale, has to be fitted. It
turns out that, the relevant observation window of the experiments falls in
the intermittent regime of the Cahn-Hilliard model, exhibiting an algebraic
scaling exponent of $1/4$. Although Méhes et al. [20] suggested an algebraic
scaling exponent of $1/3$ for the experimental data, we find a good agreement
with the Cahn-Hilliard model as well, due to the short observation span. The
fact that both models, the cellular automaton and the Cahn-Hilliard model,
both agree with the experimental data, while exhibiting different scalings for
the experimental setup, corroborates the ambiguity of scaling behavior, when
segregation processes are only observed on short time spans. Since it is
difficult to extend experimental time spans, this highlights the importance of
considering additional features of segregation beyond segregation indices.

## Results

### Cellular Automaton can reproduce in vitro experiments

We compare our cellular automaton simulations with in vitro data of Méhes et
al. [20], see Fig. 1. They measured the segregation indices and cluster sizes
in the segregation of EPC (fish keratocyte cell line) with PFK (primary
goldfish keratocytes) and HaCaT (human keratocyte cell line) with EPC over two
orders of magnitude in time. The cellular automaton has five parameters, which
are calibrated to the experimental data: Three adhesion parameters
$\text{\boldmath$\beta$}=(\beta_{00},\beta_{10},\beta_{11})^{T}$, which set
the homotypic $(\beta_{11},\beta_{00})$ and the heterotypic $(\beta_{01})$
adhesion strengths, the cell type ratio $N_{0}/N_{1}$, which reflects the
ratio of all numbers $N_{i}$ of each cell type in the segregation experiments,
$i\in\\{0,1\\}$, and the time scale of migration $\tau$, which relates to the
dimensionless time of the cellular automaton to physical time. Note that the
three adhesion parameters can be reduced to two effective parameters, the
difference of homotypic adhesion $db$ and the difference between average
homotypic and heterotypic adhesion $b^{*}$, see Materials and Methods for
details.

In the experiments, equal areas are covered by each cell type, which results
in different cell numbers due to slightly different cell sizes for each type.
We show that the ratio of cell type numbers $N_{0}/N_{1}$ is set by the ratio
of the segregation indices $\gamma_{1}(t)/\gamma_{0}(t)$ and thus can be
obtained directly from the experimental data, see Eq. (12) and Fig. 7. We
check that this ratio is consistent with the ratio of cell sizes of each type
and that the total numbers of cells of in the experiments and the simulations
are comparable, see Materials and Methods.

Fig 1: The cellular automaton simulations reproduce the biological cell
segregation experiments of Méhes et al. [20]. The segregation indices
$\gamma_{i}(t)$ for the two experiments PFK with EPC (dark red and dark green
lines) and HaCaT with EPC (dark blue and dark green lines), within the
observed time interval $50-1000$ min, match with the segregation indices
predicted by the cellular automaton (lines with corresponding brighter
colors). Within the given time interval (grey box and inset), the segregation
indices seem to decay algebraically with exponent $1/3$ (black dashed line) as
expected asymptotically for fluid segregation. For the simulation of the
segregation indices $\gamma_{i}(t)$ of PFK ($i=0$) mixed with EPC ($i=1$), we
obtain a cell type ratio of $N_{\text{PFK}}/N_{\text{EPC}}=41.2/58.8$ and fit
the adhesion parameters $\textbf{$\beta$}=(-8,-6.5,0)^{T}$ and the time scale
of migration $\tau\approx 3.24\text{ min}$. For the simulation of the
segregation indices $\gamma_{i}(t)$ of HaCaT ($i=0$) mixed with EPC ($i=1$) we
obtain a cell type ratio of $N_{\text{HaCaT}}/N_{\text{EPC}}=35.2/64.8$ and
fit the parameters $\textbf{$\beta$}=(-8,-5.5,0)^{T}$ and $\tau\approx
13.32\text{ min}$. In both cases $140^{2}$ cells are simulated, comparable to
the cells visible in the experiments, starting from a random mixture.
Snapshots of the cell mixtures at the points marked with crosses labeled A),
C) and B), D) are displayed in Fig. 3.

Simulations and experiments match well for both cell mixtures, see Fig. 1.
This match is surprising, as the cellular automaton displays a logarithmic
scaling which contradicts the proposed algebraic scaling of the data in Méhes
et al. [20]. However, over just two orders of magnitude in time a logarithmic
decay may appear as almost straight line in a log-log plot, making it
difficult to distinguish it from a power law. Therefore, we denote an increase
or decay which approximately follows a straight line in a log-log plot, but
only for a limited time span, as pseudo-algebraic scaling. The match between
the prediction of our model and the experiment in Fig. 1 demonstrates, that it
is not possible to decide in the limited observation time of experiments
whether the segregation indices decay algebraically or logarithmically.

The fact that we are able to reproduce in vitro cell segregation with a model
which uses only adhesive forces suggests that additional complex processes or
forces to segregate cells, like collective motion, are not required and their
presence may only slightly modify the observations.

The resulting time scales $\tau=3$ min for PFK/EPC and $13$ min for HaCaT/EPC
correspond to migration velocities $\Delta x/\tau=370~{}\mu$m/h and
$86~{}\mu$m/h, where $\Delta x$ is the width of a lattice node in the cellular
automaton, which is consistent with the range of reported average velocities
for PFK ($500~{}\mu$m/h), EPC ($30~{}\mu$m/h), and HaCaT ($34~{}\mu$m/h) [20].
While the fitted adhesion parameters may suggest that the homotypic adhesion
of HaCaT and PFK is weaker than EPC as well as that the homotypic adhesion of
HaCaT equal to that of PFK, this fit is not unique. In fact, due to the short
time span the fit is based on, a wide range of adhesion parameters can
reproduce the experimental observations (data not shown). In order to refine
the fit, additional data would have to be incorporated, for instance single
cell measurements of adhesion forces of each cell type.

### Cahn-Hilliard can reproduce in vitro experiments too

We also compare the segregation experiments of Méhes et al. [20] with fluid
segregation. For this, we use the Cahn-Hilliard model which well describes
fluid segregation in the diffusive regime in terms of a phase-field
formulation, see SI text for details. To fit the parameters of the spatially
continuous Cahn-Hilliard model to the experimental data, which is based on
discrete cells, we develop a mapping between the agent-based cellular
automaton and the Cahn-Hilliard model, see SI text. Due to this mapping, only
the mobility constant $D$ of the Cahn-Hilliard model has to be fitted to match
the time scale of the experiments, while the remaining parameters can be
inferred from the parameters of the cellular automaton used for Fig. 1.

Fig 2: The Cahn-Hilliard simulations reproduce the biological cell segregation
experiments of Méhes et al. [20]. The segregation indices $\gamma_{i}(t)$ for
the two experiments PFK with EPC (dark red and dark green lines) and HaCaT
with EPC (dark blue and dark green lines), within the observed time interval
$50-1000$ min, match with the segregation indices predicted by the cellular
automaton (lines with corresponding brighter colors). The segregation indices
of the Cahn-Hilliard simulation decay algebraically with exponent $1/4$ (black
dotted line) rather than $1/3$ (black dashed line), which implies that the
segregation process is in an intermittent regime of fluid segregation, see
text for details. Using a mapping from the cellular automaton model to the
Cahn-Hilliard model, see Materials and Methods. Parameters are set analogous
to the parameters used in Fig. 1 except for the mobility constant $D$, which
is fitted to $D=36\text{$\mu$m}^{2}/\text{min}$ for the mixture of PFK with
EPC and $D=18\text{$\mu$m}^{2}/\text{min}$ for the mixture of HaCaT with EPC.
Snapshots of the cell mixtures at the points marked with crosses labeled A),
C) and B), D) are displayed in Fig. 3.

Cahn-Hilliard simulations and experiments match well for both cell mixtures,
see Fig. 2. The model fits EPC and PFK very well. However, a small discrepancy
can be observed at the end of the fit of HaCaT from the Cahn-Hilliard model,
which nevertheless reproduces the data as well as the by Méhes et al. [20]
suggested $1/3$ algebraic scaling exponent. This match is surprising, as the
Cahn-Hilliard simulations rather display an algebraic decay with exponent of
$1/4$ than $1/3$, which was proposed for the data in Méhes et al. [20].
However, within just two orders of magnitude in time it is hard to distinguish
a power-law decay with exponent $1/4$ and one with exponent $1/3$.

Note that the segregation indices resulting from the Cahn-Hilliard model
follow only asymptotically $(t\rightarrow\infty)$ an algebraic scaling with
the exponent of $1/3$. This asymptotic decay is usually referred to when cell
segregation is compared to fluid segregation and the exponent does not depend
on the parameters of the Cahn-Hilliard model. However, the intermittent decay
of the segregation indices, before the asymptotic regime is reached, displays
a slower algebraic scaling, with exponents down to $1/6$ [28], and can even
exhibit logarithmic decay, see Fig. 2. This intermittent regime can last for
several orders of magnitude in time, and an uneven cell type ratio can
increase the duration of this regime [28]. Furthermore, while the mobility
constant $D$ primarily rescales the physical time in the Cahn-Hilliard model,
we observe that it can also alter the duration of the intermittent decay. For
instance, the simulation displayed in Fig. C, which is based on parameters
comparable to the ones used for Fig. 1 except that the mobility constant $D$
is several orders of magnitude bigger, already exhibits an algebraic scaling
with the exponent of $1/3$ at segregation indices $\leq 1/2$.

The determined mobility constants of $D=36\text{$\mu$m}^{2}/\text{min}$ and
$18\text{$\mu$m}^{2}/\text{min}$ for PFK with EPC and HaCaT with EPC,
respectively, are consistent with the range of the experimentally measured
mobility constants for each cell type, i.e. PFK
($132\text{$\mu$m}^{2}/\text{min}$), EPC ($1.29\text{$\mu$m}^{2}/\text{min and
}$), HaCaT ($1.61\text{$\mu$m}^{2}/\text{min}$) [20]. In particular, the
fitted mobility constant for PFK with EPC is greater than that of HaCaT with
EPC, as expected from the individual mobility constants of each cell type.

In conclusion, we observe that two fundamentally different models both match
the experimental segregation indices on the limited time span, see Fig. 1 and
Fig. 2. Since the segregation indices are not sufficient to distinguish
between both models with respect to the experimental observations, we
additionally compare the distribution of cluster sizes, the morphology of the
clusters, and the average cluster diameter at two different levels of
segregation qualitatively, see Fig. 3: In all three cases, the cell type that
is less abundant, here PFK shown in red, forms clusters surrounded by a single
contiguous domain of the more abundant cell type, here EPC shown in green. The
Cahn-Hilliard model results in a rather narrow distribution of cluster sizes
while clusters form circular shapes or slightly elongated bulges, see Fig. 3
E) and F). In contrast, the cells in the experiment of Méhes et al. [20]
display a wider distribution of cluster sizes with different shapes of
clusters, see Fig. 3 C) and D). Interestingly, the configurations of the
cellular automaton exhibit features very similar to the experiment, see Fig. 3
A) and B). Note, that we extracted the images for the experiment from Video S5
in Méhes et al. [20] according to the time points, which are coarsely labeled
in the video. Thus, the images may not exactly represent the chosen
segregation indices.

For the Cahn-Hilliard model, we find that the average cluster diameter is
inverse-proportional to the segregation indices (not shown), which is expected
for narrow cell size distribution, see Materials and Methods. For the cellular
automaton, the average cluster diameter displays a slightly steeper slope as
predicted by the inverse of the segregation indices (not shown), but still
close to $\sim t^{1/3}$. In contrast, Méhes et al. [20] measured an algebraic
increase of their average cluster size with an exponent close to $0.74$
instead of the suggested value of $1/3$, which means that their average
cluster size is not inverse-proportional to their segregation indices. They
suspect that this is a consequence of collective motion, meaning that
collective motion contributes to wider distribution of cluster size.

Fig 3: The cellular automaton reproduces morphology and size distribution of
the cell clusters in the experiments of Méhes et al. [20] of EPC (green) with
PFK (red) closer than the Cahn-Hilliard model. The snapshots of the cell
mixtures A), C), E) of the first row are taken at a segregation index of EPC
$\gamma_{\text{EPC}}=0.25$, at the start of the experimental recording, while
the pictures B), D), F) in the second row are at a segregation index of EPC
$\gamma_{EPC}=0.1$, at the end of the recording. A) and B) show the cellular
automaton, C) and D) show the experiments and are taken from video Video S5 in
Méhes et al. [20], and E) and F) show the Cahn-Hilliard model. The images A),
B), E) and F) show a detail from the simulations, such that approximately
$100^{2}$ cells are visible, to match the spatial scale of the images C) and
D) of the experiments. The time points corresponding to the images are marked
by black crosses in Fig. 1 and Fig. 2.

### Parameter influence of the cellular automaton on the segregation

We have already shown that the pseudo algebraic scaling in the segregating
experiments can be explained both by a logarithmic scaling law from the
cellular automaton and by an algebraic scaling with exponent of $1/4$ from the
Cahn-Hilliard model. Thus, despite the fact that both models only incorporate
adhesion forces, as proposed by Steinberg, the resulting segregation differs
fundamentally between both models. In addition, in the time frame of the
experiment neither model generates an algebraic scaling with an exponent
$1/3$, which is usually associated with fluid-like segregation. Firstly this
highlights, that not only an algebraic exponent of $1/3$ corresponds to fluid-
like segregation, but exponents between $1/6$ and $1/3$ may indicate it as
well. Secondly, this implies that in contrast to implicit suggestions of
previous works, an exponent differing from $1/3$ does not necessitate other
intercellular interactions or mechanical forces besides adhesion.

To relate the range of segregation dynamics displayed by the cellular
automaton to previous models, we study numerically the pseudo algebraic
scaling exponents, which can be generated by the automaton, and how they
depend on the adhesion parameters. The reparameterized effective adhesion
parameters $db$ and $\beta^{*}$ determine the kinetics of the segregation
results and therefore, by adjusting these parameters we are able to study the
impact of those on the possible exponents. The exact influence of $db$ and
$\beta^{*}$ on the scaling behavior is not trivial. However, the cellular
automaton is capable of producing a wide range of pseudo algebraic scalings,
see Fig. 4. We observe an upper bound for the exponent of the pseudo-algebraic
scaling at $1/3$, consistent with exponents observed in previous particle
models [18, 17, 19, 31]. Due to the logarithmic decay, the pseudo algebraic
scaling exponent over two orders of magnitudes increased with increasing
starting time of the observation window, i.e., it is maximal if the
segregation indices at the start of the observation are small.

Fig 4: A $1/3$ exponent forms an upper bound for the pseudo algebraic scaling
in cellular automaton. Segregation indices obtained from the simulation are
shown for a range of the effective parameters $db$ and $\beta^{*}$, but only
for the last two orders of magnitudes in time $\tilde{t}$ before the
segregation index reaches $\gamma=0.1$ in each simulation. For comparability,
the time scale of each simulation is set such that all simulations reach
$\gamma=0.1$ at $\tilde{t}=1$. For each simulation we use $100^{2}$ cells, a
cell type ratio of 50/50, periodic boundary conditions, and a random mixture
$\gamma=0.5$ as initial configuration.

However, in contrast to the parameters $db$ and $\beta^{*}$, the cell type
ratio does not influence the exponent of the pseudo-algebraic scaling. As
shown in Materials and Methods and visualized in Fig. 5 A), the cell type
ratio just increases the distance between $\gamma_{0}$ and $\gamma_{1}$, but
never the slope in the last orders of magnitudes in time, Fig. 5 B).

Fig 5: Segregation indices obtained from the simulation shown for a range of
cell type ratios. For each simulation we use $100^{2}$ cells, periodic
boundary conditions, $db=0$, $\beta^{*}=3$, and a random mixture as initial
configuration. For comparability, the time scale of migration $\tau$ of each
simulation is set such that all simulations reach segregation indices
$\gamma_{0}$ and $\gamma_{1}$ with $\gamma_{0}N_{0}=\gamma_{1}N_{1}=500$ at
dimensionless time $\tilde{t}=1$. Every color represents a specific cell type
ratio, while each cell type ratio was simulated five times. A) shows the raw
data of the simulations. The black lines correspond to an even cell type
ratio, for which both segregation indices match, while for uneven ratios the
segregation index of the more abundant cell type is below the black line and
the other above. B) shows the same data where each segregation index
$\gamma_{i}$ is rescaled to a segregation index $\tilde{\gamma}_{i}$ at an
even ratio according to $\tilde{\gamma}_{0}=2\gamma_{0}N_{0}/(N_{1}+N_{0})$.

## Materials and Methods

### Cellular automaton

For simulating cell segregation, we use a cellular automaton based on Voss-
Böhme et al. [21]. We use a 2D-quadratic lattice $S$ with $N\in\mathbb{N}$
nodes for each dimension, $S=\\{1,...,N\\}\times\\{1,...,N\\}$. We assign
exactly one cell to each node. Each cell has the area of $(\Delta x)^{2}$,
which leads to lattice lengths $N\Delta x$ for each side. Every cell is mapped
to a specific cell type $W=\\{0,1\\}$ with $\xi:S\rightarrow W$ defining a
specific configuration of cells on the lattice. Based on two possible cell
types, we define three adhesion parameters
$\text{\boldmath$\beta$}=(\beta_{11},\beta_{10},\beta_{00})^{T}$ which set the
stickiness of two directly neighboring cells depending of their type. The more
two neighboring cells stick to each other, the larger the associated
$\beta_{ij}$ parameter. $\tau$ denotes a parameter to adjust the time scale of
migration in the simulation. Further, based on these parameters, the rate
$r(\textbf{x},\textbf{y})$ of two cells at neighboring positions
$\textbf{x},\textbf{y}\in S,|\textbf{x}-\textbf{y}|=\Delta x$ swapping their
locations is given by:

$r(\textbf{x},\textbf{y},\xi)=\begin{cases}\tau^{-1}\exp\\{\beta_{\text{sum}}(\textbf{x},\textbf{y},\xi)\\}&\text{,
if }\xi(\textbf{y})\neq\xi(\textbf{x})\text{ and
}|\textbf{x}-\textbf{y}|=\Delta x\\\ 0&\text{, otherwise}\end{cases}$ (1)

where

$\beta_{\text{sum}}(\textbf{x},\textbf{y},\xi)=-\sum_{\textbf{z}:|\textbf{z}-\textbf{x}|=\Delta
x}\beta_{\xi(\textbf{x})\xi(\textbf{z})}-\sum_{\textbf{z}:|\textbf{z}-\textbf{y}|=\Delta
x}\beta_{\xi(\textbf{y})\xi(\textbf{z})}\text{.}$ (2)

The model is formulated in continuous time, see SI text for the details, on
the implementation.

Further, the cellular automaton simulates segregation and thus the segregation
indices $\gamma(\tilde{t})$ in a dimensionless time $\tilde{t}$. The time
scale of migration $\tau$ which transforms this dimensionless time $\tilde{t}$
into physical time $t$, as $t=\tau\tilde{t}$, is calibrated based on the
experimental data. By matching the physical time
$t_{\gamma=\gamma_{\text{match}}}$ at which the experimental segregation
indices first reaches a particular value $\gamma_{\text{match}}$, such that
$\gamma_{\text{exp}}(t_{\gamma=\gamma_{\text{match}}})=\gamma_{\text{match}}$
and the dimensionless time $\tilde{t}_{\gamma=\gamma_{\text{match}}}$ at which
the simulated segregation indices first reaches this value
$\gamma_{\text{sim}}(\tilde{t}_{\gamma=\gamma_{\text{match}}})=\gamma_{\text{match}}$,
and estimate
$\tau=t_{\gamma=\gamma_{\text{match}}}/\tilde{t}_{\gamma=\gamma_{\text{match}}}$.
If no value for $\tau$ is provided, it is set to $1$ dimensionless and
therefore neglected.

Voss-Böhme et al. [21] proposed an effective parameter $\beta^{*}$ for two
cell types, which determines the asymptotic sorting behavior, where

$\beta^{*}=\beta_{00}+\beta_{11}-2\beta_{10}\text{.}$ (3)

The impact of this parameter has been numerically confirmed and generalized to
an arbitrary number of cell types by Rossbach et al. [37]. We reparametrize
the adhesion parameters $\beta$ based on the effective parameter $\beta^{*}$
to better describe the impact of the parameters on the segregation behavior:

$db=\beta_{11}-\beta_{00}\text{,}$ (4)
$d=\beta_{00}+\beta_{10}+\beta_{11}\text{.}$ (5)

This leads to the following invertible transformation equation:

$\text{\boldmath$\beta$}=\begin{pmatrix}\beta_{00}\\\ \beta_{10}\\\
\beta_{11}\end{pmatrix}=\frac{\beta^{*}}{3}\begin{pmatrix}\frac{1}{2}\\\ -1\\\
\frac{1}{2}\end{pmatrix}+\frac{d}{3}\begin{pmatrix}1\\\ 1\\\
1\end{pmatrix}+\frac{db}{2}\begin{pmatrix}-1\\\ 0\\\ 1\end{pmatrix}\text{.}$
(6)

The parameter $d$ rescales the rates in a trivial way, since an increase of
$d$ by $\Delta d$ will increase all $\beta_{ij}$ by the same amount $1/3\Delta
d$ and therefore decrease all rates by a factor $\exp\\{-8/3\Delta d\\}$,
independently of $\xi,x,y$:

$\begin{split}r(\textbf{x},\textbf{y})&=\tau^{-1}\exp\left\\{-\sum_{\textbf{z}:|\textbf{z}-\textbf{x}|=\Delta
x}{\left(\beta_{\xi(\textbf{x})\xi(\textbf{z})}+\frac{1}{3}\Delta
d\right)}-\sum_{\textbf{z}:|\textbf{z}-\textbf{y}|=\Delta
x}{\left(\beta_{\xi(\textbf{y})\xi(\textbf{z})}+\frac{1}{3}\Delta
d\right)}\right\\}\\\ &=\tau^{-1}\exp\left\\{-\frac{8}{3}\Delta
d\right\\}\exp\left\\{-\sum_{\textbf{z}:|\textbf{z}-\textbf{x}|=\Delta
x}{\beta_{\xi(\textbf{x})\xi(\textbf{z})}}-\sum_{\textbf{z}:|\textbf{z}-\textbf{y}|=\Delta
x}{\beta_{\xi(\textbf{y})\xi(\textbf{z})}}\right\\}\end{split}\text{.}$ (7)

The factor $\exp\\{-8/3\Delta d\\}$ just rescales the time scale of migration
$\tau$.

The effects of the parameter $db$ on the model system are more complex and
have been examined numerically. We initialize with a random configuration
$\xi$ and measure successively, for each subsequent configuration $\xi_{t}$
the sum $\lambda_{t}$ of all heterotypic transition rates in the whole system
at this time. The value $\lambda_{t}$ sets the current average waiting time
$\Delta t_{\text{swap}}=1/\lambda_{t}$ between two cell switches, see
implementation of the cellular automaton in SI text. We find that on average
an increase of the parameter $db$ will increase $\lambda_{t}$ and therefore
decrease the average waiting time $\Delta t_{\text{swap}}$. As illustration we
show the dependency of $\lambda_{0}$ on the parameters for a random
configuration $\xi$ in Fig. 6. Further, we notice that an increase of $db$
will also increase the cpu time, i.e., the number of cell switches, required
to reach the same level of segregation.

Fig 6: The influence on the time scale is trival for the parameter $d$, and
non trival for the parameter $db$. Shown is, for a constant $\beta^{*}=3$ and
random initial conditions $\xi$ with a $50/50$ cell type ratio, the color
coded sum $\lambda_{0}$ of all heterotypic transitions rates. In direction of
$(1,1,1)^{T}$, the value of $\lambda$ decreases and therefore the simulation
time $\Delta t_{\text{swap}}$ for two neighboring heterotypic cells to change
positions increases. In direction of $(-1,0,1)^{T}$, the time dependency is
nontrivial, but symmetric to $db=0$.

### Segregation index

As in the experiment, we use type specific segregation indices $\gamma_{i}$ to
determine the degree of segregation over time in the cell-based model. For
type $i\in W$, the index $\gamma_{i}$ is the average of the amount
$n_{\neq}(k)$ of heterotypic neighbors, where the average is taken over all
positions $k$ carrying cells of type $i$, in relation to the maximum possible
numbers of neighbors, which is $4$ for a von-Neumann neighborhood,

$\gamma_{i}=\frac{1}{4}\left\langle
n_{\neq}(k)\right\rangle_{\begin{subarray}{c}k\in S\\\
\eta(k)=i\end{subarray}}=\frac{1}{4N_{i}}\sum_{\begin{subarray}{c}k\in S\\\
\eta(k)=i\end{subarray}}{n_{\neq}(k)}=\frac{1}{4N_{i}}I\text{,}$ (8)

where $N_{i}$ denotes the total number of cells of type $i$ and $I$ denotes
the interface length,

$I=\sum_{\begin{subarray}{c}k\in S\\\
\eta(k)=i\end{subarray}}{n_{\neq}(k)}\text{,}$ (9)

which is another commonly used measure of segregation. Further, if an even
cell type ratio is given ($50/50$), it applies $N_{i}=N^{2}/2$, where
$N^{2}=|S|$. The resulting prefactor $2N^{2}$ is equal to the maximum
achievable interface length in the cellular automaton, which corresponds to a
checkerboard configuration where each cell has four heterotypic neighbors,
$I_{\text{max}}=2N^{2}$. Based on this, the relative interface length $I_{r}$
can be defined as the interface length $I$ normalized by the maximal interface
length $I_{\text{max}}$,

$I_{r}=\frac{I}{I_{\text{max}}}=\frac{1}{2N^{2}}\sum_{\begin{subarray}{c}k\in
S\\\ \eta(k)=i\end{subarray}}{n_{\neq}(k)}\text{.}$ (10)

Thus, for an even cell type ratio $N_{0}=N_{1}$ the relative interface length
$I_{r}$ is equal to the segregation indices $\gamma_{0}=\gamma_{1}=I_{r}$. If
the numbers of cells of each type $N_{i}$ are not equal, it follows from Eq.
(8) that the segregation indices $\gamma_{i}$ are inverse-proportional to the
cell type ratio

$\frac{\gamma_{0}(t)}{\gamma_{1}(t)}=\frac{N_{1}}{N_{0}}\text{.}$ (11)

Therefore, the scaling exponents of $\gamma_{i}$ and $I_{r}$ are always
identical.

The previous equations Eq. (8) and Eq. (11) apply exactly for periodic
boundary conditions. For other boundary conditions, the cell type ratio still
approximates the type specific segregation indices ratio
$\gamma_{0}(t)/\gamma_{1}(t)\sim N_{1}/N_{0}$, and the segregation indices
approximate the relative interface length $I_{r}\sim\gamma_{i}$. This is due
to the fact that boundary cells at the edge and in the corners have less than
$4$ neighbors, but their contribution gets less with rising lattice size $N$,
since the boundary size scales with $O(N)$ and the lattice size scales with
$O(N^{2})$.

Fig 7: The cell type ratio for the simulation can be obtained from the
experiments of Méhes et al. [20]. Shown is a comparison of the cell type ratio
$r=\gamma_{0}(t)/\gamma_{1}(t)$ from the experiments of Méhes et al. [20]
(dotted in color) and from the cellular automaton (dashed in corresponding
color). The cell type ratio in the cellular automaton is set to the
approximate mean of the ratios observed in the experiment (dashed black
lines). The cell type ratio was calculated by the ratio of the type specific
segregation indices for each time $t$ per experiment EPC with PFK (red) and
HaCaT with EPC (blue).

Based on Eq. (11) it is possible to calculate for every pair of type specific
segregation indices $\gamma_{0}(t)$ and $\gamma_{1}(t)$ the corresponding cell
type ratio and vice versa. We define the cell type ratio as

$r:=\frac{\gamma_{0}}{\gamma_{1}}\quad\text{with }N_{0}\geq N_{1}\text{,}$
(12)

where we assume without loss of generality that $r\leq 1$.

In this sense, the cell type ratio in the experiment can be obtained from the
ratio of the corresponding segregation indices, see Fig. 7. Indeed, the ratio
of segregation indices is relatively constant over the time span of the
experiment, and we set the cell type ratio of the cellular automaton
accordingly, as indicated in Fig. 7. In their experiments, Méhes et al. [20]
choose the cell type ratio such that initially the same amount of area is
covered by each type. Since the cells of each type are similar in size, as EPC
is $300\mu m^{2}$ [38], PFK is $400\mu m^{2}$ [11] and HaCaT is $80-400\mu
m^{2}$ [39], it is reasonable to simulate the segregation with the cellular
automaton, where every cell has the same space of the grid. However, the small
differences in size imply that the number of cells of each type in the
experiment is not equal. Instead, one estimates from the cell size ratios
$0.75=A_{\text{EPC}}/A_{\text{PFK}}=A_{\text{EPC}}/A_{\text{HaCaT}}$, cell
type ratios which are consistent with the ones obtained from the ratio of
segregation indices $0.70=N_{\text{PFK}}/N_{\text{EPC}}$ and
$0.54=N_{\text{HaCaT}}/N_{\text{EPC}}$ for EPC with PFK and HaCaT with EPC,
respectively, see Fig. 7. Note that the specific cell type sizes where not
reported by Méhes et al. [20] and cells can vary in size during an experiment
as well as depending on the experimental setup.

Additionally to the segregation indices and the interface length, the average
cluster diameter is a third commonly used measure to determine order in
segregation processes. For the cellular automaton it can be shown that the
average cluster diameter $d$ is inverse-proportional to the interface length
$I$, assuming the cell type ratio equals $50/50$, and the cluster size
distribution is narrow, i.e. $\langle d_{l}\rangle^{2}\approx\langle
d_{l}^{2}\rangle$, and the total area $A_{\text{sum}}$ of all clusters is
constant. In the following $n_{c}$ denotes the number of clusters,
$A_{l},l\in\mathbb{N}|1\leq l\leq n_{c}$ the size and $U_{l}$ the scope of the
$l$-th cluster. Approximating the clusters as circles, we have

$\begin{split}&\sum_{l}^{n_{c}}{A_{l}}=A_{\text{sum}}=n_{c}\langle
A_{l}\rangle=n_{c}\frac{4}{\pi}\langle d_{l}^{2}\rangle\Leftrightarrow
n_{c}=\frac{4A_{\text{sum}}}{\pi\langle d_{l}^{2}\rangle}\\\
&\sum_{l}^{n_{c}}{U_{l}}=I=n_{c}\langle U_{l}\rangle=n_{c}\pi\langle
d_{l}\rangle\Leftrightarrow n_{c}=\frac{I}{\langle
d_{l}\rangle\pi}\text{.}\end{split}$ (13)

By combining the two expressions for $n_{c}$ in Eq. (13), we get

$\langle d_{l}\rangle=\frac{I\langle
d_{l}^{2}\rangle}{4A_{\text{sum}}}\Rightarrow\langle
d_{l}\rangle\underset{\langle d_{l}\rangle^{2}\approx\langle
d_{l}^{2}\rangle}{\approx}\frac{\langle d_{l}^{2}\rangle}{\langle
d_{l}\rangle}=\frac{4A_{\text{sum}}}{I}\text{,}$ (14)

where the last approximation is only valid for a narrow distribution of
cluster sizes. Since $4A_{\text{sum}}$ is constant, it results that the
average cluster diameter is inverse-proportional to the interface length
$\langle d_{l}\rangle\sim 1/I$. We infer from the fact that the average
cluster diameter in both, the cellular automaton and the Cahn-Hilliard model,
is inverse-proportional to the interface length, that their distribution of
cluster sizes is sufficiently narrow.

## Discussion

In conclusion, we modeled cell segregation with a cellular automaton, which
solely incorporates differential adhesion, to reproduce experimental
segregation indices measured by Méhes et al. [20]. While Méhes et al.
interpreted the decay of the segregation indices as an algebraic scaling with
the exponent of $1/3$, the cellular automaton model exclusively generates a
logarithmic decay. We attribute this contradiction to the limited time span
observable in the experiment, which is insufficient to determine the long term
scaling of the segregation indices. Thus, we refer to the seemingly algebraic
decay observed on a limited span as pseudo-algebraic scaling. The match of the
experimental results and the ones from the cellular automaton highlights the
possible ambiguity of scalings on short time spans. We quantify the range of
exponents possible with the pseudo-algebraic scaling of the cellular automaton
model and find the exponent of $1/3$ to be an upper bound. Since an algebraic
decay of the segregation indices with exponent $1/3$ is commonly considered
for fluid-like segregation, and Steinberg [15] proposed that cell segregation
is similar to that of fluids, we additionally compare the experimental results
with fluid segregation expressed by the Cahn-Hilliard model. In order to
adjust the spatial scale of the Cahn-Hilliard model to the cell segregation
experiment, we developed a mapping between the cellular automaton and the
Cahn-Hilliard model. The resulting segregation indices from the Cahn-Hilliard
model fit well the experimental ones, although they rather follow an algebraic
decay with exponent $1/4$ than $1/3$ on the relevant time interval, which is
again hard to distinguish on the short time span of the experimental data.

Note that the Cahn-Hilliard model, which well describes fluid segregation in a
diffusive regime, displays an algebraic decay of segregation indices with
exponent $1/3$ only asymptotically. There is also an intermittent regime,
which can last several orders of magnitude in time, during which exponents
down to $1/6$ are possible [28]. The smaller exponent of $1/4$ observed by us
means that for the setup of the experiment the corresponding fluid-like
segregation dynamics are in the intermittent regime. On the one hand, this
highlights that experimentally observed exponents smaller than $1/3$ do not
necessarily rule out fluid like segregation. On the other hand, this
demonstrates the importance of calibrating segregation models to actual
experimental data, as only limited time regimes of the model may be
experimentally relevant.

In this context, we present a way to fit both models to experimental data,
which can be applied to future experiments. Since the cell type ratio can
directly obtained from the segregation indices ratio and the time scale of
migration rescales the time scale by a factor, only two parameters remain to
be fitted for the cellular automaton. With respect to the mapping we developed
between the Cahn-Hilliard model and the cellular automaton model, only the
mobility constant $D$ as a single parameter has to be fitted for the Cahn-
Hilliard model.

We point out that in many experiments and simulations only two or less orders
of magnitudes in time are considered to determine the scaling behavior [19,
20, 18, 31, 34]. Our results suggest, that scaling behavior of segregation
indices on such short time spans is ambiguous, and algebraic scaling on these
time spans should be rather called pseudo-algebraic scaling, since it may also
be a misinterpreted logarithmic decay. This possible ambiguity has already
been hinted at before: Nakajima and Ishihara mentioned, that their segregation
can also be interpreted as a logarithmic scaling since the algebraic decay was
measured only in the last orders of magnitude [18]. Belmonte et al. indicated
that a logarithmic decay might be possible, if no coordinated motion of
neighbor cells is present [17]. The cellular automaton, which is solely based
on differential adhesion, can generate pseudo-algebraic decays which cover the
same range of exponents $\leq 1/3$ as models which additionally incorporate
other mechanisms like collective motion or differential velocities [17, 18,
19, 31]. Due to this ambiguity on short time scales, it is not possbile to
distinguish between specific models and thus to determine which mechanisms
govern the segregation solely based on the scaling behavior.

Many studies conclude from the scaling behavior of the segregation indices to
the impact of certain cell mechanisms like collective motion. In contrast, our
results strongly suggest to utilize additional metrics of segregation to
compare between simulations and experiments to overcome the ambiguous
interpretations of the segregation indices of experimental data on limited
time spans. Such segregation metrics could be the cluster size distribution,
morphology of the clusters, and the average cluster diameter. Note, that
despite the match of segregation indices between cellular automaton model,
Cahn-Hilliard model and experimental data, the distribution and morphology of
clusters observed in the experiment are closer reproduced by the cellular
automaton than the Cahn-Hilliard model. Note however that neither model did
reproduce the scaling of the average cluster diameter of the experiments by
Méhes et al. [20]: The average cluster diameter in the cellular automaton and
the Cahn-Hilliard model scale inverse-proportional to the segregation indices,
which is consistent with observations in CPM models [18]. In contrast, Méhes
et al. measured a steeper algebraic scaling with the exponent of $0.74$ for
their average cluster diameter [20].

Another interesting feature of cell segregation is which cell type encloses
the other. While it seems reasonable that the more abundant cell type should
enclose the other, Beatrici and Brunnet [19] have found different behavior
depending on the cell type ratios and the cells’ velocities. In addition, to
resolve the contradicting logarithmic decay found by Glazier and Graner [14]
in the CPM and the algebraic scalings found in successive studies with CPM
[18] and particle models [17, 19], Nakajima and Ishihara [18] proposed that
the number of cells considered in a simulation affects the scaling behavior.
However, we suggest that the considered time period is the crucial factor,
since many simulations display a logarithmic decay at the start [19, 17, 34,
18]. This is supported by Beatrici and Brunnet [19], which found no difference
in the scaling behavior for a wide range of cell numbers ($500$ to $8000$) in
their simulations. Likewise, we observed the same logarithmic decay for a
range of $25^{2}$ to $140^{2}$ cells per simulation, while only the
fluctuations of the segregation indices are diminished by using more cells.

## Acknowledgments

We thank Tamás Vicsek for his helpful communication about his publication.
This research is co-financed by the EU, the European Social Fund (ESF) and by
tax funds on the basis of the budget passed by the Saxon state parliament.

## Supporting Information (SI)

### Cahn-Hilliard Navier-Stokes

The Cahn-Hilliard Navier-Stokes model accurately describes the evolution of
two immiscible fluids under flow and diffusion [40, 41]. The model is also
used as a typical choice for simulating fluid segregation [27, 25, 24]. This
model is known for producing an algebraic scaling with the exponent of $1/3$
for segregation processes with small length scale [28, 26]. For larger length
scales, an exponent of 2/3 can be observed [26]. The complete model is given
by the following set of differential equations:

$\displaystyle\partial_{t}\Phi+\textbf{u}\cdot\nabla\Phi=D\Delta\mu\text{,}$
(S1) $\displaystyle\mu=-\epsilon^{2}\Delta\Phi+\Phi^{3}-\Phi\text{,}$
$\displaystyle\rho(\partial_{t}\textbf{u}+(\textbf{u}\cdot\nabla)\textbf{u})+\nabla
p=\nabla\cdot(\eta(\nabla\textbf{u}+\nabla\textbf{u}^{T}))+\tilde{\sigma}\epsilon^{-1}\mu\nabla\Phi\text{,}$
$\displaystyle\nabla\cdot\textbf{u}=0\text{.}$

At small length scales, predefined by the size of the biological cells, the
fourth order diffusion in the Cahn-Hilliard equations dominates, such that the
influence of flow can be neglected and Eq. (S1) simplifies to:

$\displaystyle\partial_{t}\Phi=D\Delta\mu\text{,}$ (S2)
$\displaystyle\mu=-\epsilon^{2}\Delta\Phi+\Phi^{3}-\Phi\text{.}$

These equations are defined for a domain $\Omega=[0,L_{x}]\times[0,L_{y}]$
where $L_{x}$ and $L_{y}$ denotes the maximal size of the domain. We define a
phase $\Phi:\Omega\rightarrow[-1,1]$ on this domain, where $\Phi\approx 1$
denotes the first fluid, like water, and $\Phi\approx-1$ denotes the second
fluid, like oil. Values of $|\Phi|<\Phi_{0}$ are defined as interface area,
e.g. $\Phi_{0}=0.9$. The width of the interface area is proportional to the
parameter $\epsilon$:

$\delta=\text{arctanh}(\Phi_{0})\sqrt{2}\epsilon\text{,}$ (S3)

see Fig. A below.

Fig A: Schematic representation of the relationship between $\Phi$ and
$\delta$ and their correlation to the interface area.

The parameter $D$ is the mobility constant and influences the time scale of
the diffusion process, which is set by $\tau\sim\frac{\epsilon^{2}}{D}$. With
the phase $\Phi$ and the parameter $\epsilon$ the chemical potential
$\mu(\Phi,\epsilon)$ can be calculated. Each fluid has its own typical
parameters like viscosity $\eta_{\Phi=-1}$ and $\eta_{\Phi=_{1}}$ and density
$\rho_{\Phi=-1}$ and $\rho_{\Phi=_{1}}$. By linear interpolation of the
viscosity and density pairs, the functions $\eta(\Phi)$ and $\rho(\Phi)$ can
be calculated. Dependent on the types of fluids there is a surface tension
$\sigma$. Which enters the Navier-Stokes equations after a rescaling as
$\tilde{\sigma}:=\frac{3}{2{\sqrt{2}}}$ [40, 41].

The implementation of this model follows a special pressure projection scheme
with incomplete pressure iterations and an explicit Euler approach described
in the paper of Adam et al. [42].

The interface length $I$ for this model is approximated by the Cahn-Hilliard
surface energy,

$I\approx\frac{3}{2\sqrt{2}}\int_{\Omega}{\frac{1}{\epsilon}W(\Phi)+\frac{\epsilon}{2}|\nabla\Phi^{2}|}\text{,}$
(S4) $W(\Phi)=\frac{1}{4}(1-\Phi^{2})^{2}\text{.}$ (S5)

### Mapping of the cellular automaton model and the Cahn-Hilliard model

We developed a mapping process, to equally start a simulation of the cellular
automaton model and the Cahn-Hilliard model and compare their segregation
behaviors. Therefore a match of the time or length scale between both models
is needed. Since the cellular automaton always segregates logarithmically and
the Cahn-Hilliard model algebraically, the time scale can never be the same.
However, the length scale can be matched. The length scale in the Cahn-
Hilliard model is set by the parameter for the width of the interface area
$\epsilon$. Since the interface area in the cellular automaton is sharp, we
define an equally wide transition to be from the middle of one cell
$\Phi\approx-0.9$ to the middle of a neighboring heterotypic cell $\Phi\approx
0.9$. The length of a side $\triangle x$ of one cell in the cellular automaton
can be calculated by the square root of the average of the specific cell areas
$A=(A_{0}+A_{1})/2\approx 350\mu m^{2}$. The transition area in the cellular
automaton refers to $2\delta$, see Eq. (S3), and is equal to one cell side
length $\triangle x\approx\sqrt{350}\mu m$, see Fig. B.

Fig B: Visualization of the mapping of the transition area between the
cellular automaton model and the Cahn-Hilliard model. Each square represents a
cell in the cellular automaton, where red denotes $i=0$ and blue $i=1$ for
$i\in W$. The length of one side is defined as $\triangle x$ and can be
calculated by the square root of the specific cell type area $A_{0},A_{1}$. We
define the middle of each cell to be equally to each phase $\Phi_{0}=0.9,-0.9$
of the Cahn-Hilliard model. Therefore can the absolute width of the transition
area $\delta$ be calculated and $\epsilon$ can be determined.

Further, the cell type ratio of the cellular automaton can be directly
integrated into the Cahn-Hilliard model, by initializing the simulation with
an equal phase ratio. The domain size $L_{x}=L_{y}$ can be obtained by the
number of cells per dimension multiplied with the average size of a cell
$N\delta=N\sqrt{350}\mu m$. Only the mobility constant $D$ needs to be fitted.
Since the interface width $\epsilon$ is set, the time scale gets only
influence by $D$. Therefore if $D$ gets doubled, the time scale halved. If the
length scales is matched, the initial interface length $I$ will be equal in
both models, if both start from a random field. The comparison of both models
can be seen in Fig. C.

Fig C: The length scale of the cellular automaton model and the Cahn-Hilliard
model can be matched. Due to diverse time scaling behavior the course of the
segregation of each model differs. The simulation of the cellular automaton
started with a random cell mixture, adhesion parameter
$\textbf{$\beta$}=(-1.56,-3.06,-1.56)^{T}$, time scale of migration $\tau=1s$,
$100^{2}$ cells, and a cell type ratio of $50/50$. The Cahn-Hilliard
simulation start from a random even distributed field with a mobility constant
$D=0.24$. Note that the Cahn-Hilliard model has always a settling process at
the start, where the phases get formed, see Fig A),B),C). Therefore, to verify
the match, the interface length at point C) of the Cahn-Hilliards model must
be compared to the beginning of the cellular automaton. Due to the continuous
interface evolution of the Cahn-Hilliard model, different methods are usually
used to measure the interface length in both models.

### Pseuodo algorithm for the cellular automaton

In order to simulate the cellular automaton in an effective way, we
implemented a version with continuous time. This improved the performance of
our simulations drastically, in comparison to a simple algorithm with discrete
time steps. The pseudo algorithm used for the cellular automaton reads as
follows:

  1. 1.

Initialise the lattice.

  2. 2.

Choose random one heterogeneous transition ($\textbf{x},\textbf{y}\in
S,|\textbf{x}-\textbf{y}|=\Delta x\wedge\xi\neq\xi^{\textbf{x},\textbf{y}}$).
Transitions with a higher rate $r(\textbf{x},\textbf{y})$, will be chosen with
a linear higher probability $P(\xi\rightarrow\xi^{\textbf{x},\textbf{y}})$.

  3. 3.

The two cells of the selected transition will swap there position on the
lattice.

  4. 4.

$\Delta t_{\text{swap}}$ is added to the time.

  5. 5.

If an end condition is reached, stop here.

  6. 6.

Else, return to step 2.

$\begin{aligned}
&P(\textbf{x},\textbf{y})=\frac{r(\textbf{x},\textbf{y})}{\sum_{B}r(e)},e\in
B,r(e):=r(\textbf{x},\textbf{y})\text{,}\\\ &B=\\{\textbf{x},\textbf{y}\in
S,|\textbf{x}-\textbf{y}|=\Delta
x\wedge\xi\neq\xi^{\textbf{x},\textbf{y}}\\}\end{aligned}\text{,}$ (S6)
$\Delta t_{\text{swap}}\sim Exp(\lambda(\xi))\text{,}$ (S7)
$\lambda(\xi)=\sum_{B}\exp(\beta_{\text{sum}}(\textbf{x},\textbf{y},\xi))\text{.}$
(S8)

Here $B$ denotes the set of all possible heterogeneity transitions.

## References

  * [1] Townes PL, Holtfreter J.  Directed movements and selective adhesion of embryonic amphibian cells.  J Exp Zool. 1955;128(1):53–120.  doi:10.1002/jez.1401280105.
  * [2] Cerchiari AE, Garbe JC, Jee NY, Todhunter ME, Broaders KE, Peehl DM, et al.  A strategy for tissue self-organization that is robust to cellular heterogeneity and plasticity.  Proc Natl Acad Sci USA. 2015;112(7):2287–2292.
  * [3] Xiong F, Tentner AR, Huang P, Gelas A, Mosaliganti KR, Souhait L, et al.  Specified neural progenitors sort to form sharp domains after noisy Shh signaling.  Cell J. 2013;153(3):550–561.  doi:10.1016/j.cell.2013.03.023.
  * [4] Kay RR, Thompson CRL.  Forming Patterns in Development without Morphogen Gradients: Scattered Differentiation and Sorting Out.  Cold Spring Harb Perspect Biol. 2009;1(6).  doi:10.1101/cshperspect.a001503.
  * [5] Canty L, Zarour E, Kashkooli L, François P, Fagotto F.  Sorting at embryonic boundaries requires high heterotypic interfacial tension.  Nat Commun. 2017;8(1):1–15.  doi:10.1038/s41467-017-00146-x.
  * [6] Rieu JP, Sawada Y.  Hydrodynamics and cell motion during the rounding of two dimensional hydra cell aggregates.  Eur Phys J B. 2002;27(1):167–172.  doi:10.1140/epjb/e20020142.
  * [7] Steinberg MS.  Differential adhesion in morphogenesis: a modern view.  Curr Opin Genet Dev. 2007;17(4):281–286.  doi:10.1016/j.gde.2007.05.002.
  * [8] Krieg M, Arboleda-Estudillo Y, Puech PH, Käfer J, Graner F, Müller DJ, et al.  Tensile forces govern germ-layer organization in zebrafish.  Nat Cell Biol. 2008;10(4):429–436.  doi:10.1038/ncb1705.
  * [9] Mombach n, Glazier n, Raphael n, Zajac n.  Quantitative comparison between differential adhesion models and cell sorting in the presence and absence of fluctuations.  Phys Rev Lett. 1995;75(11):2244–2247.  doi:10.1103/PhysRevLett.75.2244.
  * [10] Foty RA, Pfleger CM, Forgacs G, Steinberg MS.  Surface tensions of embryonic tissues predict their mutual envelopment behavior.  Development. 1996;122(5):1611–1620.
  * [11] Méhes E, Vicsek T.  Segregation mechanisms of tissue cells: from experimental data to models.  Complex Adapt Syst Model. 2013;1:4.  doi:10.1186/2194-3206-1-4.
  * [12] Krens SFG, Heisenberg CP.  Cell sorting in development.  Curr Top Dev Biol. 2011;95:189–213.  doi:10.1016/B978-0-12-385065-2.00006-2.
  * [13] Graner F, Glazier JA.  Simulation of biological cell sorting using a two-dimensional extended Potts model.  Phys Rev Lett. 1992;69(13):2013–2016.  doi:10.1103/PhysRevLett.69.2013.
  * [14] Glazier JA, Graner F.  Simulation of the differential adhesion driven rearrangement of biological cells.  Phys Rev E. 1993;47(3):2128–2154.  doi:10.1103/PhysRevE.47.2128.
  * [15] Steinberg MS.  Does differential adhesion govern self-assembly processes in histogenesis? Equilibrium configurations and the emergence of a hierarchy among populations of embryonic cells.  J Exp Zool. 1970;173(4):395–433.  doi:10.1002/jez.1401730406.
  * [16] Méhes E, Vicsek T.  Segregation mechanisms of tissue cells: from experimental data to models.  Complex Adapt Syst Model. 2013;1(1):1–13.  doi:10.1186/2194-3206-1-4.
  * [17] Belmonte JM, Thomas GL, Brunnet LG, de Almeida RMC, Chaté H.  Self-Propelled Particle Model for Cell-Sorting Phenomena.  Phys Rev Lett. 2008;100(24):248702.  doi:10.1103/PhysRevLett.100.248702.
  * [18] Nakajima A, Ishihara S.  Kinetics of the cellular Potts model revisited.  New J Phys. 2011;13(3):033035.  doi:10.1088/1367-2630/13/3/033035.
  * [19] Beatrici CP, Brunnet LG.  Cell sorting based on motility differences.  Phys Rev E. 2011;84(3):031927.  doi:10.1103/PhysRevE.84.031927.
  * [20] Méhes E, Mones E, Németh V, Vicsek T.  Collective Motion of Cells Mediates Segregation and Pattern Formation in Co-Cultures.  PLoS ONE. 2012;7(2):e31711.  doi:10.1371/journal.pone.0031711.
  * [21] Voss-Böhme A, Deutsch A.  The cellular basis of cell sorting kinetics.  J Theor Biol. 2010;263(4):419–436.  doi:10.1016/j.jtbi.2009.12.011.
  * [22] Steinberg MS.  Reconstruction of tissues by dissociated cells. Some morphogenetic tissue movements and the sorting out of embryonic cells may have a common explanation.  Science. 1963;141(3579):401–408.  doi:10.1126/science.141.3579.401.
  * [23] Lamorgese AG, Mauri R.  Diffuse-interface modeling of phase segregation in liquid mixtures.  Int J Multiph Flow. 2008;34(10):987–995.  doi:10.1016/j.ijmultiphaseflow.2008.03.003.
  * [24] Voorhees PW.  The theory of Ostwald ripening.  J Stat Phys. 1985;38(1):231–252.  doi:10.1007/BF01017860.
  * [25] Hardy SC, Voorhees PW.  Ostwald ripening in a system with a high volume fraction of coarsening phase.  Metall Trans A. 1988;19(11):2713–2721.  doi:10.1007/BF02645806.
  * [26] Naso A, Náraigh L.  A flow-pattern map for phase separation using the Navier-Stokes Cahn-Hilliard model.  Eur J Mech. 2017;72.  doi:10.1016/j.euromechflu.2018.08.002.
  * [27] Witkowski T, Backofen R, Voigt A.  The influence of membrane bound proteins on phase separation and coarsening in cell membranes.  Phys Chem Chem Phys. 2012;14(42):14509–14515.  doi:10.1039/C2CP41274H.
  * [28] Garcke H, Niethammer B, Rumpf M.  Transient Coarsening Behaviour In The Cahn-Hilliard Model.  Acta Mater. 2003;51.  doi:10.1016/S1359-6454(03)00087-9.
  * [29] Zhang Y, Thomas GL, Swat M, Shirinifard A, Glazier JA.  Computer Simulations of Cell Sorting Due to Differential Adhesion.  PLoS ONE. 2011;6(10):e24999.  doi:10.1371/journal.pone.0024999.
  * [30] Osborne JM, Fletcher AG, Pitt-Francis JM, Maini PK, Gavaghan DJ.  Comparing individual-based approaches to modelling the self-organization of multicellular tissues.  PLoS Comput Biol. 2017;13(2):e1005387.  doi:10.1371/journal.pcbi.1005387.
  * [31] Strandkvist C, Juul J, Baum B, Kabla AJ, Duke T.  A kinetic mechanism for cell sorting based on local variations in cell motility.  Interface Focus. 2014;4(6).  doi:10.1098/rsfs.2014.0013.
  * [32] Schötz EM, Burdine RD, Jülicher F, Steinberg MS, Heisenberg CP, Foty RA.  Quantitative differences in tissue surface tension influence zebrafish germ layer positioning.  HFSP J. 2008;2(1):42–56.  doi:10.2976/1.2834817.
  * [33] Beysens DA, Forgacs G, Glazier JA.  Cell sorting is analogous to phase ordering in fluids.  Proc Natl Acad Sci USA. 2000;97(17):9467–9471.
  * [34] Cochet-Escartin O, Locke TT, Shi WH, Steele RE, Collins EMS.  Physical Mechanisms Driving Cell Sorting in Hydra.  Biophys J. 2017;113(12):2827–2841.  doi:10.1016/j.bpj.2017.10.045.
  * [35] Vishwakarma M, Spatz JP, Das T.  Mechanobiology of leader–follower dynamics in epithelial cell migration.  Curr Opin Cell Biol. 2020;66:97–103.  doi:10.1016/j.ceb.2020.05.007.
  * [36] Fujimori T, Nakajima A, Shimada N, Sawai S.  Tissue self-organization based on collective cell migration by contact activation of locomotion and chemotaxis.  Proc Natl Acad Sci USA. 2019;116(10):4291–4296.  doi:10.1073/pnas.1815063116.
  * [37] Rossbach P, Böhme HJ, Lange S, Voss-Böhme A.  Model-based prediction of an effective adhesion parameter guiding multi-type cell segregation. 2021;.
  * [38] Fijan N, Sulimanović D, Bearzotti M, Muzinić D, Zwillenberg LO, Chilmonczyk S, et al.  Some properties of the Epithelioma papulosum cyprini (EPC) cell line from carp cyprinus carpio.  Annales de Virologie. 1983;134(2):207–220.  doi:10.1016/S0769-2617(83)80060-4.
  * [39] Boelsma E, Verhoeven MCH, Ponec M.  Reconstruction of a Human Skin Equivalent Using a Spontaneously Transformed Keratinocyte Cell Line (HaCaT).  J Invest Dermatol. 1999;112(4):489–498.  doi:10.1046/j.1523-1747.1999.00545.x.
  * [40] Aland S, Boden S, Hahn A, Klingbeil F, Weismann M, Weller S.  Quantitative comparison of Taylor flow simulations based on sharp-interface and diffuse-interface models.  Int J Numer Methods Fluids. 2013;73.  doi:10.1002/fld.3802.
  * [41] Aland S, Voigt A.  Benchmark computations of diffuse interface models for two-dimensional bubble dynamics.  Int J Numer Methods Fluids. 2012;69(3):747–761.  doi:10.1002/fld.2611.
  * [42] Adam N, Franke F, Aland S.  A Simple Parallel Solution Method for the Navier–Stokes Cahn–Hilliard Equations.  J Math. 2020;8(8).  doi:10.3390/math8081224.