Disease in prey population and body size of intermediate predator reduce the prevalence of chaos-conclusion drawn from Hastings–Powell model

Ecological Complexity 6 (2009) 363–374

Contents lists available at ScienceDirect

Ecological Complexity journal homepage: www.elsevier.com/locate/ecocom

Disease in prey population and body size of intermediate predator reduce the prevalence of chaos-conclusion drawn from Hastings–Powell model Krishna pada Das, Samrat Chatterjee, J. Chattopadhyay * Agricultural and Ecological Research Unit, Indian Statistical Institute, 203, B.T. Road, Kolkata 700108, West Bengal, India

A R T I C L E I N F O

A B S T R A C T

Article history: Received 7 February 2008 Received in revised form 27 January 2009 Accepted 3 March 2009 Available online 18 April 2009

In ecology the disease in the prey population plays an important role in controlling the dynamical behaviour of the system. We modify Hastings and Powell’s (HP) [Hastings, A., Powell, T., 1991. Chaos in three-species food chain. Ecology 72 (3), 896–903] model by introducing disease in the prey population. The conditions for which the modified HP model system represents extinction, permanence or impermanence of population are worked out. The modified model is analyzed to obtain different conditions for which the system exhibits stability around the biologically feasible equilibria. Through numerical simulations we display that the modified system enters into stable solutions depending upon the force of infection in prey population as well as body size of intermediate predator. Our results demonstrate that disease in prey population and body size of intermediate predator are the key parameters for controlling the chaotic dynamics observed in original HP model. ß 2009 Elsevier B.V. All rights reserved.

Keywords: Disease in prey Body size Intermediate predator Chaos Stable Permanence Impermanence

1. Introduction Disease in ecological system plays an important role and cannot be ignored. This leads to a new field of study known as ecoepidemiology. Though the term ‘eco-epidemiology’ was coined by Chattopadhyay and Arino (1999) but Anderson and May (1986) were probably the first who considered the disease factor in predator–prey dynamics. Quite a good number of articles have appeared after then, for example see Arino et al. (2004), Freedman (1990), Hadeler and Freedman (1989), Hethcote et al. (2002) and Venturino (2002). The infected prey leads to an array of interesting dynamics in ecological interaction which was studied by many authors (Dobson, 1988; Hethcote et al., 2004; Xiao and Chen, 2001). But as far as our knowledge goes a little effort has been made to understand such dynamics in eco-epidemiological systems. Recently Chatterjee et al. (2006) observed that rate of infection and the predation rate are two prime factors that govern the chaotic dynamic in eco-epidemiological system. Though the study of chaos in an eco-epidemiological system is new but the literature of chaos in ecological system is very rich.

* Corresponding author. Tel.: +91 33 25753231; fax: +91 33 25773049. E-mail addresses: [email protected] (K.p. Das), [email protected] (S. Chatterjee), [email protected] (J. Chattopadhyay). 1476-945X/$ – see front matter ß 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.ecocom.2009.03.003

Chaotic dynamics is common in tri-trophic food chain and is of common interest to both modelers and experimental ecologists. Hastings and Powell (HP) (Hastings and Powell, 1991) produced a chaotic population system in a simple tri-trophic food chain model with type-II functional response. It is to be noted that natural system seem to have no difficulty to switch from one state to another (chaos to order and order to chaos). After the work of Hastings and Powell (1991), many researchers explored their model by including various ecological factors. Rai and Sreenivasan (1993) reformulated HP model for aquacultural research. Ruxton (1994, 1996) showed that the system of linked populations or imposition of a population floor on HP model has a stabilizing effect on tri-trophic food chain. Eisenberg and Maszle (1995) revisited the HP model and observed that gradual addition of refugia provide a stabilizing influence for which the chaotic dynamics collapsed to stable limit cycles. McCann and Hastings (1997) stabilized the food web by eliminating chaotic dynamics and limit cycles behaviour of the HP model by introducing omnivore on top-predator and also by increasing the strength of omnivore. Varriale and Gomes (1998) analyzed HP model in two different approaches. Firstly they observed the asymptotic states of the system resulting from numerical integration of the equations. Secondly they applied the embedding procedure to extract the relevant dynamical exponents from a time series for only one scalar variable. By considering intraspecific density dependence (IDD), an important ecological factor, in food chain model (deterministic and

364

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

stochastic), Xu and Li (2002) showed that the deterministic tritrophic food chain model stabilized the food chain and the very famous teacup chaotic attractor disappear. By including the effect of the toxin, produced by Toxin producing phytoplankton, as extra mortality in zooplankton population, Chattopadhyay and Sarkar (2003) observed that chaotic behaviour is less likely to occur in a real food chain dynamics. Lonngren et al. (2004) also visited the HP model with different set of parameters by using techniques from modern control theory. They observed that the synchronization of two diverse food chain systems are possible. Maionchi and Reis (2005) proposed a spatial version of the HP model in which predators seeked their preys only in a finite neighborhood of their home location. They showed that local predation leads to spontaneous pattern formation that leads to qualitative changes in the global dynamics of the system. They also observed that the chaotic strange attractor present in the original HP model could be replaced by a stable fixed point or by an attractor of simpler topology depending on the size of the predator neighborhood. The above reviews clearly indicate that a lot of research has already been performed on HP model. As far as our knowledge goes none of the studies take into account the disease factor in the system. The main objective of the present paper is to observe the dynamics of original HP model in the presence of disease in the prey population. We like to see how the disease in prey population influence the dynamics of tri-trophic food chain model proposed by Hastings and Powell (1991). Another important aspect of the present work is to study the relationship between the body size of the intermediate predator population and the dynamical nature of the system. This study becomes a prominent focus in present ecology and evolutionary biology. Many authors showed all physiological and ecological aspects of an organism to bring the importance of body size as a biological variable forcibly to the attention of biologists (La Barbera, 1986; McMahon, 1973). System dynamics are regulated mainly by the physiological rate parameters of all living state variables and body size contributes significantly to control these parameters (Hrbacek, 1962; Tag, 1995). Blueweiss et al. (1978) and Lehman (1976) reported that the half saturation constant of zooplankton grazing on phytoplankton increases logarithmically with the body size. In the present work we vary the body size of the intermediate predator to find the critical value of body size for which the system do not enter into the chaotic region. 2. Model formulation The HP model (Hastings and Powell, 1991) with pairwise interactions between three species, namely, X; Y; Z, which incorporates a Holling type-II functional response in both consumer species, namely Y and Z is as follows:   dX X YX ¼ RX 1   C 1 A1 dT K B1 þ X dY dT

¼

A1

YX YZ A2 D1 Y B1 þ X B2 þ Y

dZ dT

¼

C 2 A2

(1)

YZ D2 Z: B2 þ Y

Here X is the numbers of the species at the lower level of the food chain, Y the numbers of the species that preys upon X, and Z is the numbers of the species that preys upon Y. Here T is time. The constant R is the intrinsic growth rate and the constant K is the carrying capacity of the species X. The constants C11 and C 2 are conversion rates of prey to predator for the species Y and Z respectively; D1 and D2 are constant death rates for species Y and Z respectively. The constants Ai and Bi for i ¼ 1; 2 are maximal

predation rate and half saturation constants for Y and Z respectively. Hastings and Powell (1991) demonstrated that the dynamic interaction between prey and predators in simple threespecies food chain is chaotic in a certain region of parametric space. We modify the HP model (1) by introducing disease in prey population. We make the following assumptions: 1. The prey population is divided into two classes, viz. (a) susceptible class whose population density is denoted by S and (b) infected class whose population density is denoted by I. The intermediate predator whose population density is denoted by P 1 and the density of top-predator is denoted by P 2 . 2. A part of the susceptible prey population becomes infected at a rate a, following the law of mass action. 3. Infected population is not in a state of reproduction and also does not compete for the resource. 4. Behaviour of the entire community is assumed to arise from the coupling of these interacting species, where P 1 preys on both susceptible prey and infected prey in the form of Holling type-II and type-I (mass action law) respectively. This combination of mixed selection of functional form is taken because the capturing of infected prey is easier than the susceptible prey. Top-predator preys intermediate predator in the form of Holling type-II. 5. The infected prey population dies at the rate D1 and the intermediate predator and top-predator die at the rate D2 and D3 respectively. With the above assumptions, the HP model takes the following form:   dS S P1 S ¼ RS 1   aIS  C 1 A1 dT K B1 þ S dI ¼ aSI  A2 P 1 I  D1 I dT (2) dP1 P1 S P1 P2 ¼ A1 þC 2 A2 IP 1  A3 D2 P1 B1 þ S dT B2 þ P 1 dP2 P1 P2 ¼ C 3 A3 D3 P 2 dT B2 þ P 1 System (2) has to be analyzed with the following initial conditions: Sð0Þ > 0;

Ið0Þ > 0;

P 1 ð0Þ > 0;

P 2 ð0Þ > 0

Here S; I; P 1 ; P 2 are respectively the susceptible prey population, infected prey population, the intermediate predator population and the top-predator population respectively. A1 and A2 are the maximal predation rate of intermediate predator for susceptible and infected prey respectively; A3 is the maximal predation rate of top-predator for intermediate predator; B1 and B2 are the half saturation constant for functional response of intermediate and the top-predator respectively; C 1 1 is the conversion rate of susceptible prey to intermediate predator; C 2 is the conversion rate of infected prey to intermediate predator; C 3 is the conversion rate of intermediate predator to top-predator. We nondimensionalize the system (2) with S I P1 P2 i¼ ; p1 ¼ s¼ ; ; p2 ¼ ; t ¼ RT; K K K K then the system (2) takes the form ds dt di dt d p1 dt d p2 dt

p1 s  F 1 ðs; i; p1 ; p2 Þ 1 þ cs

¼

sð1  sÞ  asi  b

¼

asi  d p1 i  ei  F 2 ðs; i; p1 ; p2 Þ

¼ ¼

p1 s p1 p2 þgi p1  h  j p1  F 3 ðs; i; p1 ; p2 Þ 1 þ cs 1 þ m p1 p1 p2 k l p2  F 4 ðs; i; p1 ; p2 Þ 1 þ m p1 f

(3)

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

System (3) has to be analyzed with the following initial conditions: sð0Þ > 0;

ið0Þ > 0;

p1 ð0Þ > 0;

g¼

aK R

;

b¼

C 2 A2 K ; R

C 1 A1 K ; RB1

h¼

A3 K ; B2 R

c¼

K ; B1

m¼

3.3. Extinction criterion of population

p2 ð0Þ > 0 Lemma 1. If 1  aiðtÞ, then lim t ! 1 sðtÞ ¼ 0. If asðtÞ  e, then lim t ! 1 iðtÞ ¼ 0. If f  c j and giðtÞ  j, then lim t ! 1 p1 ðtÞ ¼ 0. If k  ml, then lim t ! 1 p2 ðtÞ ¼ 0.

where a¼

d¼

K ; B2

A2 K ; R

e¼

D1 ; r

D2 ; R

k¼

C 3 A3 K ; B2 R

j¼

f ¼

A1 K B1 R l¼

D3 R

Proof. ds p s ¼ sð1  sÞ  asi  b 1  sð1  aiÞ: dt 1 þ cs

3. Some preliminary results ; sðtÞ  sðt 0 Þexp

3.1. Positive invariance

Z

t

 ð1  aiðrÞÞ dr:

t0 T

4

Setting X ¼ ðs; i; p1 ; p2 Þ 2 R , we may write (3) in a vector form as T

FðXÞ ¼ ½F 1 ðXÞ; F 2 ðXÞ; F 3 ðXÞ; F 4 ðXÞ ; 4

365

1

(4)

Hence lim t ! 1 sðtÞ ¼ 0 provided 1  aiðtÞ. di ¼ asi  d p1 i  ei: dt

4

where F : C þ ! R and F 2 C ðR Þ. Then Eq. (3) becomes X˙ ¼ FðXÞ

; iðtÞ ¼ iðt 0 Þ exp  iðt 0 Þ exp

with Xð0Þ ¼ X 0 2 R4þ . It is easy choosing Xð0Þ 2 R4þ such that X i

Z

t

 ðasðrÞ  d p1 ðrÞ  eÞ dr  ðasðrÞ  eÞ dr:

t0

Thus, lim t ! 1 iðtÞ ¼ 0 provided asðtÞ  e. d p1 p1 s p1 p2 þ gi p1  h ¼ f  j p1 : 1 þ cs dt 1 þ m p1

3.2. Boundedness Theorem 1. All the solutions of the system (3) which start in R4þ are uniformly bounded.

t

t0

(5)

to check in Eq. (4) that whenever ¼ 0 then F i ðXÞjX i ¼0  0, (i=1, 2, 3, 4). Now any solution of Eq. (5) with X 0 2 R4þ , say XðtÞ ¼ Xðt; X 0 Þ, is such that XðtÞ 2 R4þ for all t > 0 (Nagumo, 1942).

Z

  sðrÞ p2 ðrÞ þgiðrÞh  j dr f t0 1 þ m p 1 ðrÞ    1 þ csðrÞ Rt sðrÞ þgiðrÞ  j dr  p1 ðt 0 Þ exp t0 ð f  c jÞ 1 þ csðrÞ

; p1 ðtÞ ¼ p1 ðt 0 Þ exp

  Rt

Proof. Let

Thus if f  c j and giðtÞ  j, then lim t ! 1 p1 ðtÞ ¼ 0.

W ¼ s þ i þ p1 þ p2 :

d p2 p1 p2 ¼k  l p2 : dt 1 þ m p1

Using Eq. (3), we have dW p s p1 p2 ¼ sð1  sÞ  ðb  f Þ 1  ðd  gÞ p1 i  ðh  kÞ  ei dt 1 þ cs 1 þ m p1

; p2 ðtÞ



 j p1  l p2 : Since b  f , d  g and h  k we get following expression: dW  sð1  sÞ  ei  j p1  l p2 : dt ;

¼

dW þ LW  ðs  1Þ2 þ 1; dt

  p1 ðrÞ l dr 1 þ m p1 ðrÞ p2 ð0Þexp ðlðt  t 0 ÞÞ; p2 ð0Þexp



Rt

t0



k

by the given condition. Thus lim t ! 1 p2 ðtÞ ¼ 0 provided k  ml.

where L ¼ min f1; e; j; lg:

&

4. Existence and local stability of equilibrium points The system has seven equilibrium points. E0 ð0; 0; 0; 0Þ is the trivial equilibrium point. The axial equilibrium point is E1 ð1; 0; 0; 0Þ. Disease free planar equilibrium point is E2 ðu 1 ; 0; u 2 ; 0Þ, where

Hence, dW þ LW  1: dt Applying the theory of differential inequality we obtain, 1  exp ðLtÞ 0 < Wðs; i; p1 ; p2 Þ < L

u1 ¼

j ; f cj

u2 ¼

f f f  ðc þ 1Þ jg 2

bð f  c jÞ

:

The existence conditions of disease free planar equilibrium point E2 are

þ Wðsð0Þ; ið0Þ; p1 ð0Þ; p2 ð0ÞÞexp ðLtÞ: For t ! 1, we have 0 < W < 1=L.

f  cj>0 R4þ

and

f  ðc þ 1Þ j > 0:

are confined in the Hence all solutions of (3) that initiate in region   1 B ¼ ðs; i; p1 ; p2 Þ 2 R4þ : W ¼ þ j for j > 0 : L

  The endemic planar equilibrium point is E3 e=a; ða  eÞ=a2 ; 0; 0 ; where a  e > 0. E4 ðs¯; 0; p¯ 1 ; p¯ 2 Þ is disease free space equilibrium point where

This proves the theorem.

p¯ 1 ¼

&

l ; k  ml

p¯ 2 ¼

kð f s¯  jð1 þ csÞÞ ¯ hð1 þ cs¯Þ

366

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

and s¯ is the positive root of the equation cs2 þ ð1  cÞsþ b p¯ 1  1 ¼ 0. The disease free space equilibrium point E4 exists if k  ml > 0;

kð f s¯  jð1 þ csÞÞ ¯ > 0;

1 < c;

2

ð1  cÞ ¼ 4cðb p¯ 1  1Þ ˆ pˆ 1 ; 0Þ is the top-predator free equilibrium point where E5 ðˆs; i; sˆ ¼

e þ d pˆ 1 ; a

A p21 þ B p1 þ C ¼ 0; where 2

Hence E2 is stable if the condition given in the theorem is satisfied. The Jacobian matrix at E3 is given by 2 3 e be  0 7 6  a e ð1 þ ceÞ 6 7 6 7 6a  e 7 dða  eÞ 6 0  0 7 6 a 7 2 a J3 ¼ 6 7 6 7 6 7 fe gða  eÞ 6 0 þ 0 j 0 7 6 7 1 þ ce a2 4 5

e f þ d f pˆ 1 iˆ ¼ j  a þ ce þ cd pˆ 1

and pˆ 1 is the positive root of the equation

A ¼ cd ;

The characteristic roots of the Jacobian matrix J 2 are ku 2 =ð1 þ mu 2 Þ  l; au1  du2  e; and roots of the equation   bcu1 u2 b f u1 u2 l2   u1 l þ ¼0 3 1 þ cu 1 ð1 þ cu 1 Þ

B ¼ facd þ a2 d f  ð2cde þ ad þ a2 jcd þ a2 bÞg;

0

and

0

0

l

C ¼ aða þ ceÞ þ a2 e f  eða þ ceÞ  a2 jða þ ceÞ:

The characteristic roots of the Jacobian matrix J3 are l; fe=ð1 þ ceÞ þ ðgða  eÞÞa2  j and the roots of the equation

The top-predator free equilibrium point E5 exists if

l2 þ l þ

j>

e f þ d f pˆ 1 ; a þ ce þ cd pˆ 1

B2 ¼ 4AC;

B < 0: 





The interior equilibrium point is given by E ðs ; i ; d p1 þ e s ¼ ; a and p2 ¼

a e

b p1 ð1  s Þ  ; i ¼ a að1 þ cs Þ 

p1 ;

p2 Þ,

where

l p1 ¼ k  ml

as¯ < e þ d p¯ 1 ;



The interior equilibrium point E exists if b p1 ; 1 þ cs

The Jacobian ðs; i; p1 ; p2 Þ is 2 F 1s F 1i 6 F 2s F 2i 6 J4 F 3s F 3i 0 0

k > ml;

fe=ð1 þ ceÞþ

s¯ >

hm p¯ 1 p¯ 2 ð1 þ m p¯ 1 Þ2

þ

bc p¯ 1 s¯ ð1 þ cs¯Þ2

; (7)

bc p¯ 1 > ð1 þ cs¯Þ2

fs þ gi > j 1 þ cs

and the top-predator free equilibrium point E5 is stable if

matrix J of the system (3) at any arbitrary point given by 3 F 1 p1 0 F 2 p1 0 7 7: (6) F 3 p1 F 3 p2 5 F 4 p1 F 4 p2

Theorem 2. The trivial equilibrium point E0 is always unstable. Existence of E2 and E3 ensures the instability of the axial equilibrium point E1 . The disease free planar equilibrium point E2 is locally stable if ku2 =ð1 þ mu 2 Þ < l; au1 < du2 þ e; bcu2 < 1 þ cu1 . The endemic planar equilibrium point E3 is locally stable if fe=ð1 þ ceÞ þ ðgða  eÞÞ=a2 < j.

Proof. Since an eigenvalue associated with the Jacobian matrix at E0 is 1, so E0 is an unstable equilibrium point. From the Jacobian matrix computed around E1 it is observed that the equilibrium point E1 is unstable if a > e and f > ðc þ 1Þ j, which are the existence conditions for the equilibrium points E2 and E3 . The Jacobian matrix at E2 is given by 2 bcu2 bu 1 s au1  2 6 ð1 þ cu 1 Þ ð1 þ c u Þ 6 1 6 6 0 a u  d u  e 0 1 2 6 6 J2 ¼ 6 6 f u2 6 g u2 0 6 6 ð1 þ cu1 Þ2 6 4 0 0 0

Hence, the equilibrium point E3 is stable if ðgða  eÞÞ=a2 < j. &

Theorem 3. The disease free space equilibrium point E4 is locally stable if

 1 þ m p1 fs  þ gi  j : h 1 þ cs

1 > s þ

eða  eÞ ¼ 0: a

0

3

7 7 7 7 0 7 7 7 7 hu 2 7  ð1 þ mu2 Þ 7 7 7 5 ku2 l 1 þ mu2

k pˆ 1 < l; 1 þ m pˆ 1

2

bc pˆ 1 < ð1 þ cˆsÞ ;

bgð1 þ cˆsÞ  d f

Proof. The Jacobian matrix at E4 is given by 2 bc p¯ 1 s¯ b¯s s¯ a¯s  6 ð1 þ cs¯Þ 6 ð1 þ cs¯Þ2 6 6 0 as¯  d p¯ 1  e 0 6 6 6 J4 ¼ 6 f p¯ 1 hm p¯ 2 p¯ 2 6 g p¯ 1 6 6 ð1 þ cs¯Þ2 ð1 þ m p¯ 1 Þ2 6 6 k p¯ 2 4 0 0 ð1 þ m p¯ 1 Þ2

(8)

3 0

7 7 7 7 0 7 7 7 7 h p¯ 1 7  7 ð1 þ m p¯ 1 Þ 7 7 7 5 0

The characteristic roots of the Jacobian matrix J 4 are as¯  d p¯ 1  e; and the roots of the equation given by

l3 þ V1 l2 þ V2 l þ V3 ¼ 0; where

"

V1 ¼  V2 ¼ V3 ¼

hm p¯ 1 p¯ 2

ð1 þ m p¯ 1 Þ2

þ

bc p¯ 1 s¯

#

s¯ ð1 þ cs¯Þ2(

) hm p¯ 1 p¯ 2 bc p¯ 1 s¯ b f s¯ p¯ 1 þ ¯ s þ 2 3 ð1 þ m p¯ 1 Þ3 ( ð1 þ m p¯ 1 Þ2) ð1 þ csÞ ð1 þ csÞ ¯ ¯ hk p¯ 1 p¯ 2 bc p¯ 1 s¯ s¯ ð1 þ m p¯ 1 Þ3 ð1 þ cs¯Þ2 hk p¯ 1 p¯ 2

Hence, by Routh–Hurwitz criterion the equilibrium point E4 is stable if the condition (7) holds.

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

Using the Routh–Hurwitz criteria we observe that the system (3) is stable around the positive equilibrium point E if the conditions stated in the theorem hold.

The Jacobian matrix at E5 is given by 2 3 bc pˆ 1 sˆ bˆs ˆ s aˆ s  0 6 7 ð1 þ cˆsÞ 6 ð1 þ cˆsÞ2 7 6 7 6 7 ˆ ˆ ai 0 di 0 6 7 6 7 6 7 J5 ¼ 6 7 ˆ ˆ f p h p 1 1 6 7 ˆ g p 0  1 6 2 1 þ m pˆ 1 7 6 ð1 þ cˆsÞ 7 6 7 4 5 k pˆ 1 0 0 0 l 1 þ m pˆ 1

5. Permanence and impermanence of the system (3)

The characteristic roots of the Jacobian matrix J 5 are k pˆ 1 =ðð1 þ m pˆ 1 ÞÞ  l and the roots of the equation given by

l3 þ Q1 l2 þ Q2 l þ Q3 ¼ 0; where

"

Q1 ¼ 

Q3 ¼

 sˆ ; 2 ð1 þ cˆsÞ

Q2 ¼ gd pˆ 1 Iˆ þ a2 sˆIˆ þ

b f pˆ 1 sˆ 3

ð1 þ cˆsÞ

ˆ ˆ bc pˆ 1 sˆ ˆ  ad f sˆ pˆ 1 i þ agbˆs pˆ 1 i :  sˆ ðgd pˆ 1 IÞ 2 1 þ cˆs 1 þ cˆs ð1 þ cˆsÞ

Hence, by Routh–Hurwitz criterion the equilibrium point E5 is stable if the condition (8) holds. & Theorem 4. The interior point E ðs ; i ; p1 ; p2 Þ for the system (3) is locally asymptotically stable if the following conditions are hold as follows:

s 1 s 2  s 3 > 0;

fe gða  eÞ þ > j; 1 þ ce a2 ðiiÞ as¯ > ð1 þ c p¯ 1 Þ2 ; k pˆ 1 ðiiiÞ > l: 1 þ m pˆ 1

ðiÞ



s 1 > 0;

From biological point of view, permanence of a system means the survival of all populations of the system in future time. Mathematically, permanence of a system means that strictly positive solutions do not have omega limits points on the boundary of the non negative cone.

Theorem 5. Let f > ðc þ 1Þ j and a > e and the following conditions are satisfied

#

bc pˆ 1 sˆ

367

s 3 ðs 1 s 2  s 3 Þ  s 4 s 21 > 0;

Further if there exists a finite number of periodic solutions s ¼ fr ðtÞ; p1 ¼ cr ðtÞ;r ¼ 1; 2; . . . ; n, in the s  p1 plane, then system (3) is uniformly persistent provided for each periodic solutions of period T, Z 1 T hr ¼ e þ ðafr  dcr Þ dt > 0; T 0 r ¼ 1; 2; . . . ; n.

where s i ’s are given in the proof of the theorem. Proof. The Jacobian matrix at the interior point E ðs ; i ; p1 ; p2 Þ is 2 3 A11 A12 A13 A14 6 7 6A 7 6 21 A22 A23 A24 7 7 V¼6 6 7 6 A31 A32 A33 A34 7 4 5 A41

A42

A43

A44

where A11 A21 A31 A41

bcs p1

bs ¼ s ; A12 ¼ as ; A13 ¼  ; A14 ¼ 0 1 þ cs ð1 þ cs Þ2   ¼ ai ; A22 ¼ 0; A23 ¼ di ; A24 ¼ 0 f p1 hm p1 p2 h p1 ¼ ; A ¼ g p ; A ¼ ; A34 ¼  32 33 1 2 2  1 þ m p1 ð1 þ cs Þ ð1 þ m p1 Þ k p2 ¼ 0; A42 ¼ 0; A43 ¼ ; A44 ¼ 0 ð1 þ m p1 Þ2 



The characteristic equation of Jacobian matrix is given by

l4 þ s 1 l3 þ s 2 l2 þ s 3 l þ s 4 ¼ 0 where "

s3 ¼ s4 ¼

&

We claim that E0 2 = VðxÞ. If E0 2 VðxÞ then by the Butler– T McGehee lemma there exists a point P in VðxÞ W s ðE0 Þ where W s ðE0 Þ denotes the stable manifold of E0 . Since o ðPÞ lies in VðxÞ and W s ðE0 Þ is the i  p1  p2 space, we conclude that o ðPÞ is unbounded, which is a contradiction. Next E1 2 = VðxÞ, for otherwise, since E1 is a saddle point (which follows from the existence of E2 and E3 ), by the Butler–McGehee T lemma there exists a point P in VðxÞ W s ðE1 Þ. Now W s ðE1 Þ is the s  p2 plane implies that an unbounded orbit lies in VðxÞ, a contradiction. Next we show that E3 2 = VðxÞ. If E3 2 VðxÞ, the condition fe=ð1 þ ceÞ þ ðgða  eÞÞ=a2 > j implies that E3 is saddle point. W s ðE3 Þ is the s  i  p2 space and hence the orbits in this space emanate from either E0 or E1 or an unbounded orbit lies in VðxÞ, again a contradiction. The condition as¯ > e þ d p¯ 1 implies that E4 is a unstable point and also the condition k pˆ 1 =ð1 þ m pˆ 1 Þ > l implies that E5 is

# hm p1 p2  s þ ; ð1 þ cs Þ2 ð1 þ m p1 Þ2 ( ) bc p1 s hm p1 p2 b f p1 s kh p1 p2 þ s þg p1 di þ a2 s i  ; 3 2 2 ð1 þ m p1 Þ ð1 þ cs Þ ð1 þ m p1 Þ ð1 þ cs Þ3 " #" # bc p1 s kh p1 p2 hm p1 p2 f p1 i p1 s  g p1 di  a2 s i þads i þagbs 2 3 2 2     ð1 þ cs Þ ð1 þ cs Þ ð1 þ m p1 Þ ð1 þ m p1 Þ ð1 þ cs Þ   kh p1 p2 a2 s i ð1 þ m p1 Þ3

s1 ¼  s2 ¼

Proof. Let x be a point in the positive quadrant and oðxÞ be orbit through x and V be the omega limit set of the orbit through x. Note that VðxÞ is bounded.

bc p1 s

= VðxÞ unstable. So by similar arguments we can show that E4 2 and E5 2 = VðxÞ.

368

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

Lastly we show that no periodic orbit in the s  p1 plane or E2 2 VðxÞ. Let r i i ¼ 1; 2; . . . ; n denote the closed orbit of the periodic solution ðfr ðtÞ; cr ðtÞÞ in s  p1 plane such that r i lies inside r i1 . Let, the Jacobian matrix J given in (6) corresponding to r i is denoted by J r ðfr ðtÞ; 0; cr ðtÞ; 0Þ. Computing the fundamental matrix of the linear periodic system,

Proof. Proof is obvious and hence omitted. Before obtaining the conditions for impermanence of system (3), we briefly define the impermanence of a system. Let x ¼ ðx1 ; x2 ; x3 ; x4 Þ be the population vector, let D ¼ fx : x1 ; x2 ; x3 ; x4 > 0g, and @D is the boundary of D. rð:; :Þ is the distance in R4þ . Let us consider the system of equations is

X ¼ J r ðtÞX; Xð0Þ ¼ X 0 :

x˙ ¼ xi f i ðxÞ;

We find that its Floquet multiplier in the i direction is eh1 T . Then proceeding in an analogous manner like Kumar and Freedman (1989), we conclude that no r i lies on VðxÞ. Thus, VðxÞ lines in the positive quadrant and system (3) is persistent. Finally, since only the closed orbits and the equilibria from the omega limit set of the solutions on the boundary of R4þ and system (3) is dissipative. Now using a theorem of Butler et al. (1986), we conclude that system (3) is uniformly persistent. &

where f i : R4þ ! R and f i 2 C 1 .

0

Theorem 6. Let f > ðc þ 1Þ j and a > e and the following conditions are satisfied fe gða  eÞ þ > j; 1 þ ce a2 2 ðiiÞ a¯s > ð1  þ c p¯ 1 Þ ;  k pˆ 1 ku2 ðiiiÞ min ; >l 1 þ m pˆ 1 1 þ mu 2

ðiÞ

and if there exists no limit cycle in the s  p1 plane, then system (3) is uniformly persistent.

i ¼ 1; 2; 3; 4

The semi orbit g þ is defined by the set fxðtÞ : t > 0g where xðtÞ is the solution with initial value xð0Þ ¼ x0 . The above system is said to be impermanent (Hutson and Law, 1985) if and only if there is an x 2 D such that lim t ! 1 rðxðtÞ; @DÞ ¼ 0. Thus a community is impermanent if there is at least one semi-orbit which tends to the boundary. & Theorem 7. Let f > ðc þ 1Þ j and a > e and if the condition fe=ð1 þ ceÞ þ ðgða  eÞÞ=a2 < j holds, then the system (3) is impermanent. Proof. The conditions f > ðc þ 1Þ j and a > e are obtained from the existence of the equilibria points E2 and E3 . The given condition fe=ð1 þ ceÞ þ ðgða  eÞÞ=a2 < j implies that E3 is a saturated equilibrium point on the boundary. Hence, there exists at least one orbit in the interior that converges to the boundary (Hofbauer, 1986). Consequently the system (3) is impermanent (Hutson and Law, 1985). &

Fig. 1. The figure depicts the solution of system (3) in the absence of the disease. (a) Depicts the steady state stable distribution with b ¼ 5:0; f ¼ 5:0; h ¼ 0:1; k ¼ 0:1; m ¼ 2:0; j ¼ 0:4; l ¼ 0:01; c ¼ 0:5. (b) Shows the limit cycle oscillation for c ¼ 2:2. (c) Depicts the period doubling of the system for c ¼ 2:3. Finally, (d) shows the chaotic tea-cup for c ¼ 3.

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

6. Numerical results In this section we perform the numerical experiments to observe the dynamics of the system (3) with the following hypothetical set of parameter values and most of the values are taken from Hastings and Powell (1991). The parameter values are

369

disease. We observe that for 0 < c < 2:1 the system (3) is stable around the positive steady state (Fig. 1 a) and for 2:1 < c < 2:3 it shows limit cycle oscillations (Fig 1b). The period doubling is observed at c ¼ 2:3 (Fig. 1c). Finally, chaotic dynamics is observed for 2:3 < c  6:2 (Fig. 1d). 6.2. Effect of disease in the system (3)

b ¼ 5; m ¼ 2;

d ¼ 3; j ¼ 0:4;

e ¼ 0:5; k ¼ 0:1;

f ¼ 5;

g ¼ 2:5;

h ¼ 0:1;

l ¼ 0:01:

This parametric set is being kept fixed throughout the numerical experiments. In our numerical experiment we have varied a and c. For c ¼ 3 and a ¼ 1:3, we obtain the positive equilibrium E ð0:6731; 0:0922; 0:1250; 11:8157Þ. For the above set of parameter values we have s 1 ¼ 0:3456 > 0; s 4 ¼ 0:00079347 > 0; s 1 s 2 s 3 ¼ 0:0046 > 0; s 3 ðs 1 s 2  s 3 Þ  s 4 s 21 ¼ 0:00016281 > 0, which means that the system (3) is locally asymptotically stable around positive equilibrium E . 6.1. Dynamics of the disease free system of (3) Hastings and Powell (1991) observed stable focus, limit cycle, period doubling and chaotic dynamics of the system by changing half saturation constant c. Our aim is, first to observe the exchange of states (stability-limit cycle-period doubling-chaos) as shown in HP model for different values of c and subsequently to observe the system under the variation of infection rate a. Thus we first run the system (3) with the assumption that the system is free from

To observe the effect of disease in the system (3) we increase the force of infection a, keeping the half-saturation constant c fixed at c ¼ 3:0. For a  1:15 the system shows chaotic dynamics (Fig. 2 a). For a ¼ 1:16 the system enters into period doubling position from chaotic region and very famous tea cup chaotic attractor is disappeared (Fig. 2b). This shows that the presence of disease in the system prevents it from entering into chaotic region. The system goes into the limit cycle oscillation from period-doubling for a ¼ 1:2 (Fig. 2c). The system goes into the steady stable position from limit-cycle for a ¼ 1:3 (Fig. 2d). To make it more clear we plot the bifurcation diagram for all the populations with a as the bifurcating parameter, see Fig. 3. We observed that the system shows period-doubling in the range 1:16  a < 1:2, limit cycle in the range 1:2  a < 1:3 and finally system (3) settles down to steady state solutions, a stable situation for a  1:3 keeping c ¼ 3:0; m ¼ 2; g ¼ 2:5 and other parameters fixed. So it is clear that increasing the strength of force of infection up to a certain level reduce the prevalence of chaos, perioddoubling and limit cycle oscillations. This observations indicate that the disease factor could be used as a biological control parameter for persistence of the species.

Fig. 2. (a) Depicts the chaotic dynamics of other three species for b ¼ 5:0; f ¼ 5:0; h ¼ 0:1; k ¼ 0:1; m ¼ 2:0; j ¼ 0:4; l ¼ 0:01; c ¼ 3; d ¼ 3; e ¼ 0:5; g ¼ 2:5; a ¼ 1:15. (b) Depicts the period-doubling of the system (3) for a ¼ 1:17. (c) Depicts the co-existence of all four species in limit cycle oscillation for a ¼ 1:2. Finally, (d) Depicts the steady state stable distribution of all four species for a ¼ 1:3.

370

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

Fig. 3. Figure depicts the bifurcation for a 2 ½1:16; 1:5 and c ¼ 3; m ¼ 2; g ¼ 2:5, showing that system enters into order from chaos.

But while increasing the infection rate to control the chaos one must be careful otherwise the increase in the infection rate may bring some form of instability into the system which leads to the extinction of a population. This is clear from the bifurcation diagram shown in Fig. 4. Here we plot the populations with respect

to the infection rate a and we vary a from 1.5 to 8.5. We observe that when a crosses 3 (approx.) the system changes from stable steady state solution to limit cycle solution. If we go on increasing the infection rate the top-predator is finally washed away from the system, see Fig. 4.

Fig. 4. Figure depicts the bifurcation for a 2 ½1:5; 8:5 and c ¼ 3; m ¼ 2; g ¼ 2:5, showing that system enters into unstable situation from stable condition.

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

371

ecological allometric principles of body size by using logarithmic scale (log10) of this size for parameterizations of half-saturation constant of zooplankton. There are some authors (Jogensen et al., 2002; Ray et al., 2001a,b) who used the cell or body volume as measure of size for the scaling of allometric relationship with half-saturation constant of zooplankton and observed that normally the zooplankton comprises between the smaller 10 mm3 and larger 104 mm3 . We adopt the following relation from Mandal et al. (2007), c ¼ 1:2log ðV z Þ  0:5:

Fig. 5. Figure depicts the chaotic dynamics of system (3) with b ¼ 5:0; f ¼ 5:0; h ¼ 0:1; k ¼ 0:1; m ¼ 2:0; j ¼ 0:4; l ¼ 0:01; d ¼ 3; e ¼ 0:5; g ¼ 2:5; a ¼ 1:16; c ¼ 2:92 (corresponding intermediate predator body volume Vz2:85 .

6.3. Effect of body size of intermediate predator in system (3) In this section we observe that body size of intermediate predator controls the chaos of system (3). Here we vary the halfsaturation constant in such a way that the body size of the intermediate predator population gradually shifted from small body size to large. We proceed in the same way as deduced by Mandal et al. (2007). To make the calculation easier here we assume that (3) represents aquatic system with phytoplankton as prey, zooplankton as intermediate-predator and fish as toppredator. Mandal et al. (2007) make deductions from general

(9)

Here we denote V z as body volume of zooplankton and Vzi ¼ log ð10i mm3 Þ. We have increased the body size of zooplankton by incorporating the gradual increment of half-saturation constant of zooplankton following the Eq. (9). The system (3) shows chaotic dynamics with b ¼ 5:0; f ¼ 5:0; h ¼ 0:1; m ¼ 2:0; j ¼ 0:4; l ¼ 0:01; g ¼ 2:5; a ¼ 1:16; c ¼ 2:92 (corresponding intermediate predator body volume Vz2:85 ), see Fig. 5. To observe the effect of different body size of intermediate predator on the System (3) we again plot the bifurcation diagram of all the populations with V z as the bifurcation parameter, see Fig. 6. We observe that increase in the body size of intermediate predator also reduce the prevalence of chaotic dynamics originally observed in HP model. Thus, from the above results we observe that both the strength of force of infection a and body size of intermediate predator V z play important role in maintaining the stability of the system around the positive steady state. Now it would be interesting to observe how the force of infection a vary with the body size V z to maintain the stability of the system around the positive steady state. For this we produce a 2D-parametric figure where we plot the body size V z along the x-axis and the force of infection a along the y-axis, see Fig. 7. The idea behind generating Fig. 7 is very simple. For each value of V z we integrate the system for different values of a. We collect the minimum and maximum value of a (for a

Fig. 6. Figure shows the bifurcation for V z 2 ½2:85; 5:8 and a ¼ 1:16 and it also shows that system enters into order from chaos.

372

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

population is stable (see Fig. 8). Now we vary the parameter value a to see other dynamics. We note that oscillatory behaviour observed in original HP model stabilizes in the presence of infection in prey population and this stabilization phenomena is observed for an internal of force of infection a. From Fig. 9 we have observed that there is a critical value of force of infection a (ac ¼ 3:2) above which all four species coexist in oscillating position. A bifurcation diagram is presented in Fig. 10 to observe the global dynamics of the system. We observe from Fig. 10 that all species coexist in stable position for 2:0  a < 3:2 and all species survive in oscillating position for a  3:2. 7. Discussion

Fig. 7. Figure shows the plot between the body size V z and the minimum value of the infection rate a for which the system (3) is stable around the positive steady state.

particular V Z ) for which the system is stable around the interior equilibrium point. Finally we plot this maximum and minimum values of a in y-axis and the corresponding value of V z in the x-axis. To make the curve smooth we use cubic spline data interpolation technique. We observe that with the increase in the value of V z initially there is an increase in the area of the parameter a for which the system is stable. After certain critical value of V z this area becomes constant. This findings show that with the increase of the body size of the intermediate predator, the system can tolerate more disease in order to maintain the stability of the system around its steady state. Now we observe the dynamical behaviour of HP model with disease in prey population (i.e. our proposed model) by choosing another set of parameter values. We choose the parameter values a ¼ 2:0; b ¼ 2:5; c ¼ 5:5; d ¼ 3:5; e ¼ 0:6; f ¼ 2:2; g ¼ 2:2; h ¼ 0:15; m ¼ 4:2; j ¼ 0:3; k ¼ 0:12; l ¼ 0:01. It is interesting to note for these set of parameter values the HP model disease in prey

In the present paper we modified the HP model by introducing disease in the prey population. We worked out the conditions for which the population will go to extinction and the conditions for which the system is stable around different equilibria. We have also derived the conditions for which the system is permanent or impermanent. Hastings and Powell (1991) observed stability, limit cycle, period-doubling and finally chaos in a simple tritrophic food chain model by increasing half saturation constant c. However, the introduction of disease in the prey population changes the dynamics of the original model. From our numerical simulations we observe that by increasing the force of infection a the system enters into period-doubling from chaos, limit cycle solution from period-doubling and finally settles down to a steady state solution, a stable situation when other parameters are fixed. Thus increase in the infection rate not only prevent the system from entering into chaotic region, but also make the system stable around the coexistence steady state. But while increasing the infection rate one must keep in mind that it should not become too high, otherwise this may cause instability to the system leading to extinction of population. We also observe that body size of intermediate predator control chaotic dynamics. The system (3) settles down into stable from chaos by increasing the body size of predator. The above findings clearly indicate that moderate infection and body size of intermediate predator are prime factors for disappearance of chaotic dynamics observed in HP model.

Fig. 8. Figure shows the stable distribution of all four species for a ¼ 2:0; b ¼ 2:5; c ¼ 5:5; d ¼ 3:5; e ¼ 0:6; f ¼ 2:2; g ¼ 2:2; h ¼ 0:15; m ¼ 4:2; j ¼ 0:3; k ¼ 0:12; l ¼ 0:01.

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

373

Fig. 9. Figure shows that all populations coexist in oscillating position for a ¼ 3:2 and other parameters fixed as in Fig. 8.

To make our conclusion robust, we choose another set of parameter values. This experiment also confirms that moderate infection in prey population has a stabilizing effect in HP model. We like to mention that our findings are in well agreement with others observations in different context. For example, Chattopadhyay and Arino (1999) showed that a classical eco-epidemiological model entered into Hopf bifurcation for increasing the force of infection. Anderson and May (1986) showed that invasion of a resident predator–prey system by a new strain of parasites could cause destabilization and exhibits limit cycle oscillation. Oscillatory coexistence are particularly worthwhile to study, because prey-

predator (or host-parasitoid, plant-herbivore, consumer-resource) systems are well-known examples of inherently fluctuating populations (Turchin, 2003). Oscillations are also a concern of biological conservation, because populations could reach such small abundances that they are likely to go extinct (Rosenzweig, 1971). We may conclude that increasing the strength of force of infection and body size of intermediate predator reduce the prevalence of chaos, period-doubling and limit cycle oscillations. We also observe that with increase in the body size the area of infection rate increases under which the system is stable around the coexistence steady state.

Fig. 10. Figure shows that all four species settle down into oscillating position from stable situation for a 2 ½2; 4 and other parameter fixed as in Fig. 8.

374

K.p. Das et al. / Ecological Complexity 6 (2009) 363–374

Acknowledgement The authors are grateful to the reviewers for their critical comments and suggestions which have immensely improved the content and presentation of this version. References Anderson, R.M., May, R.M., 1986. The invasion, persistence and spread of infectious diseases within animal and plant communities. Philos. Trans. R. Soc. Lond. B 314, 533–570. Arino, O., Abdllaoui, A., Mikram, J., Chattopadhyay, J., 2004. Infection on prey population may act as a biological control in ratio-dependent predator–prey model. Nonlinearity 17, 1101–1116. Blueweiss, L., Fox, H., Kudzma, V., Nakashima, D., Peters, R.H.S., 1978. Relationships between some body size and some life history parameters. Oecologia (Berl.) 37, 257–272. Butler, G.J., Freedman, H., Waltman, P., 1986. Uniformly persistent systems. Proc. Am. Math. Soc. 96, 425–429. Chatterjee, S., Bandyopadhyay, M., Chattopadhyay, J., 2006. Proper predation makes the system disease free- conclusion drawn from an eco-epidemiological model. J. Biol. Syst. 14 (4), 599–616. Chattopadhyay, J., Sarkar, R.R., 2003. Chaos to order: preliminary experiments with a population dynamics models of three trophic levels. Ecol. Model. 163 (1–2), 45–50. Chattopadhyay, J., Arino, O., 1999. A predator–prey model with disease in the prey. Nonlinear Anal. 36, 747–766. Dobson, A.P., 1988. The population biology of parasite induced changes in host behaviour. Q. Rev. Biol. 63, 139–165. Eisenberg, J.N., Maszle, D.R., 1995. The structural stability o f a three-species food chain model. J. Theor. Biol. 176, 501–510. Freedman, H.I., 1990. A model of predator–prey dynamics modified by the action of parasite. Math. Biosci. 99, 143–155. Hadeler, K.P., Freedman, H.I., 1989. Predator–prey populations with parasitic infection. J. Math. Biol. 27, 609–631. Hastings, A., Powell, T., 1991. Chaos in three-species food chain. Ecology 72 (3), 896– 903. Hethcote, H.W., Ma, Z., Lio, S., 2002. Effects of quarantine in six endemic models for infectious diseases. Math. Biosci. 180, 141–160. Hethcote, H.W., Han, W.W.L., Zhien, Ma., 2004. A predator–prey model with infected prey. Theor. Pop. Biol. 66, 259–268. Hofbauer, J., 1986. Saturated equilibria, permanence and stability for ecological systems. In: Groos, L., Hallam, T., Levin, S. (Eds.), Mathematical Ecology, Proc. Trieste. World Scientific, Singapore. Hrbacek, J., 1962. Species composition and amount of the zooplankton in relation to fish stock. Rozpr. Cesk. Akad. Ved., Rada. Mat. Prir. Ved. 10, 1–116. Hutson, V., Law, R., 1985. Permanent coexistence in general models of three interacting species. J. Math. Biol. 21, 289–298.

Jogensen, S.E., Ray, S., Berec, L., Straskraba, M., 2002. Improved calibration of a eutrophication model by use of the size variation due to succession. Ecol. Model. 153, 269–277. Kumar, R., Freedman, H., 1989. A mathematical model of facultative mutualism with populations interacting in a food chain. Math. Biosci. 97, 235–261. La Barbera, M., 1986. The evolution and ecology of body size. In: Ranp, M., Jablonsky, D. (Eds.), Patterns and Process in the History of Life. Springer-Verlag, Berlin, Heidelberg, pp. 69–98. Lehman, J., 1976. The filter feeder as an optimal forager and the predicted shapes of feeding curves. Limnol. Oceanogr. 21, 501–516. Lonngren, K.E., Bai, E.W., Ucar, A., 2004. Dynamics and synchronization of the Hastings–Powell model of the food chain. Chaos Solitons Fract. 20, 387–393. Maionchi, D.O., Reis, S.F.d., Aguiar, M.A.M.d., 2005. Chaos and pattern formation in a spatial tritrophic food chain. Ecol. Model., doi:10.1016/j.ecolmodel.2005.04.028. Mandal, S., Ray, S., Roy, S., Jorgensen, S.E., 2007. Investigation of thermodynamic properties in an ecological model developing from ordered to chaotic states. Ecol. Model. 204, 40–46. McCann, K., Hastings, A., 1997. Re-evaluating the omnivory-stability relationship in food webs. Proc. R. Soc. Lond. B 264, 1249–1254. McMahon, T.A., 1973. Size and shape in biology. science 179, 1201–1204. Nagumo, M., 1942. Uber die Lage der Integralkurven gew onlicher Differentialgleichungen. Proc. Phys. Math. Soc. Japan 24, 551. Rai, V., Sreenivasan, R., 1993. Period-doubling bifurcations leading to chaos in a model food chain. Ecol. Model. 69 (1–2), 63–77. Ray, S., Berec, L., Straskraba, M., Jogensen, S.E., 2001a. Optimization of exergy and the implications of body sizes of phytoplankton and zooplankton in an aquatic ecosystem model. Ecol. Model. 140, 219–234. Ray, S., Berec, L., Straskraba, M., Ulanowicz, R.E., 2001b. Evaluation of system performance through optimizing ascendancy in an aquatic ecosystem model. J. Biol. Syst. 9, 1–22. Rosenzweig, M.L., 1971. Paradox of enrichment: destabilization of exploitation ecosystems in ecological time. Science 171, 385–387. Ruxton, G.D., 1994. Low levels of immigration between chaotic populations can reduce system extinctions by inducing asynchronous cycles. Proc. R. Soc. Lond. B 256, 189–193. Ruxton, G.D., 1996. Chaos in a three-species food chain with a lower bound on the bottom population. Ecology 77 (1), 317–319. Tag, E.P.Y., 1995. The allometry of algal growth rates. J. Plank. Res. 17, 1325–1335. Turchin, P., 2003. Complex Population Dynamics. A Theoretical empirical Synthesis. Princeton University Press, Princeton, NJ. Varriale, M.C., Gomes, A.A., 1998. A study of a three species food chain. Ecol. Model. 110 (2), 119–133. Venturino, E., 2002. Epidemics in predator–prey model: disease in the predators. IMA J. Math. Appl. Med. Biol. 19, 185–205. Xiao, Y., Chen, L., 2001. Modelling and analysis of a predator–prey model with disease in the prey. Math. Biosci. 171, 59–82. Xu, C., Li, Z., 2002. Influence of intraspecific density dependence on a three-species food chain with and without external stochastic disturbances. Ecol. Model. 155, 71–83.