\documentclass[twoside]{article} \pagestyle{myheadings} \markboth{\hfil Asymptotic and transient analysis \hfil}% {\hfil Thomas C. Gard \hfil} \begin{document} \setcounter{page}{51} \title{\vspace{-1in}\parbox{\linewidth}{\footnotesize\noindent Fourth Mississippi State Conference on Differential Equations and\newline Computational Simulations, Electronic Journal of Differential Equations, \newline Conference 03, 1999, pp 51-62. \newline URL: http://ejde.math.swt.edu or http://ejde.math.unt.edu \newline ftp ejde.math.swt.edu (login: ftp)} \vspace{\bigskipamount} \\ % Asymptotic and transient analysis of stochastic core ecosystem models \thanks{ {\em Mathematics Subject Classifications:} 92D25, 60H10. \hfil\break\indent {\em Key words:} stochastic plankton-fish models, stochastic differential equations, \hfil\break\indent asymptotic properties, ultimate boundedness, transient behavior, exit probabilities. \hfil\break\indent \copyright 2000 Southwest Texas State University and University of North Texas. \hfil\break\indent Published July 10, 2000. } } \date{} \author{ Thomas C. Gard } \maketitle \begin{abstract} General results on ultimate boundedness and exit probability estimates for stochastic differential equations are applied to investigate asymptotic and transient properties of models of plankton-fish dynamics in uncertain environments \end{abstract} \newtheorem{thm}{Theorem} \section{Introduction} Opposing general points of view on whether or not populations ultimately survive are succinctly expressed recently by Halley and Iwasa ([10]) and Jansen and Sigmund ([12]). Is extinction certain if random variability is taken into account? Or if parameters are restricted to realistic ranges and the mitigating effects of population communities are explicitly considered, will persistence occur possibly after some initial risk period? Answers to such questions for real systems are more than pedagogical niceties. They can, for example, lead to proposed strategies for maintaining large scale bio-physico-chemical systems such as the highly utilized natural systems constituting watershed ecosystems. Central to dynamical models for watershed ecosystems are what might be called core lake ecosystem models, such as plankton-fish models which describe the dynamics of a limiting nutrient $P$, and algae $A$, zooplankton $Z$, and small fish $F$ populations \begin{eqnarray} \frac{dP}{dt} &=& \delta(P_I(t) - P) + g_P(t, P, A, Z, F)\nonumber\\ \frac{dA}{dt} &=&g_A(t,P,A,Z)-\delta A\nonumber\\ \frac{dZ}{dt}&=&g_Z(t,A,Z,F)-\delta Z \\ \frac{dF}{dt}&=&g_F(t,Z,F)+F_I(t)\nonumber \end{eqnarray} recently discussed in the literature.(See, for example Doveri et.\,al.\,[3] where specific functional forms for the interaction portions $g_P$, $g_A$, $g_Z$, and $g_F$ of the net growth rates are given.) In (1) $P_I$ denotes the nutrient input rate, and $F_I$ the small fish recruitment rate from large fish. Simplified submodels of (1) have been discussed; the $PA$ submodel is a resource-consumer model with similar dynamics to the simple chemostat model ([18]), the $PAZ$ submodel has been discussed by Ruan ([16],[17]) and others, and the $AZF$ model is a three-species food chain. The relative novelty of (1) is the explicit inclusion of small fish dynamics - the timing and size of large annual recruitment peaks simultaneously effecting and being determined by $PAZ$ levels. Temperature and other seasonality time variations of parameters together with cyclic nonlinearities in the model can lead to chaotic regimes ([15]). All models of real biological systems account for uncertainty in parameters and structure in one way or another. There is always a variability ansatz, although in many cases such assumptions are implicit. An explicit approach on the other hand is to formulate models with well-defined stochastic features in order to account for random variability. The class of stochastic differential equation models for interacting populations is such a class which can take into account environmental randomness: the SDE model analogous to (1) has the general form \begin{eqnarray} dP&=&[\delta(P_I-P)+g_P] dt+\sigma_P \,dW_P\nonumber\\ dA&=&[g_A-\delta A] dt+\sigma_A \,dW_A\nonumber\\ dZ&=&[g_Z-\delta Z] dt+\sigma_Z \,dW_Z \\ dF&=&[g_F+F_I] dt+\sigma_F \,dW_F\,,\nonumber \end{eqnarray} where $W_P$, $W_A$, $W_Z$, and $W_F$ are standard Brownian motions and $\sigma_P$, $\sigma_A$, $\sigma_Z$, and $\sigma_F$ denote the corresponding intensities of the noise fluctuations; the $\sigma$'s may be functions of the state variables and time. A specific example of (2), motivated by a stochastic model of two competitors in a chemostat suggested by Stephanopoulos, Aris, and Fredrickson [19], is given by \begin{eqnarray} dP&=&[\delta_0(P_I-P)+g_P] dt+\delta_1(P_I-P)\, dW\nonumber\\ dA&=&[g_A-\delta_0 A] dt+\delta_1A\, dW\nonumber\\ dZ&=&[g_Z-\delta_0 Z] dt+\delta_1Z \,dW\\ dF&=&[g_F+F_I] \, dt\,.\nonumber \end{eqnarray} System (3) arises when the dilution or washout rate $\delta$ is viewed as the sum of an average value $\delta_0$ plus a random noise fluctuation with intensity $\delta_1$ about the average: \begin{equation} \delta=\delta_0+\delta_1N\,. \end{equation} In (4) $N$ represents standard white noise - in a generalized sense \begin{equation} N=\frac{dW}{dt} \end{equation} with $W$ a standard Brownian motion. In the next section we give a result which obtains asymptotic estimates for the average values of the state variables in (2) which is analogous to uniform persistence for the corresponding deterministic model (1). Application to specific $PAZF$ models is incomplete at this time and the subject of ongoing work. Transient behavior of the model may be important whether or not the former result applies. The third section contains a result which gives estimates for first exit location probabilities from certain bounded sets in the feasible region which may indicate initial survival or extinction tendencies of populations. We show that this result can be applied to models of the form (3). \section{Persistence in the mean} Permanence (uniform persistence together with dissipativity) is the most basic general qualitative feature to verify for interacting population models ([11],[20]); it is the model analog of mutual survival and non-explosion of the populations represented in the model. Permanence means that there are positive constants $K$ and $L$ such that for any component population $X(t)$ with any positive initial value $X(0)$ \begin{equation} K\leq \liminf_{t \to \infty} X(t)\leq\limsup_{t\to \infty} X(t)\leq L. \end{equation} If \begin{equation} Y(t) = \ln X(t) \end{equation} or some other transformation of the ray $(0, \infty)$ to the line $(-\infty, \infty)$, permanence of $X$ is equivalent to dissipativity or ultimate boundedness of $Y$: there is a positive constant $M$ such that for any initial value $Y(0)$ \begin{equation} \limsup_{t\to\infty} |Y(t)|\leq M. \end{equation} There are well-known theorems in differential equations which give ultimate boundedness if a Liapunov function exists. In this section we will apply an analogous theorem of Miyahara ([14]) for stochastic differential equations. It is convenient to change notation here: let \begin{eqnarray} X&=&(X_{1,}X_2,X_3,X_4)=(P,A,Z,F)\,, \\ Y&=&\ln(X)\leftrightarrow Y_i=\ln(X_i),\; i =1,2,3,4. \end{eqnarray} Applying Ito's formula to (2) yields a transformed system of the form \begin{equation} dY=H(t,Y) dt+\Gamma(t,Y) dW \end{equation} where here $W $= ($W_P$, $W_A$, $W_Z$, $W_F), H$ is a 4-d vector function, and $\Gamma$ is a $4\times 4$ diagonal matrix function. \begin{thm} {\rm (Miyahara [14])} Suppose there exists a scalar function $V(t,y)$ which is ${\cal C}^1$ in $t$ and ${\cal C}^2$ in $y$ and a number $p\geq 1$ such that for some constants $a_1$ and $a_2$ and positive constants $c_1$ and $c_2$ and all $y$ \begin{enumerate} \item $V(t, y) \geq-a_1+c_1\|y\|^p$ \item ${\cal L} V(t, y) \leq a_2-c_2 V(t, y)$ \end{enumerate} \noindent where ${\cal L} V = V_t+H\cdot\nabla V+\frac{1}{2} \mathop{\rm trace}(\Gamma\Gamma^T V_{yy})$ and $V_t$ is the partial derivative of $V$ with respect to $t,\ \nabla V$ is the $y$-gradient of $V$ and $V_{yy}$ denotes the matrix of second partial derivatives of $V$ with repect to $y$. % Then for, any solution $Y(t)$ of (11), \begin{equation} \limsup_{t\to\infty} E \|Y(t)\|^p \leq \frac{a_1}{c_1}+\frac{a_2}{c_1 c_2}% \end{equation} where $E(\cdot)$ denotes the expected value or mean. \end{thm} Equations (6) - (8) suggest that the conclusion (12) of Theorem 1 could be called persistence (or permanence) in the mean for $X$. Applying Theorem 1 requires a candidate for the Liapunov function $V$. It has been shown by the author ([8]) that functions of the form \begin{equation} V(y)=\exp \{U(y)\} \end{equation} with \begin{equation} U(y)=\sum_{i = 1}^n\alpha_i[e^{y_i}-y_i-1] \end{equation} for some positive constants $\alpha_i$ can be applied to predator-prey models. The function $U$ in (14) is Volterra's Liapunov function transformed from the positive cone $R_{+}^n$ to all of $R^n$ by (10). Utilizing the $1$-norm $$ \|y\|=\sum_{i = 1}^n|y_i| $$ we obtain \begin{equation} U(y) \geq-\beta+\alpha \|y\| \end{equation} where $$ \beta=\sum_{i = 1}^n\alpha _i{\rm \ \ and \ \ }\alpha={\rm min\ } \alpha _i. $$ From (15) and (13) it follows that $$ V(y) \geq exp\{-\beta+\alpha \|y\|\}\geq e^{-\beta}(1+\alpha \|y\|). $$ Condition 1 of Theorem 1 is verified for $p=1$ with \begin{equation} a_1=0{\rm \ \ and \ \ } c_1=\alpha e^{-\beta}. \end{equation} Obtaining Condition 2 is more difficult. Generally one expects \begin{equation} a_2=a_2(\Gamma), c_2=c_2(\Gamma)% \end{equation} if the function $V$ has negative definite derivative $\dot{V}$ along trajectories of the corresponding deterministic system. The following example shows that Condition 2 can be obtained for a single population dynamics model with a Liapunov function of the form (13) even when $\Gamma$ is not necessarily small. Consider the stochastic logistic model \begin{equation} dX=X(1-X) dt+\frac{1}{\sqrt{2}}X dW \end{equation} which transforms to \begin{equation} dY=\Bigl(\frac{3}{4}-e^Y\Bigr) dt+\frac{1}{\sqrt{2}} dW \end{equation} under $Y=\ln (X)$. The Liapunov function (13) here is \begin{equation} V(y)=\exp \{e^y-y-1\}. \end{equation} It is easy to see that $$ V(y)\geq|y| $$ and so we can actually do a little better than (16): we can take \begin{equation} a_1=0{\rm \ \ and\ \ } c_1=1. \end{equation} To get Condition 2 we need to estimate \begin{eqnarray} {\cal L} V(y)&=&\Bigl(1-\frac{1}{4}-e^y\Bigr)V'(y)+\frac{1}{2}\Bigl(% \frac{1}{2}\Bigr)V''(y)\nonumber\\ &=&\Bigl(1-\frac{1}{4}-e^y\Bigr)(e^y-1)V(y)+% \frac{1}{4}[(e^y-1)^2+e^y]V(y). \end{eqnarray} A brief calculation gives \begin{equation} {\cal L} V(y)=\frac{1}{4}[1-3(e^y-1)^2]V(y)\leq\frac{1}{2}-% \frac{1}{4}V(y), \end{equation} i.e., we have Condition 2 satisfied with \begin{equation} a_2=\frac{1}{2}{\rm \ \ and\ \ } c_2=\frac{1}{4}. \end{equation} Using (21) and (24), the conclusion of the theorem yields \begin{equation} \limsup_{t\to\infty} E \|\ln (X(t))\| \leq\frac{1/2}{1/4}=2. \end{equation} Volterra's Liapunov function does not seem to work for resource-consumer models; the counterpart to inequality (23) does not hold. So there remain some problems in applying Miyahara's result to stochastic $PAZF$ models. We conclude this section by summarizing the conditions which would have to be met in order to get a specific result here. For clarity in stating the conditions, we will consider only the following simplified version of (2) \begin{eqnarray} dP&=&[\delta(P_I-P)+Pf_P(P,A)] dt+(P_I-P)\mu_P(P) \,dW_P\nonumber\\ dA&=&[Af_A(P,A,Z)-\delta A] dt+A\mu_A(A) \,dW_A\nonumber\\ dZ&=&[Zf_Z(A,Z,F)-\delta Z] dt+Z\mu_Z(Z) \,dW_Z\\ dF&=&[Ff_F(Z,F)+F_I] dt+F\mu_F(F)\, dW_F\,.\nonumber \end{eqnarray} In particular, we are ignoring nutrient recycling and we are assuming that parameters are constants. Under the log transformation: \begin{equation} Y_1=\ln (P/P_I),Y_2=\ln (A),Y_3=\ln (Z),Y_4=\ln (F) \end{equation} the system becomes \begin{eqnarray} dY_1&=&[\delta(e^{-Y_1}-1)-\frac{\mu_P^2}{2}(e^{-Y_1}-1)^2+f_P] dt+(e^{-Y_1}-1)\mu_P \,dW_P\nonumber\\ > dY_2&=&[f_A-\delta-\frac{\mu_A^2}{2}] dt+\mu_A \,dW_A\nonumber\\ > dY_3&=&[f_Z-\delta-\frac{\mu_Z^2}{2}] dt+\mu_Z \,dW_Z\\ dY_4&=&[f_F+F_Ie^{-Y_4}-\frac{\mu_F^2}{2}] dt+\mu_F \,dW_F\,,\nonumber \end{eqnarray} where in (28) $f_P=f_P(P_Ie^{Y_1},e^{Y_2})$, $\mu_P=\mu_P(P_Ie^{Y_1})$,\dots. If we choose a ${\cal C}^2$ function $V$ which satisfies Condition 1 of Miyahara's Theorem: for some number $p\geq 1$ there is a constant $a_1$ and a positive constant $c_1$ such that \begin{equation} V(y) \geq-a_1+c_1\|y\|^p \end{equation} we need to verify Condition 2. Condition 2 in Miyahara's Theorem is, for some positive constant $c$ \begin{eqnarray} {\cal L} V+cV &=&\frac{\partial V}{\partial y_1}[\delta(e^{-y_1}-1)+f_P]+\frac{\partial V}{\partial y_2}[f_A-\delta]\nonumber\\ &&+\frac{\partial V}{\partial y_3}[f_Z-\delta]+\frac{\partial V}{\partial y_4}[f_F+F_Ie^{-y_4}]\nonumber\\ &&+\frac{1}{2}\Bigl\{\mu_P^2\Bigl(% \frac{\partial^2V}{\partial y^2_1}-\frac{\partial V}{\partial y_1}\Bigr) (e^{-y_1}-1)^2+\mu_A^2\Bigl(\frac{\partial^2V}{\partial y^2_2}-\frac{\partial V}{\partial y_2} \Bigr)\\ &&+\mu_Z^2\Bigl(\frac{\partial^2V}{\partial y^2_3}-\frac{\partial V}{\partial y_3} \Bigr)+\mu_F^2\Bigl(\frac{\partial^2V}{\partial y^2_4}-\frac{\partial V}{\partial y_4} \Bigr)\Bigr\}+cV\nonumber \end{eqnarray} is bounded. The result then is \begin{thm} Suppose there exists a ${\cal C}^2$ function $V$ defined on $R^4$ which satisfies (29) and (30). Then, for any solution $Y(t)$ of (28), \begin{equation} \limsup_{t\to \infty} E \|Y(t)\|^p \leq \frac{a_1}{c_1}+\frac{b}{c_1 c} \end{equation} where $b$ is a bound for ${\cal L} V+cV$ i. e., for any solution $(P(t)$,$A(t),Z(t),F(t))$ of (26) and any positive number $\epsilon$, \begin{equation} E \|(\ln (P(t)/P_I),\ln (A(t)),\ln (Z(t)),\ln (F(t)))\|^p \leq \frac{a_1}{c_1}+% \frac{b}{c_1 c}+ \epsilon \end{equation} for all sufficiently large $t$. \end{thm} \section { Exit probabilities} Even when some form of stability can be verified for a model, transient behavior may still be important to investigate. If trajectories enter a region of state space where one or more model components are small, features neglected in the model could lead to collapse before a predicted recovery can occur. For models which attempt to account for random effects, the situation is particularly critical. Estimating certain exit statistics is a natural first approach to deal with this problem ([6],[7],[9],[13]). Suppose \begin{equation} X=(X_{1,}X_2,\ldots,X_n) \end{equation} represents the $n$ components of a stochastic dynamical population model taking values in the usual positive cone \begin{equation} R^n_{+}=\{x=(x_{1,}x_2,\ldots,x_n):x_i>0,i=1,2,\ldots,n\} \end{equation} in $n$-dimensional space, and $ B\subseteq R_{+}^n$ is a bounded set. Then for any fixed $ x\in B$, we can consider the realization $$ X=X(t,x), t\geq 0 $$ of the model with $ X(0,x)=x\in B$, and the corresponding first exit time of $X$, \begin{equation} \tau=\tau_x(B)=inf\bigl\{t:X(t,x)\notin B\bigr\} \end{equation} from $B$. The first exit time $\tau$ or even its mean or expected value \begin{equation} u(x)=E(\tau_x) \end{equation} gives an indication of persistence of $X$ relative to the set $B$ ([13]). For example, if $\tau_x=\infty$ for all $x\in B$, then $B$ is positive invariant for $X$, and if also the boundary $\partial B$ of $B$ is contained in $R_{+}^n$, then the set $B$ is a candidate for a practical persistence estimate for the model. If it could also be shown that each realization $X$ which begins at an $x$ outside $B$ hits $B$ in a finite time before hitting the boundary of $R^n_{+}$, verification of practical persistence would be complete. If the model for $X $takes the form of a stochastic differential equation \begin{equation} dX=G(X) dt+\Lambda(X) dW \end{equation} as discussed in the previous section, then it is known that the expected exit time $u$ solves the boundary value problem \begin{eqnarray} {\cal L} u(x) = -1, \quad x\in B\\ u(x) = 0, \quad x\in\partial B\,,\nonumber \end{eqnarray} where, as in Theorem 1 above, \begin{equation} {\cal L} u =G\cdot\nabla u+\frac{1}{2} \mathop{\rm trace}(\Lambda\Lambda^T u_{xx})\,. \end{equation} For example, for the simple scalar problem, \begin{equation} dX = \sqrt{\varepsilon}\ dW, X(0) = x \in B = (0,1) \end{equation} with $\epsilon$ any positive number, the boundary value problem for $u(x)=E(\tau_x)$ is \begin{equation} -1 = {\cal L}U(x) = \frac{\varepsilon}{2} u''(x), u(0) = 0 = u(1). \end{equation} The solution is easily calculated: \begin{equation} u(x)=\frac{1}{\epsilon}(x-x^2). \end{equation} Note that $\tau_x$ is finite almost surely in this example. The unit interval $B$ here is not an estimate for practical persistence; in fact persistence fails in this example. Although the above example is not a very interesting population model, it does exhibit what has become anticipated behavior of randomly perturbed deterministic models - loss of stability. In this situation the size of $\tau_x$ or $E(\tau_x)$ still can indicate relative persistence. One can also try to determine other exit statistics such as exit point location probabilities \begin{equation} v(x)=P\{X(\tau_x)\in\partial_\eta B\} \end{equation} where $\partial_\eta B$ is some particular subset of the boundary of $B$. This can be both physically relevant and mathematically tractable if the set $B$ and the boundary portion $\partial_\eta B$ are suitably chosen. Suppose $V$ is a ${\cal C}^2$ function, and \begin{equation} B\subseteq\{x:\eta\leq V(x)\leq\gamma\} \end{equation} and \begin{equation} \partial_\eta B=\partial B\cap\{x:V(x)=\eta\}. \end{equation} We have the following result (See also [8].) \begin{thm} Let $p =v(x)$, and $q =u(x)=E(\tau_x)$. Suppose there is a constant $c\geq 0$ such that \begin{equation} {\cal L} V(x)\geq c{\rm ,\ for\ all\ } x\in B \end{equation} where ${\cal L}$ is the operator given by (39). Then \begin{equation} p\leq[\gamma-V(x)-cq]/[\gamma-\eta]. \end{equation} \end{thm} \noindent{\bf Remark 1.} Note that, if, for example, \begin{equation} V(x) = \frac{\eta + \gamma}{2} \end{equation} then (47) becomes \begin{equation} p\leq\frac{1}{2}-c q/[\gamma-\eta] \end{equation} i.\ e., the term $-c q/[\gamma-\eta]$ gives an estimate of the net bias due to the drift $G(x)$ and diffusion $\Lambda(x)$ terms in (37). \bigskip \noindent Proof of Theorem 3. By Dynkin's formula ([4]) (or by Ito's formula and taking expected values - see [5], for example) applied to the process $V(X(t))$ on the random interval $[0,\tau]$ we have \begin{equation} EV(X(\tau))-V(x)=E\int_0^\tau{\cal L} V(X(s)) ds \end{equation} and then (46) yields \begin{equation} E\int_0^\tau{\cal L} V(X(s)) ds\geq cE(\tau). \end{equation} From (50) and (51) then we have \begin{equation} EV(X(\tau))\geq V(x)+cE(\tau). \end{equation} On the other hand, taking into account (44) and (45) , we get \begin{equation} EV(X(\tau))\leq p \eta+(1-p) \gamma.\end{equation} Inequalities (52) and (53) give \begin{equation} p \eta+(1-p) \gamma\geq V(x)+cE(\tau).\end{equation} or \begin{equation} p\leq[\gamma-V(x)-cq]/[\gamma-\eta]. \end{equation} \hfill$\diamondsuit$\medskip Returning to the example (40): $dX=\sqrt{\epsilon}\ dW, X(0)=x\in B=(0,1)$ \linebreak for a simple application of Theorem 3, we take $$ \eta=0,\gamma=1,\ {\rm and\ } V(x)=x^r, $$ for any number $r$ satisfying $10$ and $\delta_1$ sufficiently large, \begin{eqnarray} {\cal L} V(x)&=&\{[\delta_0(P_I-x_1)+g_P+\frac{1}{2}(\delta_1% (P_I-x_1))^2]/x_1\nonumber\\ &&+ r_2[g_A-\delta_0 x_2+\frac{1}{2}(\delta_1x_2)^2]\\ &&+ r_3[g_Z-\delta_0 x_3+\frac{1}{2}(\delta_1x_3)^2]\nonumber\\ &&+r_4g_F\}V(x)\ \geq \ c\nonumber \end{eqnarray} for some positive constant c and for all $x\in B$, since everything in (63) is bounded, and all of the terms involving $\delta_1$ are positive. We remark finally that it should be noted that both $c$ and $q$ in (47) generally will depend on $\delta_1$. \begin{thebibliography}{00} \bibitem{c1} Y. Cao and T. C. Gard, Uniform persistence and net functions. J. Dynamics and Differential Equations 7(1995), 491-520. \bibitem{c2} Y. Cao and T. C. Gard, Practical persistence for differential delay models of population interactions. Electronic J. Differential Equations (1998) Conference 01, 1997, 41-53. \bibitem{d1} F. Doveri, M. Scheffer, S. Rinaldi, S. Muratori, and Y. Kuznetsov, Seasonality and chaos in a plankton-fish model. Theor. Pop. Biol. 43(1993), 159-183. \bibitem{d2} E. B. Dynkin, Markov Processes I, II. Springer-Verlag, New York, 1965. \bibitem{g1} T. C. Gard, Introduction to Stochastic Differential Equations. Marcel Dekker, New York, 1988. \bibitem{g2} T. C. Gard, Stochastic population models on a bounded domain: exit probabilities and persistence. Stochastics and Stochastics Reports 33(1990), 63-73. \bibitem{g3} T. C. Gard, Exit probabilities for stochastic population models: initial tendencies for extinction, explosion, or permanence. Rocky Mountain J. Math. 20(1990), 1-15. \bibitem{g4} T. C. Gard, Practical persistence in stochastic population models. Fields Institute Communications 21(1999), 205-215. \bibitem{g5} T. C. Gard, Transient effects of stochastic multi-population models. preprint. \bibitem{h1} J. M. Halley and Y. Iwasa, Extinction rate of a population under both demographic and environmental stochasticity. Theor. Pop. Biol. 53(1998), 1-15. \bibitem{h2} V. Hutson and K. Schmitt, Permanence and the dynamics of biological systems. Math. Biosci. 111(1992), 1-71. \bibitem{j1} V. A. A. Jansen and K. Sigmund, Minireview Shaken not stirred: on permanence in ecological communities. Theor. Pop. Biol. 54(1998), 195-201. \bibitem{l1} D. Ludwig, Persistence in dynamical systems under random perturbations. SIAM Rev. 17(1975), 605-640. \bibitem{m1} Y. Miyahara, Ultimate boundedness of the systems governed by stochastic differential equations. Nagoya Math. J. 47(1972), 111-144. \bibitem{r1} S. Renaldi and C. Solidoro, Chaos and peak-to-peak dynamics in a plankton-fish model. Theor. Pop. Biol. 54(1998), 62-77. \bibitem{r2} S. Ruan, A three-trophic-level model of plankton dynamics with nutrient recycling. Can. Appl. Math. Quart. 1(1993), 529-553. \bibitem{r3} S. Ruan, Persistence and coexistence in zooplankton-phytoplankton-nutrient models with instantaneous nutrient recycling. J. Math. Biol. 31(1993), 633-654. \bibitem{s1} H. Smith and P. Waltman, The Theory of the Chemostat: Dynamics of Microbial Competition. Cambridge University Press, Cambridge, 1995. \bibitem{s2} G. Stephanopoulos, R. Aris and A. G. Fredrickson, A stochastic analysis of the growth of competing microbial populations in a continuous biochemical reactor. Math. Biosci. 45(1979), 90-135. \bibitem{w1} P. Waltman, A brief survey of persistence in dynamical systems. in Delay Differential Equations and Dynamical Systems (S. Busenberg and M. Martelli, eds.), Lect. Notes in Math. 1475, Springer-Verlag, Berlin, 1991, 31-40. \end{thebibliography} \noindent{\sc Thomas C. Gard} \\ Department of Mathematics \\ The University of Georgia \\ Athens, GA 30602, USA \\ e-mail: tgard@uga.edu \end{document}