- Email: [email protected]
Persistence and global stability of a ratio-dependent predator–prey model with stage structure
Applied Mathematics and Computation 158 (2004) 729–744 www.elsevier.com/locate/amc
Persistence and global stability of a ratio-dependent predator–prey model with stage structure Rui Xu *, M.A.J. Chaplain, F.A. Davidson Department of Mathematics, University of Dundee, 23 Perth Road, Dundee DD1 4HN, UK
Abstract A ratio-dependent predator–prey model with stage structure for prey is investigated. First, sufficient conditions are derived for the uniform persistence and impermanence of the model. Next, by constructing appropriate Lyapunov functions, a set of easily verifiable sufficient conditions are obtained for the global asymptotic stability of nonnegative equilibria of the model. Numerical simulations are presented to illustrate the validity of our main results. Ó 2003 Elsevier Inc. All rights reserved. Keywords: Stage structure; Ratio-dependence; Predator–prey model; Uniform persistence; Global stability
1. Introduction The traditional Lotka–Volterra type predator–prey systems are very important in the models of multi-species population dynamics. Standard LotkaVolterra type models, on which a large body of existing predator–prey theory is built by assuming that the per capita rate of predation depends on the prey numbers only. Recently, the traditional prey-dependent predator–prey models have been challenged by several biologists (see, for example, [1–4]) based on the fact that functional and numerical response over typical ecological time scales * Corresponding author. Permanent address: Department of Mathematics, Institute of Shijiazhuang Mechanical Engineering, Shijiazhuang 050003, Hebei Province, PR China. E-mail address: [email protected] (R. Xu).
0096-3003/$ - see front matter Ó 2003 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2003.10.012
730
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
ought to depend on the densities of both prey and predator, especially when predators have to search, share or compete for food. A more suitable general predator–prey model should be based on the ‘‘ratio-dependent’’ theory. This roughly states that the per capita predator growth rate should be a function of the ratio of prey to predator abundance. Moreover, as the number of predators often changes slowly (relative to prey number), there is often competition among the predators, and the per capita rate of predation should therefore depend on the numbers of both prey and predator, most probably and simply on their ratio. These hypotheses are strongly supported by numerous field and laboratory experiment and observations [2,3,5,6]. Based on the Michaelis–Menten or Holling type-II function, Arditi and Ginzburg [5] proposed a ratio-dependent function of the form x cðx=yÞ cx P ¼ ; ¼ y m þ ðx=yÞ my þ x and the following ratio-dependent predator–prey model: x_ ¼ xða bxÞ cxy=ðmy þ xÞ; y_ ¼ yðd þ fx=ðmy þ xÞÞ:
ð1:1Þ
Here xðtÞ and yðtÞ represent the densities of the prey and the predator at time t, respectively. a=b is the carrying capacity of the prey, d > 0 is the death rate of the predator, and a; c; m and f =c are positive constants that stand for the intrinsic growth rate of the prey, capturing rate, half saturation constant and conversion rate of the predator, respectively. The ratio-dependent predator–prey models with or without time delays have been studied by many researchers recently and very rich dynamics have been observed (see, for example, [7–17] and the references cited therein). It is assumed in the standard ratio-dependent predator–prey model that each individual prey admits the same risk to be attacked by predator. This assumption is obviously unrealistic for many animals. In the natural world, there are many species whose individuals have a life history that take them through two stages, immature and mature. In this paper, we classify individuals of prey as belonging either the immature or the mature and suppose that the immature prey does not have a risk to be attacked by the predator. This seems reasonable for a number of mammals, where the immature preys concealed in the mountain cave, are raised by their parents; they do not necessarily go out for seeking food, the rate they are attacked by the predators can be ignored. Stage-structured models have received much attention in recent years (see, for example, [18–29]). In [18], a stage-structured model of single species growth consisting of immature and mature individuals was proposed and discussed. In [19], it was further assumed that the time from immaturity to maturity is itself state dependent. An equilibrium analysis and eventual lower and upper bounds
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
731
of positive solutions for the model were addressed. Recently, Wang and Chen [27], Magnusson [25], Zhang and Chen [29] proposed and investigated predator–prey models with stage structure for prey or predator to analyze the influence of a stage structure for the prey or the predator on the dynamics of predator–prey models. The main purpose of this paper is to study the effect of stage structure for prey on the dynamics of ratio-dependent predator–prey system. To do so, we study the following differential system: x_ 1 ðtÞ ¼ ax2 ðtÞ r1 x1 ðtÞ bx1 ðtÞ; a1 x2 ðtÞyðtÞ ; x_ 2 ðtÞ ¼ bx1 ðtÞ b1 x22 ðtÞ myðtÞ þ x2 ðtÞ a2 x2 ðtÞ y_ ðtÞ ¼ yðtÞ r þ ; myðtÞ þ x2 ðtÞ
ð1:2Þ
where x1 ðtÞ represents the density of immature individual preys at time t, and x2 ðtÞ denotes the density of mature individual preys at time t, yðtÞ represents the density of the predator at time t. The model is derived under the following assumptions: (A1) The immature prey population: the birth rate into the immature population is proportional to the existing mature prey population with a proportionality a > 0; the death rate is proportional to the existing immature prey population with proportionality r1 . (A2) The mature prey population: b1 is the death and overcrowding rate of the population; the transformation rate from the immature prey individuals to mature prey individuals is proportional to the existing immature prey population with proportionality b. (A3) The predator population: the predator feed only on the mature prey. r is the death rate of the predator; a1 is the capturing rate of the predator; a2 =a1 is the rate of conversion of nutrients into the reproduction of the predator; m is the half saturation rate of the predator. The initial conditions for system (1.2) take the form of x1 ð0Þ > 0;
x2 ð0Þ > 0;
yð0Þ > 0:
ð1:3Þ
It is easy to show that system (1.2) has a unique solution zðtÞ ¼ ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ satisfying initial conditions (1.3), and all solutions of system (1.2) corresponding to initial conditions (1.3) are defined on ½0; þ1Þ and remain positive for all t P 0. In this paper, the solution of system (1.2) satisfying initial conditions (1.3) is said to be positive. The organization of this paper is as follows. In the next section, sufficient conditions are derived for the uniform persistence and impermanence of system
732
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
(1.2). In Section 3, by means of suitable Lyapunov functions, sufficient conditions are obtained for the global stability of nonnegative equilibria of system (1.2). Numerical simulations are also presented respectively in Sections 2 and 3 to illustrate the feasibility of our main results. In Section 4, a short discussion is given to conclude this work.
2. Uniform persistence In this section, we will establish conditions for the persistence and impermanence of system (1.2) with initial conditions (1.3). Definition (i) System (1.2) is said to be uniformly persistent if there exists a compact region D Int R3þ such that every solution X ðtÞ of system (1.2) with initial conditions (1.3) eventually enters and remains in the region D. (ii) System (1.2) is said to be impermanent if there is a positive solution ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ of (1.2) satisfying minflim inf x1 ðtÞ; lim inf x2 ðtÞ; lim inf yðtÞg ¼ 0: t!þ1
t!þ1
t!þ1
In order to prove the uniform persistence of system (1.2), we first give a result on the boundedness of solutions of system (1.2). Lemma 2.1. Positive solutions of system (1.2) with initial conditions (1.3) are ultimately bounded. Proof. Let ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be any solution of system (1.2) with initial conditions (1.3). Let qðtÞ ¼ x1 ðtÞ þ x2 ðtÞ þ
a1 yðtÞ: a2
Calculating the derivative of qðtÞ along solutions of system (1.2), we obtain _ ¼ ax2 ðtÞ r1 x1 ðtÞ b1 x22 ðtÞ qðtÞ
a1 r yðtÞ a2
6 AqðtÞ þ ða þ AÞx2 ðtÞ b1 x22 ðtÞ 6 AqðtÞ þ where A ¼ minfr; r1 g. It follows from (2.1) that lim sup qðtÞ 6 t!1
aþA :¼ M : 4Ab1
aþA ; 4b1
ð2:1Þ
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
733
Therefore, there exist positive constants M > M and T1 > 0 such that if t P T1 , qðtÞ < M. This completes the proof. h We are now in a position to state our result on the uniform persistence of system (1.2). Theorem 2.1. System (1.2) is uniformly persistent provided that (H1) a2 > r, (H2) ab=ðr1 þ bÞ > a1 =m. Proof. Let ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be any positive solution of system (1.2) with initial conditions (1.3). It follows from the first and the second equations of system (1.2) that x_ 1 ðtÞ ¼ ax2 ðtÞ r1 x1 ðtÞ bx1 ðtÞ; a1 x_ 2 ðtÞ P bx1 ðtÞ b1 x22 ðtÞ x2 ðtÞ: m
ð2:2Þ
Consider the following auxiliary system u_ 1 ðtÞ ¼ au2 ðtÞ r1 u1 ðtÞ bu1 ðtÞ; a1 u_ 2 ðtÞ ¼ bu1 ðtÞ b1 u22 ðtÞ u2 ðtÞ: m
ð2:3Þ
System (2.3) has two equilibria: E10 ð0; 0Þ, E1 ðu1 ; u2 Þ, where au2 1 ab a1 ; u2 ¼ u1 ¼ : b1 r 1 þ b m r1 þ b It is easy to show that E10 is always a saddle point, E1 is always stable. System (2.3) can be rewritten as a ½u2 ðtÞðu1 ðtÞ u1 Þ þ u1 ðtÞðu2 ðtÞ u2 Þ; u1 b u_ 2 ðtÞ ¼ ½u1 ðtÞðu2 ðtÞ u2 Þ þ u2 ðtÞðu1 ðtÞ u1 Þ b1 u2 ðtÞðu2 ðtÞ u2 Þ: u2
u_ 1 ðtÞ ¼
ð2:4Þ Consider the following Lyapunov function u1 u2 V1 ðtÞ ¼ c1 u1 u1 u1 ln þ c2 u2 u2 u2 ln ; u1 u2 where c1 , c2 are constants to be determined.
ð2:5Þ
734
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
Calculating the derivative of V1 ðtÞ along solutions of (2.4) we derive 2 dV1 ðtÞ X u_ i ðtÞ ¼ ci ðui ðtÞ ui Þ dt ui ðtÞ i¼1 ¼
c1 aðu1 ðtÞ u1 Þ ½u2 ðtÞðu1 ðtÞ u1 Þ þ u1 ðtÞðu2 ðtÞ u2 Þ u1 u1 ðtÞ þ
c2 bðu2 ðtÞ u2 Þ ½u1 ðtÞðu2 ðtÞ u2 Þ þ u2 ðtÞðu1 ðtÞ u1 Þ u2 u2 ðtÞ 2
c2 b1 ðu2 ðtÞ u2 Þ : Set c1 ¼
ð2:6Þ
bu1 =ðau2 Þ,
c2 ¼ 1. It follows from (2.6) that sffiffiffiffiffiffiffiffiffiffi "sffiffiffiffiffiffiffiffiffiffi #2 dV1 ðtÞ b u2 ðtÞ u1 ðtÞ ¼ ðu1 ðtÞ u1 Þ þ ðu2 ðtÞ u2 Þ dt u2 u1 ðtÞ u2 ðtÞ 2
2
b1 ðu2 ðtÞ u2 Þ 6 b1 ðu2 ðtÞ u2 Þ :
ð2:7Þ
Integrating both sides of (2.7) on ½0; t yields Z t 2 V1 ðtÞ þ b1 ðu2 ðsÞ u2 Þ ds 6 V1 ð0Þ; 0 2
which implies ðu2 ðtÞ u2 Þ 2 L1 ½0; 1Þ. In addition, with a similar argument in the proof of Lemma 2.1 we can prove that u1 ðtÞ, u2 ðtÞ are bounded. This, together with (2.3), implies that u2 ðtÞ u2 and u_ 2 ðtÞ are uniformly continuous. Applying BarbalatÕs lemmas (Lemma 1.2.2 and 1.2.3, Gopalsamy [30]), we conclude that lim ðu2 ðtÞ u2 Þ2 ¼ 0;
t!1
which imply that lim u2 ðtÞ ¼ u2 ; t!1
lim u_ 2 ðtÞ ¼ 0;
t!1
lim u1 ðtÞ ¼ u1 :
t!1
Therefore, there exists a T2 P T1 such that if t P T2 , u1 ðtÞ > u1 e=2, u2 ðtÞ > u2 e=2; where e > 0 is chosen sufficiently small. A standard comparison argument shows that e e lim inf x1 ðtÞ P u1 ; lim inf x2 ðtÞ P u2 : t!þ1 t!þ1 2 2 Thus, there is T3 P T2 such that if t P T3 , x1 ðtÞ > u1 e, x2 ðtÞ > u2 e; where e > 0 is chosen sufficiently small. Using a similar argument in the proof of Theorem 2.2 in [13], we can derive from the third equation of system (1.2) that lim inf yðtÞ P t!þ1
u2 ða2 rÞ : 2mr
This completes the proof.
h
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
735
Theorem 2.2. If a2 < r or r1 > r, a1 =m a > r, then system (1.2) is impermanent. Proof. Let ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be a positive solution of system (1.2) with initial conditions (1.3). We first assume a2 < r. It follows from the third equation of system (1.2) that y_ ðtÞ < yðtÞða2 rÞ; which implies yðtÞ < yð0Þeða2 rÞt : Thus, we have limt!1 yðtÞ ¼ 0. In this case, the predator population will go to extinction. If r1 > r, a1 =m a > r, then there exists b > 0 such that a1 ¼ a þ r: mþb Let d ¼ ðx1 ð0Þ þ x2 ð0ÞÞ=yð0Þ < b. We now claim that x1 ðtÞ þ x2 ðtÞ
for all t P 0;
and
limt!þ1 ðx1 ðtÞ þ x2 ðtÞÞ ¼ 0:
Otherwise, there is a first time t1 such that x1 ðt1 Þ þ x2 ðt1 Þ ¼ b; yðt1 Þ
x1 ðtÞ þ x2 ðtÞ
for all t 2 ½0; t1 Þ:
Then for any t 2 ½0; t1 , we derive from the first and the second equations of (1.2) that d a1 x2 ðtÞyðtÞ ðx1 ðtÞ þ x2 ðtÞÞ ¼ ax2 ðtÞ r1 x1 ðtÞ b1 x22 ðtÞ dt myðtÞ þ x2 ðtÞ a1 x2 ðtÞ 6 ax2 ðtÞ r1 x1 ðtÞ m þ x2 ðtÞ=yðtÞ a1 6 r1 x1 ðtÞ þ a x2 ðtÞ 6 rðx1 ðtÞ þ x2 ðtÞÞ; mþb which yields x1 ðtÞ þ x2 ðtÞ 6 ðx1 ð0Þ þ x2 ð0ÞÞert :
ð2:8Þ
On the other hand, it follows from the third equation of system (1.2) that for all t P 0 y_ ðtÞ P ryðtÞ; which yields yðtÞ P yð0Þert :
ð2:9Þ
736
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
Thus, we derive from (2.8) and (2.9) that x1 ðtÞ þ x2 ðtÞ x1 ð0Þ þ x2 ð0Þ 6 ¼d
for all t 2 ½0; t1 ;
which leads to a contradiction. Hence, we obtain x1 ðtÞ þ x2 ðtÞ 6 x1 ð0Þ þ x2 ð0Þert
for all t P 0:
We therefore derive that lim ðx1 ðtÞ þ x2 ðtÞÞ ¼ 0;
t!þ1
which yields lim x1 ðtÞ ¼ lim x2 ðtÞ ¼ 0:
t!þ1
t!þ1
The proof is complete.
h
Theorem 2.3. If r1 > r, a1 =m a > r, then there exists a positive solution ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ of system (1.2) satisfying lim ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ ¼ ð0; 0; 0Þ:
t!þ1
Proof. Let b ¼ a1 =ða þ rÞ m. Assume that ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ is a solution of system (1.2) with ðx1 ð0Þ þ x2 ð0ÞÞ=yð0Þ < b. It follows from the proof of Theorem 2.2 that lim x1 ðtÞ ¼ lim x2 ðtÞ ¼ 0:
t!þ1
t!þ1
Since yðtÞ is positive and bounded, we have 0 6 lim inf yðtÞ 6 lim sup yðtÞ ¼ s < þ1: t!þ1
t!þ1
We now show that limt!þ1 yðtÞ ¼ 0. Otherwise, we have s > 0. Since limt!þ1 x2 ðtÞ ¼ 0, there is a T4 > 0 such that mrs x2 ðtÞ < for all t P T4 : 2a2 On the other hand, by the definition of s, there is a t1 > T4 such that s and y_ ðt1 Þ > 0: ð2:10Þ yðt1 Þ > 2 It follows from (2.10) that a2 x2 ðt1 Þyðt1 Þ ryðt1 Þ > 0; myðt1 Þ þ x2 ðt1 Þ
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
737
which yields x1 ðt1 Þ >
mr mrs yðt1 Þ > ; a2 2a2
which leads to a contradiction. Hence, we must have s ¼ 0, i.e., lim y2 ðtÞ ¼ 0:
t!þ1
The proof is complete.
We now give an example to show the feasibility of Theorem 2.3. Example 1. Consider the following system x_ 1 ðtÞ ¼ 3x2 ðtÞ x1 ðtÞ 2x1 ðtÞ; 6x2 ðtÞyðtÞ ; x_ 2 ðtÞ ¼ 2x1 ðtÞ 2x22 ðtÞ yðtÞ þ x2 ðtÞ 1:5x2 ðtÞ y_ ðtÞ ¼ yðtÞ 0:5 þ : yðtÞ þ x2 ðtÞ
ð2:11Þ
It is easy to show that r1 ¼ 1 > r ¼ 0:5, a1 =m ¼ 6 > a þ r ¼ 3:5. Thus, by Theorem 2.3 we see that there is a positive solution ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ of system (2.11) with ðx1 ð0Þ þ x2 ð0ÞÞ=yð0Þ < b ¼ 1=3 satisfying lim x1 ðtÞ ¼ lim x2 ðtÞ ¼ lim yðtÞ ¼ 0:
t!þ1
t!þ1
t!þ1
Numerical integration can be carried out using standard algorithms. We used the package ode45 in MATLAB. As shown in Fig. 1, numerical 2 x1 x2 y
solution
1.5
1
0.5
0
0.5
0
5
10
15
20
25
time t
Fig. 1. The temporal solution found by numerical integration of system (2.11) with initial conditions ðx1 ð0Þ; x2 ð0Þ; yð0ÞÞ ¼ ð0:25; 0:25; 2Þ.
738
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
simulation also confirms this fact for system (2.11) with initial conditions ðx1 ð0Þ; x2 ð0Þ; yð0ÞÞ ¼ ð0:25; 0:25; 2Þ.
3. Global stability In this section, we are concerned with the global asymptotic stability of nonnegative equilibria of system (1.2). The strategy is to construct appropriate Lyapunov functions. System (1.2) always has equilibria E0 ð0; 0; 0Þ, E1 ðx1 ; x2 ; 0Þ, where x1 ¼
a2 b b1 ðr1 þ bÞ
2
x2 ¼
;
ab : b1 ðr1 þ bÞ
It is easy to show that system (1.2) has a positive equilibrium E ðx1 ; x2 ; y Þ if and only if (H1) and the following hold: (H3) ab=ðr1 þ bÞ > a1 ða2 rÞ=ðma2 Þ,where x1 ¼
ax2 ; r1 þ b
x2 ¼
ab a1 ða2 rÞ ; b1 ðr1 þ bÞ ma2 b1
y ¼
x2 ða2 rÞ : mr
The Jacobian matrix of system (1.2) at E1 ðx1 ; x2 ; 0Þ takes the form of 0
r1 b @ b 0
a 2b1x2 0
1 0 a1 A: a2 r
It is easy to verify that if a < r, E1 is locally asymptotically stable, if a > r, it is locally unstable. The Jacobian matrix of system (1.2) at E ðx1 ; x2 ; y Þ takes the form of 0
r1 b B b @ 0
a 2 2 2b1 x2 ma1 ðy Þ =ðmy þ x2 Þ 2
ma2 ðy Þ =ðmy þ x2 Þ
2
1 0 2 2 a1 ðx2 Þ =ðmy þ x2 Þ C A: 2 ma2 x2 y =ðmy þ x2 Þ
It is easy to show that if ab=ðr1 þ bÞ > 2a1 =m, the positive equilibrium E is locally asymptotically stable. We are now able to state and prove our main result on the global asymptotic stability of the positive equilibrium E ðx1 ; x2 ; y Þ of system (1.2). Theorem 3.1. Let (H1) hold. Assume further that ab=ðr1 þ bÞ > 2a1 =m. Then the positive equilibrium E ðx1 ; x2 ; y Þ of system (1.2) is globally asymptotically stable.
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
739
Proof. Let zðtÞ ¼ ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be any positive solution of system (1.2) with initial conditions (1.3). System (1.2) can be rewritten as a ½x2 ðtÞðx1 ðtÞ x1 Þ þ x1 ðtÞðx2 ðtÞ x2 Þ; x1
x_ 1 ðtÞ ¼
b ½x1 ðtÞðx2 ðtÞ x2 Þ þ x2 ðtÞðx1 ðtÞ x1 Þ x2 þ x2 ðtÞ b1 ðx2 ðtÞ x2 Þ
x_ 2 ðtÞ ¼
ð3:1Þ
þ _ ¼ yðtÞ yðtÞ
a1 y ðx2 ðtÞ x2 Þ a1 x2 ðyðtÞ y Þ ; ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
ma2 ½y ðx2 ðtÞ x2 Þ x2 ðyðtÞ y Þ : ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
Define a Lyapunov function
V ðtÞ ¼
2 X xi y ci xi xi xi ln þ c3 y y y ln ; y xi i¼1
ð3:2Þ
where ci ði ¼ 1; 2; 3Þ are positive constants to be determined. Calculating the derivative of V ðtÞ along solutions of system (3.1), it follows that 2 dV ðtÞ X x_ i ðtÞ y_ ðtÞ ¼ þ c3 ðyðtÞ y Þ ci ðxi ðtÞ xi Þ dt x yðtÞ ðtÞ i i¼1
¼
c1 aðx1 ðtÞ x1 Þ ½x2 ðtÞðx1 ðtÞ x1 Þ þ x1 ðtÞðx2 ðtÞ x2 Þ x1 x1 ðtÞ c2 bðx2 ðtÞ x2 Þ ½x1 ðtÞðx2 ðtÞ x2 Þ þ x2 ðtÞðx1 ðtÞ x1 Þ x2 x2 ðtÞ þ c2 ðx2 ðtÞ x2 Þ b1 ðx2 ðtÞ x2 Þ þ
a1 y ðx2 ðtÞ x2 Þ a1 x2 ðyðtÞ y Þ þ ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ þ c3 ma2 ðyðtÞ y Þ
y ðx2 ðtÞ x2 Þ x2 ðyðtÞ y Þ : ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
ð3:3Þ
740
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
Set c1 ¼ b1 x1 =ða1 x2 Þ, c2 ¼ 1, c3 ¼ a1 x2 =ðma2 y Þ. We derive from (3.3) that sffiffiffiffiffiffiffiffiffiffi "sffiffiffiffiffiffiffiffiffiffi #2 dV ðtÞ b x2 ðtÞ x1 ðtÞ ¼ ðx1 ðtÞ x1 Þ þ ðx2 ðtÞ x2 Þ dt x2 x1 ðtÞ x2 ðtÞ 2
b1 ðx2 ðtÞ x2 Þ2 þ
a1 y ðx2 ðtÞ x2 Þ ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
2
2
a1 ðx2 Þ ðyðtÞ y Þ : y ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
ð3:4Þ
Choose a positive constant e such that ab=ðr1 þ bÞ 2a1 =m b1 e > 0. It follows from the proof of Theorem 2.1 that there is a T > 0 such that if t P T , x2 ðtÞ > u2 e. Thus we derive from (3.4) and Lemma 2.1 that for t P T , 2 2 dV ðtÞ a1 a1 ðx2 Þ ðyðtÞ y Þ 2 6 b1 : ðtÞ x Þ ðx 2 2 dt mðu2 eÞ y Mðm þ 1Þðmy þ x2 Þ ð3:5Þ Integrating both sides of (3.5) on ½0; t gives Z t a1 2 V ðtÞ þ b1 ðx2 ðsÞ x2 Þ ds mðu2 eÞ 0 Z t 2 a1 ðx2 Þ þ ðyðsÞ y Þ2 ds 6 V ð0Þ; y Mðm þ 1Þðmy þ x2 Þ 0 2
2
which implies ðx2 ðtÞ x2 Þ ; ðyðtÞ y Þ 2 L1 ½0; 1Þ. In addition, by Lemma 2.1, x2 ðtÞ, yðtÞ are bounded. This, together with (1.2), implies that x2 ðtÞ x2 , yðtÞ y and x_ 2 ðtÞ are uniformly continuous. Applying BarbalatÕs lemmas (Lemmas 1.2.2 and 1.2.3, Gopalsamy [30]), we conclude that lim ðx2 ðtÞ x2 Þ2 ¼ 0;
t!1
lim x_ 2 ðtÞ ¼ 0;
t!1
lim ðyðtÞ y Þ2 ¼ 0;
t!1
which imply that lim x2 ðtÞ ¼ x2 ;
t!1
lim x1 ðtÞ ¼ x1 ;
t!1
This completes the proof.
lim yðtÞ ¼ y :
t!1
h
In the following we discuss the global stability of the nonnegative equilibrium E1 ðx1 ; x2 ; 0Þ. Theorem 3.2. Let (H2) hold. If a2 < r, then the nonnegative equilibrium E1 ðx1 ; x2 ; 0Þ is globally asymptotically stable.
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
741
Proof. Let ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be a positive solution of system (1.2) with initial conditions (1.3). By the proof of Theorem 2.1, we see that there is a T3 > 0 such that if t P T3 , x2 ðtÞ > u2 =2. By the proof of Theorem 2.2, we have limt!þ1 yðtÞ ¼ 0. Thus, for any given e 2 ð0; ab=ðr1 þ bÞÞ, there is a T5 P T3 such that if t > T5 , a1 yðtÞ=ðmyðtÞ þ x2 ðtÞÞ < e. We consider the following two auxiliary systems: u_ 1 ðtÞ ¼ au2 ðtÞ r1 u1 ðtÞ bu1 ðtÞ; u_ 2 ðtÞ ¼ bu1 ðtÞ b1 u22 ðtÞ eu2 ðtÞ;
ð3:6Þ
and v_ 1 ðtÞ ¼ av2 ðtÞ r1 v1 ðtÞ bv1 ðtÞ; v_ 2 ðtÞ ¼ bv1 ðtÞ b1 v22 ðtÞ:
ð3:7Þ
Using similar arguments in the proof of Theorem 2.1, we can show that ae e lim u1 ðtÞ ¼ x1 ; lim u2 ðtÞ ¼ x2 ; t!þ1 t!þ1 b1 ðr1 þ bÞ b1 lim v1 ðtÞ ¼ x1 ; lim v2 ðtÞ ¼ x2 : t!þ1
t!þ1
A standard comparison argument shows that ae e x1 6 lim x1 ðtÞ 6 x1 ; x2 6 lim x2 ðtÞ 6 x2 : t!þ1 b1 ðr1 þ bÞ t!þ1 b1 Letting e ! 0; we finally obtain lim x1 ðtÞ ¼ x1 ;
t!þ1
lim x2 ðtÞ ¼ x2 :
t!þ1
This completes the proof.
h
We now give two examples to illustrate the validity of Theorems 2.1, 3.1 and 3.2. Example 2. As an example, we consider the following system: x_ 1 ðtÞ ¼ 4x2 ðtÞ 0:1x1 ðtÞ 2x1 ðtÞ; x2 ðtÞyðtÞ ; x_ 2 ðtÞ ¼ 2x1 ðtÞ x22 ðtÞ 2yðtÞ þ x2 ðtÞ 1:5x2 ðtÞ y_ ðtÞ ¼ yðtÞ 0:5 þ : 2yðtÞ þ x2 ðtÞ
ð3:8Þ
System (3.8) has a unique positive equilibrium E ð2920=441; 73=21; 73=21Þ. It is easy to verify that the coefficients of system (3.8) satisfy the assumptions of Theorem 2.1, 3.1. By Theorem 2.1 we see that system (3.8) is uniformly persistent. By Theorem 3.1 we see that the positive equilibrium E of system (3.8)
742
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744 7
x1 x2 y
6
solution
5
4
3
2
1
0
10
20
30
40
50 time t
60
70
80
90
100
Fig. 2. The temporal solution found by numerical integration of system (3.8) with initial conditions ðx1 ð0Þ; x2 ð0Þ; yð0ÞÞ ¼ ð1; 1; 1Þ.
is globally asymptotically stable. Numerical simulation also suggests these facts (see, Fig. 2). Example 3. Consider the following ratio-dependent predator–prey model x_ 1 ðtÞ ¼ 6x2 ðtÞ x1 ðtÞ 3x1 ðtÞ; x2 ðtÞyðtÞ ; x_ 2 ðtÞ ¼ 3x1 ðtÞ 2x22 ðtÞ yðtÞ þ x2 ðtÞ x2 ðtÞ y_ ðtÞ ¼ yðtÞ 2 þ : yðtÞ þ x2 ðtÞ
ð3:9Þ
3.5 x1 x2 y
3
solution
2.5 2 1.5 1 0.5 0
0
5
10
15
20
25 time t
30
35
40
45
50
Fig. 3. The temporal solution found by numerical integration of system (3.9) with initial conditions ðx1 ð0Þ; x2 ð0Þ; yð0ÞÞ ¼ ð1; 1; 1Þ.
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
743
System (3.9) has a nonnegative equilibrium E1 ð27=8; 9=4; 0Þ. It is easy to verify that the coefficients of system (3.9) satisfy the assumptions of Theorem 3.2. By Theorem 3.2 we see that E1 is globally asymptotically stable. Numerical simulation also confirms this fact (see, Fig. 3). 4. Discussion In this paper, we have discussed the effect of stage structure for prey on the dynamics of a ratio-dependent predator–prey model. By some comparison arguments, sufficient conditions are derived for the uniform persistence and impermanence of system (1.2). By constructing suitable Lyapunov functions, sufficient conditions are obtained for the global asymptotic stability of nonnegative equilibria of system (1.2). By Theorem 2.1, we see that system (1.2) will be uniformly persistent if the birth rate into the immature prey population, the rate of immature prey becoming mature prey, and the conversion rate and the half saturation rate of the predator are high and the capturing rate of the predator and the death rates of both the immature prey and the predator are low enough satisfying (H1)–(H2). By Theorem 2.2 we see that if the conversion rate of the predator is less than its death rate, then the predator will go to extinction, system (1.2) will not be persistent. Theorem 3.1 shows the global stability of the positive equilibrium of system (1.2) with somewhat stronger assumptions than those in Theorem 2.1 on the uniform persistence of system (1.2). We would like to mention here that we are unable to show system (1.2) admits limit cycles when E exists and is unstable. This is known to be true for system (1.1) (see, Hsu et al. [10]). We leave this for future work. References [1] R. Arditi, H. Saiah, Empirical evidence of the role of heterogeneity in ratio-dependent consumption, Ecology 73 (1992) 1544–1551. [2] R. Arditi, L.R. Ginzburg, H.R. Akcakaya, Variation in plankton densities among lakes: a case for ratio-dependent models, American Nat. 138 (1991) 1287–1296. [3] R. Arditi, N. Perrin, H. Saiah, Functional response and heterogeneities: an experiment test with cladocerans, OIKOS 60 (1991) 69–75. [4] A.P. Gutierrez, The physiological basis of ratio-dependent predator–prey theory: a methbolic pool model of NicholsonÕs blowflies as an example, Ecology 73 (1992) 1552–1563. [5] R. Arditi, L.R. Ginzburg, Coupling in predator–prey dynamics: ratio-dependence, J. Theor. Biol. 139 (1989) 311–326. [6] I. Hanski, The functional response of predator: worries bout scale, TREE 6 (1991) 141–142. [7] E. Beretta, Y. Kuang, Global analyses in some delayed ratio-dependent predator–prey systems, Nonlinear Anal. 32 (1998) 381–408. [8] F. Berezovskaya, G. Karev, R. Arditi, Parametric analysis of the ratio-dependent predator– prey model, J. Math. Biol. 43 (2001) 221–246.
744
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
[9] M. Fan, K. Wang, Periodicity in a delayed ratio-dependent predator–prey system, J. Math. Anal. Appl. 262 (2001) 179–190. [10] S.B. Hsu, T.W. Hwang, Y. Kuang, Global analysis of the Michaelis–Menten ratio-dependent predator–prey system, J. Math. Biol. 42 (2001) 489–506. [11] S.B. Hsu, T.W. Hwang, Y. Kuang, Rich dynamics of a ratio-dependent one prey two predator model, J. Math. Biol. 43 (2001) 377–396. [12] C. Jost, O. Arino, R. Arditi, About deterministic extinction in ratio-dependent predator–prey models, Bull. Math. Biol. 61 (1999) 19–32. [13] Y. Kuang, E. Beretta, Global qualitative analysis of a ratio-dependent predator–prey system, J. Math. Biol. 36 (1998) 389–406. [14] D. Xiao, S. Ruan, Global dynamics of a ratio-dependent predator–prey system, J. Math. Biol. 43 (2001) 268–290. [15] D. Xiao, S. Ruan, Stability and bifurcation in a delayed ratio-dependent predator–prey system, Proc. Edin. Math. Soc. 45 (2002) 205–220. [16] R. Xu, L. Chen, Persistence and stability for a two-species ratio-dependent predator–prey system with time delay in a two-patch environment, Comput. Math. Appl. 40 (2000) 577–588. [17] R. Xu, L. Chen, Persistence and global stability for n-species ratio-dependent predator–prey system with time delays, J. Math. Anal. Appl. 275 (2002) 27–43. [18] W.G. Aiello, H.I. Freedman, A time delay model of single-species growth with stage structure, Math. Biosci. 101 (1990) 139–153. [19] W.G. Aiello, H.I. Freedman, J. Wu, Analysis of a model representing stage-structured population growth with state-dependent time delay, SIAM J. Appl. Math. 52 (1992) 855–869. [20] F. Brauer, Z. Ma, Stability of stage-structured population models, J. Math. Anal. Appl. 126 (1987) 301–315. [21] H.I. Freedman, J. Wu, Persistence and global asymptotic stability of single species dispersal models with stage structure, Quart. Appl. Math. 49 (1991) 351–371. [22] Y. Kuang, J.W.H. So, Analysis of a delayed two-stage population with space-limited recruitment, SIAM J. Appl. Math. 55 (1995) 1675–1695. [23] S. Liu, L. Chen, Z. Liu, Extinction and permanence in nonautonomous competitive system with stage structure, J. Math. Anal. Appl. 274 (2002) 667–684. [24] S. Liu, L. Chen, R. Agarwal, Recent progress on stage-structured population dynamics, Math. Comput. Model. 36 (2002) 1319–1360. [25] K.G. Magnusson, Destabilizing effect of cannibalism on a structured predator–prey system, Math. Biosci. 155 (1999) 61–75. [26] X. Song, L. Chen, Optimal harvesting and stability for a two-species competitive system with stage structure, Math. Biosci. 170 (2001) 173–186. [27] W. Wang, L. Chen, A predator–prey system with stage structure for predator, Comput. Math. Appl. 33 (1997) 83–91. [28] W. Wang, G. Mulone, F. Salemi, V. Salone, Permanence and stability of a stage-structured predator–prey model, J. Math. Anal. Appl. 262 (2001) 499–528. [29] X. Zhang, L. Chen, The stage-structured predator–prey model and optimal harvesting policy, Math. Biosci. 168 (2000) 201–210. [30] K. Gopalsamy, Stability and Oscillations in Delay Differential Equations of Population Dynamics, Kluwer Academic, Dordrecht/ Norwell, MA, 1992.
Persistence and global stability of a ratio-dependent predator–prey model with stage structure Rui Xu *, M.A.J. Chaplain, F.A. Davidson Department of Mathematics, University of Dundee, 23 Perth Road, Dundee DD1 4HN, UK
Abstract A ratio-dependent predator–prey model with stage structure for prey is investigated. First, sufficient conditions are derived for the uniform persistence and impermanence of the model. Next, by constructing appropriate Lyapunov functions, a set of easily verifiable sufficient conditions are obtained for the global asymptotic stability of nonnegative equilibria of the model. Numerical simulations are presented to illustrate the validity of our main results. Ó 2003 Elsevier Inc. All rights reserved. Keywords: Stage structure; Ratio-dependence; Predator–prey model; Uniform persistence; Global stability
1. Introduction The traditional Lotka–Volterra type predator–prey systems are very important in the models of multi-species population dynamics. Standard LotkaVolterra type models, on which a large body of existing predator–prey theory is built by assuming that the per capita rate of predation depends on the prey numbers only. Recently, the traditional prey-dependent predator–prey models have been challenged by several biologists (see, for example, [1–4]) based on the fact that functional and numerical response over typical ecological time scales * Corresponding author. Permanent address: Department of Mathematics, Institute of Shijiazhuang Mechanical Engineering, Shijiazhuang 050003, Hebei Province, PR China. E-mail address: [email protected] (R. Xu).
0096-3003/$ - see front matter Ó 2003 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2003.10.012
730
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
ought to depend on the densities of both prey and predator, especially when predators have to search, share or compete for food. A more suitable general predator–prey model should be based on the ‘‘ratio-dependent’’ theory. This roughly states that the per capita predator growth rate should be a function of the ratio of prey to predator abundance. Moreover, as the number of predators often changes slowly (relative to prey number), there is often competition among the predators, and the per capita rate of predation should therefore depend on the numbers of both prey and predator, most probably and simply on their ratio. These hypotheses are strongly supported by numerous field and laboratory experiment and observations [2,3,5,6]. Based on the Michaelis–Menten or Holling type-II function, Arditi and Ginzburg [5] proposed a ratio-dependent function of the form x cðx=yÞ cx P ¼ ; ¼ y m þ ðx=yÞ my þ x and the following ratio-dependent predator–prey model: x_ ¼ xða bxÞ cxy=ðmy þ xÞ; y_ ¼ yðd þ fx=ðmy þ xÞÞ:
ð1:1Þ
Here xðtÞ and yðtÞ represent the densities of the prey and the predator at time t, respectively. a=b is the carrying capacity of the prey, d > 0 is the death rate of the predator, and a; c; m and f =c are positive constants that stand for the intrinsic growth rate of the prey, capturing rate, half saturation constant and conversion rate of the predator, respectively. The ratio-dependent predator–prey models with or without time delays have been studied by many researchers recently and very rich dynamics have been observed (see, for example, [7–17] and the references cited therein). It is assumed in the standard ratio-dependent predator–prey model that each individual prey admits the same risk to be attacked by predator. This assumption is obviously unrealistic for many animals. In the natural world, there are many species whose individuals have a life history that take them through two stages, immature and mature. In this paper, we classify individuals of prey as belonging either the immature or the mature and suppose that the immature prey does not have a risk to be attacked by the predator. This seems reasonable for a number of mammals, where the immature preys concealed in the mountain cave, are raised by their parents; they do not necessarily go out for seeking food, the rate they are attacked by the predators can be ignored. Stage-structured models have received much attention in recent years (see, for example, [18–29]). In [18], a stage-structured model of single species growth consisting of immature and mature individuals was proposed and discussed. In [19], it was further assumed that the time from immaturity to maturity is itself state dependent. An equilibrium analysis and eventual lower and upper bounds
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
731
of positive solutions for the model were addressed. Recently, Wang and Chen [27], Magnusson [25], Zhang and Chen [29] proposed and investigated predator–prey models with stage structure for prey or predator to analyze the influence of a stage structure for the prey or the predator on the dynamics of predator–prey models. The main purpose of this paper is to study the effect of stage structure for prey on the dynamics of ratio-dependent predator–prey system. To do so, we study the following differential system: x_ 1 ðtÞ ¼ ax2 ðtÞ r1 x1 ðtÞ bx1 ðtÞ; a1 x2 ðtÞyðtÞ ; x_ 2 ðtÞ ¼ bx1 ðtÞ b1 x22 ðtÞ myðtÞ þ x2 ðtÞ a2 x2 ðtÞ y_ ðtÞ ¼ yðtÞ r þ ; myðtÞ þ x2 ðtÞ
ð1:2Þ
where x1 ðtÞ represents the density of immature individual preys at time t, and x2 ðtÞ denotes the density of mature individual preys at time t, yðtÞ represents the density of the predator at time t. The model is derived under the following assumptions: (A1) The immature prey population: the birth rate into the immature population is proportional to the existing mature prey population with a proportionality a > 0; the death rate is proportional to the existing immature prey population with proportionality r1 . (A2) The mature prey population: b1 is the death and overcrowding rate of the population; the transformation rate from the immature prey individuals to mature prey individuals is proportional to the existing immature prey population with proportionality b. (A3) The predator population: the predator feed only on the mature prey. r is the death rate of the predator; a1 is the capturing rate of the predator; a2 =a1 is the rate of conversion of nutrients into the reproduction of the predator; m is the half saturation rate of the predator. The initial conditions for system (1.2) take the form of x1 ð0Þ > 0;
x2 ð0Þ > 0;
yð0Þ > 0:
ð1:3Þ
It is easy to show that system (1.2) has a unique solution zðtÞ ¼ ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ satisfying initial conditions (1.3), and all solutions of system (1.2) corresponding to initial conditions (1.3) are defined on ½0; þ1Þ and remain positive for all t P 0. In this paper, the solution of system (1.2) satisfying initial conditions (1.3) is said to be positive. The organization of this paper is as follows. In the next section, sufficient conditions are derived for the uniform persistence and impermanence of system
732
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
(1.2). In Section 3, by means of suitable Lyapunov functions, sufficient conditions are obtained for the global stability of nonnegative equilibria of system (1.2). Numerical simulations are also presented respectively in Sections 2 and 3 to illustrate the feasibility of our main results. In Section 4, a short discussion is given to conclude this work.
2. Uniform persistence In this section, we will establish conditions for the persistence and impermanence of system (1.2) with initial conditions (1.3). Definition (i) System (1.2) is said to be uniformly persistent if there exists a compact region D Int R3þ such that every solution X ðtÞ of system (1.2) with initial conditions (1.3) eventually enters and remains in the region D. (ii) System (1.2) is said to be impermanent if there is a positive solution ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ of (1.2) satisfying minflim inf x1 ðtÞ; lim inf x2 ðtÞ; lim inf yðtÞg ¼ 0: t!þ1
t!þ1
t!þ1
In order to prove the uniform persistence of system (1.2), we first give a result on the boundedness of solutions of system (1.2). Lemma 2.1. Positive solutions of system (1.2) with initial conditions (1.3) are ultimately bounded. Proof. Let ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be any solution of system (1.2) with initial conditions (1.3). Let qðtÞ ¼ x1 ðtÞ þ x2 ðtÞ þ
a1 yðtÞ: a2
Calculating the derivative of qðtÞ along solutions of system (1.2), we obtain _ ¼ ax2 ðtÞ r1 x1 ðtÞ b1 x22 ðtÞ qðtÞ
a1 r yðtÞ a2
6 AqðtÞ þ ða þ AÞx2 ðtÞ b1 x22 ðtÞ 6 AqðtÞ þ where A ¼ minfr; r1 g. It follows from (2.1) that lim sup qðtÞ 6 t!1
aþA :¼ M : 4Ab1
aþA ; 4b1
ð2:1Þ
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
733
Therefore, there exist positive constants M > M and T1 > 0 such that if t P T1 , qðtÞ < M. This completes the proof. h We are now in a position to state our result on the uniform persistence of system (1.2). Theorem 2.1. System (1.2) is uniformly persistent provided that (H1) a2 > r, (H2) ab=ðr1 þ bÞ > a1 =m. Proof. Let ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be any positive solution of system (1.2) with initial conditions (1.3). It follows from the first and the second equations of system (1.2) that x_ 1 ðtÞ ¼ ax2 ðtÞ r1 x1 ðtÞ bx1 ðtÞ; a1 x_ 2 ðtÞ P bx1 ðtÞ b1 x22 ðtÞ x2 ðtÞ: m
ð2:2Þ
Consider the following auxiliary system u_ 1 ðtÞ ¼ au2 ðtÞ r1 u1 ðtÞ bu1 ðtÞ; a1 u_ 2 ðtÞ ¼ bu1 ðtÞ b1 u22 ðtÞ u2 ðtÞ: m
ð2:3Þ
System (2.3) has two equilibria: E10 ð0; 0Þ, E1 ðu1 ; u2 Þ, where au2 1 ab a1 ; u2 ¼ u1 ¼ : b1 r 1 þ b m r1 þ b It is easy to show that E10 is always a saddle point, E1 is always stable. System (2.3) can be rewritten as a ½u2 ðtÞðu1 ðtÞ u1 Þ þ u1 ðtÞðu2 ðtÞ u2 Þ; u1 b u_ 2 ðtÞ ¼ ½u1 ðtÞðu2 ðtÞ u2 Þ þ u2 ðtÞðu1 ðtÞ u1 Þ b1 u2 ðtÞðu2 ðtÞ u2 Þ: u2
u_ 1 ðtÞ ¼
ð2:4Þ Consider the following Lyapunov function u1 u2 V1 ðtÞ ¼ c1 u1 u1 u1 ln þ c2 u2 u2 u2 ln ; u1 u2 where c1 , c2 are constants to be determined.
ð2:5Þ
734
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
Calculating the derivative of V1 ðtÞ along solutions of (2.4) we derive 2 dV1 ðtÞ X u_ i ðtÞ ¼ ci ðui ðtÞ ui Þ dt ui ðtÞ i¼1 ¼
c1 aðu1 ðtÞ u1 Þ ½u2 ðtÞðu1 ðtÞ u1 Þ þ u1 ðtÞðu2 ðtÞ u2 Þ u1 u1 ðtÞ þ
c2 bðu2 ðtÞ u2 Þ ½u1 ðtÞðu2 ðtÞ u2 Þ þ u2 ðtÞðu1 ðtÞ u1 Þ u2 u2 ðtÞ 2
c2 b1 ðu2 ðtÞ u2 Þ : Set c1 ¼
ð2:6Þ
bu1 =ðau2 Þ,
c2 ¼ 1. It follows from (2.6) that sffiffiffiffiffiffiffiffiffiffi "sffiffiffiffiffiffiffiffiffiffi #2 dV1 ðtÞ b u2 ðtÞ u1 ðtÞ ¼ ðu1 ðtÞ u1 Þ þ ðu2 ðtÞ u2 Þ dt u2 u1 ðtÞ u2 ðtÞ 2
2
b1 ðu2 ðtÞ u2 Þ 6 b1 ðu2 ðtÞ u2 Þ :
ð2:7Þ
Integrating both sides of (2.7) on ½0; t yields Z t 2 V1 ðtÞ þ b1 ðu2 ðsÞ u2 Þ ds 6 V1 ð0Þ; 0 2
which implies ðu2 ðtÞ u2 Þ 2 L1 ½0; 1Þ. In addition, with a similar argument in the proof of Lemma 2.1 we can prove that u1 ðtÞ, u2 ðtÞ are bounded. This, together with (2.3), implies that u2 ðtÞ u2 and u_ 2 ðtÞ are uniformly continuous. Applying BarbalatÕs lemmas (Lemma 1.2.2 and 1.2.3, Gopalsamy [30]), we conclude that lim ðu2 ðtÞ u2 Þ2 ¼ 0;
t!1
which imply that lim u2 ðtÞ ¼ u2 ; t!1
lim u_ 2 ðtÞ ¼ 0;
t!1
lim u1 ðtÞ ¼ u1 :
t!1
Therefore, there exists a T2 P T1 such that if t P T2 , u1 ðtÞ > u1 e=2, u2 ðtÞ > u2 e=2; where e > 0 is chosen sufficiently small. A standard comparison argument shows that e e lim inf x1 ðtÞ P u1 ; lim inf x2 ðtÞ P u2 : t!þ1 t!þ1 2 2 Thus, there is T3 P T2 such that if t P T3 , x1 ðtÞ > u1 e, x2 ðtÞ > u2 e; where e > 0 is chosen sufficiently small. Using a similar argument in the proof of Theorem 2.2 in [13], we can derive from the third equation of system (1.2) that lim inf yðtÞ P t!þ1
u2 ða2 rÞ : 2mr
This completes the proof.
h
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
735
Theorem 2.2. If a2 < r or r1 > r, a1 =m a > r, then system (1.2) is impermanent. Proof. Let ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be a positive solution of system (1.2) with initial conditions (1.3). We first assume a2 < r. It follows from the third equation of system (1.2) that y_ ðtÞ < yðtÞða2 rÞ; which implies yðtÞ < yð0Þeða2 rÞt : Thus, we have limt!1 yðtÞ ¼ 0. In this case, the predator population will go to extinction. If r1 > r, a1 =m a > r, then there exists b > 0 such that a1 ¼ a þ r: mþb Let d ¼ ðx1 ð0Þ þ x2 ð0ÞÞ=yð0Þ < b. We now claim that x1 ðtÞ þ x2 ðtÞ
for all t P 0;
and
limt!þ1 ðx1 ðtÞ þ x2 ðtÞÞ ¼ 0:
Otherwise, there is a first time t1 such that x1 ðt1 Þ þ x2 ðt1 Þ ¼ b; yðt1 Þ
x1 ðtÞ þ x2 ðtÞ
for all t 2 ½0; t1 Þ:
Then for any t 2 ½0; t1 , we derive from the first and the second equations of (1.2) that d a1 x2 ðtÞyðtÞ ðx1 ðtÞ þ x2 ðtÞÞ ¼ ax2 ðtÞ r1 x1 ðtÞ b1 x22 ðtÞ dt myðtÞ þ x2 ðtÞ a1 x2 ðtÞ 6 ax2 ðtÞ r1 x1 ðtÞ m þ x2 ðtÞ=yðtÞ a1 6 r1 x1 ðtÞ þ a x2 ðtÞ 6 rðx1 ðtÞ þ x2 ðtÞÞ; mþb which yields x1 ðtÞ þ x2 ðtÞ 6 ðx1 ð0Þ þ x2 ð0ÞÞert :
ð2:8Þ
On the other hand, it follows from the third equation of system (1.2) that for all t P 0 y_ ðtÞ P ryðtÞ; which yields yðtÞ P yð0Þert :
ð2:9Þ
736
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
Thus, we derive from (2.8) and (2.9) that x1 ðtÞ þ x2 ðtÞ x1 ð0Þ þ x2 ð0Þ 6 ¼d
for all t 2 ½0; t1 ;
which leads to a contradiction. Hence, we obtain x1 ðtÞ þ x2 ðtÞ 6 x1 ð0Þ þ x2 ð0Þert
for all t P 0:
We therefore derive that lim ðx1 ðtÞ þ x2 ðtÞÞ ¼ 0;
t!þ1
which yields lim x1 ðtÞ ¼ lim x2 ðtÞ ¼ 0:
t!þ1
t!þ1
The proof is complete.
h
Theorem 2.3. If r1 > r, a1 =m a > r, then there exists a positive solution ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ of system (1.2) satisfying lim ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ ¼ ð0; 0; 0Þ:
t!þ1
Proof. Let b ¼ a1 =ða þ rÞ m. Assume that ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ is a solution of system (1.2) with ðx1 ð0Þ þ x2 ð0ÞÞ=yð0Þ < b. It follows from the proof of Theorem 2.2 that lim x1 ðtÞ ¼ lim x2 ðtÞ ¼ 0:
t!þ1
t!þ1
Since yðtÞ is positive and bounded, we have 0 6 lim inf yðtÞ 6 lim sup yðtÞ ¼ s < þ1: t!þ1
t!þ1
We now show that limt!þ1 yðtÞ ¼ 0. Otherwise, we have s > 0. Since limt!þ1 x2 ðtÞ ¼ 0, there is a T4 > 0 such that mrs x2 ðtÞ < for all t P T4 : 2a2 On the other hand, by the definition of s, there is a t1 > T4 such that s and y_ ðt1 Þ > 0: ð2:10Þ yðt1 Þ > 2 It follows from (2.10) that a2 x2 ðt1 Þyðt1 Þ ryðt1 Þ > 0; myðt1 Þ þ x2 ðt1 Þ
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
737
which yields x1 ðt1 Þ >
mr mrs yðt1 Þ > ; a2 2a2
which leads to a contradiction. Hence, we must have s ¼ 0, i.e., lim y2 ðtÞ ¼ 0:
t!þ1
The proof is complete.
We now give an example to show the feasibility of Theorem 2.3. Example 1. Consider the following system x_ 1 ðtÞ ¼ 3x2 ðtÞ x1 ðtÞ 2x1 ðtÞ; 6x2 ðtÞyðtÞ ; x_ 2 ðtÞ ¼ 2x1 ðtÞ 2x22 ðtÞ yðtÞ þ x2 ðtÞ 1:5x2 ðtÞ y_ ðtÞ ¼ yðtÞ 0:5 þ : yðtÞ þ x2 ðtÞ
ð2:11Þ
It is easy to show that r1 ¼ 1 > r ¼ 0:5, a1 =m ¼ 6 > a þ r ¼ 3:5. Thus, by Theorem 2.3 we see that there is a positive solution ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ of system (2.11) with ðx1 ð0Þ þ x2 ð0ÞÞ=yð0Þ < b ¼ 1=3 satisfying lim x1 ðtÞ ¼ lim x2 ðtÞ ¼ lim yðtÞ ¼ 0:
t!þ1
t!þ1
t!þ1
Numerical integration can be carried out using standard algorithms. We used the package ode45 in MATLAB. As shown in Fig. 1, numerical 2 x1 x2 y
solution
1.5
1
0.5
0
0.5
0
5
10
15
20
25
time t
Fig. 1. The temporal solution found by numerical integration of system (2.11) with initial conditions ðx1 ð0Þ; x2 ð0Þ; yð0ÞÞ ¼ ð0:25; 0:25; 2Þ.
738
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
simulation also confirms this fact for system (2.11) with initial conditions ðx1 ð0Þ; x2 ð0Þ; yð0ÞÞ ¼ ð0:25; 0:25; 2Þ.
3. Global stability In this section, we are concerned with the global asymptotic stability of nonnegative equilibria of system (1.2). The strategy is to construct appropriate Lyapunov functions. System (1.2) always has equilibria E0 ð0; 0; 0Þ, E1 ðx1 ; x2 ; 0Þ, where x1 ¼
a2 b b1 ðr1 þ bÞ
2
x2 ¼
;
ab : b1 ðr1 þ bÞ
It is easy to show that system (1.2) has a positive equilibrium E ðx1 ; x2 ; y Þ if and only if (H1) and the following hold: (H3) ab=ðr1 þ bÞ > a1 ða2 rÞ=ðma2 Þ,where x1 ¼
ax2 ; r1 þ b
x2 ¼
ab a1 ða2 rÞ ; b1 ðr1 þ bÞ ma2 b1
y ¼
x2 ða2 rÞ : mr
The Jacobian matrix of system (1.2) at E1 ðx1 ; x2 ; 0Þ takes the form of 0
r1 b @ b 0
a 2b1x2 0
1 0 a1 A: a2 r
It is easy to verify that if a < r, E1 is locally asymptotically stable, if a > r, it is locally unstable. The Jacobian matrix of system (1.2) at E ðx1 ; x2 ; y Þ takes the form of 0
r1 b B b @ 0
a 2 2 2b1 x2 ma1 ðy Þ =ðmy þ x2 Þ 2
ma2 ðy Þ =ðmy þ x2 Þ
2
1 0 2 2 a1 ðx2 Þ =ðmy þ x2 Þ C A: 2 ma2 x2 y =ðmy þ x2 Þ
It is easy to show that if ab=ðr1 þ bÞ > 2a1 =m, the positive equilibrium E is locally asymptotically stable. We are now able to state and prove our main result on the global asymptotic stability of the positive equilibrium E ðx1 ; x2 ; y Þ of system (1.2). Theorem 3.1. Let (H1) hold. Assume further that ab=ðr1 þ bÞ > 2a1 =m. Then the positive equilibrium E ðx1 ; x2 ; y Þ of system (1.2) is globally asymptotically stable.
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
739
Proof. Let zðtÞ ¼ ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be any positive solution of system (1.2) with initial conditions (1.3). System (1.2) can be rewritten as a ½x2 ðtÞðx1 ðtÞ x1 Þ þ x1 ðtÞðx2 ðtÞ x2 Þ; x1
x_ 1 ðtÞ ¼
b ½x1 ðtÞðx2 ðtÞ x2 Þ þ x2 ðtÞðx1 ðtÞ x1 Þ x2 þ x2 ðtÞ b1 ðx2 ðtÞ x2 Þ
x_ 2 ðtÞ ¼
ð3:1Þ
þ _ ¼ yðtÞ yðtÞ
a1 y ðx2 ðtÞ x2 Þ a1 x2 ðyðtÞ y Þ ; ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
ma2 ½y ðx2 ðtÞ x2 Þ x2 ðyðtÞ y Þ : ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
Define a Lyapunov function
V ðtÞ ¼
2 X xi y ci xi xi xi ln þ c3 y y y ln ; y xi i¼1
ð3:2Þ
where ci ði ¼ 1; 2; 3Þ are positive constants to be determined. Calculating the derivative of V ðtÞ along solutions of system (3.1), it follows that 2 dV ðtÞ X x_ i ðtÞ y_ ðtÞ ¼ þ c3 ðyðtÞ y Þ ci ðxi ðtÞ xi Þ dt x yðtÞ ðtÞ i i¼1
¼
c1 aðx1 ðtÞ x1 Þ ½x2 ðtÞðx1 ðtÞ x1 Þ þ x1 ðtÞðx2 ðtÞ x2 Þ x1 x1 ðtÞ c2 bðx2 ðtÞ x2 Þ ½x1 ðtÞðx2 ðtÞ x2 Þ þ x2 ðtÞðx1 ðtÞ x1 Þ x2 x2 ðtÞ þ c2 ðx2 ðtÞ x2 Þ b1 ðx2 ðtÞ x2 Þ þ
a1 y ðx2 ðtÞ x2 Þ a1 x2 ðyðtÞ y Þ þ ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ þ c3 ma2 ðyðtÞ y Þ
y ðx2 ðtÞ x2 Þ x2 ðyðtÞ y Þ : ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
ð3:3Þ
740
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
Set c1 ¼ b1 x1 =ða1 x2 Þ, c2 ¼ 1, c3 ¼ a1 x2 =ðma2 y Þ. We derive from (3.3) that sffiffiffiffiffiffiffiffiffiffi "sffiffiffiffiffiffiffiffiffiffi #2 dV ðtÞ b x2 ðtÞ x1 ðtÞ ¼ ðx1 ðtÞ x1 Þ þ ðx2 ðtÞ x2 Þ dt x2 x1 ðtÞ x2 ðtÞ 2
b1 ðx2 ðtÞ x2 Þ2 þ
a1 y ðx2 ðtÞ x2 Þ ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
2
2
a1 ðx2 Þ ðyðtÞ y Þ : y ðmyðtÞ þ x2 ðtÞÞðmy þ x2 Þ
ð3:4Þ
Choose a positive constant e such that ab=ðr1 þ bÞ 2a1 =m b1 e > 0. It follows from the proof of Theorem 2.1 that there is a T > 0 such that if t P T , x2 ðtÞ > u2 e. Thus we derive from (3.4) and Lemma 2.1 that for t P T , 2 2 dV ðtÞ a1 a1 ðx2 Þ ðyðtÞ y Þ 2 6 b1 : ðtÞ x Þ ðx 2 2 dt mðu2 eÞ y Mðm þ 1Þðmy þ x2 Þ ð3:5Þ Integrating both sides of (3.5) on ½0; t gives Z t a1 2 V ðtÞ þ b1 ðx2 ðsÞ x2 Þ ds mðu2 eÞ 0 Z t 2 a1 ðx2 Þ þ ðyðsÞ y Þ2 ds 6 V ð0Þ; y Mðm þ 1Þðmy þ x2 Þ 0 2
2
which implies ðx2 ðtÞ x2 Þ ; ðyðtÞ y Þ 2 L1 ½0; 1Þ. In addition, by Lemma 2.1, x2 ðtÞ, yðtÞ are bounded. This, together with (1.2), implies that x2 ðtÞ x2 , yðtÞ y and x_ 2 ðtÞ are uniformly continuous. Applying BarbalatÕs lemmas (Lemmas 1.2.2 and 1.2.3, Gopalsamy [30]), we conclude that lim ðx2 ðtÞ x2 Þ2 ¼ 0;
t!1
lim x_ 2 ðtÞ ¼ 0;
t!1
lim ðyðtÞ y Þ2 ¼ 0;
t!1
which imply that lim x2 ðtÞ ¼ x2 ;
t!1
lim x1 ðtÞ ¼ x1 ;
t!1
This completes the proof.
lim yðtÞ ¼ y :
t!1
h
In the following we discuss the global stability of the nonnegative equilibrium E1 ðx1 ; x2 ; 0Þ. Theorem 3.2. Let (H2) hold. If a2 < r, then the nonnegative equilibrium E1 ðx1 ; x2 ; 0Þ is globally asymptotically stable.
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
741
Proof. Let ðx1 ðtÞ; x2 ðtÞ; yðtÞÞ be a positive solution of system (1.2) with initial conditions (1.3). By the proof of Theorem 2.1, we see that there is a T3 > 0 such that if t P T3 , x2 ðtÞ > u2 =2. By the proof of Theorem 2.2, we have limt!þ1 yðtÞ ¼ 0. Thus, for any given e 2 ð0; ab=ðr1 þ bÞÞ, there is a T5 P T3 such that if t > T5 , a1 yðtÞ=ðmyðtÞ þ x2 ðtÞÞ < e. We consider the following two auxiliary systems: u_ 1 ðtÞ ¼ au2 ðtÞ r1 u1 ðtÞ bu1 ðtÞ; u_ 2 ðtÞ ¼ bu1 ðtÞ b1 u22 ðtÞ eu2 ðtÞ;
ð3:6Þ
and v_ 1 ðtÞ ¼ av2 ðtÞ r1 v1 ðtÞ bv1 ðtÞ; v_ 2 ðtÞ ¼ bv1 ðtÞ b1 v22 ðtÞ:
ð3:7Þ
Using similar arguments in the proof of Theorem 2.1, we can show that ae e lim u1 ðtÞ ¼ x1 ; lim u2 ðtÞ ¼ x2 ; t!þ1 t!þ1 b1 ðr1 þ bÞ b1 lim v1 ðtÞ ¼ x1 ; lim v2 ðtÞ ¼ x2 : t!þ1
t!þ1
A standard comparison argument shows that ae e x1 6 lim x1 ðtÞ 6 x1 ; x2 6 lim x2 ðtÞ 6 x2 : t!þ1 b1 ðr1 þ bÞ t!þ1 b1 Letting e ! 0; we finally obtain lim x1 ðtÞ ¼ x1 ;
t!þ1
lim x2 ðtÞ ¼ x2 :
t!þ1
This completes the proof.
h
We now give two examples to illustrate the validity of Theorems 2.1, 3.1 and 3.2. Example 2. As an example, we consider the following system: x_ 1 ðtÞ ¼ 4x2 ðtÞ 0:1x1 ðtÞ 2x1 ðtÞ; x2 ðtÞyðtÞ ; x_ 2 ðtÞ ¼ 2x1 ðtÞ x22 ðtÞ 2yðtÞ þ x2 ðtÞ 1:5x2 ðtÞ y_ ðtÞ ¼ yðtÞ 0:5 þ : 2yðtÞ þ x2 ðtÞ
ð3:8Þ
System (3.8) has a unique positive equilibrium E ð2920=441; 73=21; 73=21Þ. It is easy to verify that the coefficients of system (3.8) satisfy the assumptions of Theorem 2.1, 3.1. By Theorem 2.1 we see that system (3.8) is uniformly persistent. By Theorem 3.1 we see that the positive equilibrium E of system (3.8)
742
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744 7
x1 x2 y
6
solution
5
4
3
2
1
0
10
20
30
40
50 time t
60
70
80
90
100
Fig. 2. The temporal solution found by numerical integration of system (3.8) with initial conditions ðx1 ð0Þ; x2 ð0Þ; yð0ÞÞ ¼ ð1; 1; 1Þ.
is globally asymptotically stable. Numerical simulation also suggests these facts (see, Fig. 2). Example 3. Consider the following ratio-dependent predator–prey model x_ 1 ðtÞ ¼ 6x2 ðtÞ x1 ðtÞ 3x1 ðtÞ; x2 ðtÞyðtÞ ; x_ 2 ðtÞ ¼ 3x1 ðtÞ 2x22 ðtÞ yðtÞ þ x2 ðtÞ x2 ðtÞ y_ ðtÞ ¼ yðtÞ 2 þ : yðtÞ þ x2 ðtÞ
ð3:9Þ
3.5 x1 x2 y
3
solution
2.5 2 1.5 1 0.5 0
0
5
10
15
20
25 time t
30
35
40
45
50
Fig. 3. The temporal solution found by numerical integration of system (3.9) with initial conditions ðx1 ð0Þ; x2 ð0Þ; yð0ÞÞ ¼ ð1; 1; 1Þ.
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
743
System (3.9) has a nonnegative equilibrium E1 ð27=8; 9=4; 0Þ. It is easy to verify that the coefficients of system (3.9) satisfy the assumptions of Theorem 3.2. By Theorem 3.2 we see that E1 is globally asymptotically stable. Numerical simulation also confirms this fact (see, Fig. 3). 4. Discussion In this paper, we have discussed the effect of stage structure for prey on the dynamics of a ratio-dependent predator–prey model. By some comparison arguments, sufficient conditions are derived for the uniform persistence and impermanence of system (1.2). By constructing suitable Lyapunov functions, sufficient conditions are obtained for the global asymptotic stability of nonnegative equilibria of system (1.2). By Theorem 2.1, we see that system (1.2) will be uniformly persistent if the birth rate into the immature prey population, the rate of immature prey becoming mature prey, and the conversion rate and the half saturation rate of the predator are high and the capturing rate of the predator and the death rates of both the immature prey and the predator are low enough satisfying (H1)–(H2). By Theorem 2.2 we see that if the conversion rate of the predator is less than its death rate, then the predator will go to extinction, system (1.2) will not be persistent. Theorem 3.1 shows the global stability of the positive equilibrium of system (1.2) with somewhat stronger assumptions than those in Theorem 2.1 on the uniform persistence of system (1.2). We would like to mention here that we are unable to show system (1.2) admits limit cycles when E exists and is unstable. This is known to be true for system (1.1) (see, Hsu et al. [10]). We leave this for future work. References [1] R. Arditi, H. Saiah, Empirical evidence of the role of heterogeneity in ratio-dependent consumption, Ecology 73 (1992) 1544–1551. [2] R. Arditi, L.R. Ginzburg, H.R. Akcakaya, Variation in plankton densities among lakes: a case for ratio-dependent models, American Nat. 138 (1991) 1287–1296. [3] R. Arditi, N. Perrin, H. Saiah, Functional response and heterogeneities: an experiment test with cladocerans, OIKOS 60 (1991) 69–75. [4] A.P. Gutierrez, The physiological basis of ratio-dependent predator–prey theory: a methbolic pool model of NicholsonÕs blowflies as an example, Ecology 73 (1992) 1552–1563. [5] R. Arditi, L.R. Ginzburg, Coupling in predator–prey dynamics: ratio-dependence, J. Theor. Biol. 139 (1989) 311–326. [6] I. Hanski, The functional response of predator: worries bout scale, TREE 6 (1991) 141–142. [7] E. Beretta, Y. Kuang, Global analyses in some delayed ratio-dependent predator–prey systems, Nonlinear Anal. 32 (1998) 381–408. [8] F. Berezovskaya, G. Karev, R. Arditi, Parametric analysis of the ratio-dependent predator– prey model, J. Math. Biol. 43 (2001) 221–246.
744
R. Xu et al. / Appl. Math. Comput. 158 (2004) 729–744
[9] M. Fan, K. Wang, Periodicity in a delayed ratio-dependent predator–prey system, J. Math. Anal. Appl. 262 (2001) 179–190. [10] S.B. Hsu, T.W. Hwang, Y. Kuang, Global analysis of the Michaelis–Menten ratio-dependent predator–prey system, J. Math. Biol. 42 (2001) 489–506. [11] S.B. Hsu, T.W. Hwang, Y. Kuang, Rich dynamics of a ratio-dependent one prey two predator model, J. Math. Biol. 43 (2001) 377–396. [12] C. Jost, O. Arino, R. Arditi, About deterministic extinction in ratio-dependent predator–prey models, Bull. Math. Biol. 61 (1999) 19–32. [13] Y. Kuang, E. Beretta, Global qualitative analysis of a ratio-dependent predator–prey system, J. Math. Biol. 36 (1998) 389–406. [14] D. Xiao, S. Ruan, Global dynamics of a ratio-dependent predator–prey system, J. Math. Biol. 43 (2001) 268–290. [15] D. Xiao, S. Ruan, Stability and bifurcation in a delayed ratio-dependent predator–prey system, Proc. Edin. Math. Soc. 45 (2002) 205–220. [16] R. Xu, L. Chen, Persistence and stability for a two-species ratio-dependent predator–prey system with time delay in a two-patch environment, Comput. Math. Appl. 40 (2000) 577–588. [17] R. Xu, L. Chen, Persistence and global stability for n-species ratio-dependent predator–prey system with time delays, J. Math. Anal. Appl. 275 (2002) 27–43. [18] W.G. Aiello, H.I. Freedman, A time delay model of single-species growth with stage structure, Math. Biosci. 101 (1990) 139–153. [19] W.G. Aiello, H.I. Freedman, J. Wu, Analysis of a model representing stage-structured population growth with state-dependent time delay, SIAM J. Appl. Math. 52 (1992) 855–869. [20] F. Brauer, Z. Ma, Stability of stage-structured population models, J. Math. Anal. Appl. 126 (1987) 301–315. [21] H.I. Freedman, J. Wu, Persistence and global asymptotic stability of single species dispersal models with stage structure, Quart. Appl. Math. 49 (1991) 351–371. [22] Y. Kuang, J.W.H. So, Analysis of a delayed two-stage population with space-limited recruitment, SIAM J. Appl. Math. 55 (1995) 1675–1695. [23] S. Liu, L. Chen, Z. Liu, Extinction and permanence in nonautonomous competitive system with stage structure, J. Math. Anal. Appl. 274 (2002) 667–684. [24] S. Liu, L. Chen, R. Agarwal, Recent progress on stage-structured population dynamics, Math. Comput. Model. 36 (2002) 1319–1360. [25] K.G. Magnusson, Destabilizing effect of cannibalism on a structured predator–prey system, Math. Biosci. 155 (1999) 61–75. [26] X. Song, L. Chen, Optimal harvesting and stability for a two-species competitive system with stage structure, Math. Biosci. 170 (2001) 173–186. [27] W. Wang, L. Chen, A predator–prey system with stage structure for predator, Comput. Math. Appl. 33 (1997) 83–91. [28] W. Wang, G. Mulone, F. Salemi, V. Salone, Permanence and stability of a stage-structured predator–prey model, J. Math. Anal. Appl. 262 (2001) 499–528. [29] X. Zhang, L. Chen, The stage-structured predator–prey model and optimal harvesting policy, Math. Biosci. 168 (2000) 201–210. [30] K. Gopalsamy, Stability and Oscillations in Delay Differential Equations of Population Dynamics, Kluwer Academic, Dordrecht/ Norwell, MA, 1992.