Advances in Mathematical Physics
Volume 2016 (2016), Article ID 8938970, 7 pages
http://dx.doi.org/10.1155/2016/8938970

Research Article

Unstable Modes and Order Parameters of Bistable Signaling Pathways at Saddle-Node Bifurcations: A Theoretical Study Based on Synergetics

Till D. Frank1,2

¹CESPA, Department of Psychology, University of Connecticut, 406 Babbidge Road, Storrs, CT 06269, USA
²Department of Physics, University of Connecticut, 2152 Hillside Road, Storrs, CT 06269, USA

Received 27 July 2016; Accepted 25 September 2016

Academic Editor: Giorgio Kaniadakis

Copyright © 2016 Till D. Frank. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Mathematical modeling has become an indispensable part of systems biology which is a discipline that has become increasingly popular in recent years. In this context, our understanding of bistable signaling pathways in terms of mathematical modeling is of particular importance because such bistable components perform crucial functions in living cells. Bistable signaling pathways can act as switches or memory functions and can determine cell fate. In the present study, properties of mathematical models of bistable signaling pathways are examined from the perspective of synergetics, a theory of self-organization and pattern formation founded by Hermann Haken. At the heart of synergetics is the concept of so-called unstable modes or order parameters that determine the behavior of systems as a whole close to bifurcation points. How to determine these order parameters for bistable signaling pathways at saddle-node bifurcation points is shown. The procedure is outlined in general and an explicit example is worked out in detail.

1. Introduction

Bistable signaling pathways on the cellular level play a key role [1–7]. They can act as switches in general and determine cell fate and function as memory storage in particular [1–7]. While on a conceptual level there is a good understanding of bistable signaling pathways, there is still a need for improving our methodological tools to address such bistable components from a mechanistic point of view. Within the framework of ordinary differential equation (ODE) models, bistable signaling pathways have been described by single species or multiple species models. Single species models capture the pathway activity by a single variable. A benchmark model in this regard is a model for promoter activity under the impact of a positive feedback loop. The feedback loop is established by means of a transcription factor that acts as an activator [2, 4, 8, 9]. In the case of single variable models, the bistability is given in terms of a low and high state of the species concentration under consideration or the gene expression or promoter activity of interest. In what follows, the focus will be on two-species models as described by two coupled ODEs [2, 3, 5]. The two states of the bistable element typically describe that either species is dominant over species or species is dominant over species ; that is, there are two states of interest that will in what follows be denoted by (high), (low) and (low), (high). In order to function as a switch, the bistable component must feature a control parameter and the bistable component must be driven through a saddle-node bifurcation by varying the control parameter gradually. At the saddle-node bifurcation the bistable pathway becomes monostable such that only one of the two species can dominate the signaling pathway. In doing so, bistable models can a switch from (high), (low) to (low), (high) or vice versa from (low), (high) to (high), (low). While the saddle-node bifurcation has been studied using analytical and numerical solution methods [2, 3, 5, 10], in this context little is known about the so-called unstable mode or order parameter [11] emerging at this saddle-node bifurcation. According to the theory of synergetics, a theory of pattern formation and self-organization, the unstable mode at the bifurcation point determines the order of the whole two-species system and for this reason may be referred to as order parameter [11]. In the following section, the general approach to determine the unstable mode or order parameter will be described and subsequently an application to the bistable model proposed by Markevich et al. (2004) [5] will be presented.

2. Unstable Mode and Order Parameter of Bistable Signaling Pathways Involving Two Species

2.1. General Case

Let us first consider the general case of a two-species model of a bistable signaling pathway given in terms of two coupled ODEs. Accordingly, the pathway is described by the concentrations and of the two species under consideration measured in arbitrary units. Variations and over time are assumed to satisfy

where and are appropriately defined functions and depend on the pathway being considered. Equation (1) is assumed to involve a control parameter . Furthermore, (1) is assumed to be bistable in a particular interval of . At the boundaries of that interval a stable fixed point merges with an unstable fixed point giving rise to a saddle-node bifurcation. The stable and unstable fixed points merge and disappear and the pathway becomes monostable. Let us consider the saddle-node bifurcation in more detail. To this end, we “follow” the stable fixed point coordinates and involved in the saddle-node bifurcation along variations of . The linearized version of (1) at that stable fixed point reads

with and . The variables ε and δ denote small deviations from the fixed point (st), (st). The matrix coefficients correspond to the coefficients of the Jacobian matrix at and . The eigenvalues and of the matrix are given by

where the trace and the determinant of are given by

Note that in (3) comes with a plus sign in front of the square root function, whereas exhibits the minus sign. In general, at a saddle-node bifurcation of a two-dimensional dynamical system there is one eigenvalue equal to zero while the other is negative. Using the aforementioned plus and minus sign convention of (3), at the saddle-node bifurcation, we have , , , and . From it follows that

Moreover, substituting into (3) we obtain

It is useful to summarize this intermediate result. Let denote a critical value at which the saddle-node bifurcation occurs. Consequently, at we have and and (5) holds. If we shift next away from the critical value by a small amount such that the pathway is monostable, , where is small, then there are two directions in the two-dimensional state space. In one direction the dynamics evolves relatively slowly whereas in the direction orthogonal to that direction the dynamics evolves relatively fast because for we have as compared to . In line with synergetics, the direction associated with the slow dynamics corresponds to the unstable mode at and dominates the behavior of the whole system; that is, the mode can be considered as the order parameter [11]. In particular, the mode dominates the signaling pathway even when μ is not zero but a small parameter. This holds for perturbations that drive the system temporarily away from the fixed point (st), (st) in the bistable parameter domain and for the transient dynamics in the monostable domain that describes the transition from the disappearing fixed point (st), (st) to the alternative fixed point that is still stable in the monostable domain.

Let us determine the unstable mode. In general, that is, for arbitrary eigenvalues and , the eigenvectors normalized to unity are given in components by

for . At the saddle-node bifurcation point, we have and the corresponding eigenvector represents the unstable mode or order parameter. Consequently, at a saddle-node bifurcation the order parameter of a bistable signaling pathway described by a two-species model given in terms of two coupled ODEs is defined by

2.2. The Signaling Pathway Model by Markevich et al. (2004)

2.2.1. Definition of the Model

Markevich et al. (2004) [5] discussed the following two-species signaling pathway model:

with , , , , and . The concentration levels are normalized such that and are restricted to the interval . The model and generalizations of it have also been discussed in subsequent studies [10, 12].

2.2.2. Nullcline Approach

For sake of simplicity, we introduce the new parameters:

The nullclines are given by

with and . Using the nullclines we would like to illustrate the characteristic features of the pathway model. After that, a more detailed analysis will be conducted.

In order to illustrate the characteristic features of the pathway model, we consider first the symmetric case given by . In this case, the model can exhibit a single asymptotically stable fixed point (Figure 1) or three fixed points (Figure 2). In the latter case, the pathway model is bistable. The fixed point on the diagonal is unstable. Note that in the symmetric case the number of fixed points only depends on a single parameter as we will see later. In particular, at a critical value of , there is a pitchfork bifurcation connecting the monostable and bistable parameter domains. We are interested in the switch-like functioning of bistable pathway models and the order parameter that characterizes the switching from one state to another state. With respect to the model by Markevich et al. (2004) [5], we need to consider the asymmetric case. As we will see below, if is much larger than the model is monostable with the species dominating the species . In contrast, if is much larger than then the model is monostable with the species dominating the species . Consequently, if we start initially with a large ratio then the pathway exhibits a high state for and a low state for . If we gradually decrease the ratio, then the fixed point for (high) and (low) will converge towards the unstable fixed point of the pathway model. The two fixed points will merge and disappear in a saddle-node bifurcation. This is shown in Figures 3 and 4. The pathway will switch to the alternative state with a relatively high state for and a relatively low state for .

Figure 1: Nullclines as defined by (11) (dashed line) and (12) (solid line) for the monostable case. Parameters: (as in [10]) and . All parameters are in arbitrary units.

Figure 2: It is as in Figure 1 but for the bistable case. Parameters: and .

Figure 3: Illustration of the saddle-node bifurcation (part I) for the asymmetric pathway model. The stable fixed point given by (high) and (low) is about to merge with the unstable fixed point. Parameters: , , and .

Figure 4: Illustration of the saddle-node bifurcation (part II). The pathway has become monostable. Magnification (not shown) of the top left corner of the phase space reveals that the nullclines do not intersect each other in this area. This is because the asymmetry of the pathway model has increased as compared to Figure 3. Parameters: , , and .

Figures 1, 2, 3, and 4 illustrate the two basic mechanisms leading to the saddle-node bifurcation underlying the switch-like functioning of the model. First, in the symmetric case, the model is bistable. Second, when breaking the symmetry, in the asymmetric case, the bistable domain can be destroyed via a saddle-node bifurcation resulting in a switch from the state (low), (high) to the state (high), (low) or vice versa. Let us analyze these mechanisms in more detail. To this end, we derive a fixed point equation.

2.2.3. Cubic Fixed Point Equation and Bifurcation Diagram

In order to determine the fixed point semianalytically, we put the left-hand sides of (9) equal to zero. From the two relations in (9) we can obtain the following relations:

It is useful to introduce the new variables and like

Moreover, we introduce the new parameters and like

The parameter measures the degree of asymmetry. In the symmetric case we have . From (14) we obtain the intermediate result:

and finally

In deriving (18) we used . The variable as a function of defined by (18) is a monotonically decaying function for any parameter . For we have . Importantly, from (17), it follows that for any parameter the solutions of (18) are symmetric; that is, if and is a solution, then and is a solution as well. From (14) it follows that

with

Equation (19) is a cubic equation in the variable . Therefore, there are maximal three solutions. If then the equation is symmetric in and . Therefore, in view of the symmetry property of (18) for any parameter and the symmetry property of (19) in the special case , we conclude that if then the number of solutions is either 1 or 3. From (19) it follows that in the symmetric case (which implies by ) the parameter acts as control parameter and induces a pitchfork bifurcation. For example, for () and numerical evaluation of (19) shows (see below) that for (19) exhibits three admissible solutions (i.e., solutions in the interval between 0 and 1).

For any given set of parameters , and , the stable and unstable fixed points and of the signaling pathway model (9) can be determined numerically by finding the solutions (roots) of the cubic equation defined by (19) in combination with (18) and (20). In doing so, the bifurcation diagram can be constructed when defining a suitable control parameter . In order to focus on switch-like functioning we vary the degree of asymmetry like

with . In particular, for we have the symmetric case with and . We choose the parameters and sufficiently large such that in this (symmetric) case the model is bistable. Figure 5 presents the bifurcation diagram thus obtained for a fixed value of . Figures 6, 7, and 8 depict the functions and of the implicit equation (19) for the symmetric case (Figure 7) and the critical values at which saddle-node bifurcations occur (Figures 6 and 8).

Figure 5: Bifurcation diagram of the bistable pathway model (9) acting as a switch for a particular value of when varying the degree of asymmetry. The stationary solutions are plotted versus the control parameter values . The solutions were computed from (18), (19), (20), and (21). Solid line: asymptotically stable fixed points of the state (low), (high). Dotted line: unstable fixed points in the bistable parameter domain. Dashed line: asymptotically stable fixed points of the state (high), (low). Parameters: , , and .

Figure 6: Illustration of the saddle-node bifurcation using the cubic fixed point equation (19). The functions (solid line) and (dashed line) of (19) are plotted versus at the saddle-node bifurcation point in which the states (high) and (low) emerge or disappear. Parameters: ; all other parameters are as in Figure 5.

Figure 7: The functions (solid line) and (dashed line) of (19) in the symmetric case. The signaling pathway model is bistable. Parameters: (); all other parameters are as in Figure 5.

Figure 8: Saddle-node bifurcation illustrated by means of the cubic fixed point equation (19). The functions (solid line) and (dashed line) of (19) are plotted versus at the saddle-node bifurcation point in which the state (low), (high) emerges or disappears. Parameters: ; all other parameters are as in Figure 5.

2.2.4. Order Parameters

Our next objective is to determine the order parameters or unstable modes at the saddle-node bifurcation points shown in the bifurcation diagram; see Figure 5. From (1), (2), and (9) it follows that the coefficients of the Jacobi matrix read

We computed the matrix coefficients for the asymptotically stable fixed points (low), (high) and (high), (low) shown in Figure 5 as functions of the control parameter . Using (3) we computed the eigenvalues. Figure 9 presents the eigenvalues as functions of . We found that both eigenvalues were real-valued in the considered parameter domain. Importantly, as expected, converged to zero when approaching the saddle-node bifurcation points, whereas assumed a finite negative value. The convergence of one of the eigenvalues to zero can be seen even more clearly when studying the determinant of the Jacobian matrix as function of α. We computed from (4). The result is shown in Figure 10. As expected, converged to zero at the saddle-node bifurcation points. Finally, we computed the eigenvector associated with from (7) for the asymptotically stable fixed points (low), (high) and (high), (low) shown in Figure 5 and as functions of the control parameter . Figure 11 presents the eigenvector components thus obtained. The order parameters of the signaling pathway model correspond to the eigenvectors at the two saddle-node bifurcation points. These order parameters (order parameter vectors) are highlighted in Figure 11 by circles. Figure 11 demonstrates explicitly how the general methodology outlined in Section 2.1 can be used to determine order parameters of bistable signaling pathways at saddle-node bifurcation points.

(a)

(b)

Figure 9: Eigenvalues (a) and (b) computed from (3), (4), and (22) for the state (low), (high) (solid lines) and the state (high), (low) (dashed line) shown in Figure 5 as functions of . Magnifying (not shown) the top parts of the two panels reveals that converges to zero when approaching the saddle-node bifurcation points, whereas assumes a finite negative value. Parameters: ; all other parameters are as in Figure 5.

Figure 10: Determinant of the Jacobian matrix as computed from (4) and (22) for the state (low), (high) (solid line) and the state (high), (low) (dashed line) shown in Figure 5 and as function of . Parameters: ; all other parameters are as in Figure 5.

(a)

(b)

Figure 11: Eigenvector components of the first eigenvector and order parameters at the saddle-node bifurcation points. The eigenvector components (a) and (b) as computed from (3) and (4) are shown for the state (low), (high) (solid line) and the state (high), (low) (dashed line) reported in Figure 5 and as functions of . At the saddle-node bifurcation points (circles), these components constitute the respective order parameters.

3. Discussion

Bistable signaling pathways that function as switches and can be described by means of two-variable ODE models were considered. It was demonstrated that in general the switching behavior involves so-called order parameters that identify directions in which the dynamics evolves relatively slowly. As it is known from synergetics [11], a theory of self-organization and pattern formation, these directions play a special role because they dominate in a particular sense the dynamics of the whole system. For a model of a bistable signaling pathway proposed by Markevich et al. (2004) [5], the order parameters at the two switching points described by the pathway model were determined explicitly.

In previous studies, much effort has been devoted to determine the bifurcation diagrams of bistable multispecies signaling pathways either numerically or by means of analytical methods [3, 5–7, 10, 12]. The present study adds a novel perspective to these studies and to similar studies. Accordingly, in the mathematical space spanned by the biochemical species of bistable signaling pathways exhibiting switch-like functioning there exist particular directions that determine (in the sense of synergetics) the emerging order when a switch occurs. For signaling pathways described by means of two-variable models these directions can be determined in terms of order parameter vectors (or unstable modes) as outlined in Section 2.1.

The method presented in Section 2.1 in principle can be generalized and applied to signaling pathways characterized by more than two biochemical species. For such higher dimensional problems, it might be difficult to determine analytical expressions for eigenvalues and eigenvectors. If so, it might be still possible to carry out the procedure sketched in Section 2.1 by means of numerical solution methods. Future work might be devoted to address this issue.

Competing Interests

The author declares that there is no conflict of interests regarding the publication of this paper.

References

J. E. Ferrell Jr. and W. Xiong, “Bistability in cell signaling: how to make continuous processes discontinuous, and reversible processes irreversible,” Chaos, vol. 11, no. 1, pp. 227–236, 2001.

J. E. Ferrell Jr., “Bistability, bifurcations, and Waddington's epigenetic landscape,” Current Biology, vol. 22, no. 11, pp. R458–R466, 2012.

A. Goldbeter and D. E. Koshland Jr., “An amplified sensitivity arising from covalent modification in biological systems,” Proceedings of the National Academy of Sciences of the United States of America, vol. 78, no. 11, pp. 6840–6844, 1981.

J. E. Lisman, “A mechanism for memory storage insensitive to molecular turnover: a bistable autophosphorylating kinase,” Proceedings of the National Academy of Sciences of the United States of America, vol. 82, no. 9, pp. 3055–3057, 1985.

N. I. Markevich, J. B. Hoek, and B. N. Kholodenko, “Signaling switches and bistability arising from multisite phosphorylation in protein kinase cascades,” The Journal of Cell Biology, vol. 164, no. 3, pp. 353–359, 2004.

C. Schröter, P. Rué, J. P. Mackenzie, and A. M. Arias, “FGF/MAPK signaling sets the switching threshold of a bistable circuit controlling cell fate decisions in embryonic stem cells,” Development, vol. 142, no. 24, pp. 4205–4216, 2016.

K. M. Byrne, N. Monsefi, J. C. Dawson et al., “Bistability in the Rac1, PAK, and RhoA signaling network drives actin cytoskeleton dynamics and cell motility switches,” Cell Systems, vol. 2, no. 1, pp. 38–48, 2016.

J.-R. Kim, Y. Yoon, and K.-H. Cho, “Coupled feedback loops form dynamic motifs of cellular networks,” Biophysical Journal, vol. 94, no. 2, pp. 359–365, 2008.

S. Chiangga, W. Pornkaveerat, and T. D. Frank, “On a Fitzhugh-Nagumo type model for the pulse-like jasmonate defense response in plants,” Mathematical Biosciences, vol. 273, pp. 80–90, 2016.

F. Ortega, J. L. Garcés, F. Mas, B. N. Kholodenko, and M. Cascante, “Bistability from double phosphorylation in signal transduction: kinetic and structural requirements,” The FEBS Journal, vol. 273, no. 17, pp. 3915–3926, 2006.

H. Haken, Synergetics: An Introduction, Springer, Berlin, Germany, 1987.

P. Smolen, D. A. Baxter, and J. H. Byrne, “Bistable MAP kinase activity: a plausible mechanism contributing to maintenance of late long-term potentiation,” American Journal of Physiology—Cell Physiology, vol. 294, no. 2, pp. C503–C515, 2008.