The Mean-Field Survival Model for Stripe Formation in Zebrafish Exhibits Turing Instability

TL;DR

The paper develops and analyzes a ring-based mean-field survival ODE model for zebrafish stripe formation, governed by two coupled populations $X_j$ and $M_j$ on $N$ sites with nonlocal interactions at distance $h$. Through linear stability analysis and a discrete Fourier transform, it derives a critical polynomial $F(x,y,b,d)$ and a biologically motivated region $ ext{$\\mathcal P$}$ that dictates when Turing-like patterns emerge, showing that for large $N$ and $h$ the coexistence equilibrium becomes unstable to spatially periodic perturbations with wavenumbers in a specified band; otherwise the homogeneous state is stable. The work contrasts finite-$N$ ODE behavior with PDE limits, revealing significant discrepancies and providing explicit asymptotic formulas for neutral stability curves and their slopes, which yield quantitative predictions for stripe width and projection lengths consistent with in-vivo observations. Numerical simulations corroborate the analytic predictions, showing out-of-phase melanophore and xanthophore stripes and dominant Fourier modes that match the theoretically predicted bands. Overall, the study demonstrates that distant-neighbor signaling in a discrete lattice can robustly generate zebrafish-like stripes and offers a rigorous framework for analyzing pattern formation in coupled ODE networks at physiologically relevant scales.

Abstract

Zebrafish have been used as a model organism in many areas of biology, including the study of pattern formation. The mean-field survival model is a coupled ODE system describing the expected evolution of chromatophores coordinating to form stripes in zebrafish. This paper presents analysis of the model focusing on parameters for the number of cells, length of distant-neighbor interactions, and rates related to birth and death of chromatophores. We derive the conditions on these parameters for a Turing bifurcation to occur and show that the model predicts patterns qualitatively similar to those in nature. In addition to answering questions about this particular model, this paper also serves as a case study for Turing analysis on coupled ODE systems. The qualitative behavior of such coupled ODE models may deviate significantly from continuum limit models. The ability to analyze such systems directly avoids this concern and allows for a more accurate description of the behavior at physically relevant scales.

Figures (9)

• Figure 1: Schematic of ring geometry of model \ref{['1D ODE']}. Red arrows represent short-range interactions and blue represent long-range interactions.
• Figure 2: The curves $y=T_{h}(x)$ (solid black) and $F(x,y,b,d)=0$ (dashed black) in the $xy$-plane. The set with $F<0$ is the shaded region below and to the right of the $F=0$ curve. Turing instability occurs when one or more of the points $(x_{k},T_{h}(x_{k}))$ (black markers) lie in the region with $F<0$ for some $(b,d)$. All markers are labeled with the mode(s) that they correspond to. For both plots, $h=3$ and $N=6$. (\ref{['fig: h3N6b3d10 Cheb plot']}) For $(b,d)=(3,10)$, the equilibrium is linearly stable. (\ref{['fig: h3N6b3d100 Cheb plot']}) For $(b,d)=(3,100)$, the $k=\pm1$ modes are unstable. Changing the value of $b$ has a similar effect of changing the shape of the $F=0$ curve.
• Figure 3: The curves $y=T_{h}(x)$ (solid black) and $F(x,y,b,d)=0$ (dashed black) in the $xy$-plane. The set with $F<0$ is the shaded region below and to the right of the $F=0$ curve. Turing instability occurs when one or more of the points $(x_{k},T_{h}(x_{k}))$ (black markers) lie in the region with $F<0$ for some $(b,d)$. The rightmost marker corresponds to the $k=0$ mode, the next marker corresponds to the $k=\pm 1$ modes, etc. The modes which are unstable are labeled in all plots (in the inset for Figure \ref{['fig: h11N23b3d10 Cheb plot']}). $(b,d)=(3,10)$ for all plots. (\ref{['fig: h4N8b3d10 Cheb plot']}) For $h=4$ and $N=8$, the $k=\pm 1$ modes are unstable. (\ref{['fig: h4N23b3d10 Cheb plot']}) For $h=4$ and $N=23$, the $k=\pm 2,\pm 3$ modes are unstable. Different choices of $(b,d)$ may cause only the $k=\pm 2$ or only the $k=\pm 3$ modes to be unstable. (\ref{['fig: h11N23b3d10 Cheb plot']}) For $h=11$ and $N=23$, the $k=\pm 1,\pm 3$ modes are unstable. Unlike the unstable modes in (\ref{['fig: h4N23b3d10 Cheb plot']}), $k=\pm 3$ modes being unstable implies the $k=\pm 1$ modes are unstable. Inset depicts the region $[1/2,1]\times[-1,0]$ marked by a dotted rectangle in the original image.
• Figure 4: The curves, $F(x_{k},T_{h}(x_{k}),b,d)=0$, in the $bd$-plane. The bifurcation curve (solid blue) represents the boundary where the equilibrium is neutrally stable. The neutral stability curves (dashed blue) for specific modes represent where the equilibrium is neutrally stable to the given mode, but may already be unstable to other mode(s). The number next to each curve provides the corresponding Fourier mode(s). For each neutral stability curve, the set with $F<0$ is the blue shaded region above the curve. Darker blue shaded regions indicate instability to multiple Fourier modes. Each bifurcation curve has two analytically computable asymptotes (solid red) which are described in Theorem \ref{['thm: Asymptotes']}. The solid black line represents the boundary of $\mathcal{P}$ and the gray shaded region outside of $\mathcal{P}$ are the parameter values for which the equilibrium is always stable. (\ref{['fig: h4N8 Bifurcation']}) For $h=4$ and $N=8$, only the $k=\pm 1$ modes are unstable. (\ref{['fig: h4N23 Bifurcation']}) For $h=4$ and $N=23$, four modes are unstable. When $b$ is small, the $k=\pm 2$ modes lose stability first, but when $b$ is larger it is the $k=\pm 3$ modes. (\ref{['fig: h11N23 Bifurcation']}) For $h=11$, and $N=23$, four modes are unstable. The $k=\pm 1$ modes always lose stability before $k=\pm 3$. The asymptotes are omitted from this figure because they are visually indistinguishable from the bifurcation curve itself.
• Figure 5: Scatter plots of linear instability behavior in the $bh$-plane for fixed $d$ values. Each sample point is colored blue if the homogeneous state is unstable and red if it is stable The PDE model predicts instability for all $(b,d)$ above the PDE neutral stability curve (blue line). The largest and smallest unstable values are marked with larger dots. The shaded gray region with solid black boundary is the compliment of $\mathcal{P}$ in each plot, meaning that the equilibrium is always linearly stable in this region. For panels (\ref{['fig: d4N23 bhBifurcationPlot']})-(\ref{['fig: d4LargeN bhBifurcationPlot']}), the plots range over all values of $h$ for the given $N$. (\ref{['fig: d4N23 bhBifurcationPlot']}) For $d=4$ and $N=23$, the region of instability satisfies $0< b < 5$ and $5\leq h$. (\ref{['fig: d10N23 bhBifurcationPlot']}) For $d=10$ and $N=23$, the qualitative shape is similar to the $d=4$ case, but now there is an additional interval at the $h=4$ level and all intervals extend to larger values of $b$. (\ref{['fig: d4N100 bhBifurcationPlot']}) For $d=4$ and $N=100$, the instability region again includes an interval at $h=4$ and is bounded by $b=5$, but the intervals generally extend closer to this boundary. The rich structure of the right boundary demonstrates the value of performing LSA for all values of $h$ and $N$. (\ref{['fig: d4LargeN bhBifurcationPlot']}) For $d=4$ and $N=10000$, the left boundary of these intervals are qualitatively similar to the small $N$ cases, but the right boundary maintains none of the fine structure observed in the $N=100$ case.

Theorems & Definitions (14)

• Proposition 1
• Proposition 2
• Theorem 1
• Theorem 2
• Theorem 3
• Corollary 1
• Corollary 2
• Lemma 1
• proof : Proof of Theorem \ref{['thm: Discrete Bifurcation']}
• proof : Proof of Theorem \ref{['thm: Continuous Bifurcation']}