Global Dynamics of an Epidemic Model with a Ratio-Dependent Nonlinear Incidence Rate

Abstract

We study an epidemic model with a nonlinear incidence rate which describes the psychological effect of certain serious diseases on the community when the ratio of the number of infectives to that of the susceptibles is getting larger. The model has set up a challenging issue regarding its dynamics near the origin since it is not well defined there. By carrying out a global analysis of the model and studying the stabilities of the disease-free equilibrium and the endemic equilibrium, it is shown that either the number of infective individuals tends to zero as time evolves or the disease persists. Computer simulations are presented to illustrate the results.

1. Introduction

In recent years, attempts have been made to develop realistic mathematical models for the transmission dynamics of infectious diseases. In modeling of communicable diseases, the incidence function has been considered to play a key role in ensuring that the models indeed give reasonable qualitative description of the transmission dynamics of the diseases [1, 2]. Some factors, such as media coverage, density of population, and life style, may affect the incidence rate directly or indirectly [3].

Let 𝑆(𝑡) be the number of susceptible individuals, let 𝐼(𝑡) be the number of infective individuals, and let 𝑅(𝑡) be the number of removed individuals at time 𝑡. After a study of the cholera epidemic spread in Bari in 1973, Capasso and Serio [4] introduced a saturated incidence rate 𝑔(𝐼)𝑆 into epidemic models. This is important because the number of effective contacts between infective individuals and susceptible individuals may saturate at high infective levels due to crowding of infective individuals or due to the protection measures by the susceptible individuals. If the function 𝑔(𝐼) is decreasing when 𝐼 is large, it can also be used to interpret the “psychological” effects: for a very large number of infectives the infection force may decrease as the number of infective individuals increases because in the presence of large number of infectives the population may tend to reduce the number of contacts per unit time. To incorporate the effect of the behavioral changes of the susceptible individuals, Liu et al. [5, 6] used a nonlinear incidence rate given by

𝑔(𝐼)=𝑘𝐼^𝑙1+𝛼𝐼^ℎ, (1.1) where the parameters 𝑘,𝑙, and ℎ are positive constants and 𝛼 is a nonnegative constant, 𝑘𝐼^𝑙 measures the infection force of the disease, and 1/(1+𝛼𝐼^ℎ) measures the inhibition effect from the behavioral change of the susceptible individuals when their number increases or from the crowding effect of the infective individuals. The special cases when 𝑙,ℎ, and 𝛼 take different values have been used by many authors. See, for example, Derrick and van den Driessche [7], Hethcote [8], Hethcote and Levin [9] and van den Driessche [10], Alexander and Moghadas [11], Ruan and Wang [12], Xiao and Ruan [13], and so forth.

Note that the infectious force 𝑔(𝐼) in (1.1) is only a function of infective individuals. However in the transmission of communicable diseases, it involves both infective individuals and susceptible individuals. Thus we think the infectious force ought to depend on the densities of both infective individuals and susceptible individuals and it should take the form 𝑔(𝐼,𝑆). Assuming that the infectious force is a function of the ratio of the number of the infectives to that of the susceptibles, in this paper, we consider the following infectious force function:

𝐼𝑔(𝐼,𝑆)=𝑔(𝑆)=𝑘(𝐼/𝑆)^𝑙1+𝛼(𝐼/𝑆)^ℎ. (1.2)

To describe the psychological or inhibitory effect from the behavioral change of the susceptible individuals when the number of infective individuals is very large, Xiao and Ruan [13] have considered the following SIRS model:

𝑑𝑆𝑑𝑡=𝑏−𝑑𝑆−𝑘𝐼𝑆1+𝛼𝐼²+𝛾𝑅,𝑑𝐼=𝑑𝑡𝑘𝑆𝐼1+𝛼𝐼²−(𝑑+𝜇)𝐼,𝑑𝑅𝑑𝑡=𝜇𝐼−(𝑑+𝛾)𝑅, (1.3) where the infectious force 𝑔(𝐼) takes the following form:

𝑔(𝐼)=𝑘𝐼1+𝛼𝐼², (1.4) which is the special case when 𝑙=1,ℎ=2 in (1.1). All the parameters are positive numbers. Also 𝑏 is the recruitment rate of the population, 𝑑 is the natural death rate of the population, 𝑘 is the proportionality constant, 𝛾 is the natural recovery rate of the infective individuals, 𝜇 is the rate at which recovered individuals lose immunity and return to the susceptible class, and 𝛼 is the parameter measures the psychological or inhibitory effect. Global analysis carried out by them shows that either the number of infective individuals tends to zero as time evolves or the disease persists.

Instead of the infectious force given by (1.4), in this paper, we incorporate the infectious force form (1.2) (when 𝑙=1,ℎ=2) into model (1.3). Our aim is to investigate the psychological effect of certain serious diseases on the community when the ratio of the number of infectives to that of the susceptibles is getting larger. As we can see from the following section, the model is not well defined at the origin and thus it cannot be linearized at it. We want to explore if it can produce some new dynamical behaviors different from the classical epidemic models.

This paper is organized as follows: in Section 2, we present the model and its simplified form. In Section 3, we carry out a qualitative analysis of the model. The existence of all kinds of equilibria and their stability results are derived. Finally, a brief discussion and some numerical simulations are given in Section 4.

2. The Model

Assuming that the infectious force takes the form (1.2) and 𝑙=1,ℎ=2, that is,

𝐼𝑔(𝐼,𝑆)=𝑔(𝑆)=𝑘(𝐼/𝑆)1+𝛼(𝐼/𝑆)²=𝑘𝐼𝑆𝑆²+𝐼², (2.1) consider the following SIRS model:

𝑑𝑆𝑑𝑡=𝑏−𝑑𝑆−𝑘𝐼𝑆²𝑆²+𝛼𝐼²+𝛾𝑅,𝑑𝐼=𝑑𝑡𝑘𝐼𝑆²𝑆²+𝛼𝐼²−(𝑑+𝜇)𝐼,𝑑𝑅𝑑𝑡=𝜇𝐼−(𝑑+𝛾)𝑅. (2.2) All the parameters are positive and have the similar biological meaning as in (1.3). Noting that system (2.2) is not well defined at the origin (0,0,0), we redefine that when (𝑆,𝐼,𝑅)=(0,0,0),

𝑑𝑆𝑑𝑡=𝑏,𝑑𝐼𝑑𝑡=0,𝑑𝑅𝑑𝑡=0. (2.3) With this assumption, it is easy to see that the first octant 𝑅³₊={(𝑆,𝐼,𝑅)∣𝑆,𝐼,𝑅≥0} is positively invariant for system (2.2), and that system (2.2) is continuous and satisfies the Lipschitz condition in 𝑅³₊. The following result shows that the solutions of system (2.2) are bounded, and hence, lie in a compact set and are continuable for all positive time.

Lemma 2.1. The plane 𝑆+𝐼+𝑅=𝑏/𝑑 is an invariant manifold of system (2.2), which is attracting in the first octant.

Proof. Summing up the three equations in (2.2) and denoting 𝑁(𝑡)=𝑆(𝑡)+𝐼(𝑡)+𝑅(𝑡), (2.4) we have 𝑑𝑁𝑑𝑡=𝑏−𝑑𝑁. (2.5) It is clear that 𝑁(𝑡)=𝑏/𝑑 is a solution of (2.5) and for any 𝑁(𝑡₀)≥0, the general solution of (2.5) is 1𝑁(𝑡)=𝑑[𝑏−(𝑏−𝑑𝑁(𝑡₀))𝑒^{−𝑑(𝑡−𝑡0)}]. (2.6) Thus, lim_𝑡→∞𝑏𝑁(𝑡)=𝑑, (2.7) which implies the conclusion.

In the following, we consider the existence of equilibria of system (2.2). For any values of parameters, model (2.2) always has a disease-free equilibrium 𝐸⁰₃=(𝑏/𝑑,0,0). To find the positive equilibria, set

𝑏−𝑑𝑆−𝑘𝐼𝑆²𝑆²+𝛼𝐼²+𝛾𝑅=0,𝑘𝑆²𝑆²+𝛼𝐼²−(𝑑+𝜇)=0,𝜇𝐼−(𝑑+𝛾)𝑅=0. (2.8) Define the basic reproduction number as follows:

𝑅₀=𝑘.𝑑+𝜇 (2.9) From (2.8), we have

[𝛼−(𝑅₀𝜇−1)(1+)𝑑+𝛾²]𝐼²+2𝑏𝑑(𝑅₀𝜇−1)(1+𝑏𝑑+𝛾)𝐼−²𝑑²(𝑅₀−1)=0. (2.10) It is easy to see from (2.10) that

(i)if 𝑅₀≤1, then there is no positive equilibrium; (ii)if 𝑅₀>1, then there is a unique positive equilibrium 𝐸^∗₃=(𝑆^∗,𝐼^∗,𝑅^∗), called the endemic equilibrium and given by

𝑆^∗=𝛼𝑅₀𝐼−1^∗,𝐼^∗=1𝛼/(𝑅₀𝑏−1)+(1+𝜇/(𝑑+𝛾))𝑑,𝑅^∗=𝜇𝐼𝑑+𝛾^∗. (2.11)

In the following section, we will study the properties of these equilibria and perform a global qualitative analysis of model (2.2).

3. Mathematical Analysis

It is clear that the limit set of system (2.2)-(2.3) is on the plane 𝑆+𝐼+𝑅=𝑏/𝑑. Thus, we focus on the reduced system:

𝑑𝑆𝑑𝑡=(𝑏+𝛾𝑏𝑑)−(𝑑+𝛾)𝑆−𝛾𝐼−𝑘𝐼𝑆²𝑆²+𝛼𝐼²,𝑑𝐼=𝑑𝑡𝑘𝐼𝑆²𝑆²+𝛼𝐼²−(𝑑+𝜇)𝐼,𝑑𝑆𝛾𝑑𝑡=𝑏(1+𝑑),𝑑𝐼𝑑𝑡=0,when(𝑆,𝐼)=(0,0), (3.1) confined to the set

𝔇_𝑆𝐼𝑏={(𝑆,𝐼)∣𝑆,𝐼≥0,𝑆+𝐼≤𝑑}. (3.2)

To be concise in notations, rescale (3.1) by

𝑆=(𝑏+𝛾𝑏/𝑑)𝑥𝑑+𝜇,𝐼=(𝑏+𝛾𝑏/𝑑)𝑦𝛾𝜏,𝑡=,𝑑+𝜇 (3.3) and, for simplicity, we still use variable 𝑡 instead of 𝜏. Then, we obtain

𝑑𝑥𝑑𝑡=1−𝑞𝑥−𝑦−𝑎𝑥²𝑦𝑥²+𝑝𝑦²,𝑑𝑦𝑅𝑑𝑡=[⁰𝑥²𝑥²+𝑝𝑦²−1]𝑦,𝑑𝑥𝑑𝑡=1,𝑑𝑦𝑑𝑡=0,when(𝑥,𝑦)=(0,0), (3.4) where

𝑞=𝑑+𝑟𝑑+𝜇,𝑝=𝛼(𝑑+𝜇)²𝛾²𝑘,𝑎=𝛾. (3.5)

Denote

𝔇_𝑥𝑦𝑎={(𝑥,𝑦)∣𝑥,𝑦≥0,𝑥+𝑅₀1𝑦≤𝑞}. (3.6) Obviously, 𝔇_𝑥𝑦 is a positively invariant set for system (3.4) and every solution initiated from the outside of 𝔇_𝑥𝑦 will approach or enter into and stay in it forever as 𝑡 is getting larger.

Note that the equilibrium 𝐸⁰₂=(1/𝑞,0) corresponds to the disease-free equilibrium 𝐸⁰₃ and the unique positive equilibrium 𝐸^∗₂=(𝑥^∗,𝑦^∗) corresponding to the endemic equilibrium 𝐸^∗₃ exists if and only if the condition 𝑅₀>1, where

𝑥^∗=1(1+𝑎/𝑅₀)(𝑅₀−1)/𝑝+𝑞,𝑦^∗=𝑅₀−1𝑝𝑥^∗. (3.7)

Note that we are interested in the dynamics of system (3.4) in the interior of the first quadrant. Thus, we can make a time scale change:

𝑑𝑡=(𝑥²+𝑝𝑦²)𝑑𝜏 (3.8) such that system (3.4) is equivalent to the following system in the interior of the first quadrant. For simplicity, we still use variable 𝑡 instead of 𝜏:

𝑑𝑥𝑑𝑡=𝑥²+𝑝𝑦²+Φ(𝑥,𝑦)≡𝑃(𝑥,𝑦),𝑑𝑦𝑑𝑡=Ψ(𝑥,𝑦)≡𝑄(𝑥,𝑦), (3.9) where

Φ(𝑥,𝑦)=−𝑞𝑥³−(𝑎+1)𝑥²𝑦−𝑝𝑞𝑥𝑦²−𝑝𝑦³,Ψ(𝑥,𝑦)=(𝑅₀−1)𝑥²𝑦−𝑝𝑦³. (3.10) We should point out that except for the equilibria 𝐸⁰₂ and 𝐸^∗₂, system (3.9) has another more equilibrium 𝐸⁰₂=(0,0) than (3.4) since the time scaling is (3.8).

In the following, we will study the properties of the equilibria and perform a global qualitative analysis of model (3.9).

3.1. Asymptotic Behavior of the System (3.9) at 𝐸⁰₂=(0,0)

The equilibrium 𝐸⁰₂ of system (3.9) is an isolated critical point of higher order. Obviously, system (3.9) is analytic in a neighborhood of the origin. By [14, Theorem 3.10 on page 79], any orbit of (3.9) tending to the origin must tend to it spirally or along a fixed direction, which depends on the characteristic equation of system (3.9). Here, we will show that if a solution orbit of (3.9) tends to the origin, then it must tend to it along a fixed direction. We will also determine the number of solution orbits of system (3.9) that tend to 𝐸⁰₂ along a fixed direction as 𝑡→+∞ or 𝑡→−∞ in the interior of the first quadrant by using the results in [14]. Hereafter, we refer to [14] for results and explanations of several notations involved.

First of all, we introduce the polar coordinates 𝑥=𝑟cos𝜃, 𝑦=𝑟sin𝜃 and define

𝐺(𝜃)=−sin𝜃𝑋₂(cos𝜃,sin𝜃), (3.11) where 𝑋₂(𝑥,𝑦) is homogeneous polynomial in 𝑥 and 𝑦 of degree 2 in the first equation of (3.9). Then the characteristic equation of system (3.9) takes the form

𝐺(𝜃)=−sin𝜃(cos²𝜃+𝑝sin²𝜃)=0. (3.12) Clearly, 𝐺(𝜃)=0 has two roots 𝜃₀=0 and 𝜃₁=𝜋. By the results in [14, section II.2], we know that no orbit of system (3.9) can tend to the critical point 𝐸⁰₂ spirally, and any orbit of system (3.9) approaching the origin must be along the directions 𝜃₀=0 or 𝜃₁=𝜋. Since we are only interested in the interior of the first octant, we just consider the direction 𝜃₀=0.

To determine if there exists an orbit of system (3.9) which tends to the origin along the direction 𝜃₀=0 as 𝑡 tends to +∞ or −∞, we have to compute the derivative of 𝐺(𝜃) and the function 𝐻(𝜃):

𝐺^′(𝜃)=−cos³𝜃+(2−3𝑝)sin²𝜃cos𝜃,𝐻(𝜃)=cos𝜃(cos²𝜃+𝑝sin²𝜃). (3.13)

Theorem 3.1. There exist 𝜀>0 and 𝑟₀>0 such that there exists a unique orbit of system (3.9) in the sector {(𝜃,𝑟)∶0≤𝜃<𝜀,0<𝑟<𝑟₀} that tends to 𝐸⁰₂ along 𝜃₀=0 as 𝑡→−∞.

Proof. Since 𝐺^′(𝜃₀)=−1 and 𝐻(𝜃₀)=1, by [14, Theorem 3.7 on page 70], there exist 𝜀>0 and 𝑟₀>0 such that there exists a unique orbit of system (3.9) in the sector {(𝜃,𝑟)∶0≤𝜃<𝜀,0<𝑟<𝑟₀} that tends to 𝐸⁰₂ along 𝜃₀=0 as 𝑡→−∞. The proof is thus completed.

In fact, the only orbit tending to 𝐸⁰₂ along 𝜃₀=0 as 𝑡→−∞ is 𝑦=0.

3.2. The Disease-Free Equilibrium 𝐸⁰₂=(1/𝑞,0) of System (3.9)

First,we give the following result regarding the nonexistence of periodic orbits in system (3.9), which implies the nonexistence of periodic orbits of system (2.2) by Lemma 2.1.

Theorem 3.2. System (3.9) does not have nontrivial periodic orbits.

Proof. Consider system (3.9) for 𝑥>0 and 𝑦>0. Take a Dulac function 1𝐷(𝑥,𝑦)=𝑦(𝑥²+𝑝𝑦²). (3.14) We have 𝜕(𝐷𝑃)+𝜕𝑥𝜕(𝐷𝑄)=𝜕𝑦−𝑞𝑥⁴𝑦−2𝑝𝑞𝑥²𝑦³−2𝑎𝑝𝑥𝑦⁴−𝑝²𝑞𝑦⁵𝑦²(𝑥²+𝑝𝑦²)²−2𝑅₀𝑝𝑥²𝑦(𝑥²+𝑝𝑦²)²<0. (3.15) The conclusion follows.

The following theorem shows the properties of the disease-free equilibrium 𝐸⁰₂=(1/𝑞,0).

Theorem 3.3. The disease-free equilibrium 𝐸⁰₂ of system (3.9) is (i)a global asymptotically stable node if 𝑅₀<1; (ii)a saddle if 𝑅₀>1; (iii)global asymptotically stable if 𝑅₀=1.

Proof. The Jacobian matrix of system (3.9) at 𝐸⁰₂ is −1𝐴=𝑞−𝑎+1𝑞²0𝑅₀−1𝑞², (3.16) whose determinant and trace are, respectively, det𝐴=1−𝑅₀𝑞³𝑅,tr𝐴=⁰−1−𝑞𝑞². (3.17) However 𝐸⁰₂ is locally stable (unstable) if the trace tr𝐴<0(tr𝐴>0) and det𝐴>0.
If 𝑅₀<1, we have det𝐴>0 and tr𝐴<0. Hence, 𝐸⁰₂ is locally asymptotically stable. Note that the interior equilibrium does not exist and 𝐸⁰₂ is unstable to the positive quadrant. By Theorem 3.2, the disease-free equilibrium 𝐸⁰₂ is global asymptotically stable.
If 𝑅₀>1, we have det𝐴<0. So 𝐸⁰₂ is a saddle.
If 𝑅₀=1, system (3.9) becomes 𝑑𝑥𝑑𝑡=𝑥²+𝑝𝑦²−𝑞𝑥³−(𝑎+1)𝑥²𝑦−𝑝𝑞𝑥𝑦²−𝑝𝑦³,𝑑𝑦𝑑𝑡=−𝑝𝑦³. (3.18) Translating the boundary equilibrium 𝐸⁰₂ of system (3.18) to the origin, system (3.18) can be written as 𝑑𝑥𝑥𝑑𝑡=−𝑞−𝑎+1𝑞²𝑦+Φ(𝑥,𝑦),𝑑𝑦𝑑𝑡=−𝑝𝑦³, (3.19) where Φ(𝑥,𝑦)=−2(𝑥²+(𝑎+1)𝑞²𝑥𝑦)−(𝑞𝑥³+(𝑎+1)𝑥²𝑦+𝑝𝑞𝑥𝑦²+𝑝𝑦³). (3.20) In order to study the stability of (0,0) of system (3.19), we firstly carry out a coordinate transformation 𝑢=𝑦,𝑣=(𝑥/𝑞)+((𝑎+1)/𝑞²)𝑦,𝑑𝑡=𝑞𝑑𝜏, for system (3.19), we obtain 𝑑𝑢𝑑𝜏=−𝑝𝑢³,𝑑𝑣𝑑𝜏=−𝑣+𝐹(𝑢,𝑣), (3.21) where 𝐹(𝑢,𝑣) is a polynomial in 𝑢 and 𝑣 of a degree not less than 2. Obviously, (3.21) has a center manifold 𝑣=ℎ(𝑢). Since the solution on the center manifold satisfies 𝑑𝑢𝑑𝜏=−𝑝𝑢³, (3.22) we know that the origin of (3.22) is asymptotically stable. By the theorem of center manifold, we know that the origin (0,0) of (3.21) is also asymptotically stable. That is to say, the boundary equilibrium (1/𝑞,0) of system (3.18) is asymptotically stable. Note that when 𝑅₀=1, system (3.18) has no positive equilibrium and 𝐸⁰₂ is unstable to the positive quadrant. By Theorem 3.2, we obtain that all the orbits of system (3.18) in the first quadrant tend to the equilibrium 𝐸⁰₂. The proof of Theorem 3.3 is thus completed.

3.3. The Endemic Equilibrium 𝐸^∗₂=(𝑥^∗,𝑦^∗) of System (3.9)

Theorem 3.4. Suppose 𝑅₀>1, then there is a unique endemic equilibrium (𝑥^∗,𝑦^∗) of (3.9), which is a global asymptotically stable node.

Proof. The Jacobian matrix of system (3.9) at the equilibrium (𝑥^∗,𝑦^∗) is 𝐴=2𝑥^∗−3𝑞(𝑥^∗)²−2(𝑎+1)𝑥^∗𝑦^∗−𝑝𝑞(𝑦^∗)²2𝑝𝑦^∗−(𝑎+1)(𝑥^∗)²−2𝑝𝑞𝑥^∗𝑦^∗−3𝑝(𝑦^∗)²2(𝑅₀−1)𝑥^∗𝑦^∗(𝑅₀−1)(𝑥^∗)²−3𝑝(𝑦^∗)², (3.23) whose determinant and trace are respectively, det𝐴=6𝑞𝑅₀(𝑅₀−1)(𝑥^∗)⁴+6𝑝(𝑎+𝑅₀)𝑥^∗(𝑦^∗)³+4𝑅₀(1−𝑅₀)(𝑥^∗)³={6𝑅₀(𝑅₀𝑎−1)[𝑞+(1+𝑅₀)𝑅₀−1𝑝]𝑥^∗+4𝑅₀(1−𝑅₀)}(𝑥^∗)³=2𝑅₀(𝑅₀−1)(𝑥^∗)³,tr𝐴=2𝑥^∗−3𝑞(𝑥^∗)²−2(𝑎+1)𝑥^∗𝑦^∗−𝑝𝑞(𝑦^∗)²+(𝑅₀−1)(𝑥^∗)²−3𝑝(𝑦^∗)²={[2−2𝑅₀−2𝑞−𝑞𝑅₀−2(𝑎+1)𝑅₀−1𝑝]𝑥^∗+2}𝑥^∗=[2(1−𝑅₀)−𝑞𝑅₀−2𝑎𝑅₀−1𝑝1(1−𝑅₀)](𝑥^∗)². (3.24) Since 𝑅₀>1, we have det𝐴>0 and tr𝐴<0. Hence, (𝑥^∗,𝑦^∗) is locally asymptotically stable. Note that 𝐸⁰₂ is unstable to the interior of the positive quadrant, and that the two stable manifolds of 𝐸⁰₂ are on the 𝑥-axis and one of its unstable manifolds points to the the positive quadrant. By Theorem 3.2, we know that the interior equilibrium (𝑥^∗,𝑦^∗) of system (3.9) is a global asymptotically stable node. This proves the theorem.

Summarizing Theorems 3.1–3.4, we have the following results on the dynamics of model (3.1).

Theorem 3.5. Let 𝑅₀=𝑘/(𝑑+𝜇) , then the following hold. (i)If 𝑅₀<1, then model (3.1) has a unique disease-free equilibrium (𝑏/𝑑,0), which is a global attractor in the first quadrant (see Figure 1). (ii)If 𝑅₀=1, then mode (3.1) has a unique disease-free equilibrium (𝑏/𝑑,0), which attracts all orbits in the interior of the first quadrant (see Figure 2). (iii)If 𝑅₀>1, then model (3.1) has two equilibria, a disease-free equilibrium (𝑏/𝑑,0) and an endemic equilibrium (𝑥^∗,𝑦^∗). The endemic equilibrium is a global attractor in the interior of the first quadrant (see Figure 3).

Figure 1 

Topological structure of system (3.1) at (𝑏/𝑑,0) when 𝑅0<1, where 𝑏=1,𝑑=1,𝑘=1,𝛾=1,𝛼=1,𝜇=1.

Figure 2 

Topological structure of system (3.1) at (𝑏/𝑑,0) when 𝑅0=1, where 𝑏=1,𝑑=1,𝑘=1.5,𝛾=0.5,𝛼=1,𝜇=0.5.

Figure 3 

Topological structure of system (3.1) at (𝑏/𝑑,0) and (𝑥∗,𝑦∗) when 𝑅0>1, where 𝑏=1,𝑑=1,𝑘=2,𝛾=0.5,𝛼=1,𝜇=0.5.

4. Discussion

Several nonlinear incidence rates have been proposed by researchers, see, for example, Capasso and Serio [4], Liu et al. [6], Derrick and van den Driessche [7], Hethcote and van den Driessche [10], and so forth Complex dynamics have been observed in epidemiological models with nonlinear incidence rate, such as the existence of multiple equilibria and limit cycles, various types of bifurcations including Hopf, saddle-node, homoclinic, and Bagdanov-Takens bifurcations, see also Ruan and Wang [12] and references cited therein.

In this paper, different from the classical nonlinear incident rate, we assume that the infectious force is a function of the ratio of the number of the infectives to that of the susceptibles which takes the form 𝑘𝐼𝑆/(𝑆²+𝛼𝐼²). Our aim is to investigate the psychological effect of certain serious diseases on the community when the ratio of the number of infectives to that of the susceptibles is getting larger. Note that the model with such infectious force is not well defined at the origin. By redefining system (2.2) at the origin (0,0,0), reducing it on the invariant manifold 𝑆+𝐼+𝑅=𝑏/𝑑, we get system (3.1). For system (3.1), after simplification and transformation in the time variable, we get its equivalent polynomial system (3.9). Notice that (3.9) has another more equilibrium 𝐸⁰₂=(0,0) than (3.4) since the time scaling is (3.8). Qualitative study on the the model shows that the model does not exhibit complicated dynamics as other epidemic models with other types of incidence rates reported by Liu et al. [6], Derrick and van den Driessche [7], Hethcote and Levin [9], Hethcote and van den Driessche [10], and Ruan and Wang [12], We observe that if the positive equilibrium of the system is locally asymptotically stable, then the system does not have any nontrivial positive periodic solutions. Numerical simulations illustrate our results.

Acknowledgment

Research is supported by the National Natural Science Foundation of China (10871129) and the Educational Committee Foundation of Shanghai (09YZ208).