Impact of fear in a prey-predator system with herd behaviour

: Fear of predation plays an important role in the growth of a prey species in a prey-predator system. In this work, a two-species model is formulated where the prey species move in a herd to protect themselves and so it acts as a defense strategy. The birth rate of the prey here is affected due to fear of being attacked by predators and so, is considered as a decreasing function. Moreover, there is another fear term in the death rate of the prey population to emphasize the fact that the prey may die out of fear of predator too. But, in this model, the function characterizing the fear effect in the death of prey is assumed in such a way that it is increased only up to a certain level. The results show that the system performs oscillating behavior when the fear coefficient implemented in the birth of prey is considered in a small amount but it changes its dynamics through Hopf bifurcation and becomes stable for a higher value of the coefficient. Regulating the fear terms ultimately makes an impact on the growth of the predator population as the predator is taken as a specialist predator here. The increasing value of the fear terms (either implemented in birth or death of prey) decrease the count of the predator population with time. Also, the fear implemented in the birth rate of prey makes a higher impact on the growth of the predator population than in the case of the fear-induced death rate.


Introduction
Prey-predator interactions are considered the building blocks of ecosystems. Predators evolve, compete and even travel to di erent places in search of their prey. Malthus rst proposed a single species model as an exponential growth [27]. But, later this model is re ned to overcome its shortcomings and to observe different dynamics. Predator functional response on prey species represents the amount of prey consumed by each predator per unit time [16]. A US physical chemist, Alfred James Lotka and an Italian mathematician Vito Volterra rst developed the predator-prey model independently with Holling type-I functional response for predator [23,51]. Later, in this model the functional response was replaced by Holling type-II functional response [20,33] known as Rosenzweig-MacArthur model [42]. Now, functional response has an important role in prey-predator dynamics. It mainly depicts the predator's strategy for consumption of its resources. Holling type functional responses are used in most of the prey-predator models. However, there exist some other functional responses which depend not only on the prey but also on the predator population (for example, Beddington-DeAngelis response [18], Crowley-Martin response [50] etc.). In Holling type functional responses, it is considered that the predator does not interfere in each other's activities and so there is no chance of competition at the time of consumption [16]. However, for the prey population, it is also true that predator nds it di cult to capture a particular prey when the prey species form a group while moving. Each individual in a group always tries to behave according to the majority. A population forms a group or herd for di erent reasons. For example, a prey species tend to move in a group to defend itself from a predator's attack as the chance of getting a successful attack is lesser in this case [7]. On the other hand, a large group of a predator can easily win over its prey [28,29,43]. So, a species form group either for social interaction or for protecting themselves from predators or for enhancing their foraging capabilities. There are already some research works that have been done to analyse the suitable form to depict this social behaviour [7,1]. Based on the work of Ajraldi et al., Braza proposed a functional response replacing the prey population by its square root term in Holling type-II response [1,5]. There are already some works have been published on prey-predator interaction with herd behaviour [4,24,25,30].
Except for direct consumption, the predator has some indirect e ect too on the growth of the prey population. The fear e ect is one of the non-lethal e ects that occurs in presence of a predator and makes an impact on the birth and even on the death of the prey [22]. The prey species always remain under psychological stress due to fear of possible attack and even in some cases, they die only out of fear in the absence of direct killing. The experiment on three di erent zooplankton in Lake Michigan and Lake Erie shows that the fear e ect of predatory spiny water eas (Bythotrephes longimanus) reduces the growth rate of the zooplanktons almost seven times higher than direct predation [37]. Also, according to some literature, behavioural changes in hunting and reproduction occur among a prey population due to fear [44,39,56]. It is a natural tendency of any species to migrate into those places where resources are available in su cient amount and predation risk is lower [40,53]. Several experiments have shown that some prey species spent more time for watching predator and less time for collecting food and even produced fewer eggs [17,2,56]. For example, younger backswimmers maintain a safe distance from older cannibalistic backswimmers and face less hunting [46]. The same behaviour is observed in the case of Mule deer (prey)-mountain lions (predator) [2]. Candolin (1998) made an experiment to show that the scared three-spined stickleback males adjust their reproductive strategies by assessing the fear of predation and current versus future mating probability, which means the predation risk makes an impact on their reproduction [6]. Reproductive physiology of Elk (prey) is also a ected due to fear of wolves (predator) [9]. Several experiments have led to the conclusion that this indirect e ect of predator is equally important as direct killing and even more powerful to reduce the prey population [56,41,47,52]. An experiment performed by Zanette et al. (2011) has shown that the female song sparrow (Melospiza melodia) that are exposed to the predatory sounds produce % lesser edglings due to fear of predation risk than the birds exposed to non-predatory sounds [56]. To exclude the direct killing, they surrounded the nesting with electric fencing and netting. During the entire 130-days breeding season they have broadcast previously recorded sounds of predators (calls of raccoons, hawks, owls etc.) in some places and have played non-predatory sounds (calls of seals, geese, hummingbirds etc.) in some other places. It is observed that only the fear of predators a ects the reproduction of the birds and the numbers of eggs, hatchlings, edglings etc. are reduced in the next generations. Like the terrestrial ecology, the fear of predator is observed in the aquatic ecosystem too. In the absence of largemouth bushes, trout, turtles, etc., bluegill can grow by more than % [55]. Moreover, the growth rate of juvenile Cladocera (Daphnia, Simocephalus) decreases in the presence of kairomone (bluegill sun sh) [14]. In some experiments, it is observed that the mesocarnivores like raccoons reduce their food searching time and techniques by % due to fears of large carnivores like cougars, wolves, black bears, etc. [47]. Hua et al. (2014) in their paper have concluded that the reproducing female bluebirds adjust their breeding strategies according to the vocal signals of their predators [17]. The fear e ect can be observed in snowshoe hares [45] and dugongs population [54] also. According to the experiment of Laundre et al. (2001), the moose become alert when wolves are released in Yellowstone Park [21]. It is also true that not only the birth but the death rate of prey is also a ected due to fear of predation. Baby animals, when scared of fawns, rats, alligators to birds, have their heart rate decreased due to fear. And in the tension of hiding, the heart bursts in fact in some cases. The same characteristic is observed among the ocks of small birds also. A raptor usually swoops down to snatch a bird and this makes other birds who are near the victim got frightened and in fact, fainted or died in some cases. Animals may die from the stress of being hunted when they are chased by predators etc. Some of them die immediately after crumpling to the ground while some linger for weeks and then die. It is observed that when bird trappers catch forest birds and move them into capture boxes, many of them have died immediately, terri ed by the sudden capture. As the heart and the mind are inevitably linked, so emotional stress causes sudden death.
The following sections are categorized as: a mathematical model is proposed with positive initial conditions in Section 2. The region within which the system is biologically well-posed is stated in Section 3 and the conditions under which the species die out from the system in long run is mentioned in Section 4. The feasibility of equilibrium points and corresponding local stability criterion are analysed in Section 5 and Section 6 respectively. Section 7 contains the condition for which the underlying model is permanent. The bifurcation analysis around the equilibria is analysed in Section 8. Section 9 provides the global stability conditions of the predator-free equilibrium point and the numerical simulation for the system is performed in Section 10. Lastly, the work ends up with a brief discussion and conclusion.

Mathematical Model
The proposed model, in this work, consists of prey and predator species, where the growth rate as well as the death rate of the prey are a ected by fear (felt by prey) of the predator. Formulation of an ecological system in terms of mathematical equations helps us to analyse the system dynamics and to extract the essential behaviours of the model. But the environment is in uenced by so many factors which are changeable with time. So, in order to make a simple model, we have made some assumptions here.
Let the biomass of prey and predator population at time t be x and y respectively. The cost of fear is incorporated in the birth rate as well as in the death rate of the prey species. It is considered that the prey population prefer to move in a herd, and so the consumption rate of predator is denoted by the term (Holling type-II) θ √ x +c √ x . The parameters θ and c − denote the coe cient of consumption rate of a predator on prey population and half-saturation constant respectively. It is assumed that in absence of predator, the prey species grows logistically with birth rate r, death rate d and intra-speci c competition rate b. Also, a and m respectively are the conversion e ciency and the mortality rate of predator population. Now, due to fear (felt by prey) of predator, the birth rate of prey species is reduced and the modi ed birth rate becomes rF (K , y), where F (K , y) = +K y with K as the level of fear, satisfying the following conditions: (i) F ( , y) = , i.e., in absence of fear, there is no reduction in prey species reproductive rate. (ii) F (K , ) = , i.e., in absence of predator, there is no reduction in birth rate of prey population. (iii) lim K →∞ F (K , y) = , i.e., the birth rate of prey species becomes zero for increasing anti-predator behaviour (which means for massive increase of fear of predation). (iv) lim y→∞ F (K , y) = , i.e., the birth rate of prey species becomes zero when there is a large number of predator present in the system.
< , i.e., birth rate of prey decreases due to increase of the level of fear.
On the other hand, the death rate of prey is increased due to the impact of fear and so the modi ed death rate becomes d( + F (K , y)), where F (K , y) = K y +y (with < K < ) is the additional coe cient of the death rate of prey population due to the impact of fear ful lling the following conditions: (i) F ( , y) = , i.e., in absence of fear, there is no fear-induced death rate in the system. (iv) ∂ ∂K (F (K , y)) = y + y > , i.e., death rate of prey increases due to increase of fear impact.
So, considering all the assumptions, the proposed system becomes: with positive initial conditions: x( ) > and y( ) > . Note: The dimensions of the parameters are as follows: each of r, d, m has dimension time − , K has dimension biomass − , c has dimension biomass − / , θ has dimension biomass − / time − , and K , a are dimensionless.

Positivity and Boundedness
The following two theorems in this section show that the model (1) is well-posed as the system variables are positive and bounded. Proof. Functions on right hand side of system (1) are continuous and locally Lipschitzian, so, there exists an unique solution (x(t), y(t)) of the system with positive initial conditions (x( ), y( )) on [ , τ], where < τ < +∞ [13]. First we show that, Again the second equation gives, Theorem 3.2. Solutions of system ( ) starting in R + are uniformly bounded.
Proof. First equation of (1) gives: Hence, all solutions of system (1) enter into the region:

Extinction Scenarios
This section provides the conditions for which the prey and predator population goes extinct from the system in long time. Let us adopt the following notations: Proof.
Remark. If the natural death rate exceeds the birth rate of the prey population, then the species ultimately go extinct from the system with time.
Remark. If the birth rate of prey species is reduced by the fear of predators to such a level that it is dominated by the death rate coe cient of the prey species, then the prey population dies out from the system in long run.
Remark. If the fear-induced death rate exceeds the birth rate of the overall prey population, then the prey species go extinct from the system with time.
Remark. If the predator consumes the prey species at a higher rate such that it exceeds the birth rate of prey species, then ultimately the prey population dies out from the system with time. Proof.
Remark. It can be stated that if the maximal bene t of the predator by consuming prey fails to exceed their mortality rate, then the predator population goes extinct from the system in long run.

Equilibrium Points
Solving the nullclines, it is obtained that system (1) has (i) a trivial equilibrium point E ( , ), (ii) an axial equilibrium point E ( r−d b , ) and (iii) an interior equilibrium point E * (x * , y * ) where x * = m (aθ−cm) and y * is the positive solution of the equation: In this case, there exists only one equilibrium point present in the system. This is ecologically meaningful because

Local Stability Analysis
This section contains the local stability criterion of the equilibrium points which can be determined by the eigenvalues of corresponding Jacobian matrices and applying the Routh-Hurwitz criterion. Now, the Jacobian matrix of system (1) is where a = r +K y − d Proof.
So, λ = r − d > and λ = −m < , which imply E is an unstable equilibrium point.

Remark.
It is evident that E is locally asymptotically stable when r < d, which is biologically meaningful. Proof.
Here λ is always negative and so, E is locally asymptotically stable only when λ < , i.e., when aθ < m holds. holds. Proof.
Characteristic equation for E * (x * , y * ) is where B = −a and B = −a a . So, B > when a < . By Routh-Hurwitz criterion all the roots of equation (3) have negative real parts if B i > for i = , . As B is always positive, so, only a < implies the local stability of the interior equilibrium point.

Permanence
An ecological system is said to be 'permanent' when all the species (which exist initially) in the system survive in long run irrespective of any initial population size [48].
Proof. Consider the average Lyapunov function is V(x, y) = x β y β where β i for i = , are positive. In the interior of R + , we havė Now, calculating ϕ(x, y) at the boundary equilibrium points, we obtain: If we take r > d, then ϕ( , ) > for some β i > , i = , .

Bifurcation Analysis
The local bifurcations around the equilibrium points are analysed mainly with the help of Sotomayor's theorem [38] and Hopf bifurcation theorem [33]. To verify the conditions of Sotomayor's theorem, the Jacobian matrix at the bifurcating equilibrium point has to contain a simple zero eigenvalue. Let us consider, V = (v , v ) T and W = (w , w ) T respectively as the eigenvectors of J| (eq. point) and J| T (eq. point) for zero eigenvalue at the equilibrium point.
, taking a as a bifurcation parameter. Proof.
Let, a [TC] be the value of a s.t J| E has a simple zero eigenvalue at a = a [TC] . So, at a = a [TC] : Here, λ = −(r − d) < .
After some calculations: , − T and W = ( , ) T . Therefore, By Sotomayor's Theorem, the system undergoes a transcritical bifurcation around E at a = a [TC] .
Comparing the real and imaginary parts from both sides, we have ( l )l + ( l + C )l +Ċ l = .
Solving we get,l = Hence it is proved.

Global Stability
In this section we discuss the global stability of the predator-free equilibrium point only.
Here µ > by the local stability condition of E . Thus, y(t) ≤ y( ) exp(−µt) and the system is dissipative. From above, we have, lim t→∞ y = . And in the limit form, x(t) is a positive solution of the equationẋ(t) = As, x( ) > , the theorem is proved.

Numerical Simulation
To support the analytical ndings of system (1) where fear is incorporated in both birth and death of prey species, numerical simulation is performed. The pictorial scenarios in this section, help to analyze the system dynamics. We have xed some of the system parameters in Table 1.
In this work, we have not only considered that the birth rate of prey is a ected by the fear of predator but the death rate of prey is also fear-induced. Figure 1 depicts how the fear of predation makes an impact on the birth and death of the prey population. From Figure (1.a) it is observed that the fear-induced birth rate function is decreasing with increasing fear coe cient K . On the other hand, Figure (1.b) shows that the fear-induced death rate function increases with an increasing value of K . This gure mainly ensures that the birth rate   of prey species is decreased whereas the death rate is increased due to fear (felt by prey) of predation. Now Figure 2 describe the trajectory pro le of predator population for increasing value of K ( gure (2.a)) and K ( gure (2.b)). In both cases, the ultimate growth rate of the prey population decreases which leads to the declination of the predator population. In Figure 3, the nullclines of system (1) have been plotted in positive x-y quadrant for the parametric values of Table 1. The non-trivial nullclines meet each other at interior equilibrium point E * ( . , . ) while in absence of predator we get the axial equilibrium point E ( , ). Figure 4 shows that the trajectory starting from x ≡ ( , . ) ultimately converges to the interior steady state and the predator-free equilibrium acts a saddle point in this case. From this situation, if we start to decrease the biomass conversion rate of predator (a), then it is observed that below a threshold value a [TC] there does not exist any coexisting equilibrium point and so, the system contains only the predator-free equilibrium E . Figure 5 depicts the prey nullcline for a = . along with the other parametric values mentioned in Table 1. We do not get any predator nullcline in the positive quadrant in this case. Figure 6 shows that the trajectory starting from x ≡ ( , . ) converges to the axial equilibrium point. So, E becomes a stable equilibrium point from a saddle point through a transcritical bifurcation at  Table 1.   Figure (8.a)). On the other hand, the predation rate can also regulate the system dynamics. From the situation when the system is in a steady coexisting state, if we start to decrease the consumption rate of predator (θ), then there does not exist any interior equilibrium point below a threshold value θ [TC] and so, the system tends to E . Figure 7 depicts the prey nullcline for θ = . along with the other parametric values in Table 1. Here also we do not get any predator nullcline in the positive x-y quadrant.
Taking θ as a bifurcation parameter it is observed that system (1) undergoes a transcritical bifurcation around E at θ = θ [TC] = .
(see Figure (8.b)). The points "BP" and "H" in gures (8.c) and (8.d) represent the transcritical and Hopf bifurcation thresholds respectively. It is observed that E becomes unstable when the   Table 1. parameter values cross the transcritical bifurcation thresholds and also stable branch of interior equilibrium point originate from that bifurcation points. Later on, the system exhibits oscillating behaviour around the coexisting state through supercritical Hopf bifurcations. Again from the steady coexisting state, the system can change its dynamics exhibiting oscillations through Hopf bifurcation. In this work, if we increase the value of a, then after crossing the Hopf threshold value, the population starts to oscillate. Figure 9 shows that for a = .
the system exhibits oscillating behaviour and a stable limit cycle occurs around unstable equilibrium point E * ( . , . ). Figure 10 depicts that system (1) undergoes a Hopf bifurcation around E * at a = a [H] = .
) and st Lyapunov coe cient, l = − . e − < (see gure (8.d)). In order to draw Figures 11-12, we have taken a = . . These gures imply that the fear coe cient induced in prey's birth (K ) also play an important role to control the system dynamics. From the steady-state situation, if the value of K starts to decrease, then ultimately the predator population increases. But the amount of prey available in the system is limited. So, with time, the predator population starts to oscillate which causes the oscillation of prey also. Figure 11 shows that for K = the system exhibits oscillating behaviour and a stable limit cycle occurs around unstable equilibrium point E * ( . , . ) whereas Figure 12 depicts that the trajectory starting from x ( . , . ) converges to the stable equilibrium E * ( . , . ). The calculation gives that system (1)  ) and θ (a = . ) as the bifurcation parameters. In gure (8.c) and (8.d), the blue curve is the branch of predator-free equilibrium (E ) and black curve is the branch of interior equilibrium point (E * ). Other parametric values are taken from Table 1. In Figure 13 the trajectory pro les of predator population have been plotted with time for di erent values of K and K . From the gure, it is observed that if the fear a ects either the birth or death of the prey, the overall predator population decreases with the increase of fear of predation. But the e ect of fear is observed after some time in the both cases. Also, the fear implemented in the birth rate of prey makes a bigger impact on the growth of the predator population than in the case of fear-induced death rate of prey. The di erences of the predator component for increasing value of either K or K for parametric values of Table 1 are listed in Table 2 to analyse the impact of these parameters in the system. In Figure 14 the combined e ect of K and K on the growth of the predator population is observed. This gure depicts that in presence of a larger fear-induced birth rate of prey, the predator species signi cantly decreases for an increasing value of K . It is biologically relevant as the growth of the prey species is a ected because of the fear coe cients. Also, the predator population is considered as a specialist predator which depends only on the mentioned prey and so, decreases with time.
In Figure 15, we have depicted the impact of intrinsic growth rate of prey (r) on predator population. For r = , , we get the count of predator population (y * ) as . , . and . respectively. The equilibrium points, as dependent on the parameter values, change accordingly. This gure supports the fact that the increasing birth rate of prey actually stimulate the growth of predator population and we may get the y-component greater than one for higher values of r. Predator population K = .

Discussion and Conclusion
In an ecological system, prey-predator interaction is such a biological phenomenon that balances the food web. The sustainability of predator species depends on their consumption process and the searching strategy for prey. This consumption process of predator depends on the resource population size, their availability and also on the defence mechanism of resources. There are several species, from bat, sheep, wildebeest etc.
[15] to dolphins, hyena and even horse etc. [32] who prefer to move by forming a group or pack. There are some weak animals like shes, birds, etc. (that are prey to other strong animals) that usually move forming groups to protect themselves from predators [8,26]. Also, in a group, the parents can cooperate while taking care of their young [11,36]. Sometimes, the growth of prey becomes a ected by the frequent attack of their predator. Research reveals that only the fear of predation reduces the reproduction of song sparrow [56]. Also, the fear of being consumed by large carnivores leads to a declination of the foraging time and strategy of mesocarnivores [47]. Now, the baby animals who are scared into hiding, from fawns, rats, alligators to birds, have their heart rate decreased due to fear [3]. And in the tension of hiding, the heart bursts in some cases. The same characteristic is observed among the ocks of small birds also. A raptor usually swoops down to snatch a bird and this makes other birds who are near the victim got frightened and in fact, fainted or died in some cases. So, animals may die out of the stress of the sudden attack of predators [3].
In this work, we have proposed a prey-predator model where the prey species exhibits herd behaviour as a defence strategy. Till now there are some works that have been developed where only the birth rate of prey is a ected with the fear term. But, we rst have considered the fear term not only in the birth rate but also in the death rate of prey species. The study reveals that 'the impact of fear e ect' has an important role to regulate the dynamical behaviour of the system. The fear implemented in the birth rate of prey has a stabilizing e ect because the system reaches a steady-state (through Hopf bifurcation) when K goes above a threshold value. It is also observed that both the fear terms ultimately decrease the predator population as we have considered specialist predator in our model which depends only on the mentioned prey for survival. But, the numerical studies show that the declination of predator is a ected more when the fear is increased in the birth rate than in the death of the prey species. Moreover, it is considered at the time of model formulation that the fear e ect implemented in the death rate can not be increased at a very higher level as the prey is exhibiting herd behaviour to defend themselves. Further, if the fear implemented on the birth rate of prey becomes higher (i.e., for a higher value of K ), the predator species signi cantly decreases for the higher impact of that fear which is implemented on the death rate of prey (increases the value of K ).
In ecological aspects, the 'carry-over e ect' occurs in such a condition where an individual's act depends on their previous experience in a certain situation. Research reveals that these e ects do not depend on the discrete-time scale or seasonal variation or even migration etc. [31,34]. So, it is considered that the carry-over e ect works within and across seasons, years etc. and this e ect has an impact on the reproduction of prey species. Experimental evidences show that carry-over e ects occur in a single season or over short time periods in insects [10], marine sh [19], amphibians [49] and marine invertebrates [31] etc. Researchers in ecology and evolutionary biology always try to model a system so that the mechanism makes an impact on the survival of the population. In this work, we have studied the impact of fear of predation in a prey-predator interaction and analysed how fear regulates the permanence of the population in the system. Now, the fear of predation can be considered as one of the non-lethal e ects of the predator on prey which may not a ect a single generation only but can be carried-over over next generations [35]. So, let us modify system (1) by considering both the fear of predation along with its carry-over e ects. The birth rate of prey in the modi ed where K is the carry-over e ect parameter due to fear, which is quanti ed by the parameter K . In system (7), K = (i.e., the birth of prey is not a ected by the fear of predation) implies a logistic dynamics for prey species with fear-induced death rate function. Again, for K = , the model reduces to our proposed system (1).
For this model, we get the interior equilibrium as E N (x ** , y ** ), where x ** = m (aθ−cm) and y ** is the positive solution of the equation: . Like the proposed system, here also we get only one feasible interior equilibrium when d + bx ** < r holds. The stability condition corresponding to E * N gives that the equilibrium point is locally asymptotically stable when rK K x ** y ** ( +K x ** +K y ** ) + θy ** ( + c < bx ** is satis ed. So, comparing this local stability criterion with the stability criterion corresponding to E * , it can be stated that stability of system (7) ensures the stability of system (1). In Figure 16 the fear-induced birth rate (death rate) of prey in the presence and absence of carry-over e ect have been plotted for increasing value of K (K ). Figure (16.a) shows that when prey is not that much scared of predation (i.e., K is smaller), the carry-over e ect increases the birth rate of prey but when the fear increases, K hardly makes any impact on the birth rate. Furthermore, the death rate of prey, in presence of the carry-over e ect, increases at a higher rate than the case when no carry-over e ect is considered. So, it is observed that the carry-over e ect not only escalate the birth rate but also proliferate the death rate of the prey population. Figure 17 depicts how the carry-over factor a ects the growth rate of a predator. This gure shows that the predator population increases in presence of the carry-over e ect and the inclination maintains a directly proportional relation with K . Now the equilibrium component for prey species in both the systems is independent of the parameters K i for i = , , . But if we look at the trajectories of Figure 18, then it is observed that the prey population for system (7) settle down to a steady-state earlier than system (1) in presence of a carry-over e ect. So, we can conclude that implementing the carry-over e ect in prey species makes a system stable fast in presence of the fear e ect.
As days go, more researchers are showing interest in fear e ect problems but most of the cases are dealt with fear impact in the birth of prey only. Here, we have taken a step further by incorporating the fear term in death of prey species also. It can be concluded from the analytical and numerical results that the model when prey species moves in a herd and both the birth and death rates are a ected by the fear of predation exhibits rich dynamics. Though the proposed system can be modi ed in future to analyse some di erent dynamics. The fear terms may be considered as periodic functions and as fear can not e ect instantaneously, so, one can introduce one or more delay parameters in the fear terms. Also, there may exist environmental noises which lead to a stochastic analysis of the model. As the consumption of food takes some time, so introducing a delay in the functional response can also make the system more realistic. So, in future, some models can be analysed considering all these facts. where Q s pq = ∂ Q s ∂yp ∂yq ( , ,a [H] ) and Q s pqr = ∂ Q s ∂yp ∂yq ∂yr ( , ,a [H] ) for p, q, r, s = , . If l < , then the Hopf bifurcation is supercritical implying formation of a stable limit cycle. On the other hand, the Hopf bifurcation is subcritical for l > and an unstable limit cycle can be observed in this case.