STABILITY AND BIFURCATION ANALYSIS OF A CONTAMINATED SIR MODEL WITH SATURATED TYPE INCIDENCE RATE AND HOLLING TYPE-III TREATMENT FUNCTION

The paper deals with an epidemic model which is contaminated. Logistic growth input of vulnerable is taken in our SIR model. Saturated type incidence rate alongside Holling type III treatment function and time delay in growth component has been considered. The stability analysis of equilibrium points is done and the basic reproduction number (R0) is obtained. Using (R0), Local Stability around disease-free equilibrium point has been acquired. When R0 is less than one, it is locally stable which means disease doesn’t exist in the environment anymore and when R0 is greater than one, it is unstable which means disease still exists in the environment. We have shown that at R0 = 1 the model exhibits forward bifurcation. The conditions for the local stability of endemic equilibrium point have also been obtained. Conditions for the existence of Hopf bifurcation along with its direction and stability has been determined. Numerical simulations have been performed utilizing MATLAB to confirm the logical outcomes acquired. We have seen that to diminish the effect of illness, contamination ought to be decreased. Also, suitable strategies should be chosen to expand the capacity of treatment to decrease the effect of infection.


INTRODUCTION
With constant progress and increase of population in the urban areas, a wide range of problems related to environment such as air pollution, release of toxic materials, solid waste disposals, deforestation and so on, have gained attention much greater than ever before. Moreover, the amount of pollution in the environment can also be affected by the local species; the pollution gets absorbed, digested, and eventually gets metabolized by the species. Environment is one of the most important elements present to the society. The presence of pollution in the environment is a very broad topic and many authors have incorporated the effects and factors of pollution in their study [15,5,11,19]. The major cause of pollution is the excessive use of environmental resources and chemicals, without taking into consideration the safety measures to control pollution. Hence, resulting pollution directly affects the humans' health. Many scientists have shown their concerns regarding the reduction of environmental pollution. Motivated by the same novelty here, we have developed an SIR model where we have considered the effect of contaminants on populace. Since epidemiology is a discipline which deals with infection in a population. It deals with all aspects of epidemic such as spread, control, vaccination strategy etc. Mathematical modeling is an important tool in epidemiology which analyze the spread and control of infection in the species. Based on this here, an epidemic model describing the disease dynamics is developed. In [4,1,2], the researchers in their models have considered constant recruitment rate of the susceptible with saturated incidence rate. However, logistic growth input of susceptible have been considered in many problems and is practical for infections with high mortality rate. A model was proposed by Hutchinson [7] in which he considered a logistic equation with delay in death rate term or carrying capacity which is as follows [10]: Initially, the populace growth was exponential but it oscillated at higher levels of population. In mathematical community, his formulation on logistic equation involving the delay is know as delayed logistic equation.
Since, there are many examples involving delay in this organic world. One of them being, at the time of sexual reproduction, the time from the fertilization of eggs to the birth. Another example involving delay is, time involved in the digestion of food and its conversion to biomass resources. In our paper, we have incorporated delay in growth term apart from interaction term.
Since there are numerous reasons because of which infection in caused, the doctors and health organizations came up with various treatments to have control over the spread of disease among populace. These treatments involve antibiotics, vaccination, blood transfusion etc. Researchers who studied models based on epidemics worked on the assumption that treatment rate of disease is in direct relation to number of infective populace. It is also known in general that treatment effect the recovery rate. Authors like Wang and Ruan [21] studied the SIR model in which they examined the case where T (I) = 0 for I = 0 and is equal to r which is constant for all I > 0.
Wang [20] in his further work on SIR model, incorporated piecewise continuous treatment function: T (I) = min(rI, rI 0 ), where the health care system reaches capacity at infective level I 0 , and r > 0 is constant. Zhang and Liu [24], A.Kumar and Nilam [3] used the staturated treatment function: T (I) = rI 1+αI for r > 0, α ≥ 0. Dubey et al. [13] incorporated the Holling type III treatment rate, that determines the case where recovery rate grows rapidly with the increase in infectives initially, and later on, it settles to a saturated value after developing gradually. In our paper, an SIR epidemic model has been developed with t 1 as time delay in growth component and saturated Holling type III treatment rate. We have also considered pollution in the system. The paper is organised as follows: Section (2) describes the model, and we have discussed the characteristics of the model, section(3) examines the disease-free equilibrium and its stability analysis, in section (4) we have studied the endemicity of the model, section (5) explores the existence of hopf bifurcation, in Section (6) the numerical simulations of the system are discussed and analysed.

MODEL DESCRIPTION AND BASIC PROPERTIES
We have presented the model where we have supposed that entire populace is denoted by W(t) is distributed mainly into three compartments, viz. susceptible population X(t) (who are supposed to be infected by the disease), infected population Y(t) (who are capable of spreading the disease) and recovered population Z(t) (who have recovered from the disease and no more infectious). We assumed that the susceptible populace follows logistic growth with intrinsic growth rate r and environment carrying capacity K ( [9]), where t 1 is time delay and is derived by assuming that X X i.e. the net per capita rate of change is dependent on state of system t 1 (time units) (Hutchinson [7]). Model assumes saturated incidence rate of the form: where λ is the highest disease transmission rate and the factor 1/(1 + εY) represents the inhibition effect when the number of infectives is large.
Therapies, surgeries, medications etc., are important tools of treatment of disease. Proper treatment effectively reduces the infection present in the populace. These treatments have the capacity to influence the recovery rate of the infectives. Therefore, proper treatment rate must be considered in the model. In our model, we have considered Holling type-III treatment function [2,13] which is as follows: a is considered as the cure rate of infection and b is considered as the limitation rate in the treatment of infected.
Above assumptions leads to the following set of delay differential equations which includes the treatment function: Now in order to study the behavior of population under the pollutants we have considered certain assumptions: (i) Let us suppose that the pollutants penetrates into population through the intake of food and from the outside resource, (ii) Pollutants deprives as the time passes from the organism due to some metabolic processes or any other way.
Let Q be the constant exogenous rate of input of toxicant in environment, then the composition of equation with respect to the concentration of environmental pollutant and concentration of toxicants for an organism is given by: where A(t) is the concentration of contaminants present in environment, U(t) is the concentration of pollutant present in the organism, h denotes the rate by which toxicants are lost from the resource, a 1 denotes the rate of the resource pollutants per unit mass organism, η denotes level of concentration of pollutants in the environment, the mean rate of intake of food per unit mass organism is denoted by α 2 , d 1 is the uptake rate of toxicant in the food per unit mass organism, l 1 and l 2 are rate of net ingestion and depuration of pollutant in the organism. Thus our extended model is: (1) Now since the first two equations in the above system does not have Z in it, therefore our model reduces to: (2) where (ψ 1 (ϖ), ψ 2 (ϖ)) ∈ C([−t 1 , 0], R + 2 ).
C here represents the Banach space of continuous mapping in [−t 1 , 0] into R + 2 .
2.1. Basic Properties of the Model. For system (2) to be meaningful, all of (X(t), Y(t), Z(t)) must be nonnegative and r, K, λ , ε, µ, d, σ , a, b > 0. We show that all solutions of (2) will continue to remain positive for t ≥ 0 and also we will show that system is bounded.
Proof: As the system (2) can be written as: Then the model system (2) becomes: . It can be easily verified in the model (3) that whenever T (θ ) ∈ R + is considered such that X, Y are both equals to zero, [23]], using the lemma, any solution of (3) with T (θ ) ∈ C + , say T (t) = T (t, T (θ )), is such that T (t) ∈ R + 2 for all t greater than or equals to zero. Therefore, the solution of (2) exist in R + 2 and remain positive for all t greater than zero. Following lemma will show the boundedness of the above model.
is an attracting region for the disease transmission model which is positively invariant and is given by (2) with initial The solutions of (2) along with initial conditions are bounded in the region R + 2 and enters D as t → ∞. Moreover from the above lemma, we get that the model (2) is well formulated, both epidemiologically as well as mathematically .

STABILTY OF THE DISEASE-FREE EQUILIBRIUM
Here, we deal with the stability analysis of E = K(r−α 1 U * ) r , 0 (the disease-free equilibrium).
We begin by obtaining a key parameter R 0 which is the basic reproduction number. It is basically the threshold quantity for control of infection. It represents the number of infectives who got infected because of existing single infective [14]. The equilibrium solutions of system with no delay and of the system with delay are both same. (Tipsri and Chinviriyasit [18]). Here R 0 is computed as, The Jacobian of the system with no delay, evaluated at E = K(r−α 1 U * ) r , 0 = (X 0 , 0) is as follows: then, we have the following characteristic equation as: Thus, there are two eigenvalues: is clear that ζ 1 is always negative, while ζ 2 < 0 when R 0 < 1 and ζ 2 > 0 when R 0 > 1. We know by Hartman-Grobman's theorem [12] that if ζ 2 = 0, the solutions of the system (2) and its linearization equivalent near E. Then, we have the theorem as: Theorem 2: E is a stable node when R 0 < 1, and it is a saddle when R 0 > 1.

Analysis at
Here, we will analyse the stability of the undelayed model at E(X 0 , 0) and R 0 = 1.
Theorem 3: The system (2) will undergo forward bifurcation at disease free equilibrium E when R 0 crosses one. Proof: We can notice one of the eigenvalue of the system (2) vanishes and another eigenvalue is negative when evaluated at R 0 = 1 and λ = λ * . Here, stability behavior of equilibrium point at R 0 equals to one is done by center manifold theory [Sastry [17]] as linearization is not applicable. For this purpose we define X = x and Y = y, then the system (2) becomes: Let us denote λ = λ * as bifurcation parameter and J * as Jacobian matrix at R 0 = 1. which is as follows: T be the right and u = [u 1 , u 2 ] be the left eigenvector of jacobian corresponding to 0 eigenvalue. Then we have u 1 = 0, u 2 = 1, v 1 = −λ * K r , v 2 = 1 Now we evaluate the nonzero partial derivatives associated with the functions of system (4) at Then the bifurcation constants p 1 and q 1 from Castillo-Chavez and Song's [6] theorem 4.1, are as follows: Since q 1 is positive and p 1 is negative, the forward bifurcation occurs at E when R 0 crosses unity.

ENDEMICITY OF DISEASE
In epidemiology, endemicity of the infection and existence of the host populace are the two faces of persistence of infectious disease [8]. In our model (2), the disease persists in the system when (R 0 > 1).
In this segment, we focus on the existence of E * of (2) that is when disease persist in the system.
Let endemic equilibrium point be denoted by E * = (X * , Y * ). Substitute dX dt = 0 to find the value of X * , which is: Next substitute the value of X * in the equation dY dt = 0 to find the value of Y * at t 1 = 0 and we get, Theorem 4: In the system, when R 0 is greater than unity, there is either a unique or there are three endemic equilibria which are non-negative. If all equilibria are simple roots and if R 0 is less than unity then the endemic equilibrium does not exist in the system.
Then, by Descartes' rule [22], Q(Y) can have either three positive roots or Q(Y) can have a unique root. There exists unique endemic equilibrium point, if any of the H1-H4 holds, whereas H5-H8 conditions must satisfy, for the existence of three endemic equilibria. For R 0 < 1 it is easy to check that all the coefficients are greater than zero therefore, endemic equilibrium point does not exists in the system. Also when R 0 = 1 , the coefficient g 0 is zero and other coefficients g 1 , g 2 , g 3 and g 4 all are positive, then by fundamental theorem of algebra, Q(Y) can not have real root which is non-negative. Hence for R 0 ≤ 1, system (2) will not have any endemic equilibrium point however, when R 0 > 1 endemic equilibrium exists. Hence, we have the theorem.
For our paper, we will consider the case where any of H1-H4 holds, and R 0 > 1, then the system (2) admits a unique endemic equilibrium.

Local stability of Endemic Equilibrium Point.
From the previous section, it was noted that the system has unique E * (X * , Y * ) which is greater than zero. Now the local stability of E * will be analyzed. Characteristic equation of system at E * is as follows: The second degree transcedental equation is Theorem 5: E * is locally stable for t 1 = 0, if Kλ Y * ≥ rX * (1 + εY * ) and r(1 + εY * ) ≥ Kλ are both simultaneously satisfied.
Theorem 6: Let (E * ) of the system (2) exists for t 1 greater than zero, E * is locally stable when M and N are both greater than zero, where M and N are given by (19) and (20) respectively.
The t 1 j corresponding to s 0 can be obtained from equations (9) and (10) is as follows: )| ζ =ιs 0 > 0 Thus, we have the theorem as: Theorem 7: When t 1 ∈ (0,t 0 ), positive equilibrium of the system (2) is locally stable and when t 1 > t 0 , it is unstable, that is, the system will undergo hopf bifurcation at the positive equilibrium when t 1 = t 0 .

DIRECTION AND STABILITY OF HOPF BIFURCATION
In the past segment, we got certain conditions under which the given system goes under Hopf-Bifurcation, where time delay is considered as the being the critical parameter. In this part, by taking into account the normal form theory and the center manifold theorem, we will introduce the formula for deciding the direction of Hopf-Bifurcation and will get conditions for the stability of bifurcating periodic solutions, too. Since Hopf-Bifurcation happens at t 1 of t 1 , there exists a pair of pure imaginary root±ισ (t 1 ) of the characteristic equation.
Infact, we can take where, M 1 ,M 2 have already been given above, and δ (℘) is Dirac delta function.
Next, for ϕ ∈ C 1 ([0, 1], R 2 ), the adjoint operator A * of A can be defined as, , R 2 ) a bilinear inner product, in order to normalize the eigenvalues of A and A * can be defined as follows: where η(℘) = η(℘, 0), andφ is the complex conjugate of ϕ. It can be verified that the operators A and A * are adjoint operators with respect to this bilinear form. Thus, since ±ισ t 1 are eigenvalues of A(0), they are the eigenvalues of A * as well.
Suppose that b ∈ C 1 is an eigenvector of A(0) associated with iω 0 such that b(℘) = When ℘= 0, we can use (16)and (17) to obtain Similarly, suppose that the eigenvector, b * of A * corresponding to −iω 0 takes the form b * (s) = 1 D (1, b * ) T e iω 0 s . Then, −κ 2 κ 9 + iω 0 Now, in order to ensure that < b * , b >= 1, we will determine the valueD. By (19) we have, Hence,D In the remaining part of this section, we now compute the coordinates in order to describe the center manifold C 0 at µ = 0. Let e t be the solution of (18) when µ = 0.
Next, define Now, on the center manifold C 0 , we have wherez andz are local coordinates for the center manifold C 0 in the direction of b * andb * . We note that W is real if e t is real and we will be considering the real solutions only.
From (20) we have, Therefore, we can write (21) as, g(z,z) =b * (0)N 0 (z,z) = g 20z Substituting (18) and (21) intoẆ =ė t −żb −zb, we have, On the other hand, on the center manifold C 0 near the origin,Ẇ = Wzż +W¯zż Using (5.10) to compare the coefficients, we finally arrive at the following, Since, e t (t +℘) =W (z,z,℘) + z +zb, we have, (11) (℘)zz +W (2) (02) (℘)z Comparing the coefficients of the above equations with (23), we get, 2 ) Following the procedure same as that of [16], we get, where E 1 and E 2 are the vectors such that, 2 ) and I n is the identity matrix. Therefore, it follows that, Based on above, we can express g i j in terms of the parameters in the system. Therefore we get the following terms, where, (1)μ 2 determines the direction of Hopf bifurcation, for ifμ 2 > 0, the Hopf bifurcation will be supercritical, and ifμ 2 < 0, the Hopf bifurcation will be subcritical, and the bifurcating periodic solutions exist for t 1 > t 1 or t 1 < t 1 .
(2)β 2 determines the stability of the bifurcating periodic solutions, for ifβ 2 < 0, the bifurcating periodic solutions will be stable, and ifβ 2 > 0, the bifurcating periodic solutions will be unstable.
(3)T 2 determines the period of the bifurcating periodic solutions, for ifT 2 > 0, the period increases, and ifT 2 < 0, the period decreases.

NUMERICAL SIMULATION
To analyze the dynamical behaviour of the model, we have considered hypothetical and biologically feasible parametric values as given in Table 2. Using MATLAB, numerical simulations have been performed. We have assumed that X(0) = 10, Y(0) = 10, Z(0) = 10, a = 0.1 and b = 0.007.
First, considering the above values of parameters, when U * = 0 the value of disease-free equilibrium observed is E(10.0000, 0.0000) and R 0 = 0.0056 ( figure 1(a)). Similarly, when U * = 0.3, the steady state value for disease-free equilibrium point is E(7.0000, 0.0000) and R 0 = 0.0238( figure 1(b)). From both the figures 2 and 3 it is noted that system is locally stable for E as stated in theorem (2). Also it can be seen that, if there is an increase in level of pollution the disease-free equilibrium point is achieved at a lower value.
Therefore,it can be said thet pollution effects the reproductive function as well as the overall health of the individuals, so awareness should be increased and appropriate measures should be taken accordingly. From figures (5) and (6) it can be seen that increase in a diminishes the spread of diseases. From the above table it is evident that with the incline in pollution, it results in decline of the recovered populace but if we increase a i.e. the cure rate, the recovered populace increases. Hence, in the dynamics of the model, treatment plays a vital role but this may result in over usage of the resources and is not economically feasible and hence, we need to control the pollution.

CONCLUSION
In our study, an epidemic SIR model has been developed. The model includes time delay in the growth component, saturated type incidence rate along with Holling type III treatment rate. We have also included the effect of pollution in our model. Using R 0 , it was noted that E is locally stable when R 0 < 1 which indicates that the infection is no longer present in the environment at R 0 < 1. It was also noted that E is unstable when R 0 > 1, which clearly indicates that infection exists in the environment at R 0 greater than one. Using center manifold theory, it was shown that at R 0 = 1, system (2) undergo forward bifurcation. We have investigated that for time lag t 1 ≥ 0, E * i.e. the endemic equilibrium point is locally stable if conditions stated in theorem 4 and 5 holds true. Furthermore, hopf bifurcation at E * has also been discussed.
Through numerical simulations, we have observed that with the increase in the level of pollution we can see a decline in overall populace and in the recovered individuals. Also, it can be seen that the system is stable at lower value of logistic delay when the pollution increases, which is due to decrease in population size and perhaps may be due to the reason that reproduction capabilities of the population decreases. Thus, to remove the disease completely we need to find ways to reduce the level of pollution and choose favourable techniques to expand the capacity of treatment.

CONFLICT OF INTERESTS
The author(s) declare that there is no conflict of interests.