In cboettig/bs-tipping: Bifurcation or State Tipping in a Trophic Cascade Experiment
Introduction
Ecological systems prone to abrupt changes, such as a fishery collapsing, can be challenging to manage. These abrupt dynamics are often referred to as regime shifts, and result from nonlinearities in the underlying system dynamics. Therefore, ecologists have often turned to studying the dynamics of nonlinear but simple mathematical models to better understand and manage ecosystems. There are essentially two different mechanisms by which a system exhibiting alternative stable states can experience a sudden transition from one state to another [@Scheffer2001; @Ashwin2012]: state tipping (S-tipping) or a bifurcation (B-tipping). Figure \ref{fig1} illustrates these two mechanisms. B-tipping is a sudden transition that is induced by slow changes to a parameter of the underlying dynamical system that reduce the stability and ultimately destroy the stable node, as it collides with an unstable saddle-point in a saddle-node bifurcation (aka a fold bifurcation or "blue sky catastrophe," [@Strogatz2001]). In contrast, S-tipping involves a perturbation directly to a state variable that carries the system across the separatrix into an alternative stable state, without any such change to the underlying shape or existence of either attractor. S-tipping can be driven by direct exogenous forcing or purely stochastic events, such as extreme events that are common in the biological time series of ecosystems [@Batt2017]. While real-world transitions will often be a combination of both processes, one mechanism will typically have a dominant role. For instance, near enough to the bifurcation, intrinsic noise will frequently drive the system across the tipping point before the actual bifurcation event is reached. Crucially, only one of these mechanisms, B-tipping, will exhibit the phenomenon of Critical Slowing Down (CSD), a pattern that has received great attention in both theoretical and empirical studies for its potential as an early warning signal of such a regime shift [e.g. @Wissel1984; @Scheffer2009; @Scheffer2015].
Despite this attention, convincing empirical demonstrations of critical slowing down are still rare, particularly outside carefully controlled laboratory settings [e.g. @Dai2012]. A foundational example can be found in studies of lake ecosystems [e.g. @Carpenter2008], which can exist in either a "piscivore dominated" or "planktivore dominated" state (Figure \ref{fig2}). Significant previous work has established evidence for these alternative states and the potential for sudden transitions between the states that cannot be easily reversed. Then in 2011, a landmark study of @Carpenter2011 analyzed the experimental manipulation of an initially planktivore-dominated state as a definitive example of CSD phenomena providing an early warning signal prior to a sudden transition to the piscivore dominated state. Yet was this manipulation truly an example of B-tipping, where CSD is theoretically predicted, or S-tipping, which does not create the CSD phenomenon?
\begin{figure}
\centering
\includegraphics[width=.8\linewidth]{figures/scheffer-fig.png}
\caption{S-tipping vs B-tipping. Both panels show the bifurcation diagram for a saddle-node or "fold" bifurcation. These plot the equilibrium points of a system on the vertical axis (e.g. adult bass density in our example) as a function of environmental conditions (parameters of the system). Solid lines indicate stable points, dotted lines indicate unstable points. The top panel illustrates Bifurcation Tipping: a system which begins in the desirable, "high" ecosystem state (top branch) experiences gradual degradation as the underlying environmental conditions change (increasing values along the $x$). At the F2 position, this stable node collides with the unstable saddle point (dotted line) and is annihilated, leaving the system to rapidly transition to the remaining undesirable stable state (Forward shift). To restore the system to the desirable state, it is not sufficient to merely restore the environmental conditions observed at F2. Instead, environmental conditions must improved all the way to F1, where another bifurcation occurs that creates a rapid backwards shift to the desirable state. This effect is known as hysteresis. The bottom panel introduces state tipping: the transition from the desirable state to the undesirable state is triggered without changing the environmental conditions. Instead, a large perturbation is introduced directly to the ecosystem state, sufficient to push the system from one basin-of-attraction to another. Such transitions are characterized by large, sudden shocks to the state variable, rather than slow, gradual change to the environmental conditions. These transitions do not exhibit hysteresis -- a sufficiently large perturbation could restore the system. Figure reproduced from Scheffer et al (2001) with permission.}\label{fig1}
\end{figure}
B-tipping is typically considered to be the result of slow change in a parameter of [e.g. @Scheffer2009] the system which ultimately leads to the fold bifurcation (or close enough to a fold bifurcation for stochasticity to finish the job). S-tipping involves manipulation of the state variable. @Carpenter2011 describe driving the critical transition through the gradual addition of piscivores (largemouth Bass), which drives a transition characterized by low-piscivore density into one of high piscivore density. If the state (piscivore density) is being manipulated directly, then is this not an example of S-tipping, and thus unsuitable to test for the presence of CSD? Although the simplicity of ball-in-cup (potential well) based descriptions of regime shifts makes them alluring, this question is impossible to address properly without a more mathematically precise description of the system and the manipulation involved. Surprisingly, despite the central importance the piscivore transitions have played in our understanding of ecosystem regime shifts, these details have never been clearly presented in prior publications of this classic lake ecosystem experiment [@Carpenter2008; @Carpenter2011]. Here, we provide a direct analysis of the underlying model in the piscivory alternative state system to demonstrate why these studies are indeed an example of B-tipping.
Central to resolving this question is untangling a second issue that frequently arises in the discussion of alternative stable state dynamics in ecological systems: the issue of dimension. Ecosystems are notoriously complex, and this lake ecosystem is no exception. Here, we analyze a previously published model of a trophic cascade to determine whether manipulating the top predator constitutes B-tipping or S-tipping. Common versions of this model include 5 state variables, ranging from top predators to primary producers (Figure \ref{fig2}). In sharp contrast to this complexity, both intuition and theory of critical transitions have been centered on one-dimensional interpretations. The S-tipping and B-tipping transitions illustrated in Figure \ref{fig1} are fundamentally rooted in a one-dimensional analysis: the state consists of a single continuously valued dynamical variable (i.e. population density of adult bass). Similarly, Ball-in-cup style depictions of alternative stable states and sudden transitions present only a one-dimensional depiction of a system state. While it may be tempting to dismiss such simple depictions as out-of-hand as mere cartoons, it is in fact possible to map the full five-dimensional model explicitly into a one-dimensional approximation under a few simple assumptions. Here, we present the full dimension reduction on the full 5D model analytically, and compare this approximation to numerical simulations of the original system. This exercise facilitates the interpretation of the lake experiment in the context of predicting regime shifts by using a 1D model to understand a high dimensional system.
\begin{figure}
\centering
\includegraphics[width=.8\linewidth]{figures/Slide1.png}
\caption{Food Web Heuristic. Adult bass give rise to juveniles, while
eating their own young and pumpkinseed sunfish. Sunfish compete with and
eat juvenile bass. Sunfish eat herbivorous zooplankton like
\emph{Daphnia}, which in turn eat phytoplankton.}\label{fig2}
\end{figure}
Methods
Model Description
In the model of the lake food web, largemouth bass (\emph{Micropterus salmoides}, A) are the top predator. Intuitively, juvenile bass (J) mature into A, and J are produced by A. Bass feed on other fish, including their own young (J) and other planktivorous fish (F), (i.e. pumpkinseed sunfish, \emph{Lepomis gibbosus} and similar forage fish species). These planktivores and juvenile bass (which are also planktivores) eat herbivorous zooplankton (H), which in turn eat phytoplankton (P). Because bass fry are small enough to be eaten by "planktivorous" fish, J can be consumed by F. Furthermore, F and J compete with each other for food (H). Thus, there are strong pairwise interactions between the three fish groups, and the dynamics of these groups extend into the planktonic portion of the food web. The five-dimensional model of the food web includes the aforementioned state variables consisting of three fish groups (A, F, J) and two planktonic groups (H, P):
$$\begin{array}{llr} \dot{A} &= sJ - qA - m_A A &(1) \ \dot{F} &= D_F(F_o - F) - c_{FA}FA &(2)\ \dot{J} &= fA - c_{JA}JA - \frac{c_{JF}vJF}{h+v+c_{JF}F} - sJ &(3)\ \dot{H} &= D_H(H_R-H) + \alpha c_{HP}HP - c_{HF}HF &(4)\ \dot{P} &= r_P L \gamma (I_0,P)P - mP - c_{PH}HP &(5)\ \end{array}$$
Decreases in adult bass abundance result from mortality ($m_A A$) and harvest effort ($q$), and increases result from the maturation of the fraction of $J$ that survive ($s$) mortality (Eq. 1). Stocking the lake corresponds to $q < 0$. This continuous-time formulation for harvest or stocking effort can be justified when carried out throughout the year at a fixed rate, as in the experiments considered here, rather than in large shocks. Planktivorous fish diffuse ($D_F$) between a refuge of size $F_o$ and a foraging arena, and the predation rate of $A$ on $F$ is controlled by $c_{FA}$ (Eq. 2). Growth in $J$ is controlled by the fecundity ($f$) of $A$, but are subject to losses via cannibalization at a rate $c_{JA}$, and interactions with $F$ which are controlled by the rate at which $J$ hide ($h$) in a refuge or become vulnerable ($v$) to predation at a rate $c_{JF}$ as they enter an implicit foraging arena (Eq. 3). Herbivorous zooplankton diffuse between surface waters and a deep water refuge (of size $H_R$) at a rate $D_H$, grow by preying on phytoplankton at rate $c_{HP}$ and converting phytoplankton biomass (phosphorus) by a factor $\alpha$, and decline as they are preyed upon by planktivores at a rate $c_{HF}$ (Eq. 4). Phytoplankton death is due to mortality and sinking at rate $m$ or predation by herbivores at rate $c_{HF}$, while their growth is controlled by the per capita growth rate ($r_P$), the availability of phosphorus (introduced at rate $L$), and the availability of light at the surface ($I_0$) and the density dependent growth response to phytoplankton in response to light ($\gamma$, Eq. 5).
Table 1. Model parameters, description, values, and units. Parameters are taken from [@Carpenter2008] and [@Carpenter2013]. The function $\gamma$ describes phytoplankton growth, and can be found in [@Carpenter2008].
|Symbol |Description |Quantity |Units |
|:------|:-----------------------------------------------------------|:--------|:-------------------|
|$q$ |Harvest or stocking |variable |$t^{-1}$ |
|$s$ |Survival rate of juveniles to maturation |0.6 |$t^{-1}$ |
|$m_A$ |Mortality rate of adults |0.4 |$t^{-1}$ |
|$m$ |Phytoplankton mortality |0.1 |$t^{-1}$ |
|$f$ |Fecundity of adult bass |2 |$J/A$ |
|$c_{FA}$ |Predation of planktivores by adult bass |0.3 |$t^{-1} A^{-1}$ |
|$c_{JA}$ |Predation of juvenile bass on by adult bass |0.1 |$t^{-1} A^{-1}$ |
|$c_{JF}$ |Predation of juvenile bass by planktivorous fish |0.5 |$t^{-1} F^{-1}$ |
|$c_{PH}$ |Predation of phytoplankton by herbivores |0.25 |$t^{-1} H^{-1}$ |
|$c_{HF}$ |Predation of herbivores by planktivores |0.1 |$t^{-1} F^{-1}$ |
|$F_o$ |Abundance of planktivorous fish in non-foraging arena |200 |$F$ |
|$H_R$ |Reservoir of herbivores in a deep-water refuge |4 |$H$ |
|$D_F$ |Diffusion of planktivores between refuge and foraging arena |0.09 |$t^{-1}$ |
|$D_H$ |Diffusion of between shallow and deep water |0.5 |$t^{-1}$ |
|$v$ |Rate at which J enter a foraging arena, becoming vulnerable |80 |$t^{-1}$ |
|$h$ |Rate at which J hide in a refuge |80 |$t^{-1}$ |
|$r_P$ |Growth rate of phytoplankton |3 |$m^2 mg^{-1}$ |
|$L$ |Phosphorus loading rate |0.6 |$mg m^{-2} t^{-1}$ |
|$\alpha$ |Assimilation of phytoplankton phosphorus by zooplankton |0.3 |$g g^{-1}$ |
|$\gamma$ |Density-dependent phytoplankton growth response to light |n/a |n/a |
|$I_0$ |Solar irradiance incident on the surface of the water |300 |$\mu E m^{-2} s^{-1}$ |
Dimension Reduction
Our approach to determining whether the lake experiment constituted B-tipping or S-tipping is to reduce the dimensionality of the model described in Eqs. 1-5, and see if a 1-D expression of this model can produce alternate states and the corresponding fold bifurcation. Although transitions occur throughout the food web, we focus here on the fish dynamics (Eqs. 1-3) as they are most directly related to either harvest ($q$) or direct manipulation of adult bass (A). In the model, the dynamics of herbivorous zooplankton (H) and phytoplankton (P) do not influence A, F, or J. Therefore, the dimensionality of the model can easily be reduced from 5 to 3 by ignoring Eqs. 4-5 (the plankton) without loss of generality for the dynamics of Eqs. 1-3 (the fish). In addition to parameters, the dynamics of each fish state depend on itself and at least one other fish state. We observe that the dynamics of both the Juvenile bass (J) and forage fish (F) occur on fast timescales relative to the dynamics of long-lived adult bass. Therefore, we set Eqs. 2-3 to zero and solved for the estimate of the respective state variables given equilibrium:
$$\begin{array}{llr}
\hat A&= \frac{sJ}{q+m_A} &(6)\
\hat F&=\frac{D_FF_0}{D_F + c_{FA}A} &(7)\
\hat J&=\frac{fA}{c_{JA}A+s+\frac{c_{JF}vF}{h+v+c_{JF}F}} &(8)\
\end{array}$$
The equation of motion for A (an expression of A's dynamics completely in terms of A and fixed parameters), is acquired by substituting the RHS of Eq. 7 for F in Eq. 8, then substituting the resulting expression for J in Eq. 1. After some substitution, we can write the equation of motion for the adult bass, A, completely in terms of the adult bass population and fixed parameters of the system.
$$\begin{array}{llr}
\dot{A} &= sfA\Bigg(c_{JA}A+s+c_{JF}v\bigg(\frac{(h+v)(D_F+c_{FA}A)}{D_F F_o}+c_{JF}\bigg)^{-1}\Bigg)^{-1} - qA - m_A A &(9)
\end{array}$$
Solving for the roots of this equation demonstrates alternative stable states and the potential for B-tipping with changing $q$, as illustrated graphically in Figure \ref{fig4}.
Simulations
We performed simulations in which we initialized the 5D model at equilibrium with a low but negative harvest rate ($q=-0.001$), which can be interpreted as per capita stocking of $A$. We proceeded to linearly increase $q$ at each time step over 1020 years with a time step size of 0.01 years (102000 time steps total) until $q$ reached a value of 0.001. At each time step a small amount of additive, normally distributed noise was added to Eqs. 1-3, truncated to non-zero values. R code for these simulations is included in supplementary materials (https://github.com/cboettig/bs-tipping).
Results
Plotting this equation using established ranges for ecological parameters (Table 1; @Carpenter2008; @Carpenter2013), Figure \ref{fig4} shows that the system should support two alternative stable states, and can undergo a fold bifurcation (in reverse, creating a new stable point) when "catch" $q$ is decreased sufficiently. This sudden transition is felt in other dimensions as its impact cascades down to the lower trophic levels.
\begin{figure}
\centering
\includegraphics[width=.8\linewidth]{../analysis/manuscript_report/figure-5dsimulation-1.png}
\caption{Time series of the 5 state variables, with (q) (harvest rate
of adult bass) being increased linearly over time.}\label{fig3}
\end{figure}
These same dynamics can be seen in the direct numerical simulations of the full 5-dimensional model. At a $q$ of 0, the fish in the simulation underwent a critical transition (Figure \ref{fig3}). This transition occurred during time step 510, although an abrupt change in state was not observed until step 850. The changes in A and J were similar, but were inversely related to those in F. The predation effect of A on J was small compared to the role of A's fecundity in producing J and A's role as a predator of F, which negatively impacts J as well. Overall, both bass time series decreased with increasing $q$ (abruptly when stocking converted to harvest at $q = 0$), whereas the abundance of planktivorous fish was positively correlated with $q$.
\begin{figure}
\includegraphics[width=\linewidth]{../analysis/manuscript_report/figure-Amotion-ballincup-1.png}
\caption{Top row: motion of A for different values of $q$. Bottom row: the negative integral of the top row, depicting in the so-called ball-in-cup (or potential well) diagrams. Note that the peaks and troughs (local maxima and minima) of the potential wells, bottom row, correspond to the zeros of the rate equations (top row). Troughs indicate stable nodes, peaks indicate unstable saddle points. Note that under sufficient stocking of bass, $q=-0.5$, the zero point is an unstable equilibrium, indicated by the peak in the potential well (bottom right).} \label{fig4}
\end{figure}
The fish dynamics cascaded down to lower trophic levels, although H and P did not undergo abrupt transitions at the same time as the fishes. Instead, $H$ and $P$ exhibited strong nonlinear changes much earlier in the time series after $q$ had increased only slightly. The small cumulative increase in $q$ leading up to step 20 caused a slight increase in $F$ and a dramatic decline in $H$ ($F$'s prey); after its initial decline, the abundance of $H$ was 0.001 for the remainder of the time series. When $H$ abruptly declined, $P$ abruptly increased, then continued to increased approximately linearly for the remainder of the time series. Therefore, increasing $q$ caused declines in both bass states ($A$ and $J$), an increase in planktivorous fish ($F$), a decline in herbivorous zooplankton ($H$), and an increase in phytoplankton ($P$), a pattern consistent with a trophic cascade. All state variables exhibited abrupt state changes at some point in the time series.
Discussion
These results demonstrate that the system considered can indeed be driven through fold bifurcations as the parameter $q$ is slowly tuned up or down (Figure \ref{fig5}), forcing the creation or annihilation of an alternative stable state -- confirming the presence of B-tipping in this system. Simulations of the same system using previous estimates for ecological parameters confirm that the predictions of the 1-dimensional approximation (Eq 9) correspond with the dynamics observed in the full model.
\begin{figure}
\centering
\includegraphics[width=.8\linewidth]{../analysis/manuscript_report/figure-bifurcationDiagram-1.png}
\caption{Bifurcation diagram derived from Eq 9. Solid lines are stable
equilibria, the dashed line is the unstable equilibrium or separatrix.}\label{fig5}
\end{figure}
This analysis also highlights the importance of what might appear to be a merely semantic issue that in fact corresponds to fundamentally different dynamics. In this formulation, the control parameter $q$ corresponds to harvest (positive $q$) or stocking (negative $q$) of adult bass in the lake. It is this aspect that first gives the appearance of state tipping dynamics, since it suggests the state is being somehow "directly manipulated". Yet from a mathematical standpoint, it is irrelevant if the control parameter corresponds to some environmental or policy variable (fishing regulations, fish migration into the lake), which in turn increases or decreases net mortality proportionally, or if the manipulation directly adds or removes adult bass -- the feature that distinguishes this manipulation from a S-tipping transition comes down to the scales involved. This can be seen most clearly though an analogy of considering the same model applied to a larger scale, such as marine fisheries. In fisheries models, one would routinely examine the dynamics of the fish population biomass $B$ with respect to different possibly parameters for harvests or fishing effort. Even though harvest may be measured in the same units as the total biomass, it is natural to think of the two quantities very differently. A manager can adjust the harvest effort by changing quotas, limiting permits for fishing boats etc. Slowly increasing in the amount of bass stocked in the lake (i.e. decreasing $q$, effectively running the scenario of increasing fishing harvests in reverse) gives the rest of the food web time to adjust. In contrast, a single large introduction of many bass would create an S-tipping transition, driving the system immediately into an alternate state without exhibiting CSD phenomena needed to evaluate early warning signals of a tipping point.
All though we have focused on the analysis of deterministic models for simplicity, these results can easily be interpreted in the context of stochastic models as well. (Recall that in practice, the detection of early warning signals associated with CSD requires some level of stochasticity to generate the variance needed to detect slowing down). Under a stochastic model, the distinction between S-tipping and B-tipping is slightly blurred, since a system that gets close enough to the tipping point due to changing parameters (B-tipping) will ultimately be able to jump across the separatrix (S-tipping) before the actual bifurcation event. In a stochastic system, fluctuations in the state variable act as intrinsic sources of perturbations, and though part of the system, large chance deviations are no more predictable than sudden exogenous forcing: there is no CSD prior to purely stochastic transitions [@PRSB2013], even though retrospectively a large deviation will share similar statistical signatures such as rapidly increasing variance and autocorrelation [@PRSB2012].
A second issue highlighted by this analysis is the contrast between the interpretation of the five-dimensional food web model and the classic one-dimensional description of tipping points shown in Fig \ref{fig1}. Such one-dimensional descriptions of a state are frequently presented in reviews or abstract discussion of regime shifts and early warnings with little or no discussion of how a complex system comes to be represented by a one-dimensional continuous variable [e.g. @Scheffer2001; @Scheffer2009; @Barnosky2012; @Scheffer2012; @Clements2018], and even literature focusing on specific systems such as the lake-ecosystem considered here can be vague about the dimensional reduction [@Carpenter2008; @Carpenter2011]. This can in turn create confusion as to whether such one-dimensional diagrams are meant to be interpreted as mathematically explicit approximations of higher-order dynamics, or merely as cartoons meant only to mimic the qualitative pattern but bearing no quantitative relationship with the actual process. The analysis presented here for the lake-ecosystem model demonstrates that the former does indeed hold, though only under the appropriate time-scale limits -- not only that of the "slow parameter" but also that of the other dynamic variables in the system: plankton, Juvenile Bass, Forage Fish, Fig \ref{fig2}. Our analysis also highlights that how the nonlinearity required to create a fold bifurcation arises directly from this dimension-reducing approximation. In the original five-dimensional model, Eq (1), the dynamics of the adult bass $A$ appear to be completely linear. Only after accounting for the feedback in the other dimensions do we see the nonlinearity required to create a fold bifurcation. Our numerical simulations provide an illustration of how closely the dynamics of the full five-dimensional system are captured by the one-dimensional approximation using a biologically plausible parameterization. Thus the reference to the dynamics of a one-dimensional system state, as well as the interpretation of the transition as B-tipping rather than S-tipping, appears to be well-justified in the lake-ecosystem case.
This analysis does not vindicate the use of one-dimensional models more generally. The popularity and simplicity of the fold-bifurcation model has sometimes led to more casually equating any sudden shift with fold bifurcations, overlooking any assumptions necessary to approximate a complex dynamical system with such a one-dimensional model. Influential reviews on abrupt change in ecological systems [e.g. @Scheffer2015; @Ratajczak2018] focus almost exclusively on descriptions of one-dimensional state-spaces without discussion of if, when or how complex systems such as the ecosystem food-web considered here can be appropriately described by a nonlinear equation in a single dimension. More careful modeling informed by system-specific knowledge is required before any abrupt transition can be characterized as a fold bifurcation [@Nature2013]. The fact that CSD phenomena can occur in many different models has been suggested to mean that CSD is a generic or universal early warning indicator [@Scheffer2009] of a critical transition. In practice, many higher dimensional models may not be so well approximated along a single manifold, or may correspond to other bifurcations or transitions [@TE2013]. As such, analyses of regime shifts should not be too quick to sweep these details under the rug. Mechanistic understanding will always require models.
Acknowledgements
CB acknowledges support in part from USDA National Institute of Food and Agriculture, Hatch project CA-B-INS-0162-H. The authors also wish to acknowledge Alan Hastings, whose off-hand question at a seminar talk inspired the question and analysis considered here.
\newpage
References
cboettig/bs-tipping documentation built on May 5, 2019, 7:08 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Introduction
Ecological systems prone to abrupt changes, such as a fishery collapsing, can be challenging to manage. These abrupt dynamics are often referred to as regime shifts, and result from nonlinearities in the underlying system dynamics. Therefore, ecologists have often turned to studying the dynamics of nonlinear but simple mathematical models to better understand and manage ecosystems. There are essentially two different mechanisms by which a system exhibiting alternative stable states can experience a sudden transition from one state to another [@Scheffer2001; @Ashwin2012]: state tipping (S-tipping) or a bifurcation (B-tipping). Figure \ref{fig1} illustrates these two mechanisms. B-tipping is a sudden transition that is induced by slow changes to a parameter of the underlying dynamical system that reduce the stability and ultimately destroy the stable node, as it collides with an unstable saddle-point in a saddle-node bifurcation (aka a fold bifurcation or "blue sky catastrophe," [@Strogatz2001]). In contrast, S-tipping involves a perturbation directly to a state variable that carries the system across the separatrix into an alternative stable state, without any such change to the underlying shape or existence of either attractor. S-tipping can be driven by direct exogenous forcing or purely stochastic events, such as extreme events that are common in the biological time series of ecosystems [@Batt2017]. While real-world transitions will often be a combination of both processes, one mechanism will typically have a dominant role. For instance, near enough to the bifurcation, intrinsic noise will frequently drive the system across the tipping point before the actual bifurcation event is reached. Crucially, only one of these mechanisms, B-tipping, will exhibit the phenomenon of Critical Slowing Down (CSD), a pattern that has received great attention in both theoretical and empirical studies for its potential as an early warning signal of such a regime shift [e.g. @Wissel1984; @Scheffer2009; @Scheffer2015].
Despite this attention, convincing empirical demonstrations of critical slowing down are still rare, particularly outside carefully controlled laboratory settings [e.g. @Dai2012]. A foundational example can be found in studies of lake ecosystems [e.g. @Carpenter2008], which can exist in either a "piscivore dominated" or "planktivore dominated" state (Figure \ref{fig2}). Significant previous work has established evidence for these alternative states and the potential for sudden transitions between the states that cannot be easily reversed. Then in 2011, a landmark study of @Carpenter2011 analyzed the experimental manipulation of an initially planktivore-dominated state as a definitive example of CSD phenomena providing an early warning signal prior to a sudden transition to the piscivore dominated state. Yet was this manipulation truly an example of B-tipping, where CSD is theoretically predicted, or S-tipping, which does not create the CSD phenomenon?
\begin{figure} \centering \includegraphics[width=.8\linewidth]{figures/scheffer-fig.png} \caption{S-tipping vs B-tipping. Both panels show the bifurcation diagram for a saddle-node or "fold" bifurcation. These plot the equilibrium points of a system on the vertical axis (e.g. adult bass density in our example) as a function of environmental conditions (parameters of the system). Solid lines indicate stable points, dotted lines indicate unstable points. The top panel illustrates Bifurcation Tipping: a system which begins in the desirable, "high" ecosystem state (top branch) experiences gradual degradation as the underlying environmental conditions change (increasing values along the $x$). At the F2 position, this stable node collides with the unstable saddle point (dotted line) and is annihilated, leaving the system to rapidly transition to the remaining undesirable stable state (Forward shift). To restore the system to the desirable state, it is not sufficient to merely restore the environmental conditions observed at F2. Instead, environmental conditions must improved all the way to F1, where another bifurcation occurs that creates a rapid backwards shift to the desirable state. This effect is known as hysteresis. The bottom panel introduces state tipping: the transition from the desirable state to the undesirable state is triggered without changing the environmental conditions. Instead, a large perturbation is introduced directly to the ecosystem state, sufficient to push the system from one basin-of-attraction to another. Such transitions are characterized by large, sudden shocks to the state variable, rather than slow, gradual change to the environmental conditions. These transitions do not exhibit hysteresis -- a sufficiently large perturbation could restore the system. Figure reproduced from Scheffer et al (2001) with permission.}\label{fig1} \end{figure}
B-tipping is typically considered to be the result of slow change in a parameter of [e.g. @Scheffer2009] the system which ultimately leads to the fold bifurcation (or close enough to a fold bifurcation for stochasticity to finish the job). S-tipping involves manipulation of the state variable. @Carpenter2011 describe driving the critical transition through the gradual addition of piscivores (largemouth Bass), which drives a transition characterized by low-piscivore density into one of high piscivore density. If the state (piscivore density) is being manipulated directly, then is this not an example of S-tipping, and thus unsuitable to test for the presence of CSD? Although the simplicity of ball-in-cup (potential well) based descriptions of regime shifts makes them alluring, this question is impossible to address properly without a more mathematically precise description of the system and the manipulation involved. Surprisingly, despite the central importance the piscivore transitions have played in our understanding of ecosystem regime shifts, these details have never been clearly presented in prior publications of this classic lake ecosystem experiment [@Carpenter2008; @Carpenter2011]. Here, we provide a direct analysis of the underlying model in the piscivory alternative state system to demonstrate why these studies are indeed an example of B-tipping.
Central to resolving this question is untangling a second issue that frequently arises in the discussion of alternative stable state dynamics in ecological systems: the issue of dimension. Ecosystems are notoriously complex, and this lake ecosystem is no exception. Here, we analyze a previously published model of a trophic cascade to determine whether manipulating the top predator constitutes B-tipping or S-tipping. Common versions of this model include 5 state variables, ranging from top predators to primary producers (Figure \ref{fig2}). In sharp contrast to this complexity, both intuition and theory of critical transitions have been centered on one-dimensional interpretations. The S-tipping and B-tipping transitions illustrated in Figure \ref{fig1} are fundamentally rooted in a one-dimensional analysis: the state consists of a single continuously valued dynamical variable (i.e. population density of adult bass). Similarly, Ball-in-cup style depictions of alternative stable states and sudden transitions present only a one-dimensional depiction of a system state. While it may be tempting to dismiss such simple depictions as out-of-hand as mere cartoons, it is in fact possible to map the full five-dimensional model explicitly into a one-dimensional approximation under a few simple assumptions. Here, we present the full dimension reduction on the full 5D model analytically, and compare this approximation to numerical simulations of the original system. This exercise facilitates the interpretation of the lake experiment in the context of predicting regime shifts by using a 1D model to understand a high dimensional system.
\begin{figure} \centering \includegraphics[width=.8\linewidth]{figures/Slide1.png} \caption{Food Web Heuristic. Adult bass give rise to juveniles, while eating their own young and pumpkinseed sunfish. Sunfish compete with and eat juvenile bass. Sunfish eat herbivorous zooplankton like \emph{Daphnia}, which in turn eat phytoplankton.}\label{fig2} \end{figure}
Methods
Model Description
In the model of the lake food web, largemouth bass (\emph{Micropterus salmoides}, A) are the top predator. Intuitively, juvenile bass (J) mature into A, and J are produced by A. Bass feed on other fish, including their own young (J) and other planktivorous fish (F), (i.e. pumpkinseed sunfish, \emph{Lepomis gibbosus} and similar forage fish species). These planktivores and juvenile bass (which are also planktivores) eat herbivorous zooplankton (H), which in turn eat phytoplankton (P). Because bass fry are small enough to be eaten by "planktivorous" fish, J can be consumed by F. Furthermore, F and J compete with each other for food (H). Thus, there are strong pairwise interactions between the three fish groups, and the dynamics of these groups extend into the planktonic portion of the food web. The five-dimensional model of the food web includes the aforementioned state variables consisting of three fish groups (A, F, J) and two planktonic groups (H, P):
$$\begin{array}{llr} \dot{A} &= sJ - qA - m_A A &(1) \ \dot{F} &= D_F(F_o - F) - c_{FA}FA &(2)\ \dot{J} &= fA - c_{JA}JA - \frac{c_{JF}vJF}{h+v+c_{JF}F} - sJ &(3)\ \dot{H} &= D_H(H_R-H) + \alpha c_{HP}HP - c_{HF}HF &(4)\ \dot{P} &= r_P L \gamma (I_0,P)P - mP - c_{PH}HP &(5)\ \end{array}$$
Decreases in adult bass abundance result from mortality ($m_A A$) and harvest effort ($q$), and increases result from the maturation of the fraction of $J$ that survive ($s$) mortality (Eq. 1). Stocking the lake corresponds to $q < 0$. This continuous-time formulation for harvest or stocking effort can be justified when carried out throughout the year at a fixed rate, as in the experiments considered here, rather than in large shocks. Planktivorous fish diffuse ($D_F$) between a refuge of size $F_o$ and a foraging arena, and the predation rate of $A$ on $F$ is controlled by $c_{FA}$ (Eq. 2). Growth in $J$ is controlled by the fecundity ($f$) of $A$, but are subject to losses via cannibalization at a rate $c_{JA}$, and interactions with $F$ which are controlled by the rate at which $J$ hide ($h$) in a refuge or become vulnerable ($v$) to predation at a rate $c_{JF}$ as they enter an implicit foraging arena (Eq. 3). Herbivorous zooplankton diffuse between surface waters and a deep water refuge (of size $H_R$) at a rate $D_H$, grow by preying on phytoplankton at rate $c_{HP}$ and converting phytoplankton biomass (phosphorus) by a factor $\alpha$, and decline as they are preyed upon by planktivores at a rate $c_{HF}$ (Eq. 4). Phytoplankton death is due to mortality and sinking at rate $m$ or predation by herbivores at rate $c_{HF}$, while their growth is controlled by the per capita growth rate ($r_P$), the availability of phosphorus (introduced at rate $L$), and the availability of light at the surface ($I_0$) and the density dependent growth response to phytoplankton in response to light ($\gamma$, Eq. 5).
Table 1. Model parameters, description, values, and units. Parameters are taken from [@Carpenter2008] and [@Carpenter2013]. The function $\gamma$ describes phytoplankton growth, and can be found in [@Carpenter2008].
|Symbol |Description |Quantity |Units | |:------|:-----------------------------------------------------------|:--------|:-------------------| |$q$ |Harvest or stocking |variable |$t^{-1}$ | |$s$ |Survival rate of juveniles to maturation |0.6 |$t^{-1}$ | |$m_A$ |Mortality rate of adults |0.4 |$t^{-1}$ | |$m$ |Phytoplankton mortality |0.1 |$t^{-1}$ | |$f$ |Fecundity of adult bass |2 |$J/A$ | |$c_{FA}$ |Predation of planktivores by adult bass |0.3 |$t^{-1} A^{-1}$ | |$c_{JA}$ |Predation of juvenile bass on by adult bass |0.1 |$t^{-1} A^{-1}$ | |$c_{JF}$ |Predation of juvenile bass by planktivorous fish |0.5 |$t^{-1} F^{-1}$ | |$c_{PH}$ |Predation of phytoplankton by herbivores |0.25 |$t^{-1} H^{-1}$ | |$c_{HF}$ |Predation of herbivores by planktivores |0.1 |$t^{-1} F^{-1}$ | |$F_o$ |Abundance of planktivorous fish in non-foraging arena |200 |$F$ | |$H_R$ |Reservoir of herbivores in a deep-water refuge |4 |$H$ | |$D_F$ |Diffusion of planktivores between refuge and foraging arena |0.09 |$t^{-1}$ | |$D_H$ |Diffusion of between shallow and deep water |0.5 |$t^{-1}$ | |$v$ |Rate at which J enter a foraging arena, becoming vulnerable |80 |$t^{-1}$ | |$h$ |Rate at which J hide in a refuge |80 |$t^{-1}$ | |$r_P$ |Growth rate of phytoplankton |3 |$m^2 mg^{-1}$ | |$L$ |Phosphorus loading rate |0.6 |$mg m^{-2} t^{-1}$ | |$\alpha$ |Assimilation of phytoplankton phosphorus by zooplankton |0.3 |$g g^{-1}$ | |$\gamma$ |Density-dependent phytoplankton growth response to light |n/a |n/a | |$I_0$ |Solar irradiance incident on the surface of the water |300 |$\mu E m^{-2} s^{-1}$ |
Dimension Reduction
Our approach to determining whether the lake experiment constituted B-tipping or S-tipping is to reduce the dimensionality of the model described in Eqs. 1-5, and see if a 1-D expression of this model can produce alternate states and the corresponding fold bifurcation. Although transitions occur throughout the food web, we focus here on the fish dynamics (Eqs. 1-3) as they are most directly related to either harvest ($q$) or direct manipulation of adult bass (A). In the model, the dynamics of herbivorous zooplankton (H) and phytoplankton (P) do not influence A, F, or J. Therefore, the dimensionality of the model can easily be reduced from 5 to 3 by ignoring Eqs. 4-5 (the plankton) without loss of generality for the dynamics of Eqs. 1-3 (the fish). In addition to parameters, the dynamics of each fish state depend on itself and at least one other fish state. We observe that the dynamics of both the Juvenile bass (J) and forage fish (F) occur on fast timescales relative to the dynamics of long-lived adult bass. Therefore, we set Eqs. 2-3 to zero and solved for the estimate of the respective state variables given equilibrium:
$$\begin{array}{llr} \hat A&= \frac{sJ}{q+m_A} &(6)\ \hat F&=\frac{D_FF_0}{D_F + c_{FA}A} &(7)\ \hat J&=\frac{fA}{c_{JA}A+s+\frac{c_{JF}vF}{h+v+c_{JF}F}} &(8)\ \end{array}$$
The equation of motion for A (an expression of A's dynamics completely in terms of A and fixed parameters), is acquired by substituting the RHS of Eq. 7 for F in Eq. 8, then substituting the resulting expression for J in Eq. 1. After some substitution, we can write the equation of motion for the adult bass, A, completely in terms of the adult bass population and fixed parameters of the system.
$$\begin{array}{llr} \dot{A} &= sfA\Bigg(c_{JA}A+s+c_{JF}v\bigg(\frac{(h+v)(D_F+c_{FA}A)}{D_F F_o}+c_{JF}\bigg)^{-1}\Bigg)^{-1} - qA - m_A A &(9) \end{array}$$
Solving for the roots of this equation demonstrates alternative stable states and the potential for B-tipping with changing $q$, as illustrated graphically in Figure \ref{fig4}.
Simulations
We performed simulations in which we initialized the 5D model at equilibrium with a low but negative harvest rate ($q=-0.001$), which can be interpreted as per capita stocking of $A$. We proceeded to linearly increase $q$ at each time step over 1020 years with a time step size of 0.01 years (102000 time steps total) until $q$ reached a value of 0.001. At each time step a small amount of additive, normally distributed noise was added to Eqs. 1-3, truncated to non-zero values. R code for these simulations is included in supplementary materials (https://github.com/cboettig/bs-tipping).
Results
Plotting this equation using established ranges for ecological parameters (Table 1; @Carpenter2008; @Carpenter2013), Figure \ref{fig4} shows that the system should support two alternative stable states, and can undergo a fold bifurcation (in reverse, creating a new stable point) when "catch" $q$ is decreased sufficiently. This sudden transition is felt in other dimensions as its impact cascades down to the lower trophic levels.
\begin{figure} \centering \includegraphics[width=.8\linewidth]{../analysis/manuscript_report/figure-5dsimulation-1.png} \caption{Time series of the 5 state variables, with (q) (harvest rate of adult bass) being increased linearly over time.}\label{fig3} \end{figure}
These same dynamics can be seen in the direct numerical simulations of the full 5-dimensional model. At a $q$ of 0, the fish in the simulation underwent a critical transition (Figure \ref{fig3}). This transition occurred during time step 510, although an abrupt change in state was not observed until step 850. The changes in A and J were similar, but were inversely related to those in F. The predation effect of A on J was small compared to the role of A's fecundity in producing J and A's role as a predator of F, which negatively impacts J as well. Overall, both bass time series decreased with increasing $q$ (abruptly when stocking converted to harvest at $q = 0$), whereas the abundance of planktivorous fish was positively correlated with $q$.
\begin{figure} \includegraphics[width=\linewidth]{../analysis/manuscript_report/figure-Amotion-ballincup-1.png} \caption{Top row: motion of A for different values of $q$. Bottom row: the negative integral of the top row, depicting in the so-called ball-in-cup (or potential well) diagrams. Note that the peaks and troughs (local maxima and minima) of the potential wells, bottom row, correspond to the zeros of the rate equations (top row). Troughs indicate stable nodes, peaks indicate unstable saddle points. Note that under sufficient stocking of bass, $q=-0.5$, the zero point is an unstable equilibrium, indicated by the peak in the potential well (bottom right).} \label{fig4} \end{figure}
The fish dynamics cascaded down to lower trophic levels, although H and P did not undergo abrupt transitions at the same time as the fishes. Instead, $H$ and $P$ exhibited strong nonlinear changes much earlier in the time series after $q$ had increased only slightly. The small cumulative increase in $q$ leading up to step 20 caused a slight increase in $F$ and a dramatic decline in $H$ ($F$'s prey); after its initial decline, the abundance of $H$ was 0.001 for the remainder of the time series. When $H$ abruptly declined, $P$ abruptly increased, then continued to increased approximately linearly for the remainder of the time series. Therefore, increasing $q$ caused declines in both bass states ($A$ and $J$), an increase in planktivorous fish ($F$), a decline in herbivorous zooplankton ($H$), and an increase in phytoplankton ($P$), a pattern consistent with a trophic cascade. All state variables exhibited abrupt state changes at some point in the time series.
Discussion
These results demonstrate that the system considered can indeed be driven through fold bifurcations as the parameter $q$ is slowly tuned up or down (Figure \ref{fig5}), forcing the creation or annihilation of an alternative stable state -- confirming the presence of B-tipping in this system. Simulations of the same system using previous estimates for ecological parameters confirm that the predictions of the 1-dimensional approximation (Eq 9) correspond with the dynamics observed in the full model.
\begin{figure} \centering \includegraphics[width=.8\linewidth]{../analysis/manuscript_report/figure-bifurcationDiagram-1.png} \caption{Bifurcation diagram derived from Eq 9. Solid lines are stable equilibria, the dashed line is the unstable equilibrium or separatrix.}\label{fig5} \end{figure}
This analysis also highlights the importance of what might appear to be a merely semantic issue that in fact corresponds to fundamentally different dynamics. In this formulation, the control parameter $q$ corresponds to harvest (positive $q$) or stocking (negative $q$) of adult bass in the lake. It is this aspect that first gives the appearance of state tipping dynamics, since it suggests the state is being somehow "directly manipulated". Yet from a mathematical standpoint, it is irrelevant if the control parameter corresponds to some environmental or policy variable (fishing regulations, fish migration into the lake), which in turn increases or decreases net mortality proportionally, or if the manipulation directly adds or removes adult bass -- the feature that distinguishes this manipulation from a S-tipping transition comes down to the scales involved. This can be seen most clearly though an analogy of considering the same model applied to a larger scale, such as marine fisheries. In fisheries models, one would routinely examine the dynamics of the fish population biomass $B$ with respect to different possibly parameters for harvests or fishing effort. Even though harvest may be measured in the same units as the total biomass, it is natural to think of the two quantities very differently. A manager can adjust the harvest effort by changing quotas, limiting permits for fishing boats etc. Slowly increasing in the amount of bass stocked in the lake (i.e. decreasing $q$, effectively running the scenario of increasing fishing harvests in reverse) gives the rest of the food web time to adjust. In contrast, a single large introduction of many bass would create an S-tipping transition, driving the system immediately into an alternate state without exhibiting CSD phenomena needed to evaluate early warning signals of a tipping point.
All though we have focused on the analysis of deterministic models for simplicity, these results can easily be interpreted in the context of stochastic models as well. (Recall that in practice, the detection of early warning signals associated with CSD requires some level of stochasticity to generate the variance needed to detect slowing down). Under a stochastic model, the distinction between S-tipping and B-tipping is slightly blurred, since a system that gets close enough to the tipping point due to changing parameters (B-tipping) will ultimately be able to jump across the separatrix (S-tipping) before the actual bifurcation event. In a stochastic system, fluctuations in the state variable act as intrinsic sources of perturbations, and though part of the system, large chance deviations are no more predictable than sudden exogenous forcing: there is no CSD prior to purely stochastic transitions [@PRSB2013], even though retrospectively a large deviation will share similar statistical signatures such as rapidly increasing variance and autocorrelation [@PRSB2012].
A second issue highlighted by this analysis is the contrast between the interpretation of the five-dimensional food web model and the classic one-dimensional description of tipping points shown in Fig \ref{fig1}. Such one-dimensional descriptions of a state are frequently presented in reviews or abstract discussion of regime shifts and early warnings with little or no discussion of how a complex system comes to be represented by a one-dimensional continuous variable [e.g. @Scheffer2001; @Scheffer2009; @Barnosky2012; @Scheffer2012; @Clements2018], and even literature focusing on specific systems such as the lake-ecosystem considered here can be vague about the dimensional reduction [@Carpenter2008; @Carpenter2011]. This can in turn create confusion as to whether such one-dimensional diagrams are meant to be interpreted as mathematically explicit approximations of higher-order dynamics, or merely as cartoons meant only to mimic the qualitative pattern but bearing no quantitative relationship with the actual process. The analysis presented here for the lake-ecosystem model demonstrates that the former does indeed hold, though only under the appropriate time-scale limits -- not only that of the "slow parameter" but also that of the other dynamic variables in the system: plankton, Juvenile Bass, Forage Fish, Fig \ref{fig2}. Our analysis also highlights that how the nonlinearity required to create a fold bifurcation arises directly from this dimension-reducing approximation. In the original five-dimensional model, Eq (1), the dynamics of the adult bass $A$ appear to be completely linear. Only after accounting for the feedback in the other dimensions do we see the nonlinearity required to create a fold bifurcation. Our numerical simulations provide an illustration of how closely the dynamics of the full five-dimensional system are captured by the one-dimensional approximation using a biologically plausible parameterization. Thus the reference to the dynamics of a one-dimensional system state, as well as the interpretation of the transition as B-tipping rather than S-tipping, appears to be well-justified in the lake-ecosystem case.
This analysis does not vindicate the use of one-dimensional models more generally. The popularity and simplicity of the fold-bifurcation model has sometimes led to more casually equating any sudden shift with fold bifurcations, overlooking any assumptions necessary to approximate a complex dynamical system with such a one-dimensional model. Influential reviews on abrupt change in ecological systems [e.g. @Scheffer2015; @Ratajczak2018] focus almost exclusively on descriptions of one-dimensional state-spaces without discussion of if, when or how complex systems such as the ecosystem food-web considered here can be appropriately described by a nonlinear equation in a single dimension. More careful modeling informed by system-specific knowledge is required before any abrupt transition can be characterized as a fold bifurcation [@Nature2013]. The fact that CSD phenomena can occur in many different models has been suggested to mean that CSD is a generic or universal early warning indicator [@Scheffer2009] of a critical transition. In practice, many higher dimensional models may not be so well approximated along a single manifold, or may correspond to other bifurcations or transitions [@TE2013]. As such, analyses of regime shifts should not be too quick to sweep these details under the rug. Mechanistic understanding will always require models.
Acknowledgements
CB acknowledges support in part from USDA National Institute of Food and Agriculture, Hatch project CA-B-INS-0162-H. The authors also wish to acknowledge Alan Hastings, whose off-hand question at a seminar talk inspired the question and analysis considered here.
\newpage
References
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.