Structural Causes of Pattern Formation and Loss Through Model-Independent Bifurcation Analysis

TL;DR

This paper introduces a new method to analyze how cell patterns form during development using minimal assumptions about cell communication and signaling.

Contribution

The novel contribution is a model-independent ODE framework that predicts stable cell patterns based on intercellular and intracellular signaling assumptions.

Findings

Global intercellular communication networks determine possible emergent patterns in a generic system.

Adding constraints on local intracellular signaling leads to a single stable pattern based on network features.

The framework allows inferring unknown signaling interactions by analyzing tissue-level patterns.

Abstract

During development, precise cellular patterning is essential for the formation of functional tissues and organs. These patterns arise from conserved signaling networks that regulate communication both within and between cells. Here, we develop and present a model-independent ordinary differential equation (ODE) framework for analyzing pattern formation in a homogeneous cell array. In contrast to traditional approaches that focus on specific equations, our method relies solely on general assumptions about global intercellular communication (between cells) and qualitative properties of local intracellular biochemical signaling (within cells). Prior work has shown that global intercellular communication networks alone determine the possible emergent patterns in a generic system. We build on these results by demonstrating that additional constraints on the local intracellular signaling network lead to a single stable pattern which depends on the qualitative features of the network. Our framework enables the prediction of cell fate patterns with minimal modeling assumptions, and provides a powerful tool for inferring unknown interactions within signaling networks by analyzing tissue-level patterns.

Figures

Full text

Introduction

1

Pattern formation is a hallmark of developmental biology, where cells within a tissue or organism differentiate in a precise spatial arrangement to form complex structures. This process is essential for cell fate specification, tissue development, and morphogenesis; and understanding the molecular causes of pattern formation has promising applications in regenerative medicine, tissue engineering, and the treatment of congenital malformations (Chuong and Richardson, 2009).

During development, a group of cells with the same developmental potential will receive a signal called a morphogen that prompts them to communicate and form a pattern of cell types. These local patterning events occur sequentially, with each iteration creating progressively finer cellular patterns across the developing organism. Recent ordinary differential equation (ODE) models have considered one iteration, focusing on how changes in cell-communication and chemical kinetics affect the cellular pattern (Williamson et al., 2021; O’Dea and King, 2012; de Back et al., 2013; Fisher et al., 2005; Irons et al., 2010; Giurumescu et al., 2006; Vasilopoulos and Painter, 2016; Hadjivasiliou et al., 2016). Most models investigate the Notch signaling pathway, which is a key signaling pathway involved in contact-mediated cell-communication. Although the full complexity of the Notch pathway is difficult to model, simplified systems – such as the model from Collier et al. (1996) – have offered valuable insights into how cell-communication can give rise to fine-grained patterns.

Over time, modeling changes have produced further understanding. For instance, O’Dea and King (2012) incorporated morphogens that act as bifurcation parameters, showing how an initially homogeneous steady-state can be destabilized, leading to the formation of patterns. They also derived a continuum model to account for morphogen gradients, showing how tissues with different chemical gradients form different patterns. Williamson et al. (2021) used a similar approach and demonstrated that the structure of communication networks among cells determines the types of fine-grained patterns that can emerge. Others have shown via simulations that coarse-grained patterns such as spots and stripes can arise if there is long-range signaling in the Notch pathway (Vasilopoulos and Painter, 2016; Hadjivasiliou et al., 2016).

Despite the usefulness of these models, they all rely on specific equations; therefore, there are possible dynamics that the models cannot represent due to potentially incomplete assumptions about the system, including assumptions about reaction rates and the molecular pathways involved (e.g. focusing on Notch signaling). Moreover, although model investigation has helped with understanding which parameters are important for pattern formation, the inherent complexity of any model prohibits us from determining which parameters are necessary or sufficient for a pattern. Recent advances in network theory (Golubitsky and Stewart, 2023) allow a more thorough investigation of pattern formation in many cell-communication networks. By ignoring all unnecessary information and focusing on the connectivity of cells (i.e. which cells influence each other), one can predict all possible patterns of cell fates that can emerge in a generic system (Wang and Golubitsky, 2004). Additionally, we provide simple conditions – related to qualitative features of the chemical signaling network – that are both necessary and sufficient for a pattern of cell fates to form under our assumptions.

Using our theory, we find that despite the limited chemical signaling pathways involved in development, cells can generate diverse patterns by reorganizing themselves, thus validating and extending the work of Williamson et al. (2021), Vasilopoulos and Painter (2016), and Hadjivasiliou et al. (2016). On the other hand, cells whose organization is constrained can form various patterns by using different chemical signaling mechanisms, as suggested in de Back et al. (2013). We identify conditions under which cells will fail to form stable patterns, instead remaining in a homogeneous steady-state, oscillating synchronously, or forming oscillating patterns. Importantly, our theory shows that the chemical kinetics proposed in the Collier model are not the only kinetics that can lead to pattern formation. Finally, we are able to infer possible biochemical interactions from an observed pattern, which can reveal previously unknown interactions in biochemical signaling pathways.

The paper is organized as follows. In Section 2, we state our biological assumptions and introduce our mathematical framework, including the theory of network dynamics. We highlight how bifurcations in networks can give rise to patterns corresponding to eigenvectors of the system. In Section 3, we show that qualitative features of chemical signaling determine the stability of the initially homogeneous steady-state and dictate the critical eigenvectors that lead to pattern formation. In Section 3.2, we apply our theory to various biological contexts and predict patterns using minimal modeling assumptions. We show how the same chemical signaling pathway can produce different patterns depending on the cell-communication network. Finally, in Section 3.3, we demonstrate how our framework can be used to infer properties of biochemical signaling pathways from observed patterns in tissues.

Mathematical methods

2

General model assumptions

2.1

We make the following general assumptions to guide the construction of our network models and their associated internal dynamics.

1. The cells are well-mixed compartments that can be characterized by the concentration of multiple chemical species inside the cell.
2. The concentrations of chemical species within the cells change smoothly over time.
3. The cells begin in a nearly identical state.
4. The external signals to each cell are identical.
5. The cells are in a sufficiently uniform lattice structure. (The lattice may be in 1,2, or 3 spatial dimensions).
6. The cells average signals from their neighbors.

Assumptions (1) and (2) allow us to represent our system of cells with ordinary differential equations. Assumptions (3)-(6) allow us to represent our system with a strongly connected, regular network as defined in Section 2.

Network definitions and basic properties

2.2

Our approach emphasizes the role of network structure in pattern formation, so we outline key definitions and properties that will be used throughout the paper (for a fully rigorous treatment, see (Nijholt et al., 2016; Golubitsky and Stewart, 2023)).

Definition 1 (Regular Networks). A regular network is a finite directed graph $[eqn]$ where

1. $[eqn]$ is the set of nodes
2. $[eqn]$ is the set of arrows
3. $[eqn]$ is the source map
4. $[eqn]$ is the target map.

Additionally, the number of input arrows to each node is the same (i.e. $[eqn]$ for all $[eqn]$ ). We call the number of input arrows the valence of the network.

For example, each network shown in Figure 1 (ignoring colors) is a regular network because every node has two input arrows.

There is a class of ODEs called admissible ODEs that is naturally associated with a regular network $[eqn]$ . Suppose that each node $[eqn]$ has an associated “state” $[eqn]$ – where $[eqn]$ is called the node space. Then the state of an $[eqn]$ -node network with node space $[eqn]$ can be described by $[eqn]$ .

Next, suppose that each node $[eqn]$ is only influenced by itself and nodes $[eqn]$ such that there is an arrow from $[eqn]$ to $[eqn]$ (i.e. $[eqn]$ and $[eqn]$ for some $[eqn]$ ). Let $[eqn]$ denote the set of input nodes to $[eqn]$ (i.e. nodes satisfying $[eqn]$ and $[eqn]$ .

Definition 2 (Admissible ODE for a Regular Network). An ODE is called admissible for the regular network $[eqn]$ if for some smooth function $[eqn]$ , the dynamics of each node can be written in the form

$$[eqn]$$

where the overline indicates that $[eqn]$ is symmetric in the arguments 2 through $[eqn]$ .

Definition 3 (Balanced Colorings of Regular Networks). For a regular network $[eqn]$ , assign every vertex a color $[eqn]$ via the coloring map

$$[eqn]$$

This coloring is balanced if for every $[eqn]$ with $[eqn]$ there exists a color-preserving bijection between their input nodes:

$$[eqn]$$

Intuitively, a coloring is balanced if for any two nodes $[eqn]$ with the same color, the set of nodes influencing $[eqn]$ has exactly the same colors, with the same frequency, as the set of nodes influencing $[eqn]$ . See Figure 1 for examples of balanced colorings.

Definition 4 (Polysynchrony Subspace). Let $[eqn]$ be a balanced coloring of a regular network $[eqn]$ given by the coloring map $[eqn]$ . The associated polysynchrony subspace is the set

$$[eqn]$$

See Figures 1, 2 for examples of balanced coloring and corresponding polysynchrony subspaces.

Definition 5 (Flow Invariance). Let $[eqn]$ be a smooth vector field. A subspace $[eqn]$ is flow-invariant under $[eqn]$ if $[eqn]$ .

Proposition 1 (Polysynchrony Subspaces are Flow-Invariant). Let $[eqn]$ be an admissible ODE for the network $[eqn]$ . Every polysynchrony subspace given by a balanced coloring is flow-invariant under $[eqn]$ .

The balanced colorings in Figure 1 from left to right correspond to the polysynchrony subspaces

$$[eqn]$$ $$[eqn]$$ $$[eqn]$$

Coloring every node the same is always balanced, and the associated polysynchrony subspace is called the fully synchronous subspace, which we denote $[eqn]$ . If $[eqn]$ , we say that the system is fully synchronous.

Eigenvalues depend on cell-level dynamics and network structure

2.3

Consistent with our general assumptions discussed above, we assume that any admissible ODE of a regular network with $[eqn]$ nodes takes the form

$$[eqn]$$

where $[eqn]$ is a bifurcation parameter.

Let $[eqn]$ denote the total state of the network, and suppose $[eqn]$ . If we write the arguments of $[eqn]$ as

$$[eqn]$$

then the internal and coupled dynamics of the system (respectively) are given by

$$[eqn]$$ $$[eqn]$$

where $[eqn]$ represents the differential of $[eqn]$ with respect to the variable $[eqn]$ . For example,

$$[eqn]$$

We remove the arguments $[eqn]$ of $[eqn]$ and $[eqn]$ if it is clear where they are being evaluated.

Definition 6 (Adjacency Matrix). The adjacency matrix of a regular network $[eqn]$ with n nodes is a $[eqn]$ matrix $[eqn]$ with entries $[eqn]$ defined to be the number of arrows from node $[eqn]$ to node $[eqn]$ $[eqn]$ .

Proposition 2 (Eigenvalues of System via Reduced Matrices). Suppose $[eqn]$ is a regular network with adjacency matrix $[eqn]$ and admissible $[eqn]$ . If $[eqn]$ are the eigenvalues of $[eqn]$ (not necessarily distinct) with eigenvectors $[eqn]$ , and $[eqn]$ is the Jacobian of $[eqn]$ at a synchronous equilibrium $[eqn]$ , then the eigenvalues of $[eqn]$ are the union of the eigenvalues of $[eqn]$ where $[eqn]$ and $[eqn]$ are the internal and coupled dynamics, respectively. Furthermore, the corresponding eigenvectors are $[eqn]$ where $[eqn]$ is an eigenvector of $[eqn]$ .

Proof. See (Golubitsky and Stewart, 2023).

The critical eigenvalue dictates the preferred pattern

2.4

We assume that pattern formation occurs from a synchrony-breaking bifurcation and argue the first bifurcating pattern is the “preferred” pattern of a tissue. When the homogeneous steady-state is destabilized via a bifurcation, the critical eigenvalue corresponds to a pattern (Theorem 3), which is generically the only stable bifurcating pattern (Proposition 4). We assume that the bifurcation parameter $[eqn]$ changes much slower than the state of the cells. With a quasistatically (slowly) varying $[eqn]$ , the uniform state can lose stability with a stable patterned state emerging; thus, the system will move from the homogeneous state to the pattern state and will remain there unless there is a secondary bifurcation on the pattern branch. Furthermore, as in Poston and Stewart (1978), we may assume that the only bifurcations present are those that are forced by our assumptions and observed data. Extraneous secondary bifurcations would violate Occam’s Razor – the principle that the simplest model which can explain observed phenomena is best. Even with secondary bifurcations, the preferred pattern is the only stable pattern for some region of the parameter space, and it is reasonable that unusual behaviors can arise for unrealistic parameter values.

Theorem 3 (Patterns Bifurcating from Synchrony). Let $[eqn]$ be a regular network with admissible $[eqn]$ . Let $[eqn]$ denote the fully synchronous subspace, and let $[eqn]$ be the polysynchrony subspace associated with the balanced coloring $[eqn]$ . Assume that there is a synchronous equilibrium $[eqn]$ and that $[eqn]$ has a critical eigenvalue $[eqn]$ . Denote the associated generalized eigenspace $[eqn]$ . If $[eqn]$ and $[eqn]$ , then generically a unique branch of equilibria with synchrony pattern $[eqn]$ bifurcates from $[eqn]$ at $[eqn]$ .

Proof (modified from (Golubitsky and Stewart, 2023)). Since $[eqn]$ is nonsingular. Thus, by the implicit function theorem, there exists a synchronous branch of equilibria in the neighborhood of $[eqn]$ . By a change of coordinates, we can assume that the synchronous branch is given by $[eqn]$ . Since $[eqn]$ , the kernel of $[eqn]$ is 1-dimensional, so we can use Liapunov-Schmidt reduction to find a reduced equation of the restriction $[eqn]$ given by

$$[eqn]$$

for some vectors $[eqn]$ . The zeros of $[eqn]$ are in one-to-one correspondence with the zeros of $[eqn]$ in a neighborhood of the bifurcation.

Since the spatial domain and codomain of $[eqn]$ are 1-dimensional, we can write

$$[eqn]$$

for $[eqn]$ and $[eqn]$ . Liapunov-Schmidt can be chosen to preserve the existence of the trivial solution, so we may assume $[eqn]$ .

The eigenvalue crossing condition and reduction imply that $[eqn]$ is generically nonzero. Thus, by the implicit function theorem $[eqn]$ for some function $[eqn]$ of $[eqn]$ . Since we can write $[eqn]$ as a function of $[eqn]$ in a neighborhood of the bifurcation, there exists a branch of solutions to $[eqn]$ and thus to $[eqn]$ . □

Proposition 4 (Subsequent Branches are Unstable). Suppose $[eqn]$ is a regular network with admissible $[eqn]$ . Let $[eqn]$ be the adjacency matrix with eigenvalues $[eqn]$ , and assume there exists a synchronous branch of equilibria $[eqn]$ that is stable for $[eqn]$ (for some $[eqn]$ ). Suppose that $[eqn]$ has a critical eigenvalue that is an eigenvalue of the matrix $[eqn]$ . Furthermore, suppose that for $[eqn]$ has no critical eigenvalues that are eigenvalues of $[eqn]$ . If there is a bifurcation at $[eqn]$ with critical eigenvalue from $[eqn]$ for $[eqn]$ , then generically the bifurcating branches are unstable.

Proof. We initially assume that the synchronous equilibrium is stable. If there is a bifurcation at $[eqn]$ with a critical eigenvalue of $[eqn]$ , the synchronous equilibrium will generically lose stability. In particular, there will be a positive eigenvalue $[eqn]$ that is an eigenvalue of $[eqn]$ . If there are no critical eigenvalues of $[eqn]$ for $[eqn]$ , then $[eqn]$ will remain positive at the bifurcation point $[eqn]$ . Since the eigenvalues of the Jacobian $[eqn]$ continuously depend on ( $[eqn]$ ), the bifurcating branches will be unstable in some neighborhood of the bifurcation point. □

Other preliminaries

2.5

For our results, we will use the definition of a strongly connected graph and the well-known Perron-Frobenius theorem.

Definition 7 (Strongly Connected Network). A network $[eqn]$ is called strongly connected if for every pair of nodes $[eqn]$ with $[eqn]$ , there is a path from $[eqn]$ to $[eqn]$ and from $[eqn]$ to $[eqn]$ .

Theorem 5 (Perron-Frobenius). If a network $[eqn]$ is strongly connected, then the adjacency matrix $[eqn]$ has a simple positive eigenvalue $[eqn]$ equal to the spectral radius $[eqn]$ . For all other eigenvalues $[eqn]$ *of $[eqn]$ *, $[eqn]$ .

Proof. See (Golubitsky and Stewart, 2023).

Results

3

The internal and coupled dynamics determine the preferred pattern

3.1

The adjacency matrix of the cell-communication network determines the possible patterns that can form in the tissue. Qualitative features of the internal and coupled dynamics select the preferred pattern from all possible ones.

We assume that $[eqn]$ is a strongly connected regular network on $[eqn]$ nodes without self-arrows (i.e. for any $[eqn]$ there is no $[eqn]$ with $[eqn]$ and $[eqn]$ ). Denote the adjacency matrix $[eqn]$ . We consider an arbitrary admissible ODE $[eqn]$ with a synchronous branch of equilibria $[eqn]$ that is stable for $[eqn]$ . We assume that the Jacobian of the ODE is singular at $[eqn]$ . Lastly, we take $[eqn]$ and $[eqn]$ to be the linearized internal and coupled dynamics (respectively), which are defined whenever $[eqn]$ (Section 2.3).

Definition 8 (Critical Pattern Space). Suppose that $[eqn]$ has distinct eigenvalues $[eqn]$ and corresponding eigenvectors $[eqn]$ . Let $[eqn]$ be given by

$$[eqn]$$

Suppose that $[eqn]$ has critical eigenvalues that by Proposition 2 are null eigenvalues of $[eqn]$ for some $[eqn]$ . Then define the critical pattern space to be the subspace of $[eqn]$ given by

$$[eqn]$$

Proposition 6 (Stability of Synchronous 1-Dimensional Nodes). If the node space is 1-dimensional and the synchronous equilibrium $[eqn]$ is stable, then $[eqn]$ .

Lemma 6.1. The adjacency matrix $[eqn]$ always has a strictly positive eigenvalue and an eigenvalue with strictly negative real part.

Proof of lemma. The adjacency matrix of a regular graph always has a unique largest eigenvalue equal to its valence $[eqn]$ (see Definition 1). Furthermore, since we assume our graph has no self-arrows, every element on the diagonal of $[eqn]$ is zero, and $[eqn]$ . Since the trace is the sum of the eigenvalues and there exists a positive eigenvalue, there must be some eigenvalues with negative real part. In particular, the eigenvalue with smallest real part, which we denote $[eqn]$ , has negative real part. □

Proof of proposition. Recall that the eigenvalues of the system are $[eqn]$ where $[eqn]$ are the eigenvalues of the adjacency matrix $[eqn]$ (Proposition 2). By Lemma 6.1, $[eqn]$ has a positive eigenvalue $[eqn]$ and an eigenvalue $[eqn]$ with negative real part. Let each eigenvalue be written

$$[eqn]$$

where $[eqn]$ . For the remainder of this proof, we will drop the arguments of the internal and coupled dynamics, $[eqn]$ and $[eqn]$ (respectively), writing

$$[eqn]$$ $$[eqn]$$

Since $[eqn]$ is stable, the real part of each eigenvalue $[eqn]$ for all $[eqn]$ . By Lemma 6.1, we can choose $[eqn]$ such that $[eqn]$ . Then

$$[eqn]$$ $$[eqn]$$ $$[eqn]$$

□

Proposition 7 (Stability of Synchronous 2-dimensional Nodes). Let the node space be 2-dimensional, and suppose that $[eqn]$ has distinct real eigenvalues $[eqn]$ , not counting multiplicity. If the synchronous branch $[eqn]$ is stable, then $[eqn]$ .

Proof. Since $[eqn]$ is stable,

$$[eqn]$$ $$[eqn]$$

for $[eqn]$ . By Lemma 6.1, there exists positive and negative eigenvalues of $[eqn]$ . Let $[eqn]$ be an eigenvalue with the same sign as $[eqn]$ . Then

$$[eqn]$$

□

Lemma 7.1. Suppose that $[eqn]$ is a polynomial in $[eqn]$ smoothly parameterized by $[eqn]$ . Let $[eqn]$ be a set of points in $[eqn]$ . Let $[eqn]$ , and suppose that for all $[eqn]$ we have $[eqn]$ for all $[eqn]$ . Suppose that at $[eqn]$ is linear in $[eqn]$ and there exists some $[eqn]$ such that $[eqn]$ . Then we have the following cases

$[eqn]$ for all $[eqn]$ if and only if $[eqn]$ . $[eqn]$ and $[eqn]$ if and only if $[eqn]$ . $[eqn]$ and $[eqn]$ if and only if $[eqn]$ .

Proof. For the forward direction of (1), if $[eqn]$ for all $[eqn]$ , then $[eqn]$ since it is a linear function of $[eqn]$ ; therefore, $[eqn]$ . For the converse, suppose that $[eqn]$ . Then $[eqn]$ . By hypothesis, there exists some $[eqn]$ such that $[eqn]$ . Clearly, this implies that $[eqn]$ . Thus, $[eqn]$ , so $[eqn]$ for all $[eqn]$ .

For the forward direction of (2), suppose that $[eqn]$ and $[eqn]$ for all $[eqn]$ . This implies that for any $[eqn]$ ,

$$[eqn]$$

Equivalently,

$$[eqn]$$

Since $[eqn]$ , we must have that $[eqn]$ .

For the converse, suppose that $[eqn]$ . If

$$[eqn]$$

then since $[eqn]$ ,

$$[eqn]$$

but this contradicts the continuity of $[eqn]$ with respect to $[eqn]$ since $[eqn]$ for $[eqn]$ . Hence, $[eqn]$ must be the unique point from $[eqn]$ with $[eqn]$ .

For the forward direction of (3), suppose that $[eqn]$ and $[eqn]$ for all $[eqn]$ . This implies that for any $[eqn]$ ,

$$[eqn]$$

Equivalently,

$$[eqn]$$

Since $[eqn]$ , we must have that $[eqn]$ .

For the converse, suppose that $[eqn]$ . If

$$[eqn]$$

then since $[eqn]$ ,

$$[eqn]$$

but this contradicts the continuity of $[eqn]$ with respect to $[eqn]$ since $[eqn]$ for $[eqn]$ . Hence, $[eqn]$ must be the unique point from $[eqn]$ with $[eqn]$ . □

Theorem 8 (Characterizing the First Branch in 1-Dimension). Let $[eqn]$ have eigenvalues $[eqn]$ ordered by real part, and assume that the node space is 1-dimensional. Then we have the following mutually exclusive cases

The critical pattern space is $[eqn]$ if and only if $[eqn]$ .The critical pattern space is $[eqn]$ , and the bifurcation is synchrony-breaking if and only if $[eqn]$ .The critical pattern space is $[eqn]$ (a vector of ones of length $[eqn]$ ), and the bifurcation is synchrony-preserving if and only if $[eqn]$ .

Proof. Let $[eqn]$ be the eigenvalues of $[eqn]$ with distinct real parts $[eqn]$ . Take $[eqn]$ to be the linear polynomial in $[eqn]$ given by

$$[eqn]$$

so that $[eqn]$ gives the real part of some eigenvalue of $[eqn]$ . Since we assume that the synchronous branch is stable before the bifurcation, there exists $[eqn]$ such that for all $[eqn]$ for all $[eqn]$ . Furthermore, at the bifurcation point, there exists some $[eqn]$ with $[eqn]$ . Therefore, the hypotheses of Lemma 7.1 are satisfied.

Applying Lemma 7.1, we have the following three cases.

1. $[eqn]$ for all $[eqn]$ if and only if $[eqn]$ .
2. $[eqn]$ and $[eqn]$ if and only if $[eqn]$ .
3. $[eqn]$ and $[eqn]$ if and only if $[eqn]$ .

$[eqn]$ exactly when the system has the critical eigenvalue $[eqn]$ . Using Definition 8 the main conclusion immediately follows.

Lastly, the unique eigenvalue of $[eqn]$ with largest real part $[eqn]$ is the unique eigenvalue with eigenvector $[eqn]$ (a vector of ones). By Proposition 2 and Theorem 3, the bifurcation with critical eigenvalue $[eqn]$ is synchrony-preserving, whereas bifurcations with other critical eigenvalues are generically synchrony-breaking. □

Theorem 9 (Characterization of Bifurcations in 2-dimensions with $[eqn]$ ). Suppose that the node space is 2-dimensional and the adjacency matrix $[eqn]$ has real eigenvalues $[eqn]$ , where each $[eqn]$ has algebraic multiplicity $[eqn]$ and geometric multiplicity $[eqn]$ . Suppose that $[eqn]$ , and take

$$[eqn]$$

In Table 1, we enumerate all possible critical pattern spaces, the multiplicity and type of critical eigenvalues, and the dimension of the critical eigenspaces. For the first four cases listed, we give necessary and sufficient conditions on the internal and coupled dynamics for such a bifurcation to occur. For the remaining cases, we give necessary conditions, which are sufficient if we assume the bifurcation is sufficiently degenerate.

Proof. We want to use the trace and determinant to obtain information about the eigenvalues; therefore, define the following polynomials of $[eqn]$ :

$$[eqn]$$

We will drop the argument $[eqn]$ if it is obvious where the polynomial is being evaluated. Since we assume that $[eqn]$ , we have that

$$[eqn]$$

We also assume the synchronous equilibrium is initially stable, so there exists $[eqn]$ such that for all $[eqn]$ ,

$$[eqn]$$ $$[eqn]$$

for all $[eqn]$ . This setup is illustrated in Figure 3. For a bifurcation to happen, there exists $[eqn]$ such that either $[eqn]$ or $[eqn]$ (or both). Applying lemma 3.1.2 (to $[eqn]$ and $[eqn]$ ), we have that

1. $[eqn]$ and $[eqn]$ if and only if $[eqn]$ .
2. $[eqn]$ and $[eqn]$ if and only if $[eqn]$ .
3. $[eqn]$ and $[eqn]$ if and only if $[eqn]$ .
4. $[eqn]$ and $[eqn]$ if and only if $[eqn]$ .
5. $[eqn]$ for all $[eqn]$ if and only if $[eqn]$ .
6. $[eqn]$ for all $[eqn]$ if and only if $[eqn]$ .

Now, let’s consider some of the possible cases. Notice that by equations (1) and (2) and the continuity of both polynomials with respect to $[eqn]$ , if $[eqn]$ , then we must have $[eqn]$ . Similarly, if $[eqn]$ , then we must have $[eqn]$ .

If $[eqn]$ and $[eqn]$ , there are $[eqn]$ complex conjugate pairs critical eigenvalues that are eigenvalues of $[eqn]$ ; if $[eqn]$ and $[eqn]$ , then there are $[eqn]$ real critical eigenvalues that are an eigenvalue of $[eqn]$ ; and if $[eqn]$ , there are $[eqn]$ real critical eigenvalues that are the two eigenvalues of $[eqn]$ with multiplicity $[eqn]$ .

Lastly, by the Perron-Frobenius Theorem, $[eqn]$ . We use this information to enumerate all possibilities in Table 1. □

Theorem 10 (Nondegeneracy Conditions in 2-Dimensions). Assume that $[eqn]$ has valence $[eqn]$ . Suppose that the node space is 2-dimensional and the adjacency matrix $[eqn]$ has distinct real eigenvalues $[eqn]$ and corresponding algebraic multiplicities $[eqn]$ . Suppose that $[eqn]$ and take

$$[eqn]$$

1. Suppose that at $[eqn]$ where $[eqn]$ denotes the logical “or.” If

$$[eqn]$$

then the critical eigenvalues are the pair of imaginary eigenvalues from a single matrix $[eqn]$ , with multiplicity $[eqn]$ . 2. Suppose that at $[eqn]$ . If

$$[eqn]$$

then the critical eigenvalues are a real eigenvalue from a single matrix $[eqn]$ with multiplicity $[eqn]$ .

Proof. Since we can obtain important information about the eigenvalues of a 2×2 matrix using the trace and determinant, define the following polynomials of $[eqn]$ :

$$[eqn]$$

Again, this setup is illustrated in Figure 3. We will drop the arguments $[eqn]$ if it is obvious where the polynomial is being evaluated. Since we assume that $[eqn]$ , we have that

$$[eqn]$$

Suppose that $[eqn]$ at the bifurcation. Of course, this implies that $[eqn]$ . Furthermore, by the Perron-Frobenius theorem $[eqn]$ is the spectral radius of $[eqn]$ (i.e. the largest absolute value of its eigenvalues). Therefore, for any $[eqn]$ ,

$$[eqn]$$ $$[eqn]$$

and

$$[eqn]$$

Since $[eqn]$ for all $[eqn]$ , the critical eigenvalues must correspond to pairs of imaginary eigenvalues, given by matrices $[eqn]$ such that $[eqn]$ . In the case that $[eqn]$ , there is exactly one root of the line $[eqn]$ , implying that at the bifurcation, there is a unique $[eqn]$ such that $[eqn]$ . Thus, there is a pair of imaginary critical eigenvalues that are eigenvalues of $[eqn]$ . Since $[eqn]$ is an eigenvalue of $[eqn]$ with multiplicity $[eqn]$ , the critical eigenvalue is an eigenvalue of $[eqn]$ with multiplicity $[eqn]$ . Now, notice that $[eqn]$ implies $[eqn]$ (otherwise the synchronous branch would be unstable before the bifurcation), and thus our analysis of the case that $[eqn]$ reduces to the previous case that $[eqn]$ .

For (2), suppose that $[eqn]$ at the bifurcation. Then $[eqn]$ . Furthermore, for any $[eqn]$ ,

$$[eqn]$$ $$[eqn]$$

and

$$[eqn]$$

Since $[eqn]$ for all $[eqn]$ at the bifurcation, any critical eigenvalues must be real and given by matrices $[eqn]$ with $[eqn]$ . Using the same logic as previously, if $[eqn]$ , there can be exactly one root of the line $[eqn]$ , implying that at the bifurcation there is a unique $[eqn]$ with $[eqn]$ , giving a real critical eigenvalue that is an eigenvalue of $[eqn]$ . Since $[eqn]$ is an eigenvalue of $[eqn]$ with multiplicity $[eqn]$ , the critical eigenvalue is an eigenvalue of $[eqn]$ with multiplicity $[eqn]$ . We also have that $[eqn]$ implies that $[eqn]$ (otherwise the synchronous branch would be unstable before the bifurcation), thus reducing to the prior case. In both cases, since $[eqn]$ and $[eqn]$ the critical eigenvalue is real. □

We have shown that for a 2-dimensional node space, if $[eqn]$ – as is the case for Notch signaling (Section 3.2.1) – then any stable pattern that forms generically corresponds to the critical pattern space $[eqn]$ . If $[eqn]$ , however, then a stable pattern can correspond to any critical pattern space $[eqn]$ as shown below. This indicates two possible mechanisms for the diversity of patterns we see in biological systems: cells reuse chemical signaling pathways but reorganize themselves to change communication, or cells use different chemical signaling pathways.

Theorem 11 (Existence of Admissible ODE for Arbitrary First Branch, Golubitsky and Stewart). Suppose the adjacency matrix $[eqn]$ has distinct real eigenvalues $[eqn]$ . Then for any $[eqn]$ , there exists an admissible system with 2-dimensional node space such that the critical pattern space is $[eqn]$ .

Proof. See (Golubitsky and Stewart, 2023).

Theorem 12 (Sufficient Condition for Pattern in $[eqn]$ or $[eqn]$ ). Suppose that the node space is 2-dimensional and the adjacency matrix $[eqn]$ has distinct real eigenvalues $[eqn]$ not counting multiplicity. If $[eqn]$ , then the critical pattern space is either $[eqn]$ , or $[eqn]$ .

Proof. As before, let

$$[eqn]$$

If $[eqn]$ or $[eqn]$ , then there are critical eigenvalues from $[eqn]$ for all $[eqn]$ , and the critical pattern space is $[eqn]$ .

Now suppose that $[eqn]$ and $[eqn]$ . Then, there is a critical eigenvalue that is an eigenvalue of $[eqn]$ if and only if either $[eqn]$ or $[eqn]$ . From Lemma 7.1, if $[eqn]$ then $[eqn]$ .

$[eqn]$ is a quadratic polynomial in $[eqn]$ whose graph is downward opening (Figure 3), so its derivative is linear with slope 2 $[eqn]$ . Suppose that $[eqn]$ for $[eqn]$ for contradiction. We have two cases:

1. If $[eqn]$ , then $[eqn]$ since $[eqn]$ . But this contradicts the continuity of $[eqn]$ in $[eqn]$ since $[eqn]$ for $[eqn]$ , otherwise the branch $[eqn]$ would not be stable.
2. Similarly, if $[eqn]$ , then $[eqn]$ since $[eqn]$ , which again contradicts the continuity of $[eqn]$ as a function of $[eqn]$ . □

In conclusion, since $[eqn]$ or $[eqn]$ can only occur for $[eqn]$ or $[eqn]$ , there can only be critical eigenvalues from $[eqn]$ or $[eqn]$ (or both), which gives us our conclusion by Proposition 2 and Definition 8.

Predicting patterns from qualitative features of cell-communication and chemical kinetics

3.2

We use theorems in Section 3.1 to predict patterns of cell fates given qualitative features of Notch signaling and the cell-communication network.

Critical pattern space of Notch signaling is Pμ1

3.2.1

Following the model of Collier et al. (1996), the feedback between Delta and Notch in coupled cells can be represented as in Figure 4 (top) where Notch inhibits Delta within a cell and Delta activates Notch in neighboring cells. It is natural to assume that $[eqn]$ activates $[eqn]$ (represented with $[eqn]$ ) means that

$$[eqn]$$

and $[eqn]$ inhibits $[eqn]$ (represented with $[eqn]$ ) means that

$$[eqn]$$

Furthermore, we assume that Delta and Notch both decay over time, meaning that

$$[eqn]$$ $$[eqn]$$

Writing the concentrations of the biochemical species as a vector $[eqn]$ , the internal and coupled dynamics of the Delta-Notch signaling pathway take the form

$$[eqn]$$

where + denotes a positive term, and – denotes a negative term. Thus, we have $[eqn]$ ; if $[eqn]$ has real eigenvalues, we can apply Theorem 9.

We also have that

$$[eqn]$$

so from Table 1, the critical pattern space is $[eqn]$ .

Preferred pattern in mutated C. elegans vulval precursor cells (VPCs)

3.2.2

The C. elegans vulval precursor cells (VPCs) form a pattern that is mediated by Notch signaling. With a mutation of let-23, all VPCs receive approximately equal external signals (Aroian and Sternberg, 1991), so the system can be represented with the regular network in Figure 4 (bottom). The adjacency matrix is

$$[eqn]$$

Using Matlab, we find that the eigenvalues are real, and the smallest eigenvalue and its corresponding eigenvector are

$$[eqn]$$

implying that the critical pattern space $[eqn]$ has geometric multiplicity 1, so the critical eigenspace is 1-dimensional (Table 1), and the critical pattern space $[eqn]$ corresponds to the balanced coloring $[eqn]$ shown in Figure 4 (bottom), so Theorem 3 implies that there is a bifurcating branch of solutions with the pattern $[eqn]$ . Furthermore, Proposition 4 implies that this will generically be the only stable branch, so we expect the C. elegans vulva, with a mutation of let-23, to exhibit an alternating pattern of cell fates – and this is exactly what is observed experimentally.

Preferred patterns in square arrays of cells developing according to Notch signaling

3.2.3

We will compute the critical pattern space of two square arrays of cells that develop according to Notch signaling (Figures 5,6), showing that cells can change their communication to form a different pattern despite using the same signaling mechanism.

We have shown that the critical pattern space for the Notch pathway is $[eqn]$ , so we must compute $[eqn]$ to determine the preferred pattern of the tissue. $[eqn]$ depends on the cell-communication network. In Figure 5, we consider a 16 × 16 array of cells, where each cell is coupled to its neighbors as depicted on the left (for clarity, we omit arrows pointing towards the central cell), and there are periodic boundary conditions. Dotted lines represent a connection strength of 1, and solid lines represent a connection strength of 3. Since there are periodic boundary conditions, the adjacency matrix $[eqn]$ will be real and symmetric, so $[eqn]$ will have real eigenvalues. Using Matlab, we find that the smallest eigenvalue $[eqn]$ gives the checkerboard pattern on the top-right, which matches our simulations.

In Figure 6, we consider a 50 × 50 array of cells with couplings shown in the top panel and periodic boundary conditions. Using Matlab, the minimal eigenvalue $[eqn]$ has geometric multiplicity 2, and the critical pattern space is the 2-dimensional space spanned by the vectors depicted on the left of Figure 6. In Table 1, $[eqn]$ , and we do not attempt to show that $[eqn]$ for some coloring $[eqn]$ , so the hypotheses of Theorem 3 may fail. Thus, we cannot ensure that there is a solution with the pattern of the critical pattern space; however, our numerical simulations validate that the predicted pattern is realized for some parameter range, as shown on the right of Figure 6.

Inferring biochemical interactions from observed patterns

3.3

Finally, we show how the theory can be used to infer properties of the chemical signaling pathway from an observed pattern.

For a 16 × 16 array of cells communicating as in Figure 5 (left), if we observe the striped pattern in Figure 7 (left), then the critical pattern is neither $[eqn]$ , nor $[eqn]$ . Barring the extremely degenerate case that the critical pattern space is $[eqn]$ , we must have that $[eqn]$ (Theorem 12), so the development is due to a chemical signaling network with the connections given on the right of Figure 7 (top or bottom).

Now, suppose there is some regular array of uniform cells. If we have partial information about chemical signaling in a tissue, we can use the observed tissue pattern to obtain additional information about chemical signaling by using the following steps.

1. Use the partial information about chemical interactions to create matrices representing the internal dynamics $[eqn]$ and coupled dynamics $[eqn]$ with variables as entries.
2. Use Theorems 6 and 7 to determine the signs of several variables in $[eqn]$ and $[eqn]$ .
3. Check $[eqn]$ . Then reference the observed pattern against Table 1 to determine necessary conditions on $[eqn]$ and $[eqn]$ .
4. Use newfound information about $[eqn]$ and $[eqn]$ to determine the possible signs of additional variables in $[eqn]$ and $[eqn]$ .
5. Lastly, convert the information about the sign of entries in $[eqn]$ and $[eqn]$ to information about activation and inhibition of biochemical species.

For the first example, consider when one chemical $[eqn]$ in a cell influences the same chemical in a neighboring cell, and suppose that the chemical is independent of other chemical species. Then the node space is 1-dimensional, and we can represent the internal and coupled dynamics as in Figure 8 (top). This tells us that

$$[eqn]$$ $$[eqn]$$

For the uniform or synchronous state to be stable, we must have that $[eqn]$ (Proposition 6), implying that $[eqn]$ decays. Assume that $[eqn]$ decays at the same rate for all time. Then for a pattern to form, there must be a synchrony-breaking bifurcation, which by Theorem 8 occurs if and only if $[eqn]$ , meaning that $[eqn]$ in one cell must inhibit $[eqn]$ in its neighboring cells, and the inhibition must increase for a synchrony-breaking bifurcation to occur.

Now, suppose we know that the concentration of $[eqn]$ is affected by the concentration of another chemical $[eqn]$ within the cell (Figure 8, middle). Then writing the chemical species as a vector $[eqn]$ ,

$$[eqn]$$ $$[eqn]$$

for some real numbers $[eqn]$ .

If the cells are initially uniform, then from Proposition 7,

$$[eqn]$$

This is satisfied if both chemical species decay – as is standard. Furthermore, notice that

$$[eqn]$$

so we can apply Theorem 9. If the cells oscillate in sync, then we may assume that there was a synchrony-preserving Hopf bifurcation, implying that $[eqn]$ (Table 1, row 4), so $[eqn]$ . If they oscillate in a pattern, then we may assume there was a synchrony-breaking Hopf bifurcation, implying that $[eqn]$ (Table 1, row 2) and thus $[eqn]$ . If we see that a steady-state pattern has formed among the cells, then we can infer that a synchrony-breaking steady-state bifurcation occurred (Table 1, row 1) meaning that

$$[eqn]$$

Assuming that $[eqn]$ decays, $[eqn]$ , implying that $[eqn]$ and the strength of decay of $[eqn]$ or the strength of inhibition must increase for a pattern to form.

Lastly, assume that $[eqn]$ in cell $[eqn]$ influences $[eqn]$ in neighboring cells $[eqn]$ as shown by the general diagram in Figure 8 (bottom). Writing the chemical species as a vector $[eqn]$ , we have that the internal and coupled dynamics are given by

$$[eqn]$$ $$[eqn]$$

for real numbers $[eqn]$ .

If the cells are initially uniform, then from Proposition 7,

$$[eqn]$$

Again, this is satisfied if both chemical species decay. Furthermore, notice that

$$[eqn]$$

so applying Theorem 9, we are in the first or third row of Table 1. This means there can only be a steady-state bifurcation (i.e. no oscillatory behavior). And

$$[eqn]$$

If $[eqn]$ , then synchrony cannot be broken since $[eqn]$ (see Table 1). If we observe a pattern, however, then we expect that $[eqn]$ , meaning that we either have

1. $[eqn]$ and $[eqn]$
2. or $[eqn]$ and $[eqn]$ .

Discussion

4

This framework provides powerful tools for molecular biologists who are interested in uncovering the mechanisms of pattern formation. It provides a systematic approach to identifying the molecular causes of pattern failures, which can help guide protein knockdown experiments in model organisms. Additionally, our theory can help uncover unknown molecular interactions in chemical signaling networks through the analysis of patterns that emerge in tissues (Figure 8).

This framework was created by recognizing that a developing tissue can be modeled as a system of ODEs on a regular network. When the linearization has null eigenvalues, a pattern can form through a bifurcation (Theorem 3). Using the network framework, we can determine the eigenvalues by examining smaller matrices $[eqn]$ , where $[eqn]$ represents an eigenvalue of the adjacency matrix $[eqn]$ , while $[eqn]$ and $[eqn]$ represent the linearized internal and coupled dynamics, respectively (Proposition 2). Thus, the pattern is determined by both the global cell-communication structure $[eqn]$ and the local cell-level dynamics ( $[eqn]$ and $[eqn]$ ).

We use our assumptions about a developing tissue to gain information about the pattern and chemical kinetics. Since the synchronous state is initially stable, all eigenvalues are negative, which provides crucial information about the chemical kinetics described by $[eqn]$ and $[eqn]$ (Propositions 6, 7). Then for a pattern to form, $[eqn]$ must have a null eigenvalue for some $[eqn]$ . We found that in nondegenerate systems with either one or two biochemicals and one coupling between cells, the resulting steady-state pattern always corresponds to the smallest eigenvalue of $[eqn]$ (denoted $[eqn]$ ) and is given by the critical pattern space $[eqn]$ (Definition 8). The formation of the pattern, however, requires specific conditions on $[eqn]$ and $[eqn]$ (Theorems 8, 9).

Altogether, our theory suggests that (1) if the chemical kinetics $[eqn]$ and $[eqn]$ satisfy certain conditions, then a pattern $[eqn]$ will form that is dictated by the cell-communication structure $[eqn]$ (Section 3.2); and (2) if a pattern $[eqn]$ forms in an array with known communication structure $[eqn]$ , then the chemical kinetics $[eqn]$ and $[eqn]$ must satisfy certain properties (Section 3.3).

Our framework also extends beyond basic research with several potential applications. In tissue engineering, our findings suggest that controlling the range of cell-communication can guide pattern formation (Section 3.2.3). Alternatively, maintaining the same communication structure while adjusting chemical signaling can alter the tissue pattern (Figure 7). In medicine, these insights have important implications for understanding disease mechanisms, as disruptions in pattern formation are associated with congenital malformations. We provide a theoretical framework for understanding the molecular changes underlying these failures. In particular, our theory illustrates that multiple molecular factors can lead to the same pattern breakdown. This helps explain why drugs can treat a disease for some individuals but not others and suggests novel drug targets.

Although our framework relies on simplified assumptions about biological systems, future research will refine and expand its applicability. A key next step will be extending our analysis from two-dimensional signaling networks to more complex, realistic biochemical systems. This expansion will improve our ability to identify the molecular causes of pattern formation. Additionally, we plan to validate our theoretical predictions through simulations that incorporate spatial features including morphogen gradients, parameter noise, and network perturbations, enabling us to apply our theory to a broader range of biological contexts.