Effect of Time-delay on a Ratio-dependent Food Chain Model

This paper aims to study the effect of time-delay on a tritrophic food chain model with Michaelis-Menten type ratio-dependent functional responses. Boundedness of the time-delayed system is established. A simple criterion for deterministic extinction is derived. It has been shown that the time-delay may introduce instability in the system through Hopf bifurcation. Computer simulations are carried out to explain the analytical findings. It is discussed how these ideas illuminate some of the observed properties of real populations in the field, and explores practical implications.


Introduction
The use of mathematics in the modelling of biological phenomena has become more prevalent in recent years.The gradual invasion of the field of biology by mathematicians is already yielding a number of benefits for both the fields.A number of realistic mathematical models now exists but the analytical results about the dynamical behaviour of such models is still largely lacking.This is specially true in case of three-species models.In fact, after getting started with the pioneering works of Lotka [1] and Volterra [2], two-species models have so far dominated the field of Biomathematics.But the limited "caricature" of ecological systems by two species can account for only a small number of the phenomena that are commonly exhibited in nature.To the best of our knowledge, it was only in the late seventies that some interest in the mathematics of three-species food chain models (composed of prey, predator and superpredator) emerged (Freedman and Waltman [3], Gilpin [4], Gard [5]).Subsequently, some golden rules are derived in the the theory of three-species models (for detailed references, see Maiti and Samanta [6]).
Actually, the most crucial element in modelling of tritrophic systems is the choice of functional responses.Most of the three-species models are designed with Lotka-Volterra or Holling type prey-dependent functional responses but such functional responses have faced a number of questions since late eighties.A large number of papers are now available to get into the details of those heated phases of transition.Beside this, in many field studies, certain three-species communities have become the focus of considerable attention.In the following, we mention some of them.It will indicate the relevance of theoretical studies of three-species systems.
It is well known that approximately one third of world food production is lost due to pests.The farmers have to suffer from heavy economic losses due to these pests.Usually, pesticides are used to control pests but the long list of side-effects of pesticides forces the scientists to find alternative methods for pest control.For this reason, the method of bio-control of pests is getting much more importance in almost all developing countries.In this method, a pest is controlled (destroyed) by introducing a parasitoid or predator, which is natural enemy to the concerned pest.For example, the most important pests like Spodoptera sp. of maize (Zea mays) is controlled by using parasitoids like Cotesia marginiventris (see Ashley [7], Jalali et al. [8], Turlings and Fritzsche [9], and Fritzsche Hoballah and Turlings [10] for experimental evidences).To reduce the indiscriminate use of pesticides, recently tea scientists are using predators or pathogens to control the pests of tea (Das and Barua [11], Das et al. [12], Kabir [13], Mochizuki [14]).Similar studies can also be found in Gomez and Zamora [15] and Van Loon et al. [16].The literature abounds with many such evidences.Some other three-species systems have also got the attention of the scientists.For example, in waste treatment process, the bacteria lives on the waste (or nutrient) while other organisms as ciliates feed on the bacteria (Li and Kuang [17]).Thus three-species systems like nutrient-bacteria-ciliate, plant-herbivoreparasitoid, plant-pest-predator, et cetera are emerging in different branches of biology in their own right.These are the examples of three-species "food chains".That is, behaviour of the entire community arise from the coupling of these interacting species, where the species at the third trophic level prey on the species at the second trophic level and the species at the second trophic level feed on the species at the first trophic level (see Fig. 1).A distinct feature of these food chains is the so called domino effect: if one species dies out, all the species at higher trophic level die out as well.Of course, there are many other types of three-species systems.In this paper, we are not interested in those systems.

Prey
↓ Predator ↓ Top-predator Fig. 1.The feeding relationship in a tritrophic food chain.
We have already mentioned that the most important element in predator-prey interactions is the functional response or trophic function.In recent times, the ratio-dependent functional response has become the focus of considerable attention in ecological literature.A number of research papers have analyzed the merits of ratio-dependent functional response in comparison with others (see Arditi et al. [18], Hsu et al. [19], and references cited there in).Further, predator-prey models with such ratio-dependent functional response are strongly supported by numerous field and laboratory experiments (Arditi and Ginzburg [20], Arditi and Berryman [21], Arditi et al. [22,23], Hanski [24], Arditi and Saiah [25], Gutierrez [26], Blaine and DeAngelis [27], Poggiale et al. [28], Bernstein et al. [29], Cosner et al. [30], Arditi et al. [31]).Also, there is a growing evidences that in some situations, especially when predators have to search for food (and therefore have to compete for food), a more suitable functional response should be ratio-dependent (Arditi et al. [31], Cosner et al. [30], Hsu et al. [32,33], Xiao and Ruan [34]).
We now make a general comment on the realistic modelling of ecosystems.Usually, a model in theoretical ecology is considered under the framework of a set of ordinary differential equations.It is a fact that ordinary differential equations have long played an important role in the history of theoretical ecology.However, they are generally the first approximations of the considered real system.More realistic models should include some of the past states of these systems; that is, ideally, a real system should be modelled by differential equations with time-delays.Now it is well understood that many of the processes, both natural and manmade, in biology, medicine, et cetera, involve time-delays.Time-delays occur so often, in almost every situation, that to ignore them is to ignore reality.Kuang [35] mentioned that animals must take time to digest their food before further activities and responses take place and hence any model of species dynamics without delays is an approximation at best.Now it is beyond doubt that in an improved analysis, the effect of time-delay due to the time required in going from egg stage to the adult stage, gestation period, et cetera, has to be taken into account.Detailed arguments on importance and usefulness of time-delays in realistic models may be found in the classical books of Macdonald [36], Gopalsamy [37] and Kuang [35].
In recent years, theoretical ecologists have laid down immense importance on the effect of time-delay on realistic models.Still, there is no denying that a lack of analytic results on time-delayed models persists.In fact, after the recent developments in the field of computers, it has become a strong tendency of the theoretical ecologists to study the dynamical bahaviours of models numerically, without any analytic guide.No doubt, numerical simulations provide a basic understanding of these systems, even when analytic results are unavailable.However, analytic results are important because they can show the dynamical behaviour of a class of models, rather than a particular model.This is particularly important in ecology, where confidence in exact form of any particular model is weak.Not that we neglect the increased computation capacity and speed of present day computers, far from it, but it is beyond doubt that a better analytic understanding of the models would make the use of numerics even more useful.
The present paper is designed to study the effect of discrete time-delay on a tritrophic food chain model with ratio-dependent functional responses.The rest of the paper is structured as follows.In Section 2, we have introduced the basic deterministic model.Section 3 is equipped with the dynamical behaviours of the system.In Section 4, our analytical results are validated through numerical simulation.Section 5 contains the general discussions of the paper and biological implications of our mathematical findings.

The model
In the present paper we are concerned with a tritrophic food chain model with ratiodependent functional responses.Our basic model is actually formulated by Hsu et al. [19].The model equations are where X, Y, Z denote the population densities of prey, predator and top predator, respectively.Here r > 0 and K > 0 are the prey intrinsic growth rate and carrying capacity, respectively.For i = 1, 2; m i , a i , d i , η i are positive constants denoting maximal predator growth rates, half saturation constants, predator death rates, and yield constants, respectively.
The model we have just specified has ten parameters, which makes analysis difficult.To reduce the number of parameters and to determine which combinations of parameters control the behavior of the system, we non-dimensionalize the system (1).We choose Then the system (1) takes the form (after some simplification) where b = m1 η1a1r , m = m1 r , p = d1 r , c = m2 η2a2r , e = m2 r , q = d2 r .We have already discussed the usefulness of time-delay in realistic modelling of ecosystems.Here we consider the following modification of the model ( 2) incorporating discrete time-delay in it.
, z(0) > 0. ( The delay τ in (3) may be regarded as a gestation period or reaction time of the toppredator Z.We think it worth mentioning that the necessity of such models has already been emerged from the field experiment also.Reeve [38] conducted experiments to determine the nature of the functional response of clerid beetle Thanasimus dubius on bark beetle Dendroctonus frontalis during attack of the host tree.Analyzing the experimental findings, he has suggested that a model including ratio-dependence and time-delay for T. dubius and D. frontalis interaction may determine the net stabilizing or destabilizing effect on D. frontalis dynamics.
3 Dynamical behaviours 3.1 Dynamics when τ = 0 Hsu et al. [19] have already studied some important dynamical behaviours of the system (3) when τ = 0. We mention here some of their results.
Theorem 1.All solutions of the system (3) that start in R 3 + are uniformly bounded when τ = 0. Theorem 2. The system (3) always possesses the trivial equilibrium E 0 (0, 0, 0) and equilibrium E 1 (1, 0, 0).The third boundary equilibrium point E 2 (x, ŷ, 0) exists if and only if m > p and m − bm + bp > 0. When this condition is satisfied, x, ŷ are given by Theorem 3. The interior equilibrium E * (x * , y * , z * ) of the system (3) exists if and only if the following conditions are satisfied: Furthermore x * , y * , z * are given by x * , and z * = (e − q)y * q , where A = m c{(e−q)/e}+p .
On deterministic extinction, we have the following theorems : Theorem 6.If e > q and 0 < A ≤ 1, then lim t→∞ y(t) = 0 and lim t→∞ z(t Remark 1.The above theorems give conditions for the extinction of the top predator, the extinction of both middle and top predators (but not the prey species x, and total extinction of all the three species.The Theorem 6 states that if the middle predator is a high capacity and aggressive consumer (characterized by large values of b) and there is a shortage of prey to begin with, then all three species will go extinct.The last part of Theorem 6 suggests that if middle predator is a low capacity consumer, then prey species will persist.Now we study the stability of the most important equilibrium point E * (x * , y * , z * ).The variational matrix of ( 3) at E * is given by The characteristic equation is where Then we have the following theorem on local stability of E * .

Theorem 7. E * is locally asymptotically stable if and only if
The theorem directly follows from Routh Hurwitz criterion (as A 3 > 0).A sufficient condition for local asymptotic stability of E * is the following: Hence, the theorem follows from Routh Hurwitz criterion.

Dynamics when τ = 0
When τ = 0, the system (3) has the same equilibria as in the previous case.The main purpose of this section is to study the stability behavior of E * (x * , y * , z * ) in the presence of discrete delay (τ = 0).Before doing that we present some important results on boundedness and deterministic extinction.Theorem 9.All solutions of the time-delayed system (3) that start in R 3 + are uniformly bounded.
Proof.Let (x(t), y(t), z(t)) be any solution of the system (3) that start in R 3 + .Since Applying a theorem on differential ineqalities (Birkhoff and Rota, [39]), we obtain and for t → ∞, Therefore, it is possible to find two positive numbers κ and t ′ such that y(t) ≤ κ for t > t ′ .If e ≤ q, then dz dt ≤ z(e − q) = −kz, where k = q − e.
A standard comparison argument shows that lim sup t→∞ z(t) ≤ (e − q)κ q e (e−q)τ .
Hence the theorem.
The following theorem gives a criterion for deterministic extinction of the timedelayed system (3).So, for t > τ , y ≤ y(0)e −lt .Therefore, It is already shown in the proof of Theorem 9 that if e < q then lim t→∞ z = 0. Hence the theorem is established.

Remark 2.
The above theorem indicates that is the death rates exceed the corresponding maximal growth rates (i.e.d i > m i , i = 1, 2), then the predator and top-predator both will die out but the prey persists.It is worth mentioned that for the time-delayed system, the result of Theorem 5 hold as well.From this, we may roughly say that if the maximal predator growth rate is high enough, then there might be a case of total extinction.
We now study the stability behaviour of E * (x * , y * , z * ) for the system (3).We linearize the system (3) by using the following transformations : Then the linear system is given by where , and and all other a ij = 0, and other b ij = 0.
A necessary condition for a stability change of E * is that the characteristic equation ( 5) has purely imaginary solutions.Hence to obtain the stability criterion, we set ξ = 0 in ( 6) and ( 7) to obtain and Eliminating τ by squaring and adding ( 8) and ( 9), we get the equation for determining η as where Substituting η 2 = σ in (10), we get a cubic equation given by We notice that F is continuous everywhere with F (0) < 0 and lim σ→∞ F (σ) = ∞.Therefore, the cubic (11) always has at least one positive root.Consequently, the stability criteria of the system for τ = 0 will not necessarily ensure the stability of the system for τ = 0. Without loss of generality, we assume that equation ( 11) has three positive roots, denoted by σ 1 , σ 2 and σ 3 , respectively.Then equation ( 10) has three positive roots, say Then ±iη k is a pair of purely imaginary roots of the equation ( 5) with τ = τ

Thus, we can define
and η 0 = η k0 .Now we state a lemma which was proved by Ruan and Wei [40].
Then we have the following theorem on the stability of E * .Theorem 11.Suppose that E * exists with A 1 > 0 and ∆ > 0. Then E * is locally asymptotically stable for all τ ∈ [0, τ * ).
Proof.Since E * exists with A 1 > 0 and ∆ > 0, therefore all the roots of the equation ( 5) have negative real parts when τ = 0. Now it is easy to see that when τ = τ (j) k , k = 1, 2, 3; j = 0, 1, 2 . .., equation ( 5) has no root with zero real part and τ * is the minimum value of τ , so that equation ( 5) has purely imaginary roots.Hence applying the above lemma, we conclude that all roots of the equation ( 5) have negative real parts when τ ∈ [0, τ * ).This proves the theorem.Now we show the existence of Hopf bifurcation near E * by taking τ as bifurcation parameter.

Numerical simulation
In this section we present computer simulation of some solutions of the system (3) using MATLAB.The numerical solutions of the ordinary differential equation are obtained by using the fourth order Runge-Kutta method.Beside verification of our analytical findings, these numerical solutions are very important from practical point of view.
Thus using the time-delay as control, it is possible to break the stable behaviour of the system and drive it to an unstable state.Also it is possible to keep the populations at a desired level using the above control.

Concluding remarks
In this paper, we have studied the effect of discrete time-delay on a ratio-dependent tritrophic food chain model.It is shown that the time-delayed system is uniformly bounded, which, in turn, implies that the system is biologically well behaved.It has long been recognized that most of the studies of continuous time deterministic models reveal two basic patterns: approach to an equilibrium or to a limit cycle.The basic rationale behind such type of analysis was the implicit assumption that most food chains we observe in nature correspond to stable equilibria of the model.From this viewpoint, we have presented the stability and bifurcation analysis of the most important equilibrium point E * .The nonlinear differential equation ( 3) may be looked upon as the mathematical model for tea plant (Camellia sinensis L.)-pest (e.g.Looper caterpillar)-beneficial predator (natural enemy of the pest) (e.g.Bacillus thuringiensis or Sarcophaga sp.) (Kabir,[13]).Then we observe that the size of the tea plant (x) in the absence and presence of beneficial predator (z) are x and x * , respectively so that x * − x > 0. Also y * − ŷ < 0 where ŷ and y * are respectively the size of the pest population in the absence and presence of beneficial predators.This implies that predator attack of pests enhance fitness of tea plants and cause depression for the pests.This gives a strong theoretical support to the approach of bio-control of pests.A host of experiments had been laid out by tea scientists to study the possibility of implementing the technique of bio-control of pests in tea.In an experiment conducted by Das and Barua [11], Bacillus thuringiensis was found to be very effective against looper caterpillars.Mochizuki [14] found that at a few tea plantations in central Japan, populations of Amblyseius womersleyi were suppressing Kanzawa spider mite Tetranychus kanzawai even after they were sprayed with synthetic pyrethroid pesticides.Also the above result x * − x > 0 is in good agreement with the experimental findings of Gomez and Zamora [15] (on Hormathophylla spinosa-Ceutorhynchus sp.chalcid parasitoid interactions), Van Loon et al. [16] (on Arabidopsis thaliana-Pieris rapae-Cotesia rubecula interactions) and Fritzsche-Haballash and Turlings [10] (on Zea mays-Spodoptera littoralis-Cotesia marginiventris interactions).
It is mentioned by several researchers that the effect of time-delay must be taken into account in order to have a biologically useful mathematical model (Macdonald [36], Gopalsamy [37], Kuang [35]).In the model (3) the delay may be looked upon as the gestation period or reaction time of the top-predator.If the system (3) is the model for plant-pest-parasitoid interaction, then Theorem 10 implies that the pest will be completely eradicated when the death rates exceed the corresponding maximal growth rates.Then a rigorous analysis leads us to the Theorem 11 and 12 which mention that the stability criteria in absence of delay is no longer enough to guarantee the stability in presence of delay, rather there is a critical value of the delay such that the system is stable if the delay lies below the critical value and become unstable when the delay exceeds it.Such regulatory impact of the time-delay is also illustrated through computer simulation.
It is well known that natural populations of plants and animals neither increase indefinitely to blanket the world nor become extinct (except in some rare cases due to some rare reasons).Hence, in practice, we often want to reduce the predator to an acceptable level in finite time.In order to accomplish this we strongly suggest that in realistic field situations (where effect of time-delays can never be violated), the parameters of the system should be regulated in such a way that the conditions of Theorem 11 are satisfied.

Theorem 8 .
If m > b and e > c, then E * is locally asymptotically stable.Proof.We notice that (i) m > b and e > c ⇒ A 1 > 0, (ii) A 3 > 0 for all values of the parameters, (iii) m > b and e > c ⇒ ∆ = A 1 A 2 − A 3 > 0.

Fig. 4 .
Fig. 4.Here all other parameter values are same as in Fig. 3 except τ = 2.4 > τ * .(a) Phase portrait of the system (3) showing a limit cycle which grows out of E * ; (b)-(d) Oscillations of x, y, z populations respectively in finite time.