Effects of allochthonous inputs in the control of infectious disease of prey.

Allochthonous inputs are important sources of productivity in many food webs and their influences on food chain model demand further investigations. In this paper, assuming the existence of allochthonous inputs for intermediate predator, a food chain model is formulated with disease in the prey. The stability and persistence conditions of the equilibrium points are determined. Extinction criterion for infected prey population is obtained. It is shown that suitable amount of allochthonous inputs to intermediate predator can control infectious disease of prey population, provided initial intermediate predator population is above a critical value. This critical intermediate population size increases monotonically with the increase of infection rate. It is also shown that control of infectious disease of prey is possible in some cases of seasonally varying contact rate. Dynamical behaviours of the model are investigated numerically through one and two parameter bifurcation analysis using MATCONT 2.5.1 package. The occurrence of Hopf and its continuation curves are noted with the variation of infection rate and allochthonous food availability. The continuation curves of limit point cycle and Neimark Sacker bifurcation are drawn by varying the rate of infection and allochthonous inputs. This study introduces a novel natural non-toxic method for controlling infectious disease of prey in a food chain model.

Introduction

Mathematical epidemiology is an important branch of science as it is essential to determine the patterns, causes, and effects of health and disease conditions in defined populations. Populations are prone to various bacterial and viral diseases in every ecosystem. These diseases have a significant influence on the dynamical evolution of the concerned population. This has prompted the need to study the impact of epidemiological parameters from mathematical as well as ecological point of view. There were significant developments in the theory and applications of epidemiological models of predator–prey systems [1], [2], [3], [4], [5], [6], [7] in the last two decades. Mathematical models have a long history in infectious disease ecology starting from Bernoulli’s [8] modelling of smallpox and including Ross’s [9] analysis of malaria. Pioneer works on mathematical modelling of infectious disease was done by Kermack and McKendrick [5], [6], [7]. These papers had a major influence on the development of mathematical models forinfectious disease and are still relevant in many epidemic situations. There are many famous mathematical models which are related with important human infectious diseases e.g., HIV [10], malaria [11], [12], SARS-coronavirus [13] and influenza [14]. Mathematical models are also being used to explore wildlife disease dynamics [15], [16] and possible routes of zoonotic disease emergence. A predator–prey system with infected prey in polluted environment was investigated by Sinha et al. [17]. Anderson and May [18] were probably the first who considered the disease factor in a predator–prey dynamics and observed that pathogen tends to destabilize the predator–prey interaction. Hadeler and Freedman [19] determined a threshold of infection rate above which an infected equilibrium or an infected periodic solution appear in Rosenzweig prey-predator model. Many eco-epidemiological models are investigated by several researchers considering various types of interactions between the populations [20], [21]. Chattopadhyay and Arino [22] formulated three species eco-epidemiological model and they studied local stability of equilibrium points, extinction criteria of species and determined the conditions for Hopf bifurcation in an equivalent two-dimensional model. Haque and Chattopadhyay [23] investigated the role of transmissible diseases in a predator–prey system with infected prey. Bhattacharyya et al. [24] proposed and investigated an epidemiological model considering nonlinearity in infection incidence. Liu et al. [25] had done the bifurcation analysis of a predator–prey model with disease in prey and they reported various bifurcation scenarios of the system. Yongzhen et al. [26] investigated effect of delay on a predator prey model with parasitic infection. Recently, various disease control strategies are introduced to control infectious disease in food chain model [27], [28]. Seasonal variation of infection rate is an important factor in disease modelling. Modelling and analysis of the transmission dynamics of tuberculosis without and with seasonality was investigated by Bowong et al. [29].

Disease control has relied heavily on the application of chemical fungicides, bactericides and soil fumigants over the last 50 years. Now there are many problems associated with their continued deployment forces us to reduce use of chemicals in the environment in general. Consequently, applications of alternative non-chemical methods of disease control continues to gain importance. It is well known that most predators not only feed a single prey rather move for alternative foods also. Polis and Strong [30] reported that the alternative resources and pathways plays central role in food web dynamics. Unexpected consequences of predator removal was observed experimentally by Sih et al. [31] and Cote & Sutherland [32]. They found that predator removal reduced prey populations in many cases. Several experimental studies [33], [34] had shown that predators can easily capture diseased prey. If predators eliminate the most infectious individuals from the prey population then that will have an outcome equivalent to quarantine, whereby infectious individuals are removed from the ecosystem and thereby prevented from disease spreading. Therefore, predator control has significant effect on food chain model with infected prey. Allochthonous inputs to predators help to increase predator population in a food chain. Some investigations [35], [36], [37], [38], [39], [40], [41] were done either providing alternative food or considering allochthonous inputs to predator in a predator–prey system. Recently, Sahoo and Poria [40], [41] reported the consequences of alternative food for predator in a diseased predator–prey system. They showed that an diseased system may becomes disease free in presence of additional food up to some threshold infection rate. But disease controlling aspects of allochthonous inputs to intermediate predator in a three species food chain model is unexplored. Starting from a well studied tritrophic food chain model [37], we investigate the disease controllability aspects of allochthonous inputs in this paper. The bifurcation phenomena of the model are analysed with respect to important biological parameters for finding periodic and other behaviours [42] of the model.

In this paper, we consider Huxel and McCann [37] three species food chain model with infected prey population and investigate the effects of supplying allochthonous foods to intermediate predator. The main aim of this paper is to analyse the role of allochthonous inputs to intermediate predator for controlling disease in prey population. The section-wise split of this paper are as follows. In next Section, the model formulation is done. In Section 3 the conditions for positivity and boundedness of solutions are derived. Stability conditions of all equilibrium points are derived and extinction criterion of infected prey population are calculated in the Section 3. Persistence condition of the system is derived in Section 4. In Section 5 the numerical simulation results of the model are presented and the effects of seasonal variation of contact rate on the model is also investigated. Numerically, one and two parameter bifurcation results using MATCONT 2.5.1. software [43], [44], [45] are presented in Section 6. Finally, conclusion is written in Section 7.

Model

A tritrophic food chain model composed of a logistic prey (x), a Holling type II intermediate predator (y), and a Holling type II top-predator (z), with allochthonous inputs () into the intermediate predator population and allochthonous inputs () to the top-predator population was proposed by Huxel and McCann [37]. The model is the following (see Appendix A for a full derivation and discussion of the underlying biological assumptions): Here, and are the parameters describing the feeding preference for the allochthonous inputs by the intermediate predator and top-predator population respectively. The analysis is performed over a range of feeding preference from 0 (feeding only on autochthonous/classical food web sources) to 1 (feeding only on the allochthonous sources). The supply of allochthonous inputs remain constant. The parameters can be deemed biologically plausible as they represent realistic predator–prey ratios in body size found in surveys [46], [47].

Polis and Hurd [48] suggested that for an island system most of the allochthonous inputs are available to intermediate predator population only. There are many other systems which are driven by allochthonous inputs to the intermediate predator population include marine filter-feeding communities in unidirectional currents or invective areas [49], soil communities [50], [30] and headwater streams which received leaf litter inputs [51]. Therefore, it is biologically meaningful to assume supply of allochthonous inputs to intermediate predator population only in the model (1). Setting and , the model (1) reduces to:

Now, suppose that prey population is infected by some infectious disease. To modify the model (2) in presence of infectious disease we make the following assumptions:

The prey population is divided into two classes, viz. (i) susceptible class whose population density is denoted by s and (ii) infected class whose population density is denoted by i. The intermediate predator whose population density is denoted by and the density of top-predator is denoted by .

A part of the susceptible prey population becomes infected at a rate , following the law of mass action and d be the natural death rate of infected population.

Infected population is not in a state of reproduction and also does not compete for the resources.

Behaviour of the entire community is assumed to arise from the coupling of these interacting species, where preys on both susceptible prey and infected prey in the form of Holling type-II and Holling type-I respectively. This different combinations of functional forms are taken because it is easier to capture infected prey than the susceptible prey [20]. Top-predator captures intermediate predator in the form of Holling type-II interaction. Infection does not transfer to intermediate predator from prey population.

Therefore, under the above assumptions the modified model (2) in presence of infected prey become the following:The schematic diagram of our model is shown in Fig. 1 . The system (3) should be analysed with initial conditions since it is a population model.

Fig. 1

Schematic diagram of the model (3).

Preliminaries

Basic properties

In this subsection we investigate some preliminary but relevant properties of the system (3).

The solutions of the system (3) exist in the closed positive octant, .

Observe that the right-hand side of system (3) is a Lipschitz continuous function of the variables . Therefore, existence and uniqueness properties of solutions hold in .  □

If remains non-negative, then the possible solutions of the system (3) are positive.

The proof is obvious.  □

The set is positively invariant for the system’s trajectories.

This result follows from Proposition 1, Proposition 2.  □

All the solutions of the system (3) which start from points of are bounded.

From the first equation of the system (3) we can write, ,

Therefore, , as .

Now we define, .

The time derivative of the above equation along the solutions of (3) is

.

Using equations of (3), we havewhere = min, provided .

Applying the theory of differential inequality we obtainNow, as , we have .  □

Hence all the solutions of (3) that initiate in are confined in the region.

, for any  . This proves the theorem.

Criterion for the extinction of infected prey

If the initial population size is sufficiently large so that the initial growth rate of the infected population would be negative (where the susceptible population density at ) i.e., if holds, then the infected prey population will extinct, i.e., .

From the system (3), we haveHence the theorem follows.  □

From the above theorem it is clear that infected prey population will extinct if the initial intermediate predator population is greater than some positive value (threshold value) depending on infection rate of prey () as well as allochthonous inputs (A). The threshold population size increases monotonically with the increase of infection rate.

Stability analysis of equilibrium points

We now discuss the stability behaviour of the system (3). The system (3) possesses following five biological equilibrium states.

The trivial equilibrium state . An eigenvalue associated with the Jacobian matrix at is 1, therefore is an unstable equilibrium point.

The axial equilibrium state . The Jacobian matrix at is (see Appendix B). The axial equilibrium point is stable if , i.e., if (where the basic reproductive ratio , the average number of new infections in a susceptible population over the (average) duration of infections) and , i.e., . The equilibrium point is unstable if either or or both holds (see Appendix B for details).

The predator free equilibrium state is , where and . The equilibrium point exists if .

The predator free equilibrium point is stable if and (see Appendix B for details).

The disease free equilibrium state is , where and is the positive root of the equation , where and . The disease free equilibrium exists if and .

The disease free equilibrium point of the system (3) is locally asymptotically stable if and the conditions and hold (see Appendix B for details).

The interior equilibrium state is , where .The positive interior equilibrium point exists if and .

Using the Routh–Hurwitz criterion (see Appendix C), we observe that the interior equilibrium point of the system (3) is locally asymptotically stable if the conditions ,   hold (see Appendix B for Jacobian matrix and associate eigenvalues at the equilibrium point ).

Persistence of the system (3)

Persistence of a system means long term survival of all populations of the system, no matter what the initial populations are. Mathematically, persistence of a system means that strictly positive solutions do not have omega() limit points on the boundary of the non-negative cone.

Let the following conditions are satisfied and if there exists a finite number (say, n) of periodic solutions in the space, then system (3) is uniformly persistent provided for each periodic solutions of period T,

growth rate of infected prey in absence of intermediate predator is positive, i.e., ,

allochthonous input (A) is greater than the ratio of half saturation constant of prey-intermediate predator and feeding preference for the allochthonous input by the intermediate predator, i.e., ,

intermediate predators growth rate in absence of top-predator is positive, i.e., ,

.

Let be a point in the positive quadrant and be orbit through X and be the omega limit set of the orbit through X. Note that is bounded.

We claim that . If then by the Butler-McGehee lemma [52] there exist a point P in where denotes the stable manifold of . Since lies in and is the space, we conclude that is unbounded, which is contrary to the boundedness of the system.

Next , for otherwise, since is a saddle point which follows from the conditions and , by the Butler-McGehee lemma [52] there exist a point P in . Now is the s-axis implies that an unbounded orbit lies in , a contradiction.

The condition implies that is unstable and therefore .

Lastly, we show that no periodic orbits in space or . Let denote the closed orbit of the periodic solution in space such that lies inside . Let the Jacobian matrix corresponding to is given by

Here . Computing the fundamental matrix of the linear periodic system , we find that its Floquet multiplier in the i-direction is . Then proceeding in an analogous manner like Kumar and Freedman [53], we conclude that no lies on . Thus, lies in the positive quadrant and system (3) is persistent. Finally, since only the closed orbits and the equilibria from the omega limit set of the solutions on boundary of and the system (3) is dissipative. Now using a theorem of Butler et al. [54], we conclude that system (3) is uniformly persistent.

Numerical results

All numerical simulations of the system (3) are done taking the parameters values of Table 1 . Effects of variation of allochthonous inputs A as well as infection rate are investigated numerically. Numerical simulations are done using MATLAB 7.8.0 (R2009a).

Table 1

Parameters values for the system (3).

Parameters	Description	Parameter value
k	Carrying capacity	1 [37]
a_1	Mass-specific metabolic rate of intermediate predator p_1	0.4 [55], [56]
a_2	Measure of the ingestion rate per unit metabolic rate of p_1	2.009 [55], [56]
a_3	Mass-specific metabolic rate of top-predator p_2	0.08 [55], [56]
a_4	Measure of the ingestion rate per unit metabolic rate of p_2	5 [55], [56]
b_1	Half saturation constant of intermediate predator	0.16129 [55], [56]
b_2	Half saturation constant of top-predator	0.5 [55], [56]
c_1	Preference parameter	0.3 (0–1) [37]
c_2	Preference parameter	0.2 (0–1) [37]
m	Attack rate of infected prey	0.3 (0–0.5) [57]
d	death rate of infected prey	0.1 (0.05–0.16) [58]

We first numerically calculate the initial intermediate predator population (using Theorem 3.1) for extinction of infected prey using the set of parameter values of Table 1 with the variation of allochthonous inputs A for different infection rate . From Fig. 2 , we observe that the system (3) becomes disease free above the critical and disease free region is decreasing with the increase of infection rate. We observe that there is a threshold value of initial intermediate predator population for eradication of infected prey. Intermediate population size should be greater than the threshold value for extinction of infected prey. The threshold population size increases monotonically with the increase of infection rate.

Fig. 2

Graph of with respect to allochthonous inputs (A) for various infection rate (a) , (b) , (c) and (d) . The figure represents the density of initial intermediate predator population () for various infection rate.

The phase trajectories of the system (3) around predator free equilibrium point are shown in Fig. 3 with initial populations (black line), (green line), (red line), (magenta line) and (blue line). The system (3) is stable at equilibrium points for (Fig. 3(a)), for (Fig. 3(b)) and for (Fig. 3(c)). The local stability of disease free equilibrium point with initial population (1.3,0,0.5,0.1) is presented in Fig. 4 . From Fig. 4, it is evident that the system (3) has a stable limit cycle around disease free equilibrium point (0.9102, 0, 0.1562, 0.6862) for (Fig. 4(a)) and there is a stable limit cycle around the equilibrium point (0.9071, 0, 0.1562, 0.6816) for (Fig. 4(b)).

Fig. 3

Phase trajectories of the system (3) corresponding to different initial conditions (black line), (green line), (red line), (magenta line) and (blue line) for (a) , (b) and (c) are presented. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 4

Phase trajectories of the system (3) around disease free equilibrium point with initial condition and for (a) and (b) (other parameters are taken from Table 1) are shown.

The numerical bifurcation results are presented here for random initial conditions, which are lie on positive quadrant. One parameter bifurcation diagrams of the system (3) with respect to infection rate are presented in Fig. 5 in absence of allochthonous inputs (). It is clear from Fig. 5 that the species extinction is possible at higher infection rate in absence of allochthonous inputs. The system become disease free naturally at low infection rate but with the increases of infection rate (after ) infected prey population increases in absence of allochthonous inputs. Bifurcation diagrams of the system (3) with respect to in presence of allochthonous inputs (A) are presented in Fig. 6 . It is evident from Fig. 6 that in presence of allochthonous inputs () the system becomes disease free for , disease present for . If we increase the supply of allochthonous input slightly and set then the system becomes disease free almost everywhere within (Fig. 6(b)). Next, bifurcation diagram with respect to A of the infected prey population are drawn for different infection rate in Fig. 7 . Fixing at different levels (at ), we observe that the system (3) becomes disease free when allochthonous input (Fig. 7) on the other hand whenever A lies within , the infection exists in the system.

Fig. 5

Bifurcation diagrams of the system (3) for susceptible prey , infected prey , intermediate predator and top-predator population with respect to infection rate in absence of allochthonous inputs (i.e., ) (The parameters are taken from Table 1) are plotted. It shows that infected prey extincts for and survived for .

Fig. 6

Bifurcation diagrams of the system (3) for infected prey population with respect to infection rate in presence of allochthonous inputs are plotted for for in (a) and for in (b) (other parameters are taken from Table 1) are depicted. For , the system become infection free when , disease persists for .

Fig. 7

Bifurcation diagrams of the system (3) for infected prey with respect to allochthonous inputs (A) are depicted for infection rates in (a), for in (b), for in (c) and for in (d) (other parameters are same as Table 1).

From the above bifurcation diagrams, we observe an internal relationship between infection rate () and allochthonous inputs (A). So, we have numerically calculated the infected and uninfected regions in the parameter plane in Fig. 8 . From Fig. 8, it is clear that for low allochthonous inputs, the system may have disease for higher infection rate. On the other hand for very high infection rate, the system become disease free in presence of high allochthonous inputs. Notice that in Fig. 8 the boundaries between infected and uninfected regions are not perfectly distinct, because, there is some degree of sensitivity to small changes in parameter values resulting in sharp transitions between different dynamical outcomes.

Fig. 8

Behaviour of the system (3) in the ) plane taking other parameters from Table 1 is plotted. The -plane is divided into infected and uninfected regions. We observe that for higher infection rate, high level of allochthonous inputs is required to make the system disease free.

Effects of seasonally varying contact rate

It is well known that seasonal variations in temperature, rainfall and resource availability strongly influence infectious disease dynamics due to variation of contact rate. London and Yorke [59] and Yorke et al. [60] have computed seasonal contact rates for various diseases. Both London and Yorke [59] and Yorke et al. [60] showed that for monthly data, the contact rate appears to be smooth and periodic with a period of about one year. London and Yorke [59] demonstrated that seasonality is necessary for perpetuating a recurrent epidemic, so we will concentrate on examining the effects of variation in the model (3). In order to investigate the qualitative effects of seasonality, we assume the infection-rate is set by the parameter [61], [62], [63] and incorporates time-dependent external factors affecting the contact rate between members of the population and the infectivity of the disease in the model (3). We assume [60], [61], [62], where is the average contact rate-constant, represents the strength of the seasonal forcing between zero and unity and t has units of years. Note that the sinusoidal forcing may strongly influence the long term behaviour [64]. Here, we investigate the effects of seasonal variation of contact rate in the disease dynamics of prey population in presence of allochthonous inputs.

The bifurcation diagrams with respect to of infected prey are presented in Fig. 9 with random positive initial condition using MATLAB software (see Appendix G for algorithms of numerical simulations). From Fig. 9, Fig. 9(c), we observe that the infection of prey population exists for seasonally varying in absence of allochthonous inputs to intermediate predator. On the other hand, in presence of allochthonous inputs the system becomes disease free (Fig. 9, Fig. 9(d)). Therefore, allochthonous inputs to intermediate predator may be very useful to remove disease in prey population in case of seasonally varying contact rate also.

Fig. 9

Bifurcation diagrams of the system (3) for infected prey with respect to strength of the seasonal variation are presented for in (a), for in (b), for in (c) and for in (d) (other parameters are taken from Table 1).

Bifurcation of equilibrium and Hopf continuation

In this section, we have investigated detail bifurcation nature of the system (3) with respect to the parameters and A numerically with the help of parameter values of Table 1. We study the existence of periodic orbits bifurcating from equilibrium points of the system (3). This phenomenon is know as a generalized Hopf bifurcation [65]. This is actually done by studying the change in the eigenvalue of the Jacobian matrix and also following the continuation algorithm. We consider a set of fixed point initial solution, and , corresponding to a parameter set of values as in Table 1. The characteristics of Hopf point, the limit cycle and the general bifurcations are explored using the software package MATCONT2.5.1. [43], [44], [45]. In this package we use prediction-correction continuation algorithm based on the Moore–Penrose matrix pseudo inverse for computing the curves of equilibria, limit point (LP), along with fold bifurcation points of limit cycles etc (see Appendix G for algorithms of numerical simulations).

In Fig. 10 (a) the continuation curves from the equilibrium point with respect to as the free parameter for are shown. In the Fig. 10(a) existence of two Hopf points (H), two branch point (BP) and two neutral saddle (H) of s with the variation of are observed for different A (fixed). We now discuss about the equilibrium curve of the system for fixed . For , we observe two Hopf point, one branch point in the positive quadrant. The first Hopf point () is located at (0.269218, 0.861629, 0.156250, 1.239693, 0.675168) with eigenvalues . For this Hopf point the first Lyapunov coefficient turns out to be −0.008568927, indicating a supercritical Hopf bifurcation. Indeed, there are two complex eigenvalues of the equilibrium with real at this parameter. First Lyapunov coefficient is negative implies that a stable limit cycle appears when equilibrium point looses stability (see, Appendix D for general procedure of existence of Hopf bifurcation). As the parameter is increasing, second Hopf point ) is situated at . For this second Hopf point the first Lyapunov coefficient turns out to be −0.05326441. The branch point (BP) occurs at (0.915670, 0.000000, 0.156250, 0.545396, 0.198508) with eigenvalue values . The remaining points (H) indicates the neutral saddle equilibrium. Similar bifurcation analysis are carried out with respect to for the variables and and the results are displayed in Fig. 10, Fig. 11, Fig. 11 respectively.

Fig. 10

Continuation curves of equilibrium with the variation of the parameter for (blue line), for (green line), for (magenta line) of susceptible prey and infected prey are presented. The labels are: -Hopf point, BP-branch point and neutral saddle (H). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 11

Continuation curves of equilibrium with the variation of for (blue line), for (green line), for (magenta line) of intermediate predator and top-predator are plotted. The labels are: -Hopf point, BP-branch point and neutral saddle (H). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

To proceed further we start from the first Hopf point () which is predicted in Fig. 10, Fig. 11 as the initial point at , and using the backward continuation process we get a family of stable limit cycles bifurcating from this Hopf point. This phenomenon is shown in Fig. 12 , where again the infection rate in the system is the only free parameter. One observes that at , we have a LPC point with period 17.20543 and normal form coefficient = −0.01373077. At this situation two cycles collide and disappears. The critical cycle has a double multiplier equal to 1. From this it follows that a stable branch occurs after the LPC point. At there is another LPC point with period 20.078810 and the normal form coefficient = 0.007475564. At that point two of the multiplier is equal to 1. We now start from the second Hopf point () in Fig. 12 as the initial point at , and proceeding in the opposite direction, get a LPC point with period 42.85038 and the normal form coefficient is = −0.05096683. A list of normal form coefficients and eigenvalues are given for codimension 1 bifurcation points in Table 2 . We list periods of limit cycle oscillations together with normal form coefficients (Fig. 12) in Table 3 .

Fig. 12

Family of limit cycles bifurcating from the Hopf points () of infected prey at (a)    and (b)    for . The cycles undergoing the limit point cycle bifurcations are labelled by and .

Table 2

Parameter values of and A at the bifurcation points in Fig. 10, Fig. 11, together with normal form coefficients and eigenvalues [66].

Label	δ	A	First Lyapunov coefficients	Eigenvalues
H_1	0.589465	0.8	l_1=-0.008568927	-0.0473243±i0.173118,±i0.365186
H_2	1.350987	0.8	l_1=-0.05326441	-0.00955927±i0.389868,±i0.146631
BP	0.198508	0.8	-	-0.856901,0,0.0244703±i0.131973

Table 3

Limit cycle bifurcation with respect to in Fig. 12, together with first Lyapunov coefficients [66].

Label	δ	A	Period	Normal form coefficients
LPC_1	0.675167	0.8	17.20543	−0.01373077
LPC_2	0.6645317	0.8	20.7881	0.007475564
LPC_3	1.350987	0.8	42.85038	−0.05096683

Now, we investigate codimension 2 bifurcation sequences with and A as bifurcation parameters. We now consider our starting point to be the first Hopf point () occurring at , in Fig. 10. Applying forward and backward continuation techniques, we observe (GH) points and (HH) point considering and A as free parameters and this phenomena is presented in Fig. 13 (see, Appendix F for existing conditions of different types of Hopf points). At the (GH) points, where the first Lyapunov coefficient vanishes indicating that all GH points are non-degenerate, since the second Lyapunov coefficients are non-zero. First generalized Hopf point () occurs at (0.290153, 0.980010, 0.156250, 1.277808, 0.589465, 1.002103) with second lyapunov coefficient . Using forward and backward continuation, we observe another generalized Hopf point at with second lyapunov coefficient . The double Hopf (HH) is located at . We list codimension 2 bifurcation points together with normal coefficients is presented in Table 4 . As theory predicts [66], each GH point emanates an LPC curve. Actually, there is just one LPC curve connecting the GH points; while approaching these points, the critical non-hyperbolic cycle shrinks to the equilibrium point. We now consider our starting point to be the point occurring at as shown in Fig. 12 and the corresponding limit point cycle (LPC) scenarios are observed in Fig. 14 (a). Staring from another point at , the bifurcation points LPC, PD, NS are presented in Fig. 14(b) (see, Appendix E for the study of limit point cycle). We list the bifurcation points together with period and normal form of coefficients in Table 5 (Fig. 14).

Fig. 13

Hopf point continuation curves of the system (3) are drawn in the plane : GH-generalized Hopf points, HH-double Hopf point. The diagrams (b) and (c) depict the real and imiginary part of eigenvalues corresponding to these bifurcation points w.r.to infection rate .

Table 4

Parameter values of and A at the bifurcation points in Fig. 13, together with normal form coefficients and eigenvalues [66].

Label	δ	A	Normal form coefficients	Eigenvalues
GH_1	0.589465	1.002103	l_1=0,l_2=0.01359786	-0.594702±i0.184733,±i0.184733
GH_2	0.847191	1.309871	l_1=0,l_2=0.806695	-0.0296194±i0.340657,±i0.176157
HH	2.197830	0.538660	p_11p_22=1,	±i0.117916,±i0.413422
			θ=-2.236111,
			δ=1.369807,
			Θ=25.95972,
			Δ=4327.669

Fig. 14

Limit cycle continuation curve of the system (3) is plotted starting from generalized Hopf (GH) points with the variation of the parameter : LPC-Limit point cycle, PD- Period doubling, NS-Neimark-Sacker.

Table 5

Limit point cycle continuation starting from point with respect to in Fig. 14, together with period and normal form of coefficients [66].

Label	δ	A	Period	Normal form coefficients
LPC_1	0.5347408	0.8	21.63391	0.4319366
LPC_2	0.8471917	0.8	35.67527	−0.0000396739
LPC_3	0.5894640	0.8	18.95807	−0.000005725
LPC_4	1.0889	0.8	39.34845	0.2451057
PD_1	1.47621	0.8	48.93902	−0.00094017
PD_2	1.960606	0.8	59.42899	0.000071263
NS_1(Neutralsaddle)	0.9242725	0.8	55.26618	–
NS_2	1.715313	0.8	54.70501	−2.110497

We now start taking the HH point as initial point. With the variation of and A, the Neimark-Sacker (torus) bifurcation curve is presented in Fig. 15 . The shape of the NS curve is more complicated: there are five more codimension 2 points on it. First Fold-Neimark Sacker bifurcation () is located at with period 53.28529. There are 1:3 strong resonance (R3) point at with period 54.73327 and 1:4 resonance (R4) point at with period 56.85741. A generic unfolding of this singularity has a period-3 saddle cycle that does not bifurcate for nearby parameter values, although it merges with the primary fixed point as the parameter approach R3. There is 1:1 strong resonance points (R1), where the fold bifurcation curve and NS curves meet tangentially (Fig. 16 ). Second Fold-Neimark Sacker bifurcation () is situated at (). Periods of Neimark Sacker bifurcation curve starting from double Hopf point (HH) of the system (3) with the variation of the parameters and A is listed in Table 6 and is depicted in Fig. 17 .

Fig. 15

Neimark Sacker bifurcation starting curves from double Hopf point (HH) of the system (3) is presented with the variation of the parameters and A: LPNS-Fold-Neimark Sacker bifurcation, R1-Resonance 1:1, R3-Resonance 1:3, R4-Resonance 1:4.

Fig. 16

Neimark Sacker bifurcation curves starting from double Hopf point (HH) of the system (3) is presented with the variation of the parameters and A: LPNS-Fold-Neimark Sacker bifurcation, R1-Resonance 1:1, R3-Resonance 1:3, R4-Resonance 1:4.

Table 6

NS curve bifurcation point starting from HH point with the variation of and A in Fig. 15.

Label	δ	A	Period
LPNS_1	2.197839	0.5386601	53.28529
LPNS_2	2.096277	0.792391	57.74066
R1	2.096918	0.7923527	57.74986
R3	2.028442	0.6312537	54.73327
R4	1.984677	0.7096657	56.85741

Fig. 17

Periods of Neimark Sacker bifurcation curves of the system (3) are plotted starting from double Hopf point (HH) with the variation of the parameters and A: LPNS-Fold-Neimark Sacker bifurcation, R1-Resonance 1:1, R3-Resonance 1:3, R4-Resonance 1:4.

Conclusions

We formulate a three species food chain model with infected prey and investigate the effects of supplying allochthonous inputs to intermediate predator in the model. Equilibrium points of the model are calculated and their stability criteria are determined. Conditions for extinction of infected prey population are derived analytically. The persistence conditions of the system are derived in presence of allochthonous inputs. One parameter bifurcation analysis is done numerically with respect to infection rate as well as allochthonous inputs parameter using a set of biologically meaningful data. We have also simulated two parameter bifurcation scenarios by varying infection rate and allochthonous inputs together.

Numerical simulation results show that system (3) becomes disease free in absence of allochthonous inputs for lower infection rate (). This is happening because infected prey is weaker than susceptible prey and is easily captured by intermediate predator. This will reduce the number of infected prey and ultimately the system may becomes disease free at certain stage for low infection rate without allochthonous inputs. But for higher infection rate (), the disease free system can not be obtained in absence of allochthonous inputs. However, the bifurcation analysis with respect to infection rate proves that the system becomes disease free in presence of suitable allochthonous inputs for higher infection rate (Fig. 6, Fig. 7). If intermediate predator eliminate the most infected individuals from the prey population then it has an outcome equivalent to quarantine and prevent disease spreading among prey species. Allochthonous inputs to intermediate predator helps to increase its population size which helps to eat infected prey at higher rate. In this way allochthonous inputs to intermediate predator is very useful for controlling disease of prey population in a food chain model. We also consider seasonally varying contact rate and show that the system becomes disease free there in presence of suitable allochthonous inputs. Therefore, we may conclude that allochthonous inputs can control infectious disease of prey in a food chain model. We draw infected and uninfected region in the parameter plane as a function of infection rate and allochthonous inputs.

We have demonstrated the pattern of bifurcation in the proposed epidemic model to observe dynamical behaviour. We have employed numerical continuation methods to compute curves of codimension 1 and 2 bifurcation using MATCONT 2.5.1 software package. Existence of Hopf (H), double-Hopf (HH), generalized Hopf (GH) bifurcation with the variation of infection rate and allochthonous inputs are noticed. We have shown that the system also undergoes branching, period doubling, and Neimark-Sacker bifurcations. The first Lyapunov coefficient is computed to determine the nature of Hopf bifurcations in the model. The local representation of a Neimark-Sacker (NS) bifurcation curve is determined. We observe some critical bifurcation points like: 1:1 Resonance (R1), 1:3 Resonance (R3) and 1:4 resonance (R4) in the model with the variation of infection rate and allochthonous input. Therefore, many qualitatively different dynamical behaviours are observed in the model with the variation of infection rate and allochthonous input.

The key observation is that allochthonous inputs to intermediate predator can control infectious disease of prey in a food chain model if the death rate due to predation is lower than disease induced mortality of prey. If predators eliminate the most infectious individuals from the prey population then that will have an outcome equivalent to quarantine, which can prevent disease spreading among prey. Actually, the supply of allochthonous input increase the population of intermediate predator and then it can captures infected prey at higher rate and as a result, after certain time the system becomes disease free. Considering a parasite infected food chain model Packer et al. [67] showed that predation will increase prey population size when the force of predation is lower than prey mortality. On the other hand if the death rate due to predation is higher than prey mortality then negative effects of allochthonous inputs will be expected. Notice that allochthonous inputs never increase the disease of prey but may not control the infection of the prey at high infection rate (for example see Fig. 8). Results presented in this paper are important to understand the role of allochthonous inputs as a natural disease controller of infected prey in a food chain model. The possible interesting aspects for future developments are to investigate the role of allochthonous inputs as a disease controller in food chain models with other type of functional responses between intermediate predator and infected prey and also with different types of infection spreading rule. The role of allochthonous inputs to control infectious disease of prey in a network of infected food chain models are also an important area for future developments. This novel natural non-toxic method of disease control may be very useful for biological conservation of prey species and also for developing new strategies for controlling infectious disease of prey in a real world biological systems.