\SweaveOpts{ eval=TRUE, prefix=FALSE, width=4.5, height=4.5, eps=FALSE, include=FALSE} \chapter{Multiple Basins of Attraction} <>= setwd("~/Documents/Projects/Book/draft/Chap08") rm(list=ls()) badfonts = FALSE .detach.stuff = function() { s1 = grep("package|Autoloads",search()) nattach = length(search()) xattach = search()[c(-1,-s1)] for (i in xattach) eval(substitute(detach(i),list(i=i))) } .detach.stuff() ## ps.options(family="NimbusSan") ps.options(colormodel="cmyk") options(width=75, SweaveHooks=list(fig=function() par(mar=c(5,4,.5,1.5))) ) ## palette(gray((0:8)/8)) library(primer) trellis.par.set(canonical.theme(color=FALSE)) #source("DMBpplane.R") require(deSolve); @ %% [Insert Figure at start] \section{Introduction} Historically, simple models may have helped to lull some ecologists into thinking either that (i) models are useless because they do not reflect the natural world, or (ii) the natural world is highly predictable.\footnote{One notable exception was May's work revealing that chaos erupts from a very simple model of discrete logistic growth.} Here we investigate how simple models can creaate unpredictable outcomes, in models of Lotka--Volterra competition, resource competition, and intraguild predation. In all cases, we get totally different outcomes, or alternate stable states, depending on different, stochastic, initial conditions. \subsection{Alternate stable states} \emph{Alternate stable states}, or \index{alternate stable equilibria}alternate stable equilibria (ASE), \emph{are a set of two or more possible stable attractors that can occur given a single set of external environmental conditions}. For a single list of species, there may exist more than one set of abundances that provide stable equilibria. One key to assessing alternate stable states is that \emph{the stable states must occur with the same external environmental conditions}. If the external conditions differ, then the system is merely governed by different conditions.\footnote{A very important complication is that if an abiotic factor is coupled dynamically to the biotic community, then it becomes by definition an internal part of the system and no longer external to it.} If stable attractors exist, how does the system shift from one attractor to another? The system can be shifted in a variety of ways, but the key is that the system (i.e., the set of species abundances) gets shifted into a different part of the state space, and then is attracted toward another state. System shifts may occur due to demographic stochasticity, the random variation in births and deaths that may be especially important when populations become small. System shifts may also occur due to different assembly sequences. For instance the outcome of succession may depend upon which species arrive first, second, third, etc. System shifts might also arise via a physical disturbance that causes mortality. Different abundances may also arise from the gradual change and return of an environmental factor, and the resulting set of alternate equilibria is known as \emph{\index{hysteresis}hysteresis}, and below we examine a case related to resource competition. The term \emph{\index{priority effects}priority effects} refers to the situation in which initial abundances favor one species or group of species over others. We use ``priority'' to imply that if a species is given an early advantage (``priority'') then it can suppress other species, whereas without that head start, it would not do as well. In contrast, for most stable models we considered previously, the location of the attractor was independent of any early advantage. We use the terms \emph{effects of initial conditions} and priority effects to refer to situations in which starting conditions influence the outcome. So, in a 4 billion year old world, what are ``starting conditions''? Starting conditions can refer to any point resulting from a system shift, as we described above (e.g., a disturbance). Part of our thinking about ASEs relies on two assumptions that some find untenable, or at least not very useful \cite{Hastings2004a}. First, the concept assumes that fixed stable states exist and that, if perturbed, populations can return to these states. Second, the concept assumes that communities move to these states with all due haste, and achieve these states over observable time periods. If these conditions are met, then ASEs can arise and be observed. We might recognize that ASEs might be a fuzzier concept, in which there exists more than one attractor, and where complex dynamics and long return times make the attractors difficult to observe. Further, we can imagine that short-term \emph{evolution of species' traits} cause attractors to shift through time. Last, we might also want to acknowledge that saddles (i.e., attractor--repellers) also influence dynamics and species abundances. If we embrace all of these complications, then we might prefer term other than ASE for this complex dynamical landscape --- \emph{multiple basins of attraction}. \subsection{Multiple basins of attraction} \emph{\index{Multiple basins of attraction}Multiple basins of attraction} (MBA) is a phrase describing an entire landscape of ``tendencies'' in community dynamics. Imagine for a moment a two- or three-dimensional coordinate system in which each axis is the abundance of one species. An attractor is merely a place in the coordinate system (i.e., a particular set of species abundances) that exerts a pull on, or attracts, the dynamics of the populations. A stable equilibrium is an example of a global attractor; each local minimum or maximum in a stable limit cycle is also an attractor. The unstable equilibrium we observed in a Lotka--Volterra competition model is another example --- it is a saddle, or an attractor--repellor, because it attracts from one direction, but repels in another. <>= z <- 2 * volcano x <- 10 * (1:nrow(z)) y <- 10 * (1:ncol(z)) par(mar=c(0,0,0,0)) persp(x, y, -z, theta = 120, phi = 30, col="green",shade=.75, scale = FALSE, axes = FALSE) <>= par(mar=c(0.5,0.5,0.5,0.5)) contour(x,y,z) @ \begin{figure}[ht] \centering \includegraphics[width=.48\textwidth]{valleymba} \qquad \includegraphics[width=.48\textwidth]{contourmba} \caption{Perspective and contour plots of a single complex dynamical landscape, containing multiple basins of attraction. Imagine putting a ball into this landscape, and jiggling the landscape around. The landscape represents the possible states of the community, the ball represents the actual community structure at any one point in time, and the jiggling represents stochasticity, either demographic, or environmental. \label{fig:mba}} \end{figure} \index{MBA|see{multiple basins of attraction}}MBAs can be visualized as a topographic landscape, or mountain range. We see lots of little valleys (attractors) and lots of peaks (repellers) (Fig. \ref{fig:mba}). If we think more broadly about multiple basins of attraction, then we begin to see that a strict definition of ASS lies at one end of a continuum, and it is matched at the other end by a system with one global stable equilibrium. The entire continuum, and much more besides, can be conceived of in terms of a landscape with some number of both basins of attraction (attractors) and peaks (repellors) in the landscape; each of these may act with some unique force or strength. The potential for MBAs to occur has been noted for a few decades \cite{Holling1973, Lewontin1969, May:1977yq, Noy-Meir:1975fj}, and the hunt for their signature in nature followed \cite{Chase2003b, Connell1983, Long2002, Messier:1994jb, Vandermeer2004a, Vandermeer1999, Scheffer:2003vf, Schmitz:2004rt, Schroder:2005qe}. There has been discussion about what causes more than one basin of attraction, including space, stochasticity, and predation. Below we examine two mechanisms that may reveal priority effects: strong interference competition \cite{Reynolds1993}, and \emph{intraguild predation}, an indirect interaction which has elements of both competition and predation \cite{Polis:1989ly}. We first investigate interference competition using a three-species Lotka--Volterra model, and then with a model of resource competition. \section{Lotka--Volterra Competition and MBA} You have already seen in lecture the case of two-species \index{Lotka--Volterra!multiple basins of attraction}\index{Lotka--Volterra!three-species competition}Lotka--Volterra competition where an attractor--repeller, or saddle, arises when interspecific competition is greater than intraspecific competition. When this is the case, each species can suppress the other, if it gets a head start. This is a priority effect. Having a larger negative effect on your competitor than on yourself may not be too unusual. Examples that come immediately to mind are cases where species compete preemptively for space (territories, or substrate surface) or for a resource with a unidirectional flow (drifting prey or particles in a stream, or light coming down on a forest canopy), then one species may gain the upper hand by acquiring a large proportion of resources first, that is, preempting resources. In \index{human economic systems}human economic systems, businesses can have direct negative effects on each other through anti-competitive business practices. For instance, a larger company can temporarily flood a local market with below-cost merchandise, and drive out smaller competitors, prior to raising prices again. Note that it requires the raised prices for its long term equilibrium, but uses temporary below-market prices to eliminate competitors. In contrast, economies of scale can provide a very different kind of mechanism by which larger businesses can outcompete smaller businesses at equilibrium. Here we explore how MBA can arise in a simple competitive system. We use a three-species Lotka--Volterra model to illustrate how strong interference competition may reveal \index{priority effects}priority effects. We use a slightly different representation of the Lotka--Volterra competition model. Note that the sign of each $\alpha$ must be negative, because we are adding them --- this is similar to the notation in Chapter 7, but differs from previous treatment of Lotka--Volterra competition (Chapters 3, 5). \begin{align*} \frac{\D N_1}{\D t} &= r_1N_1\left( 1+\alpha_{11}N_1 + \alpha_{12}N_2 +\alpha_{13}N_3 \right)\\ \frac{\D N_2}{\D t} &= r_2N_2\left( 1+\alpha_{21}N_1 + \alpha_{22}N_2 +\alpha_{23}N_3 \right)\\ \frac{\D N_3}{\D t} &= r_3N_3\left( 1+\alpha_{31}N_1 + \alpha_{32}N_2 +\alpha_{33}N_3 \right) \end{align*} Note two aspects of the above equations. First note that within the parentheses, the $N_i$ are all arranged in the same order --- $N_1,\,N_2,\,N_3$. This reflects their relative positions in a food web matrix, and reveals the row-column relevance of the subscripts of the $\alpha$s. Second, note that $\alpha_{ii}=1/K_i$, and in some sense $K_i$ \emph{results from} a particular $\alpha_{ii}$.\footnote{We might also note is that it differs from notation typically used in textbooks for two-species models.} We can represent these equations as \begin{equation} \label{eq:compact3} \frac{\D N_{i}}{\D t}=N_{i}\left(r_{i} + \sum_{j=1}^{3}r_i\alpha_{ij}N_{j}\right) \end{equation} in a fashion similar to what we saw in Chapter 7. Another aspect of three species \index{Lotka--Volterra!equilibrium and $r$}Lotka--Volterra models that we will state, without demonstration, is that \emph{the equilibria can depend on $r_i$}. This contrasts with what we learned about the two-species model, which depended solely on the $\alpha_{ij}$. Indeed, our example uses parameter values in which variation in $r$ contributes to the equilibria. @ \begin{figure}[ht] \centering \label{fig:LVP} \includegraphics[width=11cm]{LVpriority} \caption{Interaction between strong interference competition and initial abundances with three competing species (solid line - sp. 1; dashed line - sp. 2; dotted line - sp. 3). The species with the highest initial abundance is indicated in vertical orientation at the beginning of each trajectory (``Equal'' indicates all species started at $N=200$). The eventual winner is indicated in horizontal orientation at the top of each graph. See text for more details.} \end{figure} Consider a particular set of simulations (Fig. \ref{fig:LVP}). \begin{compactitem} \item Species differ in their intrinsic rates of increase ($r_1=0.6$, $r_2=1$, $r_3=2$). \item Species with higher $r$ also have slightly greater negative impacts on themselves ($a_{11}=0.001$, $a_{22}=0.00101$, $a_{33}=0.00102$); this constitutes a tradeoff between maximum growth rate and carrying capacity. \item All species have twice as big a negative effect on each other as they do on themselves ($a_{ij, i\neq j}=0.002$); this creates unstable equilibria, which happen to be attractor--repellers or saddles. \item Initial abundances are random numbers drawn from a normal distribution with a mean of 200, and a standard deviation of 10; in addition we also start them at precisely equal abundances ($N_i=200$). \item Species 1 has the highest carrying capacity (smallest $\alpha_{ii}$), and would therefore often considered the best competitor; in the absence of priority effects, we would otherwise think that it could always displace the others. \end{compactitem} Now ponder the results (Fig. \ref{fig:LVP}). Does each species win in at least one simulation? Which species wins most often? Does winning depend on having the greatest initial abundance? Make a $3 \times 3$ table for all species pair combinations (column = highest initial abundance, row=winner), and see if there are any combinations that never happen. Note that when they start at equal abundances (Fig. \ref{fig:LVP}), the species with the intermediate carrying capacity and intermediate $r$ displaces the other two species. However, note also that (i) with a little stochasticity in initial conditions (slight variation in $N_i$ at $t=0$), this simple model generates unpredictable outcomes, and (ii) initial abundance does not determine everything. It is critical to realize that this is occurring in part because species have larger negative competitive effects on others than they have on themselves. In this case the effects are direct, because the model is Lotka--Volterra competition. The effects may also be indirect, when species compete for more than one limiting resource. MacArthur \cite{MacArthur1972} showed, for instance, that when generalists and specialists compete so that not all resources are available to all species, alternate stable states occur \cite{Long2002}. After we work through the code for the Lotka--Volterra example we just discussed, we take up MBA in the context of resource competition. \subsection{Working through Lotka--Volterra MBA} Here we create a function for multi-species Lotka--Volterra competition, taking advantage of matrix operations. Note that we can represent the three species as we would one, $\bm{\dot{N}} = \bm{rN}\left(1-\bm{\alpha N} \right)$, because the $\bm{ \alpha N}$ actually becomes a matrix operation, \texttt{a \%*\% N}.\footnote{Recall that \%*\% is matrix multiplication in \R~because by default, \R~multiplies matrices and vectors element-wise.} <<>>= lvcompg <- function(t, n, parms) { r <- parms[[1]]; a <- parms[[2]] dns.dt <- r * n * (1 - (a%*%n)) return( list(c(dns.dt)) ) } @ We are going to use one set of parameters, but let initial abundances vary stochastically around the unstable equilibrium point, and examine the results. Next we decide on the values of the parameters. We will create intrinsic rates of increase, $r$s, and intraspecific competition coefficients, $\alpha_{ii}$, that correspond roughly to an \index{tradeoffs!$r$ \emph{vs.} $K$}$r-K$ tradeoff, that is, between maximum relative growth rate ($r$) and carrying capacity ($1/\alpha_{ii}$). Species 1 has the lowest maximum relative growth rate and the weakest intraspecific density dependence. Following these ecological guidelines, we create a vector of $r$s, and a matrix of $\alpha$s. We then put them together in a \emph{list},\footnote{A list is a specific type of \R~object.} and show ourselves the result. <<>>= r <- c(r1=.6, r2=1, r3=2) a <- matrix(c(a11=.001, a12=.002, a13=.002, a21=.002, a22=.00101, a23=.002, a31=.002, a32=.002, a33=.00102), nrow=3, ncol=3) parms <- list(r,a) parms @ Next we get ready to simulate the populations 24 times. We set the time, $t$, and the \emph{mean} and \emph{standard deviation} of the initial population sizes. We then create a matrix of initial population sizes, with one set of three species' $n_0$ for each simulation. This will create a $3 \times 24$ matrix, where we have one row for each species, and each column is one of the initial sets of population sizes. <<>>= t=seq(0,40, by=.1) ni <- 200; std=10 N0 <- sapply(1:30, function(i) rnorm(3, mean=ni,sd=std) ) @ Now let's replace the first set of initial abundances to see what would happen if they start out at precisely the same initial abundances. We can use that as a benchmark.\footnote{\R's recycling rule tells it to use the single value of \texttt{ni} for all three values in the first column of \texttt{N0}.} <<>>= N0[,1] <- ni @ When we actually do the simulation, we get ready by first creating a graphics device (and adjust the margins of each graph). Next we tell \R~to create a graph layout to look like a 6 $\times$ 4 matrix of little graphs. Finally, we run the simulation, calling one column of our initial population sizes at a time, integrate with the ODE solver, and plot the result, 24 times. As we plot the result, we also record which species has the greatest initial abundance. <>= par(mar=c(2,2,1,1)) layout(matrix(1:30, ncol=5)) for(i in 1:30) {lvout <- ode(N0[,i], t, lvcompg, parms) matplot(t, lvout[,2:4], type="l", lwd=1.5, col=1) if(all(N0[,i]==200)) { text(3, 500, "Equal", srt=90) }else { text(3, 500, paste("Sp.", which.max(N0[,i])), srt=90)} lastN <- lvout[nrow(lvout),2:4] text(3, max(lastN), paste("Sp.", which.max(lastN)), adj=c(0,1)) } @ \section{Resource Competition and MBA} Above, we explored how simple Lotka--Volterra competition could result in unstable equilibria, causing saddles, and multiple basins of attraction. Here we take a look an example of how resource competition can do the same thing. Recall that \emph{resource competition is an indirect interaction}, where species interact through shared resources. This particular example results in a type of MBA scenario, \emph{\index{hysteresis}hysteresis}, where gradual changes in the external environment result in abrupt and sometimes catastrophic changes in the biological system (Fig. \ref{fig:hyst}). Scheffer and colleagues \cite{Scheffer:2003vf} provide evidence that anthropogenically enriched (eutrophic) lakes can shift away from dominance by submerged \index{macrophytes}macrophytes\footnote{``Macrophyte'' is a term often used for aquatic rooted vascular plants.} rooted in substrate, into systems completely dominated by floating plants such as \index{floating plants}duckweed (\emph{Lemna} spp.) and water fern (\emph{Azolla} spp.). Submerged macrophytes can extract nutrients out of both sediments and the water column. At low, typically unpolluted, nutrient levels, submerged plants can draw down water nitrogen levels to a very low level, below levels tolerated by duckweed and water fern. At high nutrient levels, floating plants are no longer limited by water column nitrogen levels, and can create deep shade that kills submerged vegetation. Aside from killing these wonderful submerged macrophytes, the loss of this structure typically alters the rest of the lake food web. \begin{figure}[ht] \centering\label{fig:scheffer1} \subfloat[Submerged~Plants]{\includegraphics[width=6cm]{Scheffer6.pdf}\label{fig:sch4}} \subfloat[Floating~Plants]{\includegraphics[width=6cm]{Scheffer5.pdf}\label{fig:sch3}} \caption{The outcome of competition depends on both nutrient loading and priority effects. The lines provide equilibrial biomasses that depend on nitrogen supply rate for \subref{fig:sch4} submerged plants, and \subref{fig:sch3} floating plants. Thus both lines in \subref{fig:sch4} refer to possible equilibria for submerged plants, and both lines in \subref{fig:sch3} refer to possible equilibria for floating plants. (S = submerged plants, F = floating plants). } \end{figure} Uncertainties arise in this scenario at intermediate levels of eutrophication (Fig. \ref{fig:scheffer1}). As stated above, submerged plants dominate at low nitrogen supply rates, and floating plants dominate at high nitrogen supply rates. However, at intermediate supply rates ($\sim 1.5$--$2.5$\,mg\,L$^{-1}$), the outcome depends on priority effects. When submerged plants initially dominate a lake prior to eutrophication, they continue exclude floating plants as nitrogen supply rates increase, up until about 1.5\,mg\,L$^{-1}$ (Fig. \ref{fig:sch4}). If supply rates go higher ($\sim 1.5$--$2.5$\,mg\,L$^{-1}$), submerged plants continue to dominate, because their high abundance can draw nitrogen levels down quite low. Above this level, however, floating plants reach sufficiently high abundance to shade out the submerged plants which then become entirely excluded. Thus, at $2.5$\,mg\,L$^{-1}$ we see a catastrophic shift in the community. Once the system is eutrophic, and dominated by floating plants, reducing the nitrogen supply rate does not return the system to its original state in a way that we might expect (Fig. \ref{fig:sch4}). Once the floating plants have created sufficient shade, they can suppress submerged plants even at low nitrogen levels. If changes in watershed management bring nitrogen supply rates back down, submerged plants cannot establish until supply rate falls to 1\,mg\,L$^{-1}$. This is well below the level at which submerged plants could dominate, if they were abundant to begin with. This is an excellent example of a system prone to \emph{hysteresis} (Fig. \ref{fig:hyst}). Gradual changes in an external driver (e.g., nitrogen runoff) cause catastrophic shifts in the internal system, because at intermediate levels of the driver, the internal system is influenced by priority effects. Hysteresis is often defined as a context-dependent response to an environmental driver. As a different example, imagine that the state variable in Figure \ref{fig:hyst} is annual precipitation, and the driver is average annual temperature. Imagine that over many years, the regional average temperature increases. At first, increased temperature has little effect on precipitation. Once precipitation reaches a particular threshold, precipitation drops dramatically, and additional increased temperature has little effect. When the temperature starts to come back down, however, we find that the high levels of precipitation do not return at the same threshold where we lost it. Rather, it does not return until we bring temperature way back down to original levels. At intermediate temperatures (grey region, Fig. \ref{fig:hyst}), \emph{precipitation depends on what the temperature used to be}. This history dependence, or context dependence is the hallmark of hysteresis. <>= par(mar=c(3,3,1,1), mgp=c(1,0,0), tck=0) x <- seq(-2, 5, by=0.01) y.ex <- expression(x^3 - 4*x^2) y <- eval(y.ex) plot(y, -x, type='n', axes=F, ylab="State Variable\n(e.g. Precipitation)", xlab="Environmental Variable\n(e.g. Temperature)") polygon(c(-10, 0, 0, -10), c(-4.5, -4.5, 2, 2), col="grey", border=NA) lines(y, -x) lines(y[x>=0 & x <= 8/3], -x[x>=0 & x <= 8/3], col=0, lty=2) box() xs <- x[c(1, 60, 701, 660)] x <- xs ys <- eval(y.ex) arrows(ys[1], -xs[1]-.5, ys[2], -xs[2]-.5, length=.1) arrows(ys[3], -xs[3]+.5, ys[4], -xs[4]+.5, length=.1) #segments(-9, -4.5, -1, -4.5, lwd=10, col="grey") @ \begin{figure}[ht] \centering \includegraphics[width=8cm]{Hysteresis.pdf} \caption{Hysteresis. When the system changes from left to right, the threshold value of the predictor of catastrophic change is greater than when the system moves from right to left. At intermediate levels of the driver, the value of the response variable depends on the history of the system. Each value, either high or low, represents an alternate, stable, basin of attraction. The arrows in this figure represent the direction of change in the environmental driver.} \label{fig:hyst} \end{figure} It is important to understand that, in principle, this is \emph{not} the result of a time lag. It is \emph{not the case} that this pattern is due to a lagged or slow response by the state variable. Rather, these alternate basins (Fig. \ref{fig:hyst}, solid lines in the grey area) represent permanent stable states from which the response variable cannot ever emerge, without some external force, such as very large changes in the environmental driver. Time lags may be important in other, different, circumstances, but are not, in principle, related to hysteresis. \subsection{Working through resource competition} Scheffer and colleagues represented submerged and floating plant interactions in the following manner, where $F$ and $S$ are floating and submerged plants, respectively \cite{Scheffer:2003vf}. \begin{align} \label{eq:mba3} \frac{\D F}{\D t} &=r_fF\frac{n}{n+h_f}\frac{1}{1+a_fF}-l_fF\\ \frac{\D S}{\D t} &=r_sS\frac{n}{n+h_s}\frac{1}{1+a_sS+bF+W}-l_sS\\ \end{align} As usual, $r$ represents the maximum per capita rate of increase for $F$ and $S$ respectively. Thus $rF$ is exponential growth of floating plants. This exponential growth is modified by terms for nitrogen ($n$) limited growth, and light limited growth. This modification results in type II functional responses. Here we discuss these modifications. \begin{description} \item[{\bf Nitrogen limitation}] The first factor to modify exponential growth above is $n/\left(n+h_x\right)$, water column nitrogen limitation. It varies from 0--1 and is a function of water column nitrogen concentration, $n$. The half saturation constant $h$ controls the shape of the relation: when $h=0$ there is no nitrogen limitation, that is, the term equals 1.0 even at very low $n$. If $h>0$, then the growth rate is a Michaelis-Menten type saturating function where the fraction approaches zero as $n \rightarrow 0$, but increases toward 1.0 (no limitation) at high nutrient levels. For submerged plants, $h_s=0$ because they are never limited by water column nitrogen levels because they derive most of the nitrogen from the substrate. \item[{\bf Nitrogen concentration}] Nitrogen concentration, $n$, is a simple saturating function that declines with increasing plant biomass which achieves a maximum $N$ in the absence of any plants. \begin{equation} \label{eq:4} n=\frac{N}{1+q_sS+q_fF} \end{equation} The nutrient concentration, $n$, depends not only on the maximum $N$ nutrient concentration, but also on the effect of submerged and floating plants which take up nitrogen out of the water at rates $1/\left(1+q_SS\right)$ and $1/\left(1+q_fF\right)$ respectively; $1/q$ is the half-saturation constant. \item[{\bf Light limitation}] The second factor to modify exponential growth is light limitation, $1/(1 + aF)$; $1/a_f$ and $1/a_s$ are \emph{half-saturation constants} --- they determine the plant densities at which the growth rates are half of the maxima; $b$ represents the shade cast by a single floating plant, and $W$ represents the light intercepted by the water column. \item[{\bf Loss}] The second terms in the above expressions $l_fF$ and $l_sS$ are simply density independent loss rates due to respiration or mortality. \end{description} Scheffer and colleagues very nicely provide units for their parameters and state variables (Table \ref{tab:Scheffer}). \begin{table}[ht] \centering \caption{Parameter and variable units and base values. Plant mass (g) is dry weight. \label{tab:Scheffer}} \begin{tabular}[]{lcl} \hline \hline {\bf Parameter/Variable} & {\bf Value} & {\bf Units}\\ \hline $F$, $S$ & (varies) & g\,m$^-2$\\ $N$, $n$ & (varies) & mg\,L$^-1$\\ $a_s$, $a_f$ & 0.01 & (g\,m$^-2$)$^-1$\\ $b$ & 0.02 & (g\,m$^-2$)$^-1$\\ $q_s$, $q_f$ & 0.075, 0.005 & (g\,m$^-2$)$^-1$\\ $h_s$, $h_f$ & 0.0, 0.2 & mg\,L$^-1$\\ $l_s$, $l_f$ & 0.05 & g\,g$^-1$\,day$^-1$\\ $r_s$, $r_f$ & 0.5 & g\,g$^-1$\,day$^-1$\\ \hline \end{tabular} \end{table} Consider the meanings of the parameters (Table \ref{tab:Scheffer}). What is $a$, and why is $1/a_s=1/a_f$? Recall that they are half-saturation constants of light limited growth; their identical values indicate that both plants become self-shading at the same biomasses. What is $q$? It is the per capita rate at which plants pull nitrogen out of the water. Why is $q_s > q_f$? This indicates that a gram of submerged plants can pull more nitrogen out of the water than a gram of floating plant. Last, why is $h_s=0$? Because submerged plants grow independently of the nitrogen content in the water column. \subsubsection{Code for Scheffer \emph{et al.}} Now we are set to model this in \R, using a built-in function for the ODEs, \texttt{scheffer}. First, let's set the parameters (Table \ref{tab:Scheffer}), time, initial abundances, and see what we have. <<>>= p <- c(N=1, as=0.01, af=0.01, b=0.02, qs=0.075, qf=0.005, hs=0, hf=0.2, ls=0.05, lf=0.05, rs=0.5, rf=0.5, W=0) t <- 1:200 Initial <- c(F=10, S=10) @ We then run the solver, and plot the result. <>= S.out1 <- ode(Initial, t, scheffer, p) matplot(t, S.out1[,-1], type='l') legend('right', c("F", "S"), lty=1:2, bty='n') @ From this run, at these nutrient levels, we observe the competitive dominance of the submerged vegetation (Fig. \ref{fig:sch1a}). Let's increase nitrogen and see what happens. <>= p['N'] <- 4 S.out2 <- ode(Initial, t, scheffer, p) matplot(t, S.out2[,-1], type='l') @ \begin{figure}[ht] \centering \subfloat[Low~Nitrogen]{\includegraphics[width=6cm]{Scheffer1.pdf}\label{fig:sch1a}} \subfloat[High~Nitrogen]{\includegraphics[width=6cm]{Scheffer2.pdf}\label{fig:sch2a}} \caption[Sudden shifts in ecosystems]{The outcome of competition depends on the nutrient loading.} \end{figure} Ah-ha! At high nutrient levels, floating vegetation triumphs (Fig. \ref{fig:sch2a}). So where are the cool multiple basins of attraction? We investigate that next. Let's mimic nature by letting the effect of \emph{Homo sapiens} increasing gradually, with gradually increasing over-exploitation of the environment. We will vary $N$, increasing it slowly, and hang on to only the final, asymptotic abundances, at the final time point. <<>>= N.s <- seq(.5,4, by=0.1) t <- 1:1000 S.s <- t( sapply(N.s, function(x) {p['N'] <- x ode(Initial, t, scheffer, p)[length(t),2:3] }) ) @ Now we plot, not the time series, but rather the asymptotic abundances \emph{vs.} the nitrogen levels (Fig. \ref{fig:sch1b}). <>= matplot(N.s, S.s, type='l') legend('right', c('F','S'), lty=1:2, bty='n') @ Now we can see this catastrophic shift at around 2.7 mg\,N\,L$^-1$ (Fig. \ref{fig:sch1b}). As nitrogen increases, we first see a gradual shift in community composition, but then, wham!, all of a sudden, a small additional increase at $\approx 2.7$ causes dominance by floating plants, and the loss of our submerged plants. Now let's try to fix the situation by reducing nitrogen levels. We might implement this by asking upstream farmers to use no-till practices, for instance \cite{Vanni2001}. This is equivalent to starting at high floating plant abundances, low submerged plant abundances, and then see what happens at different nitrogen levels. <<>>= Initial.Eutrophic <- c(F=600, S=10) S.s.E <- t( sapply(N.s, function(x) {p['N'] <- x ode(Initial.Eutrophic, c(1,1000), scheffer, p)[2,2:3] }) ) @ Now we plot, not the time series, but rather the asymptotic abundances \emph{vs.} the nitrogen levels (Fig. \ref{fig:sch2b}). <>= matplot(N.s, S.s.E, type='l') @ \begin{figure}[ht] \centering \subfloat[Low~Initial~N]{\includegraphics[width=6cm]{Scheffer3.pdf}\label{fig:sch1b}} \subfloat[High~Initial~N]{\includegraphics[width=6cm]{Scheffer4.pdf}\label{fig:sch2b}} \caption[Sudden shifts in ecosystems]{The outcome of competition depends on the history of nutrient loading.} \end{figure} What a second! If we start at high floating plant biomass (Fig. \ref{fig:sch2b}), the catastrophic shift takes place at a much lower nitrogen level. This is telling us that from around $N=$0.9--2.7, the system has two stable basins of attraction, or alternate stable states. It might be dominated either by floating plants or by submerged plants. This is often described as \emph{hysteresis}, where there is more than one value of the response for each value of the predictor variable. Now let's represent these changes for the two state variables in the aquatic plant model. First we represent the floating plants. Here we plot the low abundance state for the floating plants, adjusting the figure margins to accommodate all abundances, and then add in the high abundance data (Fig. \ref{fig:sch3}). <>= quartz(,4,4) <<>>= plot(N.s[1:23], S.s[1:23,1], type='l', lwd=2, xlim=c(0,4), ylim=c(0,900), main="Floating Plants", ylab=expression("Biomass (g m"^-2*")"), xlab="Nitrogen Supply Rate") lines(N.s[-(1:5)], S.s.E[-(1:5),1], lwd=2) @ Here we reinforce the concepts of multiple basins and hysteresis, by showing where the attractors are. I will use arrows to indicate these basins. At either high nitrogen or very low nitrogen, there is a single, globally stable attractor. At low nutrients, only submerged plants exist regardless of starting conditions. At high nutrients, only floating plants persist. Let's put in those arrows (Fig. \ref{fig:sch3}). <<>>= arrows(3, 10, 3, 620, length=.1); arrows(3, 820, 3, 720, length=.1) arrows(0.5, 620, 0.5, 50, length=.1) @ Next we want arrows to indicate the alternate basins of attraction at intermediate nitrogen supply rates. Floating plants might be at kept at low abundance at intermediate nitrogen supply rates if submerged plants are abundant (Fig. \ref{fig:sch2b}). Let's indicate that with a pair of arrows. <<>>= arrows(2.5, -10, 2.5, 60, length=.1); arrows(2.5, 200, 2.5, 100, length=.1) text(2.5, 100, "Coexisting\nwith S", adj=c(1.1, 0)) @ Alternatively, if submerged plants were at low abundance, floating plants would get the upper hand by lowering light levels, which would exclude submerged plants altogether (Fig. \ref{fig:sch3}). Let's put those arrows in (Fig. \ref{fig:sch3}). <<>>= arrows(2, 480, 2, 580, length=.1); arrows(2, 750, 2, 650, length=.1) text(2, 700, "Monoculture", adj=c(1.1,0) ) <>= dev.print(pdf, "Scheffer5.pdf") @ Now let's repeat the exercise with the submerged plants. First we plot the high abundance state, and then add the low abundance state (Fig. \ref{fig:sch4}). <>= quartz(,4,4) <<>>= plot(N.s[1:23], S.s[1:23,2], type='l', lwd=2, xlim=c(0,4), ylim=c(0,900), main="Submerged Plants", ylab=expression("Biomass (g m"^-2*")"), xlab="Nitrogen Supply Rate") lines(N.s[-(1:5)], S.s.E[-(1:5),2], lwd=2) @ Now we highlight the global attractors that occur at very low or very high nitrogen supply rates (Fig. \ref{fig:sch4}). <<>>= arrows(.7, 30, .7, 830, length=.1) arrows(3.8, 830, 3.8, 30, length=.1) @ Next we highlight the local, alternate stable equilibria that occur at intermediate nitrogen supply rates; either the submerged plants are dominating due to nitrogen competition, and achieving high abundance, <<>>= arrows(2.3, 650, 2.3, 750, length=.1); arrows(2.3, 900, 2.3, 800, length=.1) text(2.4, 900, "Coexisting\nwith F", adj=c(0, 1) ) @ or they are excluded entirely, due to light competition (Fig. \ref{fig:sch4}). <<>>= arrows(2, 130, 2, 30, length=.1) text(2, 140, "Excluded\nDue to Light Comp.", adj=c(.5,-.3) ) <>= dev.print(pdf, "Scheffer6.pdf") @ Once again, we see what underlies these alternate states, or basins (Fig. \ref{fig:scheffer1}). One population gains a numerical advantage that results in an inordinately large negative effect on the loser, and this competitive effect comes at little cost to the dominant species. At intermediate nitrogen supply rates, the submerged vegetation can reduce ambient nitrogen levels in the water column to undetectable levels because it gets most of its nitrogen from sediments. On the other hand, if floating plants can ever achieve high densities (perhaps due to a temporary nutrient pulse), then the shade they cast at intermediate supply rates prevents lush growth of the submerged plants. As a consequence, the submerged plants can never grow enough to draw nitrogen levels down to reduce the abundance of the floating plants. \section{Intraguild Predation} Intraguild predation (\index{IGP|see{intraguild predation}}IGP) differs from omnivory only in degree (Chapter 7). In omnivory, a predator shares a resource with one or more of its prey (Fig. \ref{fig:igp}). Thus the top predator feeds relatively high on the food chain, getting most of its energy or resources by eating its competitor ($a>0.5$ in Fig. \ref{fig:igp}). An extension of this is the case of \emph{\index{intraguild predation}intraguild predation}, in which a species preys upon one or more of its competitors (Fig. \ref{fig:igp}). Intraguild predation is thus refers to the case in which the top predator gets most of its energy or resources from the more basal resource, eating lower on the food chain ($a<0.5$ in Fig. \ref{fig:igp}). The distinction is not qualitative, but rather quantitative. If both consumer species prey upon each other, then we could make the argument that the name we ascribe to it depends entirely upon one's perspective. In such a case, however, we generally refer to the relations as intraguild predation. \begin{figure}[ht] \centering \label{fig:igp} \includegraphics[width=.33\textwidth]{IGP} \caption{We typically use ``omnivory'' when $a>0.5$, and ``intraguild predation'' when $a<0.5$. If we remove $A$ from this model, then the species represent those of Holt and Polis \cite{Holt1997a}.} \end{figure} \subsection{The simplest Lotka--Volterra model of IGP} We can extend our good ol' Lotka--Volterra competition model to describe intraguild predation. All we do is add a term onto each competitor. For the competitor that gets eaten (the ``IG-prey''), we subtract mass action predation, with a constant attack rate, $a$. For the top predator (the ``IG-predator''), we add this same term, plus an conversion efficiency, $b\ll 1$. \begin{align} \label{eq:1} \frac{\D N_1}{\D t} &= r_1N_1\left(1-\alpha_{11}N_1 - \alpha_{12}N_2\right) + ba N_1 N_2\\ \frac{\D N_2}{\D t} &= r_2N_2\left(1-\alpha_{21}N_1 - \alpha_{22}N_2\right) - a N_1 N_2 \end{align} Here $a$ is attack rate of the IG-predator, $N_1$, on the IG-prey, $N_2$; $b$ is the conversion efficiency of the prey into predator growth. Recall that this is the classic type I predator functional response of Lotka--Volterra predation. Let's work through a little logic. \begin{itemize} \item In the absence of the other species, each species will achieve its usual carrying capacity, $1/\alpha_{ii}$. \item If we could have stable coexistence without IG-predation, then adding predation will increase the risk of extinction for the prey, and increase the abundance (if only temporarily) for the predator. \item If the poorer competitor is able to feed on the better competitor, this has the potential to even the scales. \item If the poor competitor is also the prey, then --- forget about it --- the chances of persistence by the IG-prey are slim indeed. \end{itemize} Now let's move on to a model of intraguild predation with resource competition. \subsection{Lotka--Volterra model of IGP with resource competition} Here we introduce a simple IGP model where the competition between consumers is explicit resource competition \cite{Holt1997a}, rather than direct competition as in the \index{Lotka--Volterra!intraguild predation}Lotka--Volterra model above. The resource for which they compete is a logistic population. \begin{align} \label{eq:1} \frac{\D P}{\D t}&=\beta_{PB}\alpha_{BP}PB + \beta_{PN}\alpha_{NP}PN - m_PP\\ \frac{\D N}{\D t}&= \beta_{NB}\alpha_{BN}BN - m_NN - \alpha_{NP}PN\\ \frac{\D B}{\D t}&=rB \left(1-\alpha_{BB}B \right) - \alpha_{BN}BN - \alpha_{BP}PB \end{align} Recall that the units for attack rate, $\alpha$, are number of prey (killed) per individual of prey per individual of predator; the units for \index{conversion efficiency}conversion efficiency, $\beta$, are number of predators (born) per number of prey (killed, and presumably eaten and assimilated). The consumers in this model have a type I functional response (mass action). The basal resource species exhibits logistic population growth. Holt and Polis show analytically that five equilibria are present \cite{Holt1997a}. \begin{compactenum} \item All species have zero density. \item Only the resource, $B$, is present, at $B=K$. \item Only the resource, $B$, and IG-prey, $N$, are present. \item Only the resource, $B$, and IG-predator, $P$, are present. \item All species present. \end{compactenum} We will explore how initial conditions influence the outcomes of this simple IGP model. We will focus on the last three equilibria, with two or three species present. Are there lessons we can apply from the previous competition models? In those cases, we tended to get MBAs when the negative effects each species on its competitors was greater than its negative effects on itself. How can we apply that to predation? \begin{figure}[ht] \centering \subfloat[Initial densities]{\includegraphics[width=.48\textwidth]{IGPAltRA0}\label{IGPRA0}} \subfloat[Final densities]{\includegraphics[width=.48\textwidth]{IGPAltRA500}\label{IGPRA500}} \caption{Initial abundance \subref{IGPRA0} determines whether the IG-predator, $P$, or the IG-prey, $N$, win \subref{IGPRA500}. Parameters are those set with the vector \texttt{params1} (see below).} \label{fig:IGPAlt} \end{figure} Let us think about net effects of the IG-predator and IG-prey, both competition and consumption. Recall that the IG-prey must be the superior competitor --- this means that, given similar attack rates on the basal resource $B$ ($\alpha_{BN}=\alpha_{BP}$), the IG-prey must have a greater conversion efficiency ($\beta_{NB}>\beta_{PB}$).\footnote{Think of conversion efficiency as the effect of the prey $i$ on the predator $j$, $\beta_{ij}$.} In the absence of IG-predation, the superior competitor would always exclude the inferior competitor. However, if we add predation to suppress the superior competitor, then the IG-predator could win. If the relationship is such that each species has a larger net effect on each other than on themselves, we see that initial abundance (Fig. \ref{IGPRA0}) determines the winners (Fig. \ref{IGPRA500}). This allows either species to win, but not coexist. \begin{figure}[ht] \centering \subfloat[Low $\beta_{PN}$, high $\alpha_{NP}$]{\includegraphics[width=.48\linewidth]{IGPtrial1.pdf}\label{fig:AltTrial1} } \subfloat[High $\beta_{PN}$, Low $\alpha_{NP}$]{\includegraphics[width=.48\linewidth]{IGPtrial2.pdf}\label{fig:AltTrial2} } \caption[Avoiding ASEs]{Conversion efficiencies and attack rates control coexistence. \subref{fig:AltTrial1} With low conversion effciency, and high attack rates, species do not coexist. \subref{fig:AltTrial2} By reducing the attack rate of the predator on the prey ($\alpha_{NP}=10^{-4} \to 10^{-7}$), and increasing the direct benefit of prey to the predator ($\beta_{PN}=10^{-5} \to 0.5$), we get coexistence. Parameters are otherwise the same as in Fig. \ref{fig:IGPAltRA} (see \texttt{params2} below).} \end{figure} How might we get coexistence between IG-predator and IG-prey? We have already provided an example where the IG-prey is the better resource competitor (Fig. \ref{fig:IGPAlt}). To reduce the negative effect of the IG-predator on the IG-prey, we can reduce attack rate. However, when we do that, the predator cannot increase when it is rare (Fig. \ref{fig:AltTrial1}). If we further allow the predator to benefit substantially from each prey item (increasing conversion efficiency), then we see that the IG-predator can increase when rare, but eliminate the prey. Indeed, these are the essential components suggested by Holt and Polis: species coexist when their negative effects on each other are weaker (of smaller magnitude) than their negative effects on themselves. In a consumer--resource, predator-prey context, this can translate to reduced attack rates, and greater efficiency of resource use. \subsection{Working through an example of intraguild predation} To play with IBP in \R, we start by examining an \R~function for the above Lotka--Volterra intraguild predation model, \texttt{igp}. <<>>= igp @ This code uses three-letter abbreviations ($\alpha_{NP}$ = \texttt{anp}). The first letter, \texttt{a} or \texttt{b}, stands for $\alpha$ and $\beta$. The next two lower case letters correspond to one of the populations, $B$, $N$, and $P$. Next, we create a vector to hold all those parameters. <<>>= params1 <- c(bpb= 0.032, abp=10^-8, bpn=10^-5, anp=10^-4, mp=1, bnb=0.04, abn=10^-8, mn=1, r=1, abb=10^-9.5) @ Here we get ready to actually do the simulations or numerical integration with \texttt{ode}. We set the time, and then we set four different sets (rows) of initial population sizes, label them, and look at them. <<>>= t=seq(0, 60, by=.1) N.init <- cbind(B = rep(10^9, 4), N = 10^c(2, 5, 3, 4), P = 10^c(5,2,3, 4) ) @ Now we integrate the population dynamics and look at the results. Here we first set up a graphics device with a layout of four figures and fiddle with the margins. We then use a for-loop to integrate and plot four times. Then we add a legend. <<>>= quartz(,4,4); layout(matrix(1:4, nr=2)); par(mar=c(4, 4, 1, 1)) for(i in 1:4) { igp.out <- ode(N.init[i,1:3], t, igp, params1) matplot(t, log10(igp.out[,2:4]+1), type="l", lwd=2, ylab="log(Abundance)") } <>= dev.print(pdf, file="IGPDyn.pdf") @ \begin{figure}[ht] \centering \includegraphics[width=3in]{IGPDyn} \caption{Dynamics of Lotka Volterra intraguild predation, with differing initial abundances. See code for parameter values. Solid line - basal resource, dashed line - IG-prey, dotted line - IG-predator. } \label{fig:IGPDyn} \end{figure} Clearly, initial abundances affect which species can coexist (Fig. \ref{fig:IGPDyn}). If either consumer begins with a big advantage, it excludes the other consumer. In addition, if they both start at low abundances, the IGP prey, $N$, excludes the predator; if they start at moderate abundances, the IGP predator, $P$, wins. Now we need to get more thorough and systematic. The above code and its results show us the dynamics (through time) of particular scenarios. This is good, because we need to see how the populations change through time, just to see if unexpected things happen, because sometimes unexpected dynamics happen. A complementary way to analyze this model is to vary initial conditions more systematically and more thoroughly, and then simply examine the end points, rather than each entire trajectory over time. It is a tradeoff --- if we want to look at a lot of different initial conditions, we can't also look at the dynamics. In the next sections, we examine the effects of \emph{relative} abundance of the two consumers, and then of their \emph{absolute} abundances. \subsection{Effects of relative abundance} First we will vary the relative abundances of the IGP prey and predator, $N$ and $P$. We create a slightly more complete set of initial abundances, with $B$ constant, and $N$ increases as $P$ decreases. <<>>= logNP <- seq(2, 5, by=.1) N.inits <- cbind(B=rep(10^9, length(logNP)), N=10^logNP, P=10^rev(logNP) ) @ We see (scatterplot matrix not shown) that we do have negative covariation in the starting abundances in the two consumer species, the IG-prey and IG-predator. Next, we need to perform all% \footnote{Recall that \texttt{sapply} and related functions ``apply'' a function (in this case a simulation) to each element of the first argument (in this case each row number of the initial abundance matrix).} the simulations, and hold on to all the endpoints.\footnote{We transpose the output matrix (\texttt{t()}) merely to keep the populations in columns. We also use the \texttt{hmax} argument in \texttt{ode} to make sure the ODE solver doesn't try to take steps that are too big.} We do it over a long time span to see the (hopefully) truly asymptotic outcomes. We use a little manipulation to hang on to the initial abundances, at $t=50$ and the final abundances at $t=500$, putting them each in their own column and hanging on to it. <<>>= t1 <- 1:500 MBAs <- t( sapply(1:nrow(N.inits), function(i) { tmp <- ode(N.inits[i,], t1, igp, params1, hmax=0.1) cbind(tmp[1,3:4],tmp[50,3:4],tmp[500,3:4]) } ) ) colnames(MBAs) <- c('N1','P1','N50', 'P50', 'N500','P500') @ Now we need to show our results. We are interested in how the relative initial abundances of the two consumers influence the emergence of MBA. Therefore, let's put the ratio of those two populations (actually the logarithm of the ratio, $log[N/P]$\footnote{logarithms of ratios frequently have much nicer properties than the ratios themselves. Compare \texttt{hist(log(runif(100)/runif(100)))} vs. \texttt{hist(runif(100)/runif(100))}.}) on an X-axis, and graph the abundances of those two species on the Y-axis. Finally, we plot side by side the different time points, so we can see the initial abundances, the transient abundances, and (perhaps) something close to the asymptotic abundances. <>= matplot(log10(N.inits[,"N"]/N.inits[,"P"]), log10(MBAs[,1:2]+1), type="l", col=1, lty=2:3, lwd=2, ylab="log(Abundance+1)", xlab="log[N/P]") legend("right", c("N", "P"), lty=2:3, col=1, bty="n") <>= matplot(log10(N.inits[,"N"]/N.inits[,"P"]), log10(MBAs[,3:4]+1), type="l", col=1, lty=2:3, lwd=2, ylab="log(Abundance+1)", xlab="log[N/P]") <>= matplot(log10(N.inits[,"N"]/N.inits[,"P"]), log10(MBAs[,5:6]+1), type="l", col=1, lty=2:3, lwd=2, ylab="log(Abundance+1)", xlab="log[N/P]") @ It is still amazing to me that different initial abundances can have such a dramatic effect (Fig. \ref{fig:IGPAltRA}). It is also interesting that they take so long to play out. It all really just makes you wonder about the world we live in. \subsection{Effects of absolute abundance} Now let's hold relative abundance constant and equal, and vary absolute abundance. Recall that in our first explorations, we found different outcomes, depending on different total abundances. Now instead of varying $N$ and $P$ in opposite order, we have them covary positively. <<>>= logAbs <- seq(2, 7, by=.2) N.abs.inits <- cbind(B=rep(10^9, length(logAbs)), N=10^logAbs, P=10^logAbs ) @ Now we simulate\footnote{Unfortunately 'simulate' may mean 'integrate,' as it does here, or any other kind of made up scenario.} the model, using the same basic approach as above, setting the time, and hanging on to three different time points. <<>>= t1 <- 1:500 MBA.abs <- t( sapply(1:nrow(N.abs.inits), function(i) { tmp <- ode(N.abs.inits[i,], t1, igp, params1, hmax=0.1) cbind(tmp[1,3:4],tmp[50,3:4],tmp[500,3:4]) } ) ) colnames(MBAs) <- c('N1','P1','N50', 'P50', 'N500','P500') @ We plot it as above, except that now simply use $log_{10}$-abundances on the x-axis, rather than the ratio of the differing abundances. <>= quartz(,7,3) <<>>= layout(matrix(1:3, nr=1) ) matplot(log10(N.abs.inits[,"N"]), log10(MBA.abs[,1:2]+1), type="l", main="Initial Abundances (t=1)", col=2:3, lty=2:3, lwd=2, ylab="log(Abundance+1)", xlab=expression(log[10]("N"))) legend("right", c("N", "P"), lty=2:3, col=2:3, bty="n") matplot(log10(N.abs.inits[,"N"]), log10(MBA.abs[,3:4]+1), type="l", main="At time = 50)", col=2:3, lty=2:3, lwd=2, ylab="log(Abundance+1)", xlab=expression(log[10]("N"))) matplot(log10(N.abs.inits[,"N"]), log10(MBA.abs[,5:6]+1), type="l", main="At time = 500", col=2:3, lty=2:3, lwd=2, ylab="log(Abundance+1)", xlab=expression(log[10]("N"))) <>= dev.print(pdf, file="IGPAltAbs.pdf") @ \begin{figure}[ht] \centering \includegraphics[width=\textwidth]{IGPAltRA} \caption{Initial, transient, and near-asymptotic abundances of the intraguild prey, $N$, and predator, $P$, of Lotka--Volterra intraguild predation, with differing initial abundances.} \label{fig:IGPAltRA} \end{figure} \subsection{Explanation} Now, \ldots we have to explain it! Let's begin with what we think we know from Lotka--Volterra competition --- each species has a bigger effect on the others than on itself. How do we apply that here. First let's look at the per capita direct effects, the parameters for each interaction. \begin{align} \label{eq:5} \begin {array}{lccc} & B & N & P\\ B~~ & r-r\alpha_{BB}B & -\alpha_{BN} & -\alpha_{BP}\\ \noalign{\medskip} N & \beta_{NB}\alpha_{BN} & 0 & -\alpha_{NP}\\ \noalign{\medskip} P & \beta_{PB}\alpha_{BP} & \beta_{PN}\alpha_{NP} & 0 \\ \end {array} \end{align} Then we calculate the values for these and ask if competitors have larger effects on each other than they do on themselves. <<>>= params1 with(as.list(params1), { rbind(B=c(r-r*abb*10^9, -abn, -abp), N=c( bnb*abn, 0, -anp), P=c(bpb*abp, bpn*anp, 0)) }) with(as.list(params1), { rbind(B=c(r-r*abb*10^9, -abn, -abp), N=c( bnb, 0, -anp), P=c(bpb, bpn, 0)) }) @ So, from this we are reminded that the per capita direct effects on $B$, the basal resource, by both consumers are the same. $N$, the IG-P prey, however, benefits more per capita, and so can attain a higher population size, and therefore could persist, and also exclude $P$. Thus it has a larger indirect negative effect on $P$ than on itself. $P$, on the other hand, could have a huge direct negative effect on $N$. To achieve this effect, however, $P$ has to have a sufficiently large population size. That is exactly why we get the results we do. If $N$ starts out as relatively abundant, it reduces $B$ and probably excludes $P$. If, on the other hand, $P$ is abundant, they can have a large direct negative effect on $N$, and exclude $N$. Holt and Polis suggest that coexistence is more likely when (i) the IGP prey is the better competitor (as we have above) and (ii) \emph{the IG-predator benefits substantially from feeding on the IGP prey}, that is, when the conversion efficiency of prey into predators, $\beta_{PN}$, is relatively large. Let's try increasing $\beta_{PN}$ to test this idea. Let's focus on the ASS where the predator is excluded, when both species start out at low abundances (Fig. \ref{fig:IGPDyn}, upper right panel). We can focus on {\em the invasion criterion} by starting $N$ and $B$ at high abundance and then test whether the predator invade from low abundance. We can really ramp up $\beta_{PN}$ to be 100 times $\beta_{PB}$. This might make sense if $N$ nutrient value is greater than $B$. In real food chains this seems plausible, because the C:N and C:P ratios of body tissue tend to decline as one moves up the food chain \cite{Sterner:2002mi}; this makes animals more nutritious, in some ways, than plants. <<>>= params2 <- params1 params2['anp'] <- 10^-7 @ Now we numerically integrate the model, and plot it. <>= t <- 1:500 N.init.1 <- c(B=10^9, N=10^7, P=1) trial1 <- ode(N.init.1, t, igp, params2) matplot(t, log10(trial1[,-1]+1), type='l', col=1, ylab=quote(log[10]*("Density+1"))) legend('bottomright', c('B','N','P'), lty=1:3, bty='n') @ Whoa, dude! Fig. \ref{fig:AltTrial1} reveals very different results from those in Fig. \ref{fig:IGPAltRA}, but makes sense, right? We make the predator benefit a lot more from each prey item, then the predator doesn't need to be a good competitor, and can persist even if the IG prey reduces the basal resource level. Our next step is to rein the predator back in. One way to do this is to reduce the \emph{attack rate}, so that the predator has a smaller per capita direct effect on the prey. It still benefits from eating prey, but has a harder time catching them. Let's change $\alpha_{NP}$ from $10^-4$ to $10^-7$ and see what happens. <>= params2['bpn'] <- .5 trial2 <- ode(N.init.1, t, igp, params2) matplot(t, log10(trial2[,-1]+1), type='l', col=1 , ylab=quote(log[10]*("Density+1")) ) @ Now that we have allowed the predator to benefit more from individual prey (Fig. \ref{fig:AltTrial2}), but also made it less likely to attack and kill prey, we get coexistence regardless of the predator starting at low abundance. Additional exploration would be nice, but we have made headway. In particular, it turns out that this model of intraguild predation yields some very interesting cases, and the outcomes depend heavily on the productivity of the system (i.e. the carrying capacity of $B$). Nonetheless, we have explored the conditions that help facilitate coexistence of the consumers --- the IG-prey is the superior exploitative competitor, and the IGP predator benefits substantively from the prey. \section{Summary} A few points are worth remembering: \begin{compactitem} \item Alternate stable equibria, alternate stable states, and multiple basins of attraction are defined generally as mulitple attractors in a system defined by a particular suite of species in a environment in which external fluxes are constant. \item Hysteresis is typically an example of alternate stable states that are revealed through gradual changes in an external driver, such as temperature, or nutrient supply rate. \item Alternate stable equilibria will have low invasibility by definition; this lack of invasibility might come about through large direct negative effects (high attack rate, or aggression), or through a preempted resource (e.g. light interception, an exclusive substitutable resource, or allelopathy). There could also be life history variation, where long lived adults prevent colonization by less competitive juveniles, or juveniles vulnerable to predation \cite{Polis:1989ly}. \item Alternate stable equilibria seem to be more common when species have relatively larger negative effects on each other and weaker negative effects on themselves. \end{compactitem} \section*{Problems} \addcontentsline{toc}{section}{Problems} \begin{prob} \textbf{General}\\ Compare and contrast the terms ``alternate stable states'' and ``multiple basins of attraction.'' Define each and explain how the terms are similar and how they differ. \end{prob} \begin{prob} \textbf{Lotka--Volterra competition}\\ (a) Explain what we learn from Figure \ref{fig:LVP} regarding how growth rate, initial abundance and intraspecific density dependence (or carrying capacity) influence outcomes. Specifically, which of these best predicted the final outcome of competition? Which was worst? Explain.\\ (b) Explain in non-mathematical terms why strong interference allows for priority effects.\\ (c) Create a simulation to more rigorously test the conclusions you drew in part (a) above. \end{prob} \begin{prob} \textbf{Resource competition}\\ (a) Explain hysteresis.\\ (b) Alter the equation for submerged plants to represent the hypothetical situation in which submerged plants get most or all of their resources from the water column. Explain your rationale, and experiment with some simulations. What would you predict regarding (i) coexistence and (ii) hysteresis? How could you test your predictions? \end{prob} \begin{prob} \textbf{Intraguild Predation}\\ (a) Use Figure \ref{fig:IGPAlt} to explain how initial abundances influence outcomes. Are there initial conditions that seem to result in all species coexisting? Are there things we should do to check this?\\ (b) Explain how high attack rates and low conversion efficiencies by the top predator create alternate stable states. \end{prob}