Journal of Applied Mathematics

Volume 2013 (2013), Article ID 625391, 13 pages

http://dx.doi.org/10.1155/2013/625391

Research Article

Coexistence in a One-Predator, Two-Prey System with Indirect Effects

Renato Colucci

Pontificia Universidad Javeriana, Carrera 7 No. 43-82, Bogotá, Colombia

Received 26 October 2012; Accepted 27 February 2013

Academic Editor: Erik Van Vleck

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

Abstract.

We study the dynamics of a one-predator, two-prey system in which the predator has an indirect effect on the preys. We show that, in presence of the indirect effect term, the system admits coexistence of the three populations while, if we disregard it, at least one of the populations goes to extinction.

1. Introduction

The importance of indirect effects is well established in Biology (see [1–5]), for example, in the case of predation (see [6]), the predator can alter the morphology (see [7]) or the behavior of the preys. The preys, in order to reduce the possibility of contacts with the predators, could modify their normal conduct by reducing their activity or by hiding themselves for long time. There are many types of indirect effects (see [3] for a detailed discussion); another interesting case is the refuge indirect effect (see [8], e.g.); anyway, it is of great interest trying to describe the indirect interactions in population dynamics. Recently in [9] a model including indirect effects was proposed, modeling the effects of predator Daphnia over two groups of Phytoplankton of different morphology (see [10]), having Phosphorous as a resource (see [11] or [12]). A general model describing indirect effects of predation can be written in the following form:

where is the density of population of the predator while and are the densities of population of the preys. The positive function represents the indirect effects of the predator on the prey generated by predation over the prey . In [9], a system in which the functions are linear and the function is quadratic was proposed. The system takes the following form:

where , , , , , are positive parameters and where .

In the above system represent the density of population of a predator (Daphnia-Zooplankton) that predates the preys (Phytoplankton) and that are of different size, in particular is of a smaller size than . The variable represents the amount of resources (Phosphorous) for the preys and . We will show below that the system can be put in the form of system (1).

In absence of preys, the population of predators will extinguish exponentially while the presence of them brings a positive growth rate. In absence of predator, the populations of preys will increase depending on the amount of available resources, while the presence of predators affects negatively their growth rate. Moreover, the population of predator has an effect on the prey which is described by the term which produce a negative growth rate. This effect produces an indirect effect that is positive (the term ) for the population of prey . The growth rate of the amount of resources of the system is effected by the presence of the preys while, for its renewable character, the presence of the predator lets it regenerate (the term ).

We first observe that the system is closed; in fact if we sum (2), we obtain

and as a consequence we have the following proposition.

Proposition 1. The function defined as

is a constant of motion.

By using the first integral we can reduce the degree of freedom of the problem; in fact if we fix a value of the the first integral, , then the system can be rewritten in the following way:

Since we are interested only in positive solution and using condition (4), we can limit our analysis to the following region:

First of all we show that the dynamics develops in the region (see Figure 1).

Figure 1: The positive invariant set .

Theorem 2. The set is positively invariant for the solutions of the system (5).

Proof. In order to prove that is positively invariant, it is sufficient to compute the vector field of the system on the boundary of that is made up by 4 triangular regions.

First of all we consider the vertexes of , we set

The vertexes are all fixed points except the point on which the vector field is .

We study the vector field on each triangular region separately.

The Face . In this case , then, in the interior of the triangle, the vector field points toward the interior of . We check the sides of the triangle : if we have and ; moreover, as . If we have that and . If we consider the side on which then we have that . Finally we conclude that on the vector field is tangent or point inward .

The Face . In this region we have that , that is, the plane is invariant and as a consequence we have only to check the vector field on the boundary of . If we have that the vector field points inward since , while in the side where we have that and so it is positively invariant. Finally, on the side where we have that and that as .

The Face . In this region we have that , then the plane is invariant and as a consequence we have only to check the vector field on the boundary of . On the side where we have that and as a consequence , for all . Finally on the side where we have and and as , while on the side where we have that , , and when .

The Face . In this case we have that the vector field points inward . In fact, if we introduce the function that is the sum of the three variables

and we sum the equations of the system, we obtain the following differential equation depending on :

From (9) we get that

this concludes the proof since if the solutions start on then they enter the set or they remain on .

The rest of the paper is organized as follows: in Section 2 we analyze the stability of the fixed points while in Section 3 we study the dynamics on the boundary of . In Section 4 we analyze the dynamics of the system in absence of indirect effects terms (that is ) while in Section 5 we study the asymptotic behavior of the solutions of the system for , in particular we prove the coexistence of the three species under a proper choice of the parameters. In Section 6 we present some numerical simulations while in Section 7 we give some comments and remarks for future investigation.

2. Stability of the Fixed Points

The first step in studying the dynamics of the system (5) consists in finding all the fixed points in and analyzing their stability character. In Section 1 we have pointed out the existence of three fixed points among the vertexes of :

In order to find all possible further fixed points in the faces or in the interior of , we consider the intersections between the nullclines (see [13]). Since it is a standard argument, we do not present the detailed computation. If , there exists a fixed point

in the interior of the triangle , while if , then .

Moreover, there exists a segment of fixed points that we call :

The segment is one of the sides of the face (on the plane ) and its extremes are the fixed points and .

We also have an interior fixed point if the following two conditions are fulfilled:

In details, the coordinate of the fixed point are

where . As a summary of the previous analysis we can state the following theorem.

Theorem 3. The points , , and the points of the segment are equilibria of the system (5) for any (positive) value of the parameters. The system admits that the fixed point if

and the interior fixed point if (14) and (15) hold.

Remark 4. From the above analysis we notice that is a parameter of bifurcation for the fixed points. In particular we have that fixed points and collapse to a unique fixed point if , and the fixed point if .

Then we pass to study the stability character of the fixed points. We consider the functional Jacobian of the vector field:

where .

2.1. The Fixed Point

We start with the analysis of the origin . The matrix of the linearized system at is

The point is a saddle for any values of the parameters since we have one negative eigenvalue and two positive eigenvalues and . The -axis is the stable space of the linearized system at while the instable space is the plane .

It is possible to distinguish two cases: three distinct eigenvalues if (that is ) and only two distinct eigenvalues if (that is ).

2.2. The Fixed Point

We analyze the point , and the matrix of the linearized system at is

The eigenvalues of are

The eigenvalues are always negative, while if we have

then the fixed point is instable. In the case in which , that is, , we would need to study the system on the center manifold (see [14]) of the point in order to study its stability.

2.3. The Fixed Point

The matrix of the linearized system at is

The eigenvalues of the matrix are

In this case we have at least one negative eigenvalue, that is, , along the axes , then the stability depends on the first two eigenvalues. In the case in which

we can conclude that is instable. If we have to study the system on the center manifold related to the eigenvalue (and if it is 0), in order to establish the stability of .

2.4. The Fixed Point

We analyze the point in the case in which that is if . The matrix of the linearized system at is

where we have set . The eigenvalues of the matrix are

where . Then if

then are complex conjugate with negative real part; otherwise they are both negative (coincident if the inequality holds in (28)). If the eigenvalue is positive, that is,

we have that is instable.

2.5. The Segment of Fixed Point

On the fixed points that belong to the segment we have that and , then the matrix of the linearized system at any point of is of the following form:

where . The eigenvalues of the matrix are

In this case we have that the eigenvalue is always negative since . In the case in which

the eigenvalue and then the fixed points of the segments are all instable. We note that if , then satisfies

on the contrary satisfies the opposite inequality.

2.6. The Fixed Point

The matrix of the linearized system at is

where we have set . We consider the trace of the matrix :

Then, if we have that

and then at least one of the eigenvalues of is negative.

We summarize the results of the previous analysis in the following statements.

Theorem 5. Suppose that the parameters of the system satisfy the following inequalities:

Then the fixed points , , , and the fixed points that belong to the segment are all instable.

Remark 6. We note that the instability hypothesis for the point , that is, , implies the presence of the point . Then it would be interesting to prove if it is a necessary condition for the instability of , that is, if the presence of the fixed point is necessary for the instability of all fixed points in .

Remark 7. We note that the hypothesis of the previous theorem are not compatible with the inequality

that is, the condition for which .

In conclusion, the instability of the boundary fixed points implies the existence of the interior fixed point . The presence of the interior fixed point represents an obvious case of coexistence of the three species; however, if is instable then it is not useful to describe real cases of coexistence.

Since we are not able to face the general problem of stability of point , we will follow another strategy in order to prove coexistence of the system (see Section 6 below).

Remark 8. It would be interesting to find if the interior fixed point is instable for some choice of the parameters that are compatible with the instability hypothesis for the boundary fixed points. In that case it would be interesting to look for the existence of a limit cycle surrounding or of homoclinic cycle and chaotic attractors (see e.g., [15] or [16] in which the authors proved the chaotic behavior of one-predator, two prey systems without indirect effects).

3. The Dynamics on the Boundary of

In this section we analyze the dynamics on the boundary of , in order to do that we first study the dynamics on the coordinate axes (see Figure 2) and then on the faces of .

Figure 2: The dynamic on the axes.

3.1. Dynamics on the Axes

The three axes are invariant for the dynamics; in particular, on the -axes we have

Then the solution is and we have that as . On the axes and we have logistic growth and in particular as .

In fact on the -axes we have

while on the -axes we have

For both axes we have that as and the solutions are of the form

where and are the initial conditions.

3.2. The Dynamics on the Faces

In this subsection we analyze the dynamics on the faces of .

The Face is on the invariant plane . We have a system of competitive Volterra-Lotka equation:

This is a degenerate case; in fact we have a segment (see definition (13)) of fixed points that connects and . The -limit of every orbit starting inside (that is ) is a point on ; hence, the preys do not vanish and there are not periodic orbits. The behavior of the system on the face is represented in Figure 3 above.

Figure 3: The dynamics on the face .

The Face is on the invariant plane . The system reduces to a system of Volterra-Lotka equations with intraspecific competitions between the preys:

We recall that on we have the fixed point if . In order to understand the behavior of the solution we analyze the nullclines. The -nullclines are and and the -nullclines are and . We distinguish two cases: with or without the fixed point .

In the case in which or (that is ), the -nullcline is above (or pass through) the fixed point , as a consequence inside . The -nullcline divide in two regions: below it we have while above it we have . Then, all the solutions with tend to the fixed point according to the saddle character of . The behavior of the system is represented in Figure 4.

Figure 4: The dynamic on the face without the point .

In the case in which (that is, ) the nullclines intersect at the point . The direction of the vector field is represented in Figure 5 below. Hence, the solutions wind in the clockwise direction around . To be more specific we consider the linearization of (45) at . The eigenvalues of the matrix of the linearized system are

Then the eigenvalues are both negative and in the complex case they have negative real part. We conclude that is asymptotically stable.

Figure 5: The dynamic on the face in the presence of the point .

In both cases there are no periodic orbit; in fact two dimensional Volterra-Lotka equations do not admit isolated periodic orbit (see [17]), and they admit a continuum of periodic orbits if and only if the eigenvalues at the interior fixed point (in this case ) are purely imaginary (if and only if the trace of the matrix of the linearized system is zero, that is, ).

The Face is on the plane , and we have already shown (see equation (9)) that if the solutions start on it then they enter while if they remain on .

The Face is easy to analyze, if the solutions do not start on the axes (that is if ) we have that , as a consequence the solutions leave the plane and go inside .

Then we have proved the following.

Theorem 9. The system (5) admits neither limit cycles nor periodic orbits in the boundary of .

4. Analysis for

In this section we study the behavior of the system in absence of indirect effects that is in the case .

In this case we have the presence of a further boundary fixed points:

in the case in which . If , then , while if we have that . Then we have the following.

Theorem 10. The points , , and the points of the segment are equilibria of the system (5) with , for any value of the other parameters. The system admits a further fixed point if

and a fixed point if

The dynamics on the boundary is similar to that of the case . In particular on the faces and there are no differences between the two cases and .

To analyze the dynamics on the face we can use the same method of the analysis of the case and it is easy to see that we obtain the same conclusions.

The remarkable difference is that, in the case , the plane is invariant. The dynamics on it (the faces ) is similar to that on in the case in which (see Figures 6 and 7).

Figure 6: The dynamic on the face in the case in which and in which .

Figure 7: The dynamic on the face in the case in which and in which .

We can distinguish two cases: with or without the fixed point . The matrix of the linearized system at is the following:

This matrix has at least one real eigenvalue, that is, , while the other two eigenvalues cannot be purely imaginary since the trace, of the matrix formed by the first two columns and rows of , is negative. Then, as in the case , we can conclude the following.

Theorem 11. If the system (5) admits neither limit cycles nor periodic orbits on the boundary of the set .

In the case the asymptotic behavior of the solution can be studied by using a well-known result (see [17, Theorem , page 43.]) that states that if in the interior of the forward invariant set there are not fixed points then the interior of does not contain -limit or -limit sets of the solutions. In particular if there are not fixed points in the interior of then there are not periodic orbits too (see [18], e.g., of system without periodic orbit). Then we can state the following.

Theorem 12. If , the system (5) admits neither limit cycles nor periodic orbit in the set ; that is, for each solution starting in at least one of the species goes to extinction.

Remark 13. It is interesting to note that three-dimensional predator-prey systems (without indirect effects) admit nontrivial cases of coexistence; see, for example, the paper [19] where the authors studied oscillations for a class of singularly perturbed three-dimensional predator-prey systems.

5. Asymptotic Behavior for

In Section 4 we have shown that the asymptotic dynamics on does not admit coexistence of populations if . In this section we analyze the case and show that the system admits coexistence of the three species by using persistence theory and in particular an acyclicity approach (see [20]).

First of all we give a remark on the evolution of volume elements and then we pass to the proof that the system admits coexistence of the three populations under a proper choice of the parameters.

5.1. Volume Evolution

We consider the evolution of volumes under the flow of the system (5). Let be a region of with regular boundary and define , that is, the image of under the flow at time . Let be the volume of , then by Liouville's theorem we have that the evolution of is described by the following differential equation:

where is the vector field of the system (5). In details we have that

The above expression is quite difficult to study and suggests a complex behavior. However, there exist values of the parameter for which is negative. In fact at least in the case in which

we have that so the 3-dimensional volume elements contract in .

In this case, if an attractor (or limit cycle) exists, we expect that it would be contained in a small subset of . On the contrary the attractor could be contained in a bigger subset of (see the numerical simulations below).

Remark 14. From condition (53)-(55), we have that the volume element contracts in the whole region ; however, these conditions are incompatible with condition (37). Since (37) could be not necessary conditions, it would be possible to have both volume decay and instability of the fixed points.

5.2. Coexistence of the Three Populations

The biological problem of the coexistence of the three species can be put in mathematical terms by looking for the conditions which prevent the positive solutions starting in the interior of from converging to as . These ideas can be made rigorous in the context of persistence theory (see [20]). There are many definitions of persistence; the most useful for biological applications (see [21]) is uniform persistence.

Definition 15. The system (5) is uniformly -persistent if there exists such that

with

and where

In order to prove uniform persistence we use an acyclicity approach (see [20]), then (see Theorem 8.17 page 188, hypothesis (H) page 185 and Theorem 5.2 page 126 in [20]) the following conditions are needed.

(C1) There exists a compact attractor of bounded set.(C2) The invariant sets of are weakly -repelling.(C3) The invariant sets of are acyclic.

In the following lines we give the rigorous definitions of the above ideas and prove conditions (C1)–(C3).

Theorem 16. The system (5) admits a unique compact attractor that attracts all bounded sets of .

Proof. By Theorem 2.33 page 43 in [20], we get the existence of a compact attractor of bounded sets if the flow associated with the system (5) is point dissipative, asymptotically smooth, and eventually bounded on every bounded set in .

In order to prove the thesis we consider Theorem 2; in particular we write again the system (5) in function of , that is,

If with and , then

and we conclude that there exists such that for . This proves the point dissipative property of the flow.

We consider the following family of sets:

where and . Every bounded set of is contained in a set for sufficiently large. The family of sets is forward invariant (see Theorem 2 for the set ), in fact, again we have

This proves that bounded sets have bounded orbits. Finally we get the asymptotically smoothness of the flow by noting that it is asymptotically compact on every forward invariant bounded closed set (see remark 2.26, page 39 in [20]).

Before proving condition (C2) we consider the following definition (see [20, chapter 8]).

Definition 17. Let be the solution of (5) starting at . A set in is called weakly -repelling if there is no such that and as .

In order to prove weakly -repelling property it is sufficient to show that the stable manifolds of the invariant sets of are contained in . From the analysis of the stability of fixed points and of the dynamics on we have that there are no periodic orbits on and the only invariant set are , , , , and .

Theorem 18. Suppose that hypothesis of Theorem 5 holds. Then the stable manifolds of all the fixed points of are contained in .

Proof. Following the results of Section 2 we have that the stable manifold of the fixed point is the axes: , while (see also Section 3) the stable manifold of the fixed point satisfies .

The points belonging to the segment have a stable manifold inside the face on the invariant plane . In particular the extremes and of the segment satisfy and . Finally all the points belonging to the segment have the same center manifold

and this concludes the proof.

In order to prove condition (C3) we first consider the following definitions (see [20, chapter 8]).

Definition 19. Let . is chained to in , written , if there exists a total trajectory in with , and as .

Definition 20. A finite collection of subsets of is called cyclic if after possibly renumbering or in for some . Otherwise it is called acyclic.

Theorem 21. Suppose that hypothesis of Theorem 5 holds. Then there set , of invariant sets of , is acyclic.

Proof. By the analysis of stability of the fixed points and of the dynamics on the boundary of it follows that there are no periodic orbits in . Moreover, the only invariant sets , , , , and are all contained in . By the stability analysis of fixed points we get that only the following chains are admissible:

Then the set is acyclic.

From Theorems 16, 18, and 21, we get the final result.

Theorem 22. Suppose that hypothesis of Theorem 5 holds. Then the system (5) is uniformly persistent.

From the previous theorem we obtain coexistence of the three populations under the hypothesis of Theorem 5 without analyzing the stability of the fixed point .

6. Numerical Simulations

The strategy used for proving coexistence of the three species does not exclude quasi-periodic or chaotic behavior of the system. The uniform persistence only lets us conclude the existence of a bounded attractor contained in and with positive distance from . In this section we provide several numerical experiments in order to describe the structure of the attractor.

6.1. Experiment 1

We consider the following values of the parameters:

In a first numerical experiment (represented in Figures 8, 12, 13, and 14 below) we considered the value for which the attractor is the fixed point . In a second numerical experiment (represented in Figures 9, 15, 16, and 17) we considered the value . In this case the attractor appears to be bidimensional.

Figure 8: The solution of the system with the choice of parameters as in (65) and .

Figure 9: The solution of the system with the choice of parameters as in (65) and .

6.2. Experiment 2

We consider the following values of the parameters:

Again, we perform two numerical experiments with different values of . In a first simulation we considered the value (the solution is represented in Figure 10) and the attractor is the fixed point . In a second simulation we considered the value (the solution is represented in Figure 11) and the attractor appears to be a limit cycle.

Figure 10: The solution of the system with the choice of parameters as in (66) and .

Figure 11: The solution of the system with the choice of parameters as in (66) and .

Figure 12: The function in the case of (65) and .

Figure 13: The function in the case of (65) and .

Figure 14: The function in the case of (65) and .

Figure 15: The function in the case of (65) and .

Figure 16: The function in the case of (65) and .

Figure 17: The function in the case of (65) and .

6.3. Remarks on the Numerical Experiments

These experiments suggested that the parameter could be considered as a bifurcation parameter for the attractor. Moreover, we have effectuated many simulations with different choices of the function with different values of parameters and in any cases a limit cycle arises. This suggests that the attractor may have some structural stability properties.

7. Conclusions

We have shown that in absence (that is ) of the terms that describe indirect effects, the system (5) does not admit coexistence of the three populations. The present work suggests the importance of indirect effects in describing cases of coexistence (with or without an interior fixed point); in particular the parameter is a bifurcation parameter for coexistence. It has been already pointed out (see for instance [22]) that a presence of limit cycle, oscillations, or chaotic fluctuations allows the coexistence of many species and could be beneficial for the functioning of the ecosystem (see [23]).

It would be interesting to consider indirect effects in different contexts such as spatial model (see [24]) or evolution problems (see [25]). In fact it is an important ecological problem (see [26]) to determine the consequences of indirect effects in evolution.

In the context of population dynamics, it would be interesting to consider a different form of the function ; in particular if we sum the second and the third equations we get

then the indirect effect terms have no effects on the whole population of preys; they effect only the proportion between them. In order to consider a different situation it would be sufficient to put different constants and in the function for the two preys (see Figure 18).

Figure 18: The solution with (in the range of the admissible values) and (out of the range).

Moreover, it would be interesting to consider the choice , or at least . This would be an improvement in description of indirect effects; in fact in the present work we have considered that the indirect effect of the predator over the prey has the same order of magnitude of the direct effect (compare the terms and ). In the case considered, if , the form of the third equation is the most effected while the second equation can be rewritten in the following way:

that is simply a change of the coefficient of . By considering a term

also the form of the second equation would be drastically effected. Moreover, it is an interesting question (see [3]) to establish if direct and indirect effects are of the same magnitude. Then it would be useful to consider different values of in (69) in order to describe cases in which direct effects are bigger (resp., less) than indirect effects.

Another possible improvement would be to consider the nonautonomous case for which is a function of time; for example, in the case of seasonal indirect effects, we could make the following choice:

where is a positive constant (see Figures 19 and 20). Some analysis of that data (see [3]) suggests that indirect effects could take a long time to become apparent since they occur at a slower time scale than any direct effect (see [4]). For this reason it could be interesting to consider a delayed version (see [27], e.g., of existence of periodic solutions for a predator-prey model with delay) of the system (5) where, for example, in the third equation we put

Moreover, for more general forms of the function , it could be interesting to look for chaotic behavior (see [28]) of the solutions (see [29], for example, in population dynamics and [30] for a discrete time case). In conclusion we can note that, beside the system (2) which is a simplification of the description of indirect effects, from the numerical experiments, it seems that it shares the main features with more sophisticated models.

Figure 19: The solution with replaced by and .

Figure 20: The solution with replaced by and .

Acknowledgments

The author would like to thank Professor Piero Negrini (University Sapienza of Rome) and Professor Juan José Nieto (University of Santiago de Compostela) for their helpful suggestions and remarks. This work was partially supported by project PPTA 4476 of Pontificia Universidad Javeriana, Bogotá.

References

B. Bolker, M. Holyoak, V. Křivan, L. Rowe, and O. Schmitz, “Connecting theoretical and empirical studies of trait-mediated interactions,” Ecology, vol. 84, no. 5, pp. 1101–1114, 2003.

D. Cariveau, R. E. Irwin, A. K. Brody, L. S. Garcia-Mayeya, and A. Von der Ohe, “Direct and indirect effects of pollinators and seed predators to selection on plant and floral traits,” Oikos, vol. 104, no. 1, pp. 15–26, 2004.

B. A. Menge, “Indirect effects in marine rocky intertidal interaction webs: patterns and importance,” Ecological Monographs, vol. 65, no. 1, pp. 21–74, 1995.

W. E. Snyder and A. R. Ives, “Generalist predators disrupt biological control by a specialist parasitoid,” Ecology, vol. 82, no. 3, pp. 705–716, 2001.

J. T. Wootton, “Indirect effects, prey susceptibility, and habitat selection: impacts of birds on limpets and algae,” Ecology, vol. 73, no. 3, pp. 981–991, 1992.

J. Estes, K. Crooks, and R. Holt, Ecological Role of Predators. Enciclopedia of Biodiversity, vol. 4, Academic Press, New York, NY, USA, 2001.

V. Lundgren and E. Granéli, “Grazer-induced defense in Phaeocystis globosa (prymnesiophyceae): Influence of different nutrient conditions,” Limnology and Oceanography, vol. 55, no. 5, pp. 1965–1976, 2010.

M. F. Caruselaa, F. R. Momoa, and L. Romanellia, “Competition, predation and coexistence in a three trophic system,” Ecological Modelling, vol. 220, no. 10, pp. 2349–2352, 2009.

“Indirect effects affects ecosystem dynamics,” 2011, http://www.ictp-saifr.org/

R. Margalef, “Life forms of Phytoplanktos as survival alternative in an unstable environment,” Oceanologica Acta, vol. 134, no. 1, pp. 493–509, 1978.

O. Sarnelle, “Daphnia as keystone predators: effects on phytoplankton diversity and grazing resistance,” Journal of Plankton Research, vol. 27, no. 12, pp. 1229–1238, 2005.

D. O. Hessen, T. Andersen, P. Brettum, and B. A. Faafeng, “Phytoplankton contribution to sestonic mass and elemental ratios in lakes: implications for zooplankton nutrition,” Limnology and Oceanography, vol. 48, no. 3, pp. 1289–1296, 2003.

L. R. Devaney, M. W. Hirsh, and S. Smale, Differential Equations, Dynamical Systems, and an Introduction to Chaos, Elvesier, New York, NY, USA, 2004.

S. Wiggins, Introduction to Applied Nonlinear Dynamical Systems and Chaos, Springer, New York, NY, USA, 2nd edition, 2010.

A. Klebanoff and A. Hastings, “Chaos in one-predator, two-prey models: general results from bifurcation theory,” Mathematical Biosciences, vol. 122, no. 2, pp. 221–233, 1994.

S. Gakkhar and R. K. Naji, “Existence of chaos in two-prey, one-predator system,” Chaos, Solitons and Fractals, vol. 17, no. 4, pp. 639–649, 2003.

J. Hofbauer and K. Sigmund, Evolutionary Games and Population Dynamics, Cambridge University Press, Cambridge, UK, 1998.

P. van den Driessche and M. L. Zeeman, “Three-dimensional competitive Lotka-Volterra systems with no periodic orbits,” SIAM Journal on Applied Mathematics, vol. 58, no. 1, pp. 227–234, 1998.

W. Liu, D. Xiao, and Y. Yi, “Relaxation oscillations in a class of predator-prey systems,” Journal of Differential Equations, vol. 188, no. 1, pp. 306–331, 2003.

H. L. Smith and H. R. Thieme, Dynamical Systems and Population Persistence, American Mathematical Society, Providence, RI, USA, 2011.

G. Butler, H. I. Freedman, and P. Waltman, “Uniformly persistent systems,” Proceedings of the American Mathematical Society, vol. 96, no. 3, pp. 425–430, 1986.

J. Huisman and F. J. Weissing, “Biodiversity of plankton by species oscillations and chaos,” Nature, vol. 402, no. 6760, pp. 407–410, 1999.

D. Li, J. Li, and Z. Zheng, “Measuring nonequilibrium stability and resilience in an n-competitor system,” Nonlinear Analysis. Real World Applications, vol. 11, no. 3, pp. 2016–2022, 2010.

H. Nie and J. Wu, “Coexistence of an unstirred chemostat model with Beddington-DeAngelis functional response and inhibitor,” Nonlinear Analysis. Real World Applications, vol. 11, no. 5, pp. 3639–3652, 2010.

R. A. de Assis, S. Bonelli, M. Witek et al., “A model for the evolution of parasite-host interactions based on the Maculinea-Myrmica system: numerical simulations and multiple host behavior,” Nonlinear Analysis. Real World Applications, vol. 13, no. 4, pp. 1507–1524, 2012.

M. R. Walsh and D. N. Reznick, “Interactions between the direct and indirect effects of predators determine life history evolution in a killifish,” Pnas, vol. 105, no. 2, pp. 594–599, 2008.

X. Liu and L. Huang, “Periodic solutions for impulsive semi-ratio-dependent predator-prey systems,” Nonlinear Analysis. Real World Applications, vol. 10, no. 5, pp. 3266–3274, 2009.

R. L. Devaney, An Introduction to Chaotic Dynamical Systems, Addison-Wesley, Reading, Mass, USA, 2nd edition, 1989.

J.-M. Ginoux, B. Rossetto, and J.-L. Jamet, “Chaos in a three-dimensional Volterra-Gause model of predator-prey type,” International Journal of Bifurcation and Chaos in Applied Sciences and Engineering, vol. 15, no. 5, pp. 1689–1708, 2005.

Z. He and X. Lai, “Bifurcation and chaotic behavior of a discrete-time predator-prey system,” Nonlinear Analysis. Real World Applications, vol. 12, no. 1, pp. 403–417, 2011.