On the Dynamical Behavior of Toxic-Phytoplankton-Zooplankton Model with Delay

A toxin producing phytoplankton-zooplankton model with inhibitory exponential substrate and time delay has been formulated and analyzed. Since the liberation of toxic substances by phytoplankton species is not an instantaneous process but is mediated by some time lag required for maturity of the species and the zooplankton mortality due to the toxic phytoplankton bloom occurs after some time laps of the bloom of toxic phytoplankton, we induced a discrete time delay to both of the consume response function and distribution of toxic substance term. Furthermore, based on the fact that the predation rate decreases at large toxicphytoplankton density, the system is modelled via a Tissiet type functional response. We study the dynamical behaviour and investigate the conditions to guarantee the coexistence of two species. Analytical methods and numerical simulations are used to obtain information about the qualitative behaviour of the models.


Introduction
Phytoplankton are one of the most important components of the marine ecosystem.They not only form a basis for all aquatic food chains but also perform a very useful service by producing a huge amount of oxygen for human and other living animals after absorbing carbon-dioxide from surrounding environments [1].Zooplankton are microscopic animals that eat other plankton and serve as a most favorable food source for fish and other aquatic animals.During the recent years, many authors have studied the system between zooplankton and phytoplankton.Authors in [2] have dealt with a nutrient-plankton model in an aquatic environment in the context of phytoplankton bloom.In [3] the effect of seasonality and periodicity on plankton dynamics is investigated.In [4], two plankton ecosystem models with explicit representation of viruses and virally infected phytoplankton are presented.
The most common features of the phytoplankton population is rapid increase of biomass due to rapid cell proliferation and almost equally rapid decrease in populations, separated by some fixed time period.This type of rapid change in phytoplankton population density is known as "bloom" [5].
Due to the accumulation of high biomass or to the presence of toxicity, some of these blooms, more adequately called "harmful algal blooms" [6], are noxious to marine ecosystems or to human health and can produce great socioeconomic damage.There has been a global increase in harmful plankton blooms in last two decades [7][8][9].
Because of the difficulty in measuring plankton biomass, mathematical modeling of plankton population is an important alternative method of improving our knowledge of the physical and biological processes relating to plankton ecology [10][11][12][13][14][15].In [16], nutrient-plankton-zooplankton interaction models with a toxic substance which inhibits either the growth of phytoplankton, zooplankton, or both trophic levels are proposed and studied.In [17], the authors have constructed a mathematical model for describing the interaction between a nontoxic and toxic phytoplankton with a single nutrient.
Based on the fact that (1) the release of toxin from phytoplankton species is not an instantaneous process but is mediated by some time lag required for maturity of the species and (2) the zooplankton may die after some time lapse of the bloom of toxic phytoplankton (see http://www.mote.org/,http://www.mdsg.umd.edu/),models incorporating time delay in diverse biological models are extensively reviewed by Beretta and Kuang [18], Gopalsamy [19], Cooke and Grossman [20], and Cushing [21].The discrete time delay has potential to change the qualitative behavior of the dynamical systems [22][23][24][25][26][27][28][29].Chattopadhyay et al. proposed a delay model incorporating time lag in toxin liberation by phytoplankton to avoid predation by zooplankton [28].In [28], the authors introduced distribution delay and discrete to toxin liberation term.Due to discrete time delay in toxin liberation, the local existence of periodic solution through Hopf bifurcation has been obtained in [28].
In [5,10,[28][29][30][31], the following plausible toxic-phytoplankton-zooplankton system has been studied: where  = () is the density of the toxin-producing phytoplankton (TPP) population and  = () is the density of zooplankton population at any instant of time .In model (1),  is the intrinsic growth rate and  is the environmental carrying capacity of TPP population.The term () describes the functional response for the grazing of phytoplankton by zooplankton, and  is the maximum uptake rate for zooplankton. denotes the ratio of biomass conversion and  is the natural death rate of zooplankton.The function () represents the distribution of toxic substance which ultimately contributes to the death of zooplankton populations, and the parameter  denotes the rate of toxic substances produced by per unit biomass of phytoplankton.Model (1) has been studied for the following cases: (1) when () is linear but () is Holling type II [10,28,29] or Holling type III [10]; (2) when () and () are both Holling type II [5,10] or Holling type III [10]; (3) when () is Holling type II while () is Holling type III [10].
In above cases, () and () are both increasing functions of  over the entire interval 0 <  < ∞.However, in some cases, very high substrate concentrations in the lakes actually inhabit the growth of phytoplankton cells.Moreover, with the substrate concentrations increasing unlimitedly, some kind of microorganism will die eventually [32].To describe the above phenomenon accurately, we consider () from a different point of view.Adopting the idea used in [32], we assume that there exists a constant 0 <  * < ∞ such that () is increasing over the interval 0 ≤  <  * and is decreasing on the interval  * <  < ∞.More precisely, we use the so-called Tissiet functional response of the form of () =  −/ /( + ) (see, e.g., [32]).This type of functional response takes care of the fact that the predation rate decreases at large toxic-phytoplankton density.
The remainder of this paper is organized as follows.The model is described in Section 2. In Section 3, we state and prove the positivity and boundedness of the solutions.Then, in Section 4, equilibria, their existence, and local asymptotic stability are considered.Under the aids of numerical simulation method, we further analyse the model and determine if there is a parameter range for the delay parameter , where more complicated dynamics occur.In Section 5, the permanence of the system is discussed by some analytic techniques on limit sets of differential dynamical systems.Finally, a brief discussion is presented in Section 6.

State of the Model
In a real ecological context, the interaction between phytoplankton and zooplankton will not be essentially instantaneous.Instead, the response of zooplankton to contacts with phytoplankton is likely to be delayed due to a gestation period.Another fact, during the interaction between phytoplankton and zooplankton, is that the liberation of toxic substances by phytoplankton must be mediated by some time lag which is required for the maturity of toxic-phytoplankton.Let  1 stand for the time delay in conversion of food to viable biomass for the species, and  2 is the discrete time period required for the maturity of phytoplankton cells to liberate toxic substances.Based on model 3 in [5], we intend to study a model system with the assumption that () and () are described by same type of function, namely, Tissiet functional response.Then we arrive at the following model: where the parameter  is the intrinsic growth rate and  is the environmental carrying capacity of TPP population.The constant  is the maximum per capita grazing rate,  denotes the ratio of biomass conversion,  denotes the rate of toxic substances produced by per unit biomass,  is the natural death rate of the zooplankton.All the parameters in system (2) are positive constants with their usual ecological meanings.
For the sake of simplicity of mathematical analysis, in this paper, we consider model (2) for the special case: the time delay in conversion of food to viable biomass for the species equal to the discrete time period required for the maturity of phytoplankton cells to liberate toxic substances; that is,  1 =  2 = .
By performing the following scaling for model (2): we obtain a dimensionless system in the state variables  * ,  * , which can be written, by removing the stars, in the form with initial conditions: where ((), ()) ∈ ([−, 0],  2 + ).

Positive and Boundedness
In this section, we consider the positive and boundedness of the solutions of model ( 4) with initial condition (5).One has the following theorem.
It is easy to prove that if there is a  2 > 0 such that ( 2 ) =  0 , then for all  >  2 one has () >  0 .Which implies that () ≥  0 for all  ≥ 0. This completes the proof of conclusion (c) in Theorem 1.
Remark.The conclusions (b) and (c) in Theorem 1 indicate that if the ratio of biomass conversion is less than certain value, then phytoplankton population will be persistent.
Theorem 2. The following statements hold.

Local Stability and Hopf Bifurcation.
Let  * ( * ,  * ) denote any one of the equilibrium points of system (4).Linearizing (4) about  * we obtain In what follows, existence of the interior equilibria, dynamical properties of the zooplankton-free equilibrium, and the interior equilibrium are investigated.

Permanence
In this section, we will use the same methods as [35] to prove the permanence of system (4).
Proof.From Theorem 1, we see that it is enough to consider the solution ((), ()) with initial condition  = (, ) ∈ Γ. Theorem 1 implies that the omega limit set () of ((), ()) ( ≥ 0) is nonempty, compact, and invariant and () ⊂ Γ.It follows from definition of permanence and Theorem 1, that we only need to show lim inf Here  is some positive constant which dose not depend on the initial function .

Discussion
In a real ecological context, the interaction between phytoplankton and zooplankton will not be essentially instantaneous.Instead, the response of zooplankton to contact with phytoplankton is likely to be delayed due to a gestation period.Another fact, during the interaction between phytoplankton and zooplankton, is that the liberation of toxic substances by phytoplankton must be mediated by some time lag which is required for the maturity of toxic-phytoplankton.And in some cases, very high substrate concentrations in the lakes actually inhabit the growth of phytoplankton cells.Moreover, with the substrate concentrations increasing unlimitedly, some kind of microorganism will die eventually [32].Based on the above fact, in this paper, we introduce time delay to a phytoplankton-zooplankton interaction model with exponential substrate uptake and exponential distribution of toxic substance term.The model (4) accounts for some natural phenomenon.By using comparison principle for functional differential equations and traditional analysis technique for transcendental equations [34], we give a detailed analysis on boundedness of solutions of system (4) and local asymptotic stability of the equilibria of system (4).Our results show that time delay is factually harmless for the local asymptotic stability of the zooplankton free equilibrium of ( 4), but it is not always harmless for the positive equilibrium; that is to say, because of the time delay the positive equilibrium becomes unstable (Theorem 4).Based on some known techniques on limit sets of differential dynamical systems, we show that, for any time delay, the phytoplanktonzooplankton interaction model is permanent if and only if one positive equilibrium exists.
Under the aids of the numerical simulation we further investigate the delayed model system (4).Figures 1-6 illustrate that the delay  plays crucial role in determining the asymptotic behavior of solutions of model (4).The stability or oscillatory coexistence depends upon not only the parametric restriction but also upon the gestation delay (liberation delay).Our system exhibits stable behavior when  <  0 , where  0 is the threshold value of the parameter .This threshold value of  0 is determined by our numerical simulations.This value of  0 is 3.5 for case (i) of Theorem 4 and is 5.1 for case (ii).These stable solutions are shown in Figures 2(a), 2(b), and 4. When  crosses  0 , there is a delay induced instability demonstrated in Figures 2(c) and 5 for  >  0 .Stable oscillations appear when  =  0 and the Hopf bifurcation periodic solution occurs (see Figures 3 and 5).Thus, there is a range of gestation delay (liberation delay) which initially imparts stability, then induces instability, and ultimately leads to periodic behavior.