Feedback, receptor clustering, and receptor restriction to single cells yield large Turing spaces for ligand-receptor-based Turing models

Turing mechanisms can yield a large variety of patterns from noisy, homogenous initial conditions and have been proposed as patterning mechanism for many developmental processes. However, the molecular components that give rise to Turing patterns have remained elusive, and the small size of the parameter space that permits Turing patterns to emerge makes it difficult to explain how Turing patterns could evolve. We have recently shown that Turing patterns can be obtained with a single ligand if the ligand-receptor interaction is taken into account. Here we show that the general properties of ligand-receptor systems result in very large Turing spaces. Thus, the restriction of receptors to single cells, negative feedbacks, regulatory interactions among different ligand-receptor systems, and the clustering of receptors on the cell surface all greatly enlarge the Turing space. We further show that the feedbacks that occur in the FGF10-SHH network that controls lung branching morphogenesis are sufficient to result in large Turing spaces. We conclude that the cellular restriction of receptors provides a mechanism to sufficiently increase the size of the Turing space to make the evolution of Turing patterns likely. Additional feedbacks may then have further enlarged the Turing space. Given their robustness and flexibility, we propose that receptor-ligand-based Turing mechanisms present a general mechanism for patterning in biology.

Figure 1

Ligand-receptor interactions can give rise to Turing patterns. (a) Spatial patterns via a Turing mechanism can result from cooperative receptor-ligand interactions, where $m$ receptors ($R$) and $n$ ligand molecules ($L$) form an active complex that upregulates the receptor concentration by increasing its expression, limiting its turnover or similar. Importantly, the highest receptor and ligand concentrations are observed in different places. (b) In case of the standard network [panel (a)], Turing patterns emerge only for a small subset of the parameter range of the receptor and ligand production rates, $a$ and $b$. $a_{max}$ denotes the maximal value of the receptor production rate, while $b_{min}$ and $b_{max}$ denote the minimal and maximal ligand production rates. (c) Additional feedbacks (solid, red and dashed, blue arrows) can be mediated by the ligand-receptor complex, $R^{2}L$; $\leftrightarrow$ indicates receptor-ligand interactions, $⊣$ inhibitory interactions, and $-•$ up-regulating interactions. (d) The negative feedbacks in panel (c) (network U5 in Fig. S1 [23]) result in a larger Turing space when the response threshold $p$ is lowered from $p=1$ (blue shaded area) to $p=0.1$ (solid, yellow area). (e) The size of the Turing space for the network in panel (c) (network U5 in Fig. S1 [23]) increases as the response threshold $p$ is lowered. As a measure for the size of the Turing space, we record the maximum of the receptor production rate, $a_{max}$, and the ratio of the maximal and minimal ligand production rates $b_{\mathrm{ratio}}=\frac{b_{max}}{b_{min}}$, for which Turing patterns can emerge. $a=0$ is part of the Turing space and negative values of $a$ have no physiological interpretation.

Figure 2

Negative feedbacks by receptor-ligand complexes result in Turing patterns with large Turing spaces. (a) The simulated network architecture. Two receptors $R$ interact with one dimeric ligand $L$ to form a receptor-ligand complex $R^{2}L$ (black arrows, $\leftrightarrow$). The receptor-ligand complex upregulates the presence of receptor ($-•$). In addition to these core interactions that can result in a Turing mechanism, we considered negative feedbacks ($⊣$) on the ligand production (red, solid arrow) and/or the receptor production (blue dashed arrow). (b) A negative feedback on the receptor production rate [blue dashed arrow in panel (a)] increases the Turing parameter space for the receptor production rate, $a$ [blue squares in panel (b)], compared to the standard network [black part of the network in panel (a) and black star in panel (b)]. A negative feedback on ligand production [red, solid arrow in panel (a)] enlarges the Turing parameter space for the ligand production rate, $b$ [red circles in (b)]. In the presence of both feedbacks the Turing parameter space is enlarged along both axes [green triangles in panel (b)]. The feedback effects are stronger the lower the feedback threshold, $p$ ($p=0.01$, 0.1, 1, 10, 100). The gray arrow indicates the direction in which the feedback threshold, $p$, decreases.

Figure 3

Coupling of several receptor-ligand-based Turing modules further enlarges the Turing space. (a) The simulated network architecture. Two receptor-ligand-based Turing modules, as analyzed in Fig. 2 (black arrows, $\leftrightarrow$, $-•$), are coupled via additional negative feedbacks ($⊣$) on the ligand production rates (red solid arrows) and/or the receptor production rates (blue dashed arrows). (b) A negative feedback on the receptor production rate [dashed blue line in panel (a)] increases the Turing parameter space for the receptor production rate, $a$ [blue squares in panel (b)] compared to the standard network [black part of the network in Fig. 2 and black star in panel (b)]. A negative feedback on ligand production [red solid arrow in panel (a)] enlarges the Turing parameter space for the ligand production rate, $b$ [red circles in (b)]. In the presence of both feedbacks the Turing parameter space is enlarged along both axes [green triangles in panel (b)]. The feedback effects are stronger the lower the feedback threshold, $p$ ($p=0.01$, 0.1, 1, 10, 100). The gray arrow indicates the direction in which the feedback threshold, $p$, decreases.

Figure 4

Negative feedbacks enlarge the Turing space by limiting the effective production rates. The plot of the (a) effective receptor production rate $a_{\mathrm{eff}}=\frac{a}{max\left(R^{2}L\right)/p+1}$ versus the receptor production rate $a$, and (b) the plot of the effective ligand production rate $b_{\mathrm{eff}}=\frac{b}{max\left(R^{2}L\right)/p+1}$ versus $b$ show that, as a result of the negative feedbacks, the effective production rates remain in a narrow range, even as $a$ and $b$ are greatly changed. The calculation was carried out for the symmetrically coupled Turing system, shown in green in Fig. 3.

Figure 5

The restriction of receptors to single cells enlarges the Turing space. (a) Cartoon of the computational domain: diffusion of receptors is restricted to single cells, while ligand can diffuse over the entire computational domain. (b) Solution of the receptor-ligand model on a 1D, 2D and 3D (left to right) cellularized computational domain. The ligand (upper row), receptors (middle row), and ligand-receptor complexes (bottom row) pattern the domain. We provide the concentration levels (in arbitrary units) on the vertical axis for the 1D domain (left column), and intensities as color code (blue(dark)- low; red(light)- high) on the 2D and 3D domains. To distinguish cell boundaries on the 1D domain we alternate black and gray lines. (c) The size of the Turing space increases as the domain of fixed size is split into more cells, $N$. Triangles show the results for $N=10$ and $N=100$ cells. The black star reports the Turing space for the standard model, $N=1$. (d) Patterns of receptor-ligand complexes that extend over several cells can be obtained with a diffusive component, $T$, that is produced in response to the formation of receptor-ligand complexes, and that enhances the abundance of receptors on neighboring cells. The gray arrow indicates the direction, in which the feedback threshold, $p$, decreases.

Figure 6

Receptor clustering enlarges the Turing space. (a) The simulated network architecture. Clusters of $2n$ receptors $R$ interact with $n$ dimeric ligands $L$ to form a receptor-ligand complex ${\left(R^{2}L\right)}^n$ (black arrows, $\leftrightarrow$). The receptor-ligand complex upregulates the presence of receptor (black interaction, $-•$). In addition to these core interactions that can result in a Turing mechanism, we considered negative feedbacks on the ligand production (red solid arrow, $⊣$) and/or the receptor production (blue dashed arrow, $⊣$). (b) Higher cooperativity, $n>1$, as may result from larger receptor-ligand clusters further increases the size of the Turing space. The $n$-dependent increase was calculated for $p=0.01$, 0.1, 1, 10, 100 for case U5 in Fig. S1 [23]. The gray arrow indicates the direction, in which the feedback threshold, $p$, decreases.

Figure 7

Substantially enlarged Turing spaces for physiological networks. (a) The SHH-FGF10 network in the control of lung branching morphogenesis. For details see the text. (b) Schematic representation of the regulatory network for lung branching morphogenesis in panel (a). (c) The Turing space of such a physiological model is huge and further increases as the feedback threshold, $p$, is lowered. The red triangles represent the Turing spaces for $p_1=q=0.1$ (positive feedback on ligand and receptor, respectively) and $p_2$ = 0.01, 0.1, 1, 10, 100 (negative feedback); the black star represents the size of the Turing space of the standard network in Fig. 2 (black part). The gray arrow indicates the direction in which the feedback threshold, $p$, decreases. $\leftrightarrow$ indicates binding interactions, $⊣$ indicates inhibitory interactions, and $-•$ indicates up-regulating interactions.

Figure 8

Comparison of the Turing spaces calculated numerically and those derived analytically. [(a) and (b)] The shaded regions of the parameter space indicates the area where the linear stability analysis identifies a Turing instability (yellow, light shading) or other instabilities (navy, dark shading) for Eqs. (1) and (2) with zero-flux boundary conditions. The symbols indicate the points in the parameter space where the numerical solution of Eqs. (1) and (2) with zero-flux boundary conditions yielded either pattern formation (+) or not (0). $\gamma$ was chosen sufficiently large that Turing patterns could emerge on the 1D domain. Panels (a) and (b) differ in the relative diffusion coefficient $d$, with (a) $d=100$ and (b) $d=10$.

Figure 9

Convergence of the numerical solution. (a) Typical pattern of receptor-ligand complexes ($R^{2}L$) on a domain comprising two subdomains. Ligand is produced in the upper domain but free to diffuse on the entire domains. Receptor is produced in the lower domain and its diffusion is restricted to the lower domain. (b) The maximum deviation of the receptor-ligand complex ($R^{2}L$) as computed with an FEM mesh with element size equal to 0.01 from that computed at other mesh sizes.