Noise-Induced Transitions in a Nonsmooth Producer-Grazer Model with Stoichiometric Constraints.

Stoichiometric producer-grazer models are nonsmooth due to the Liebig's Law of Minimum and can generate new dynamics such as bistability for producer-grazer interactions. Environmental noises can be extremely important and change dynamical behaviors of a stoichiometric producer-grazer model. In this paper, we consider a stochastically forced producer-grazer model and study the phenomena of noise-induced state switching between two stochastic attractors in the bistable zone. Namely, there is a frequent random hopping of phase trajectories between attracting basins of the attractors. In addition, by applying the stochastic sensitivity function technique, we construct the confidence ellipse and confidence band to find the configurational arrangement of equilibria and a limit cycle, respectively.

Introduction

Ecological stoichiometry is the study of the balance of energy (carbon or C) and multipn>le nutrients (such as phosphorus or P, nitrogen or N) in ecological interactions (Sterner and Elser 2002). Organisms are composed of key elements, including C, N, P, whose balance affects organismal reproduction and growth, nutrient cycling, and trophic interactions. Plants can be easily limited by nutrient, and herbivores are more nutrient-rich organisms than plants (Elser et al. 2000). The plant nutrient quality can dramatically affect the growth of herbivorous grazers and may even lead to their extinction (Urabe et al. 2002). Historic predator-prey models only consider energy (carbon) flow in the form of population or density. In reality, producer–grazer interaction models should incorporate both food quantity and quality.

In many ecosystems, trophic interactions can be regulated by excessively enriching or limiting some key resources for reproduction and growth (Grover 2002; Peace 2015). To better understand how nutrient availability affects the population reproduction and growth, a series of newly emerged stoichiometric population models have been proposed in the past two decades for studying the underlying laws of ecological stoichiometry (Andersen 1997; Hessen and Bjerkeng 1997; Kuang et al. 2004; Wang et al. 2007, 2008, 2012; Stecha et al. 2012; Jiang et al. 2019; Song et al. 2019; Zhao et al. 2020), and one of the well received stoichiometric producer–grazer models that track both the quantity and the quality of producer is formulated in (Loladze et al. 2000):where x, y are the densities of producer (phytoplankton/algae) and grazer (daphnia) (mg C/l), respectively; b is the intrinsic growth rate of producer (/day); K is the carrying capacity of producer, which is positively related to light intensity; e is the maximal production efficiency of grazer (no unit); d is the specific loss rate of grazer that includes metabolic losses and death (/day); q is the minimal phosphorus/carbon ratio in producer (mg P/mg C); is the constant phosphorus/carbon ratio in grazer (mg P/mg C); P the total mass of phosphorus in the entire system (mg P/l); f(x) is the consumption rate of grazer, which is usually one of Holling-type functional responses. By applying the Liebig’s Law of Minimum, the producer’s growth rate is limited by both light and nutrient, and the grazer’s growth rate is limited by both food quantity and food quality. Hence, the growth terms have minimum functions and .

Model (1.1) has complex dynamics such as multiple positive equilibria and bistability (Loladze et al. 2000). Later, Li et al. (2011) provided a rigorous mathematical analysis for global stability results of all equilibria and the existence of limit cycles with Holling-type functional responses and fixed parameters except K, and Xie et al. (2018) presented complete global and bifurcation analyses for model (1.1) with Holling type II functional response with all flexible parameters. They found that the model has four types of bistability: between an internal equilibrium and a limit cycle, between an internal equilibrium and a boundary equilibrium, between two internal equilibria, and between a boundary equilibrium and a limit cycle.

In nature, deterministic systems are inevitably affected by various environmental noises which can be important or even dominant in controlling dynamics of trophic interactions. Environmental noises can change the qualitative behavior of a deterministic model (Zhao et al. 2015; Xu and Yuan 2016; Yu et al. 2018, 2019a, b; Zhao and Liu 2019; Wang and Liu 2019; Yu and Yuan 2020). Over the past several decades, many deterministic models with noise-induced transitions have been extensively studied (Xu et al. 2016; Bashkirtseva et al. 2010; Xu et al. 2018; Wu et al. 2019). In Xu et al. (2016), noise perturbations change the coexistence state to extinction for a chemostat model. A new method based on the stochastic sensitivity functions (SSF) technique has been proposed in Bashkirtseva et al. (2010) to construct the analytical description of randomly forced equilibria and cycles of discrete-time models. The presence of coexisting attractors under random perturbations can generate new dynamic regimes, which have no analogues in the deterministic case (Anishchenko et al. 2007; pan class="Chemical">Pisarchik and Feudel 2014). Effects of perturbations in limit cycles were studied in Kurrer and Schulten (1991); Baras (1997). Nonlinear dynamical models show various new phenomena, such as stochastic resonance (Gammaitoni et al. 2009; McDonnell et al. 2010), noise-induced transitions (Horsthemke and Lefever 1984), noise-induced order (Matsumoto and Tsuda 1983; Gassmann 1997), noise-induced chaos (Gao et al. 1999), and noise-induced complexity (Zaks et al. 2005). Analysis of the noise effects on dynamical systems with multipn>le stable states attracts the attention of many researchers (Kim et al. 1998). Multistable systems exhibit complex dynamics with noise-induced hopping between coexisting attractors and their basins of attraction (Kraut and Feudel 2002; de Souza et al. 2007; Dykman et al. 1990, 1994). The sensitivity analysis of randomly forced oscillations is pivotal for investigating these transitions.

The aim of this paper is to study the phenomena of noise-induced transitions for model (1.1) with Holling-type II functional response by using the SSF technique. The rest of this paper is organized as follows. In Sect. 2, we review the deterministic producer–grazer model with stoichiometric constraints and propose its stochastic version. The analysis of noise-induced transitions and the construction of confidence ellipses for this model will be presented in Sect. 3. In Sect. 4, we provide the construction of both confidence band and confidence ellipse and show that noise-induced transitions occur when confidence domains are intersected. Finally, we conclude and discuss the paper in Sect. 5.

Model Formulation and Main Results

In this section, in order to formulate our stochastic model, we first recall the main results of model (1.1) with Holling-type II functional response from Xie et al. (2018). As in Xie et al. (2018), we always assume that in this paperwhere . Then, model (1.1) is simplified as:DenoteThen, is a positively invariant set for model (2.2) [see Li et al. (2011) for the detailed proof]. Obviously, is an open trapezoid due to the limitation of (2.1). We refer the readers to Fig. 2 in Xie et al. (2018), where the grazer nullcline (red curves) in the forward invariant region is the positive x-axis and a polygonal line consisted of two line segments: (denoted by ) and (denoted by ), where ; the producer nullcline (blue curves) in the forward invariant region is the positive y-axis and a parabola arc: , . Though the vector field defined by system (2.2) is not in , it is locally Lipschitz-continuous, which guarantees the existence and uniqueness of solutions of system (2.2).

Fig. 2

Vector field of the deterministic model and the equation of the separatrix of two subsystems: (Color figure online)

System (2.2) always has two boundary equilibria: and . It may also have none or one to three coexistence/interior equilibria. The dynamics of system (2.2) is completely determined by some critical values about K listed below [see Xie et al. (2018) for more details]: The authors in Xie et al. (2018) have provided a complete global analysis for system (2.2) without fixing any parameter. Their analysis shows that the model has far richer dynamics than those found in the previous paper (Li et al. 2011). For example, system (2.2) may have three interior equilibria , , in the forward invariant region and four types of bistability may appear: (i) between and , (ii) between the limit cycle and , (iii) between and , (iv) between the limit cycle and . The readers are referred to Xie et al. (2018) for more details. In this paper, we are only concerned with the first two types and explore the impact of noises existed in the environment on bistability between two interior attractors: one is for the stable equilibrium or the unique stable limit cycle surrounding the unstable equilibrium , the other is for the stable internal equilibrium .

;

;

;

is the value of K at which the line segment is tangent to the parabola (not always exists);

.

We remark that a necessary condition for system (2.2) to have three internal equilibria is that exists. Here, we just mention the following two results from Xie et al. (2018) with small modifications, which will be used in the sequel analysis.

Lemma 2.1

(Theorem 3.2 (3) of Xie et al. (2018)) Assume that . If , there exist three internal equilibria: and are two stable equilibria, is a saddle, and the model has no limit cycle, see Fig. 1a. Hence, bistability occurs: orbits on the left of the stable manifold of will eventually tend to , orbits on the right of the stable manifold of will eventually tend to , and the separatrix is the stable manifold of . Moreover, the model has three heteroclinic orbits from to , from to , and from to , respectively.

Fig. 1

Bistability of model (2.2) occur: a two internal equilibria; b an equilibrium and a limit cycle. Here, hollow dots denote unstable equilibria, solid dots are stable equilibria; blue closed circle denotes stable limit cycle (Color figure online)

Lemma 2.2

(Theorem 3.2 (4) of Xie et al. (2018)) Assume that . If , there also exist three internal equilibria: is a saddle, is stable, but becomes unstable. Meanwhile, the model has a unique stable limit cycle surrounding the equilibrium (cyan dashed-dotted line) for , where is the broken value of K for the limit cycle, see Fig. 1b. In this case bistability occurs: solutions of the model on the right of the stable manifold of tend to the equilibrium , solutions of the model on the left of the stable manifold of tend to the limit cycle, and the separatrix is the stable manifold of . Moreover, the model has two heteroclinic orbits from to and from to , respectively.

We are now in a position to propose our stochastic model. From the biological point of view, the sources of environmental noise affecting the producer and the grazer are different. Following the approach used in Refs. Imhof and Walcher (2005) we consider randomness into deterministic system (2.2) and then obtain the following stochastic differential equation:where and are two standard one-dimensional independent Brownian motions, are the noise intensities.

For the simplicity of discussion, we assume . Notice that system (2.3) consists of two subsystems. Therefore, we assume that when and when . The following theorem is about the existence and uniqueness of the global positive solution of model (2.3), whose proof is provided in “Appendix 1”.

Theorem 2.1

Assume that . Then, for any given positive initial value , stochastic system (2.3) admits a unique positive solution (x(t), y(t)) for and the solution will remain in with probability one, in other words, for all almost surely (a.s.).

In the following, we study the phenomenon of noise-induced transitions between stochastic attractors for stochastic model (2.3).

Analysis of Noise-Induced Transitions Between Two Internal Equilibria

In order to analyze the influence of noise on model (2.2), we take the following realistic parameter values from (Li et al. 2011):We further take , it is easy to compute from Example 2 in Xie et al. (2018) thatObviously, . Now taking the above parameter values and in model (2.2) (i.e., stochastic model (2.3) with ), there exist a washout equilibrium , a mono-culture equilibrium , and three coexistence/interior equilibria , , and . It is easy to see that . From Example 2 in Xie et al. (2018) and Lemma 2.1, we know that equilibria and are unstable, is a saddle point, coexistence equilibria and are locally asymptotically stable. Using the command “DEtools [phaseportrait]” of Maple, the vector field of the deterministic model with given parameters is drawn as in Fig. 2, in which the red dash-dotted line is the separatrix of two attraction domains and the blue dash-dotted line is the separatrix of two subsystems.

For the deterministic model, the trajectory with its initial point inside the separatrix converges to the coexistence equilibrium , and the trajectory with its initial value outside the separatrix converges to the coexistence equilibrium . However, the dynamics of the stochastic model sometimes can be difficult to predict accurately. In order to illustrate the impact of environmental noise on dynamical behaviors of model (2.2), we take equilibrium as an example. For equilibrium , we can discuss similarly. For a weak noise, the stochastic trajectory with its initial value near the deterministic coexistence equilibrium will fluctuate around this equilibrium, the densities of grazer and producer stay close to their equilibrium values (see Fig. 3). However, as the noise intensity increases and becomes sufficiently large, both grazer and producer can go extinct (see Fig. 4). Hence, there is a critical noise intensity , when , both grazer and producer persist. Furthermore, in this interval we find interesting phenomena of noise-induced transitions between two stochastic attractors. Here, Figs. 3 and 4 are drawn with Matlab software by using the Euler-Maruyama method developed in Ref. Higham (2001). The same method is also used in the numerical simulations of time series and random trajectories in the later figures of the paper.

Fig. 3

a Time series of x(t) and y(t); b Phase trajectory for stochastic model (2.3) with the initial value (0.1603, 0.4415) and the noise intensity (Color figure online)

Fig. 4

a Time series of x(t) and y(t); b Phase trajectory for stochastic model (2.3) with the initial value (0.1603, 0.4415) and the noise intensity (Color figure online)

In the following, by using the SSF method (see “Appendix 2”), we construct confidence pan class="Species">ellipses for stochastic model (2.3) to characterize the configurational arrangement of stochastic coexistence equilibria and then further estimate the threshold value of the noise intensity of state switching.

DefinewhereandThe stochastic sensitivity matrixsatisfies the following equations:It then follows from (B.3) that the confidence pan class="Species">ellipse equation iswhere is a interior equilibrium of deterministic model (2.2), and pan class="Chemical">P are respectively the noise intensity and a fiducial probability.

Taking the same parameters as in Fig. 2, for the interior equilibrium we haverespectively. Then, from Eq. (3.2) the confidence pan class="Species">ellipse equation of isFor fixing fiducial probability , we take the noise intensities , and 0.003, respectively, resulting in the corresponding confidence pan class="Species">ellipses shown in Fig. 5a. Obviously, as the noise intensity increases, the confidence ellipse starts to expand and after crossing the separatrix, it enters the attraction domain of the coexistence equilibrium . The value corresponding to the tangency of the confidence ellipse can be used as an estimation for the threshold noise intensity of the onset of noise-induced transitions. Here, . Figure 5b illustrates the confidence ellipse with , one can see that the random states of the stochastic model are distributed around the corresponding deterministic coexistence equilibrium, and they belong to the interior of the confidence ellipse with probability 0.95.

Fig. 5

a Separatrix (dashed-dotted) and confidence ellipses (solid) for (small), (middle), (large). b Random state (blue) of stochastic model (2.3) and confidence ellipse (green) for

a Separatrix (dashed-dotted) and confidence pan class="Species">ellipses (solid) for (small), (middle), (large). b Random state (blue) of stochastic model (2.3) and confidence pan class="Species">ellipse (green) for

Similarly, the stochastic sensitivity matrix corresponding to equilibrium isthen the confidence pan class="Species">ellipse equation of isBy taking , we consider two confidence pan class="Species">ellipses together. For a weak noise (), two confidence ellipses of the coexistence (green) and (blue) are distinctly separated by the separatrix of two attraction basins (see the left panel of Fig. 6a), and solutions starting from one side will eventually approach the coexistence equilibrium (red) or (green) on that side (see the middle panel of Fig. 6a). Stochastic trajectories leaving the unforced deterministic attractors concentrate in their small neighborhoods (see the right panel of Fig. 6a). Here, the dynamics of the stochastic model is almost regular and small noises have little impact on the densities of grazer and producer.

Fig. 6

Confidence ellipses (left), time series (middle) and random trajectories (right) of stochastic model (2.3) for: a, ; b, ; c, ; d, (Color figure online)

As the noise intensity increases, the confidence pan class="Species">ellipse expn>ands. We keep and increase to 0.008, then the confidence pan class="Species">ellipse (blue) of the coexistence equilibrium crosses the separatrix (see the left panel of Fig. 6b). Solutions starting from the attraction basin of will eventually approach the coexistence equilibrium on the other side with high probability (see the middle panel of Fig. 6b). The phenomenon of noise-induced transition occurs (see the right panel of Fig. 6b). We also have a symmetric result in Fig. 6c (for ).

The above three situations are all related to the initial point. However, when both confidence pan class="Species">ellipses expn>and and cross the separatrix and intersect each other (see the left panel of Fig. 6d for ), we obtain frequent random hopn>ping of phase trajectories between attraction basins of the equilibria and . On the phase plane, a place of this intersection marks a location of the “transition bridge” between basins of attraction where noise-induced transitions are most likely to occur. These results obtained by the confidence domain method are in agreement with the direct numerical simulation of time series (middle panel) and random trajectories (right panel) in Fig. 6d. The difference is that the initial value can be evaluated anywhere in the invariant set.

Confidence pan class="Species">ellipses (left), time series (middle) and random trajectories (right) of stochastic model (2.3) for: a, ; b, ; c, ; d, (Color figure online)

Analysis of Noise-Induced State Switching Between Confidence Ellipse and Confidence Band

When , and other parameter values are the same as in (3.1), Lemma 2.2 implies that the coexistence equilibrium of model (2.2) is stable, but becomes unstable. Meanwhile, the model has a unique stable limit cycle surrounding the equilibrium , where both population levels fluctuate around a coexistence equilibrium. We next analyze noise-induced state switching by constructing confidence band and confidence pan class="Species">ellipse.

In the stochastic case, the limit cycle generally disappears, but the trajectories will remain in a small neighborhood of the deterministic limit cycle for a small noise. For characterizing the configurational arrangement of this neighborhood, in what follows, we will construct the confidence band for stochastic model (2.3).

For convenience, letand denote the deterministic limit cycle by , where T is the period. Then, we can write matrices F(t), G(t) and S(t) as follows:whereandFrom (B.4), we know that the stochastic sensitivity function satisfies the following boundary problem:whereHere,are elements of a vector function orthogonal to vector . It follows from (B.5) that the boundaries of the confidence band have the following explicit parametrical form:Here, the parameter k is connected with the fiducial probability by the formula , where is the error function.

In Fig. 7a, the red line is the deterministic limit cycle, the blue points are the random states on different time and the two green lines are the boundaries of the confidence band. Obviously, the random states are distributed around the deterministic limit cycle, and they belong to the interior of the confidence band with probability 0.95.

Fig. 7

Confidence bands for stochastic model (2.3). a Random states (blue) around (red) and band (green) for , ; b and (green), (blue); c and (green), (blue)

Figure 7b, c illustrate the effects of the noise intensity and fiducial probability on the size of confidence band. It is easy to see from Fig. 7b, c that the configurational arrangement of confidence band begins to expand as the noise intensity or fiducial probability increases. This result can be deduced from the expressions of and .

Following the same logic as in Sect. 3, we can obtain the stochastic sensitivity matrix of as follows:Then, from (B.3) the confidence ellipse equation of isBy taking , we investigate the evolution process between confidence band and confidence ellipse. For a weak noise, confidence band and confidence ellipse are distinctly separated, the trajectory with its initial point in either of two attraction basins will eventually approach confidence band or confidence ellipse in that attraction basin. Stochastic trajectories leaving the unforced deterministic attractors concentrate in their small neighborhoods (see Fig. 8a for ). With , the left panel of Fig. 8b illustrates that confidence band and confidence ellipse expand and intersect each other. This leads to frequent random hopping of phase trajectories between attraction basins of the limit cycle and a coexistence equilibrium (see the middle and right panels of Fig. 8b).

Fig. 8

Confidence domains (left), time series (middle) and random trajectories (right) of stochastic model (2.3) for: a, ; b, (Color figure online)

Discussion

We study noise-induced transitions from one coexistence to another coexistence in a constant interval for a stochastically forced producer–grazer model with stoichiometric constraints. The corresponding deterministic model is nonsmooth and its global stability and bifurcation are fully analyzed in Xie et al. (2018). When appropriate parameters are chosen, bistability can occur. The attractors (stable equilibria , , or the limit cycle ) correspond to the stable coexistence of grazer and producer. Two attraction basins are distinctly separated by the separatrix, which is the stable manifold of saddle point , and solutions from either side of separatrix will eventually approach the equilibrium or the limit cycle on that side. However, when there is noise disturbance, the boundary may be damaged by noise, and the solutions in the attractive basin of one stable coexistence equilibrium can eventually approach the other stable coexistence equilibrium with a high probability.

Based on the technique of SSF, we construct confidence ellipses to characterize the phenomenon of noise-induced state switching between two stochastic coexistence equilibria. It provides us the general location of the equilibria and the distribution of random states in the stochastic model. For a weak noise, the confidence ellipses are completely contained in the attraction basins of the corresponding coexistence equilibria, and the random trajectories will not leave the confidence ellipses with a high probability. As the noise intensities increase, the confidence ellipses begin to expand and when only one ellipse crosses the separatrix, we study the phenomenon of noise-induced transitions from a coexistence to another coexistence. Furthermore, when both of confidence ellipses cross the separatrix of two attraction basins and intersect each other, we obtain that frequent random hopping of phase trajectories between the attraction basins of the equilibria.

In addition, we have established the confidence band for the limit cycle of model (2.3), which provides us the general location of the stochastic cycle and the distribution of random states around the deterministic limit cycle. Following the same logic, we analyze the phenomenon of noise-induced transitions via confidence band and confidence pan class="Species">ellipse.

Theoretically, the solution of the stochastic model with multi-stability runs long enough to cause state switching. When one confidence domain intersects another, the probability of switching is very high (it can be done in a very short time).

In fact, the phenomenon of noise-induced transitions we get is independent of the way the noise is introduced, we’re just taking one of these cases. For instance, we can introduce randomness into deterministic model (2.2) by perturbing the parameter and obtain the following stochastic differential equation:The results of numerical simulation of time series and random trajectories are similar to those presented in this paper.

We discuss the generality of model (2.2), where noised-induced state switching might not be as common in a real-world situation. There are different noise intensities and ranges for different species. Our work expounded how the confidence domain method is used to understand the qualitative changes in stochastic dynamics at which noise-induced state switching occurs with high probabilities. This method is applicable to nonsmooth competition models and more complex higher dimensional models in aquatic or terrestrial ecosystems. Our results enrich the study of asymptotic behaviors in stoichiometric producer–grazer models and help better understand the stoichiometric producer–grazer dynamics in the stochastic perspective.