# Stability and Hopf Bifurcation of a Vector-Borne Disease Model with Saturated Infection Rate and Reinfection

Authors: Zhixing Hu; Shanshan Yin; Hui Wang
Journal: Computational and Mathematical Methods in Medicine
 (2019)

---

## Abstract

This paper established a delayed vector-borne disease model with saturated infection rate and cure rate. First of all, according to the basic reproductive numberR0, we determined the disease-free equilibrium E0 and the endemic equilibrium E1. Through the analysis of the characteristic equation, we consider the stability of two equilibriums. Furthermore, the effect on the stability of the endemic equilibrium E1 by delay was studied, the existence of Hopf bifurcations of this system in E1 was analyzed, and the length of delay to preserve stability was estimated. The direction and stability of the Hopf bifurcation were also been determined. Finally, we performed some numerical simulation to illustrate our main results.


---

## Body

## 1. Introduction

Malaria is a vector-borne infectious disease [1], caused by parasites. It is popular in 102 countries and regions, especially in some countries in Africa, southeast Asia, and South America. In the 30s of this century, malaria spread throughout the country. Clinical symptoms and signs of this disease, such as typical periodic onset of malaria, secondary anemia, and spleen, can cause serious consequences, including dangerous malaria, malarial kidney disease, and black urine fever.The main way of transmission of malaria is the bite of an infected female anopheline mosquito. The mosquitoes would also be infected when uninfected mosquitoes bite infected people, and this transmission process has an incubation period [2]. The important feature of malaria is that the recovered immune system may establish immune memory for such antigens. It is this characteristic that greatly reduces the spread of malaria [3, 4]. Immune process is slow and, however, takes years or even decades [5]. As time goes by, the immune system gradually weakens, and at this time, reinfection likely occurs; therefore, considering the function of delay and immune system is necessary in the study of malaria.For the vector-borne diseases such as malaria, a large number of mathematical models have been created [2, 6, 7, 8, 9], most of which consider the local immunity and delay of the spread of malaria in the crowd. Different time delay has been used to describe the latent period in the course of disease transmission [7, 8, 9]. Local stability conditions for the equilibrium of a model with two time delays have been considered by Wan and Cui [8]. The global stability of the equilibrium has been studied for a vector-borne disease model with distributed delay by Cai et al. [10].Based on the above model, this paper considers a delayed vector-borne model with saturated infection rate and partial immunity to reinfection. We prove that the stability of this system can be changed by time delay and produce Hopf bifurcation, calculating the length of delay to preserve stability. Using the center manifold theorem [11] and norm theory, we determine the stability and bifurcation direction.


## 2. Model Formulation

N
1
t represented as the host population at time t is divided into three subclasses: the susceptible St, the infected It, and the recovered Rt. N2t represented as the vector population at time t is divided into two subclasses: the susceptible Tt and the infected Vt. The Hopf bifurcation was determined in a model with direct infection and delay by Wei et al. [9]. The mathematical formulation still needs improvements. We consider an improved model as follows:(1)dSdt=Λ1−bβ1StVt1+αVt−μ1St,dIdt=bβ1StVt1+αVt+σbβ1RtVt−μ1+γIt,dRdt=γIt−σbβ1RtVt−μ1Rt,dTdt=Λ2−bβ2Tt−τIt−τ−μ2Tt,dVdt=bβ2Tt−τIt−τ−μ2Vt,where Λ1 and Λ2 represent the recruitment rate of the host population and vector population, respectively. b represents the average number of bites per mosquito per day. The incidence rate bβ1StVt/1+αVt is the number of infections of the susceptible host caused by the infected vector, and α is the inhibitory effect rate caused by the infected vector. μ1 and μ2 represent the death rates of the host population and vector population, respectively. β1 is the infection rate from vector to human. σ0≤σ≤1 represents the degree of partial protection for recovered people given by a primary infection, where σ=0 represents complete protection and σ=1 represents no protection. γ is the per capita recovery rate of the infected host population. β2 represents the infection rate from human to vector. τ is the time delay, representing the incubation period in the vector population; that is to say, a susceptible vector that bites an infective host at time t−τ will become infective at time t.The model (1) meets the initial conditions:(2)Sθ=φ1θ,Iθ=φ2θ,Rθ=φ3θ,Tθ=φ4θ,Vθ=φ5θ,φiθ≥0,i=1,2,3,4,5,−τ≤θ≤0,where φ1θ,φ2θ,φ3θ,φ4θ,φ5θ∈C−τ,0,R+5 is the Banach space of continuous functions mapping the interval −τ,0 into R+5 with the topology of uniform convergence. The norm is defined as follows:(3)φ=sup−τ≤θ≤0φ1θ,φ2θ,φ3θ,φ4θ,φ5θ.Based on the fundamental theory of functional differential equations [12], it is easy to show that the solution of the model (1) with the initial condition (2) is unique and is nonnegative for all t≥0.By (1), we know that(4)Nk′t=Λk−μkNkt,k=1,2,and can solve it by using the integrating factor:(5)Nkt=Nk0e−μkt+Λkμk1−e−μkt.That is,(6)limt⟶∞Nkt=Λkμk,k=1,2.Form the limiting theory of differential equation [13], we can draw that model (1) is the equivalent of the following equation:(7)dSdt=Λ1−bβ1StVt1+αVt−μ1St,dIdt=bβ1StVt1+αVt+σbβ1Λ1μ1−St−ItVt−μ1+γIt,dVdt=bβ2Λ2μ2−Vt−τIt−τ−μ2Vt.Next, the model (7) can be studied in the invariant set:(8)Ω=S,I,V∈R+3∣0≤S+I≤Λ1μ1,0≤V≤Λ2μ2,S≥0,I≥0.Now, let us consider the existence of equilibrium.First, it is easy to show that system (7) always has a disease-free equilibrium E0=S0,I0,V0=Λ1/μ1,0,0. The endemic equilibrium E1=S1,I1,V1 satisfies the following equation:(9)Λ1−bβ1S1V11+αV1−μ1S1=0,bβ1S1V11+αV1+σbβ1Λ1μ1−S1−I1V1−μ1+γI1=0,bβ2Λ2μ2−V1I1−μ2V1=0.Form (9), we have S1=Λ11+αV1/bβ1V1+μ11+αV1 and I1=μ2V1/bβ2Λ2/μ2−V1, where V1 satisfies the following equation:(10)P2V2+P1V+P0=0,where(11)P2=σbβ1Λ1b2β1β2+μ1μ2bβ1+αμ1,P1=b2β1β2Λ1μ1−σbβ1Λ1b2β1β2Λ2μ2+σbβ1μ12μ2+μ1μ2μ1+γbβ1+αμ1,P0=μ12μ2μ1+γ1−R0,R0=b2β1β2Λ1Λ2μ1μ22μ1+γ,where V1=−P1+P12−4P0P2/2P2 and V2=−P1−P12−4P0P2/2P2 are the two roots of (10) since σ∈0,1, and we have(12)P1≥b2β1β2Λ1μ1+σbβ1μ12μ2+μ1μ2μ1+γbβ1+αμ1−bβ1Λ1b2β1β2Λ2μ2=b2β1β2Λ1μ1+σbβ1αμ12μ2+bβ1μ1μ2μ1+γ1−R0.Obviously,P2>0,P1>0 and P0>0 if R0<1, and P2>0 and P0<0 if R0>1. From the relationship between roots and coefficients, we know that V1 and V2 are both negative if R0<1 and V1 is positive if R0>1. According to the above discussion, we can obtain the following theorems.Theorem 2.1.
System (7) has the disease-free equilibrium E0 if R0<1. System (7) has the disease-free equilibrium E0 and an endemic equilibrium E1 if R0>1.


## 3. Stability of Equilibrium and Hopf Bifurcation

In this section, we study the stability of equilibrium and the existence of Hopf bifurcation of system (7).The characteristic equation of the linear approximate equation of the system (7) at equilibrium E=S,I,V is(13)λ+bβ1V1+αV+μ10bβ1S1+αV2σbβ1V−bβ1V1+αVλ+σbβ1V+μ1+γ−bβ1S1+αV2+Q0−bβ2Λ2μ2−Ve−λτλ+μ2+bβ2Ie−λτ=0,where Q=−σbβ1Λ1/μ1−S−I.


### 3.1. The Local and Global Stability of the Disease-Free Equilibrium

At the disease-free equilibriumE0, equation (13) can be expressed as follows:(14)λ+μ1λ2+μ1+γ+μ2λ+μ2μ1+γ−b2β1β2Λ1Λ2μ1μ2e−λτ=0.Obviously, equation (14) has a negative real root λ1=−μ1. To discuss the rest of the characteristic roots of (14), we consider the following equation:(15)λ2+μ1+γ+μ2λ+μ2μ1+γ−b2β1β2Λ1Λ2μ1μ2e−λτ=0.Whenτ=0, equation (15) is equivalent to(16)λ2+μ1+γ+μ2λ+μ2μ1+γ1−R0=0.By the using Routh–Hurwitz criterion, (16) has two eigenvalues with negative real parts if R0<1.Whenτ>0, then the roots of (15) can enter the right-half plane in the complex plane by crossing the imaginary axis as the delay τ increases.Letλ=ωiω>0 be a purely imaginary root of equation (15), then separating the real and imaginary parts yields(17)−ω2+μ2μ1+γ=b2β1β2Λ1Λ2μ1μ2cosωτ,μ2+μ1+γω=−b2β1β2Λ1Λ2μ1μ2sinωτ.Squaring and taking the sum of (17) yields(18)ω4+μ22+μ1+γ2ω2+μ2μ1+γ1−R0×μ2μ1+γ+b2β1β2Λ1Λ2μ1μ2=0.Equation (18) has no roots if R0<1. Therefore, we conclude that all eigenvalues of equation (14) have negative real parts.IfR0>1, let(19)fλ=λ2+μ1+γ+μ2λ+μ2μ1+γ−b2β1β2Λ1Λ2μ1μ2e−λτ,which implies that(20)f0=μ2μ1+γ1−R0<0,limλ⟶+∞fλ=+∞.By the continuity offλ and zero point theorem, fλ=0 has at least one positive root. So, the disease-free equilibrium E0 is unstable. Based on the results, we can draw the conclusion.Theorem 3.1.
For anyτ, the virus-free equilibrium E0 of the system (7) is locally asymptotically stable if R0<1, and it is unstable if R0>1.In fact, using a similar approach to the literature [14], we can know that E0 is globally asymptotically stable if R0<1. A detailed proof is given below.For a continuous and bounded functionft, we define(21)f∞≜limt⟶∞supft,f∞≜limt⟶∞infft.For system (7), any solution with the initial conditions is St,It,Vt, and we have(22)0≤S∞≤S∞≤∞,0≤I∞≤I∞≤∞,0≤V∞≤V∞≤∞.By the fluctuation lemma [15], we know that there is a sequence tn; when tn⟶∞, we have Stn⟶S∞ and S′tn⟶0n⟶∞. Substituting tn into the first equation of (7) yields(23)S′tn=Λ1−bβ1StnVtn1+αVtn−μ1Stn≤Λ1−μ1Stn.Let us take the limits on both sides:(24)limn⟶∞S′tn≤Λ1−μ1limn⟶∞Stn,i.e.,μ1S∞≤Λ1.Similarly,(25)μ1+γI∞≤1−σbβ1S∞V∞+σbβ1Λ1μ1V∞,μ2V∞≤bβ2Λ2μ2I∞.Combining (24) and (25), we know that(26)μ2V∞≤b2β1β2Λ1Λ2μ1μ2μ1+γV∞.SinceV∞ is the supremum of the function Vt, V∞≥0. If V∞>0, by using (26), μ2≤b2β1β2Λ1Λ2/μ1μ2μ1+γ, which contradicts R0<1. That is to say, V∞=0, which implies that limt⟶∞Vt=0. In the same way by using (25), we have limt⟶∞It=0. According to the limit theorem [13], we have limt⟶∞St=Λ1/μ1. Combined with the local asymptotic stability of E0, we can get the following theorem:Theorem 3.2.
For anyτ, the virus-free equilibrium E0 of system (7) is globally asymptotically stable if R0<1.


### 3.2. The Local Stability of the Endemic Equilibrium

From (13), the characteristic equation of linear approximate equation of the system (7) at the endemic equilibrium E1 is(27)λ3+A2λ2+A1λ+A0+B2λ2+B1λ+B0e−λτ=0,where(28)A2=bβ1V11+αV1+2μ1+σbβ1V1+γ+μ2,A1=bβ1V11+αV1+μ1σbβ1V1+μ1+γ+μ2+σbβ1V1+μ1+γμ2,A0=bβ1V11+αV1+μ1σbβ1V1+μ1+γμ2,B2=bβ2I1,B1=bβ2I1bβ1V11+αV1+σbβ1V1+2μ1+γ−bβ1S11+αV111+αV1−1μ2V1I1−μ2μ1+γ,B0=bβ2I1bβ1V11+αV1+μ1σbβ1V1+μ1+γ−bβ1S11+αV1bβ1V11+αV1+μ111+αV1−1×μ2V1I1−bβ1V11+αV1+μ1μ2μ1+γ+bβ1V11+αV1−σbβ1V1bβ1S11+αV12μ2V1I1.Whenτ=0, equation (27) is equivalent to(29)λ3+H2λ2+H1λ+H0=0,where(30)H2=A2+B2=bβ1V11+αV1+bβ2I1+2μ1+μ2+σbβ1V1+γ>0,H1=A1+B1=bβ1V11+αV1+μ1σbβ1V1+μ1+γ+μ2+σbβ1V1μ2+bβ2I1bβ1V11+αV1+σbβ1V1+2μ1+γ+bβ1S11+αV11−11+αV1μ2V1I1>0,H0=A0+B0=bβ1V11+αV1+μ1σbβ1V1+μ1+γbβ2I1+σbβ1V1μ2+1−σbβ1V1bβ1S11+αV12μ2V1I1+μ1bβ1S11+αV11−11+αV1μ2V1I1>0.Notice that(31)H2H1−H0=bβ1V11+αV1+bβ2I1+2μ1+μ2+σbβ1V1+γbβ1V11+αV1+μ1×σbβ1V1+μ1+γ+μ2+σbβ1V1μ2+bβ2I1bβ1V11+αV1+σbβ1V1+2μ1+γ+bβ1S11+αV11−11+αV1μ2V1I1−bβ1V11+αV1+μ1σbβ1V1+μ1+γbβ2I1+σbβ1V1μ2−1−σbβ1V1bβ1S11+αV12μ2V1I1−μ1bβ1S11+αV11−11+αV1μ2V1I1≥bβ1V11+αV1+μ1+μ2bβ1V11+αV1γ+μ1+σbβ1V1μ2+bβ2I1×σbβ1V1+μ1+γ+bβ1S11+αV11−11+αV1μ2V1I1bβ1V11+αV1+μ1−bβ1V11+αV1+μ1σbβ1V1+μ1+γbβ2I1+σbβ1V1μ2−1−σbβ1V1bβ1S11+αV12μ2V1I1−μ1bβ1S11+αV11−11+αV1μ2V1I1≥bβ1V11+αV1μ2γ+μ1−bβ1S11+αV12μ2V1I1bβ1V11+αV1+σbβ1V1bβ1S11+αV12μ2V1I1,(32)V1I1=μ1+γbβ1S1/1+αV1+σbβ1Λ1/μ1−S1−I1≤μ1+γbβ1S1/1+αV1=μ1+γ1+αV1bβ1S1.It follows that(33)H2H1−H0≥bβ1V11+αV1μ2γ+μ1−bβ1S11+αV12bβ1V11+αV1μ2μ1+γ1+αV1bβ1S1=bβ1V11+αV1μ2γ+μ1−bβ1V11+αV12γ+μ1μ2=bβ1V11+αV1μ2γ+μ11−11+αV1>0.By using the Routh–Hurwitz criterion, equation (29) only has eigenvalues with negative real parts if R0>1. We can obtain the following theorem:Theorem 3.3.
Forτ=0, the endemic equilibrium E1 of system (7) is locally asymptotically stable if R0>1.


### 3.3. Hopf Bifurcation

In this subsection, we devote to investigating the stability of the endemic equilibrium and the existence of Hopf bifurcation.Letλ=ωiω>0 be the root of equation (27), substituting it into equation (27) and separating the real and imaginary parts; we can obtain the following equation:(34)A2ω2−A0=B0−B2ω2cosωτ+B1ωsinωτ,ω3−A1ω=B1ωcosωτ−B0−B2ω2sinωτ.Squaring and taking the sum of (34) yields(35)ω6+p2ω4+p1ω2+p0=0,where(36)p2=A22−2A1−B22,p1=A12−2A0A2+2B0B2−B12,p0=A02−B02.Letx=ω2, then equation (35) is equivalent to(37)fx≜x3+p2x2+p1x+p0=0,then f′x=3x2+2p2x+p1. The two roots of equation 3x2+2p2x+p1=0 are(38)x∗=−p2+p22−3p13,x∗∗=−p2−p22−3p13.According to [16], the condition that equation (37) has positive roots is as follows:Lemma 3.4.
For equation (37),(i)
Ifp0<0, then equation (37) has at least one positive root(ii)
Ifp0≥0 and p22≤3p1, then equation (37) has no positive root(iii)
Ifp0≥0 and p22>3p1, then equtaion (37) has positive roots if and only if x∗>0 and fx∗≤0Based on Lemma 3.4, we concluded that if (ii) is set up, then the stability ofE1 will not change when τ increases. If equation (37) has a positive root, then the stability of E1 may change with the change in τ.Suppose that equation (37) has three positive roots, written as x1,x2, and x3. Then, equation (35) has positive roots ωk=xk k=1,2,3. By using (34),(39)cosωτ=B0−B2ω2A2ω2−A0+B1ωω3−A1ωB0−B2ω22+B1ω2.Define(40)τkj=1ωkarccosB0−B2ωk2A2ωk2−A0+B1ωkωk3−A1ωkB0−B2ωk22+B1ωk2+2jπ,where k=1,2,3,j=0,1,…. Obviously, ±ωki is a pair of pure virtual root of equation (27).Let(41)τ∗=τk00=mink∈1,2,3τk0,ω∗=ωk0.It follows thatλτ=ξ0τ+iωτ is the root of equation (35) satisfying ξ0τ∗=0 and ωτ∗=ω∗.Next, we verify the transversal condition. Differentiating the two sides of equation (35) with respect to τ, we have(42)dλdτ−1=−3λ2+2A2λ+A1λλ3+A2λ2+A1λ+A0+2λB2+B0λB2λ2+B1λ+B0−τλ,then(43)Redλdτ−1τ=τ∗=3ω∗4+2A22−2A1−B22ω∗2+A12−2A0A2+2B0B2−B12B1ω∗2+B0−B2ω∗22=f′ω∗2B1ω∗2+B0−B2ω∗22.Thus,(44)signdReλdττ=τ∗=signRedλdτ−1τ=τ∗=signf′ω∗2.Iff′ω∗2>0, the transversal condition is satisfied. Therefore, according to the above discussion and the Hopf bifurcation theorem of the differential equations [12], we can get the following result:Theorem 3.5.
(i)
Ifp0≥0 and Δ=p22−3p1≤0, then the endemic equilibrium E1 is locally asymptotically stable for all τ>0.
(ii)
Ifp0<0 or p0≥0,Δ>0,x∗>0, and fx∗≤0, then the endemic equilibrium E1 is locally asymptotically stable for 0<τ<τ∗, and if f′ω∗2≠0, the system (7) undergoes a Hopf bifurcation at E1 when τ=τ∗, where
(45)
τ
∗
=
1
ω
∗
arccos
B
0
−
B
2
ω
∗
2
A
2
ω
∗
2
−
A
0
+
B
1
ω
∗
ω
∗
3
−
A
1
ω
∗
B
0
−
B
2
ω
∗
2
2
+
B
1
ω
∗
2
.


## 3.1. The Local and Global Stability of the Disease-Free Equilibrium

At the disease-free equilibriumE0, equation (13) can be expressed as follows:(14)λ+μ1λ2+μ1+γ+μ2λ+μ2μ1+γ−b2β1β2Λ1Λ2μ1μ2e−λτ=0.Obviously, equation (14) has a negative real root λ1=−μ1. To discuss the rest of the characteristic roots of (14), we consider the following equation:(15)λ2+μ1+γ+μ2λ+μ2μ1+γ−b2β1β2Λ1Λ2μ1μ2e−λτ=0.Whenτ=0, equation (15) is equivalent to(16)λ2+μ1+γ+μ2λ+μ2μ1+γ1−R0=0.By the using Routh–Hurwitz criterion, (16) has two eigenvalues with negative real parts if R0<1.Whenτ>0, then the roots of (15) can enter the right-half plane in the complex plane by crossing the imaginary axis as the delay τ increases.Letλ=ωiω>0 be a purely imaginary root of equation (15), then separating the real and imaginary parts yields(17)−ω2+μ2μ1+γ=b2β1β2Λ1Λ2μ1μ2cosωτ,μ2+μ1+γω=−b2β1β2Λ1Λ2μ1μ2sinωτ.Squaring and taking the sum of (17) yields(18)ω4+μ22+μ1+γ2ω2+μ2μ1+γ1−R0×μ2μ1+γ+b2β1β2Λ1Λ2μ1μ2=0.Equation (18) has no roots if R0<1. Therefore, we conclude that all eigenvalues of equation (14) have negative real parts.IfR0>1, let(19)fλ=λ2+μ1+γ+μ2λ+μ2μ1+γ−b2β1β2Λ1Λ2μ1μ2e−λτ,which implies that(20)f0=μ2μ1+γ1−R0<0,limλ⟶+∞fλ=+∞.By the continuity offλ and zero point theorem, fλ=0 has at least one positive root. So, the disease-free equilibrium E0 is unstable. Based on the results, we can draw the conclusion.Theorem 3.1.
For anyτ, the virus-free equilibrium E0 of the system (7) is locally asymptotically stable if R0<1, and it is unstable if R0>1.In fact, using a similar approach to the literature [14], we can know that E0 is globally asymptotically stable if R0<1. A detailed proof is given below.For a continuous and bounded functionft, we define(21)f∞≜limt⟶∞supft,f∞≜limt⟶∞infft.For system (7), any solution with the initial conditions is St,It,Vt, and we have(22)0≤S∞≤S∞≤∞,0≤I∞≤I∞≤∞,0≤V∞≤V∞≤∞.By the fluctuation lemma [15], we know that there is a sequence tn; when tn⟶∞, we have Stn⟶S∞ and S′tn⟶0n⟶∞. Substituting tn into the first equation of (7) yields(23)S′tn=Λ1−bβ1StnVtn1+αVtn−μ1Stn≤Λ1−μ1Stn.Let us take the limits on both sides:(24)limn⟶∞S′tn≤Λ1−μ1limn⟶∞Stn,i.e.,μ1S∞≤Λ1.Similarly,(25)μ1+γI∞≤1−σbβ1S∞V∞+σbβ1Λ1μ1V∞,μ2V∞≤bβ2Λ2μ2I∞.Combining (24) and (25), we know that(26)μ2V∞≤b2β1β2Λ1Λ2μ1μ2μ1+γV∞.SinceV∞ is the supremum of the function Vt, V∞≥0. If V∞>0, by using (26), μ2≤b2β1β2Λ1Λ2/μ1μ2μ1+γ, which contradicts R0<1. That is to say, V∞=0, which implies that limt⟶∞Vt=0. In the same way by using (25), we have limt⟶∞It=0. According to the limit theorem [13], we have limt⟶∞St=Λ1/μ1. Combined with the local asymptotic stability of E0, we can get the following theorem:Theorem 3.2.
For anyτ, the virus-free equilibrium E0 of system (7) is globally asymptotically stable if R0<1.


## 3.2. The Local Stability of the Endemic Equilibrium

From (13), the characteristic equation of linear approximate equation of the system (7) at the endemic equilibrium E1 is(27)λ3+A2λ2+A1λ+A0+B2λ2+B1λ+B0e−λτ=0,where(28)A2=bβ1V11+αV1+2μ1+σbβ1V1+γ+μ2,A1=bβ1V11+αV1+μ1σbβ1V1+μ1+γ+μ2+σbβ1V1+μ1+γμ2,A0=bβ1V11+αV1+μ1σbβ1V1+μ1+γμ2,B2=bβ2I1,B1=bβ2I1bβ1V11+αV1+σbβ1V1+2μ1+γ−bβ1S11+αV111+αV1−1μ2V1I1−μ2μ1+γ,B0=bβ2I1bβ1V11+αV1+μ1σbβ1V1+μ1+γ−bβ1S11+αV1bβ1V11+αV1+μ111+αV1−1×μ2V1I1−bβ1V11+αV1+μ1μ2μ1+γ+bβ1V11+αV1−σbβ1V1bβ1S11+αV12μ2V1I1.Whenτ=0, equation (27) is equivalent to(29)λ3+H2λ2+H1λ+H0=0,where(30)H2=A2+B2=bβ1V11+αV1+bβ2I1+2μ1+μ2+σbβ1V1+γ>0,H1=A1+B1=bβ1V11+αV1+μ1σbβ1V1+μ1+γ+μ2+σbβ1V1μ2+bβ2I1bβ1V11+αV1+σbβ1V1+2μ1+γ+bβ1S11+αV11−11+αV1μ2V1I1>0,H0=A0+B0=bβ1V11+αV1+μ1σbβ1V1+μ1+γbβ2I1+σbβ1V1μ2+1−σbβ1V1bβ1S11+αV12μ2V1I1+μ1bβ1S11+αV11−11+αV1μ2V1I1>0.Notice that(31)H2H1−H0=bβ1V11+αV1+bβ2I1+2μ1+μ2+σbβ1V1+γbβ1V11+αV1+μ1×σbβ1V1+μ1+γ+μ2+σbβ1V1μ2+bβ2I1bβ1V11+αV1+σbβ1V1+2μ1+γ+bβ1S11+αV11−11+αV1μ2V1I1−bβ1V11+αV1+μ1σbβ1V1+μ1+γbβ2I1+σbβ1V1μ2−1−σbβ1V1bβ1S11+αV12μ2V1I1−μ1bβ1S11+αV11−11+αV1μ2V1I1≥bβ1V11+αV1+μ1+μ2bβ1V11+αV1γ+μ1+σbβ1V1μ2+bβ2I1×σbβ1V1+μ1+γ+bβ1S11+αV11−11+αV1μ2V1I1bβ1V11+αV1+μ1−bβ1V11+αV1+μ1σbβ1V1+μ1+γbβ2I1+σbβ1V1μ2−1−σbβ1V1bβ1S11+αV12μ2V1I1−μ1bβ1S11+αV11−11+αV1μ2V1I1≥bβ1V11+αV1μ2γ+μ1−bβ1S11+αV12μ2V1I1bβ1V11+αV1+σbβ1V1bβ1S11+αV12μ2V1I1,(32)V1I1=μ1+γbβ1S1/1+αV1+σbβ1Λ1/μ1−S1−I1≤μ1+γbβ1S1/1+αV1=μ1+γ1+αV1bβ1S1.It follows that(33)H2H1−H0≥bβ1V11+αV1μ2γ+μ1−bβ1S11+αV12bβ1V11+αV1μ2μ1+γ1+αV1bβ1S1=bβ1V11+αV1μ2γ+μ1−bβ1V11+αV12γ+μ1μ2=bβ1V11+αV1μ2γ+μ11−11+αV1>0.By using the Routh–Hurwitz criterion, equation (29) only has eigenvalues with negative real parts if R0>1. We can obtain the following theorem:Theorem 3.3.
Forτ=0, the endemic equilibrium E1 of system (7) is locally asymptotically stable if R0>1.


## 3.3. Hopf Bifurcation

In this subsection, we devote to investigating the stability of the endemic equilibrium and the existence of Hopf bifurcation.Letλ=ωiω>0 be the root of equation (27), substituting it into equation (27) and separating the real and imaginary parts; we can obtain the following equation:(34)A2ω2−A0=B0−B2ω2cosωτ+B1ωsinωτ,ω3−A1ω=B1ωcosωτ−B0−B2ω2sinωτ.Squaring and taking the sum of (34) yields(35)ω6+p2ω4+p1ω2+p0=0,where(36)p2=A22−2A1−B22,p1=A12−2A0A2+2B0B2−B12,p0=A02−B02.Letx=ω2, then equation (35) is equivalent to(37)fx≜x3+p2x2+p1x+p0=0,then f′x=3x2+2p2x+p1. The two roots of equation 3x2+2p2x+p1=0 are(38)x∗=−p2+p22−3p13,x∗∗=−p2−p22−3p13.According to [16], the condition that equation (37) has positive roots is as follows:Lemma 3.4.
For equation (37),(i)
Ifp0<0, then equation (37) has at least one positive root(ii)
Ifp0≥0 and p22≤3p1, then equation (37) has no positive root(iii)
Ifp0≥0 and p22>3p1, then equtaion (37) has positive roots if and only if x∗>0 and fx∗≤0Based on Lemma 3.4, we concluded that if (ii) is set up, then the stability ofE1 will not change when τ increases. If equation (37) has a positive root, then the stability of E1 may change with the change in τ.Suppose that equation (37) has three positive roots, written as x1,x2, and x3. Then, equation (35) has positive roots ωk=xk k=1,2,3. By using (34),(39)cosωτ=B0−B2ω2A2ω2−A0+B1ωω3−A1ωB0−B2ω22+B1ω2.Define(40)τkj=1ωkarccosB0−B2ωk2A2ωk2−A0+B1ωkωk3−A1ωkB0−B2ωk22+B1ωk2+2jπ,where k=1,2,3,j=0,1,…. Obviously, ±ωki is a pair of pure virtual root of equation (27).Let(41)τ∗=τk00=mink∈1,2,3τk0,ω∗=ωk0.It follows thatλτ=ξ0τ+iωτ is the root of equation (35) satisfying ξ0τ∗=0 and ωτ∗=ω∗.Next, we verify the transversal condition. Differentiating the two sides of equation (35) with respect to τ, we have(42)dλdτ−1=−3λ2+2A2λ+A1λλ3+A2λ2+A1λ+A0+2λB2+B0λB2λ2+B1λ+B0−τλ,then(43)Redλdτ−1τ=τ∗=3ω∗4+2A22−2A1−B22ω∗2+A12−2A0A2+2B0B2−B12B1ω∗2+B0−B2ω∗22=f′ω∗2B1ω∗2+B0−B2ω∗22.Thus,(44)signdReλdττ=τ∗=signRedλdτ−1τ=τ∗=signf′ω∗2.Iff′ω∗2>0, the transversal condition is satisfied. Therefore, according to the above discussion and the Hopf bifurcation theorem of the differential equations [12], we can get the following result:Theorem 3.5.
(i)
Ifp0≥0 and Δ=p22−3p1≤0, then the endemic equilibrium E1 is locally asymptotically stable for all τ>0.
(ii)
Ifp0<0 or p0≥0,Δ>0,x∗>0, and fx∗≤0, then the endemic equilibrium E1 is locally asymptotically stable for 0<τ<τ∗, and if f′ω∗2≠0, the system (7) undergoes a Hopf bifurcation at E1 when τ=τ∗, where
(45)
τ
∗
=
1
ω
∗
arccos
B
0
−
B
2
ω
∗
2
A
2
ω
∗
2
−
A
0
+
B
1
ω
∗
ω
∗
3
−
A
1
ω
∗
B
0
−
B
2
ω
∗
2
2
+
B
1
ω
∗
2
.


## 4. Estimation of the Length of Delay to Preserve Stability

In this section, we use a Nyquist criterion [17] to calculate the length of delay to preserve stability.Consider the system (7) and the space of the real continuous functions that is defined in −τ,+∞ and satisfied the initial conditions (2) in the interval −τ,0. Define(46)St=S1+Xt,It=I1+Yt,Vt=V1+Zt.Linearization system (7) at the endemic equilibrium E1 is expressed as follows:(47)dXtdt=−bβ1V11+αV1+μ1Xt−bβ1S11+αV12Zt,dYtdt=bβ1V11+αV1−σbβ1V1Xt−σbβ1V1+μ1+γYt+bβ1S11+αV12+σbβ1Λ1μ1−S1−I1Zt,dZtdt=−μ2Zt+bβ2Λ2μ2−V1Yt−τ−bβ2I1Zt−τ.By taking the Laplace transformation for (47), we can obtain(48)SLX−X0=−bβ1V11+αV1+μ1LX−bβ1S11+αV12LZ,SLY−Y0=bβ1V11+αV1−σbβ1V1LX−σbβ1V1+μ1+γLY+bβ1S11+αV12+σbβ1Λ1μ1−S1−I1LZ,SLZ−Z0=−μ2LZ+bβ2Λ2μ2−V1LYτ−bβ2I1LZτ,where(49)LYτ=∫0∞e−stYt−τdt=∫0τe−stYt−τdt+∫τ∞e−stYt−τdt.Lett=t1+τ, then(50)LYτ=∫−τ0e−st1+τYt1dt1+∫0∞e−st1+τYt1dt1=M1e−sτ+e−sτLY,where M1=∫−τ0e−stYtdt.Similarly,(51)LYτ=∫−τ0e−st1+τYt1dt1+∫0∞e−st1+τYt1dt1=M1e−sτ+e−sτLY,where M2=∫−τ0e−stZtdt.Thus, (48) can be written as(52)A−SILXLYLZ=B,where(53)A=−bβ1V11+αV1+μ10−bβ1S11+αV12bβ1V11+αV1−σbβ1V1−σbβ1V1+μ1+γbβ1S11+αV12+σbβ1Λ1μ1−S1−I10bβ2Λ2μ2−V1e−sτ−μ2−bβ2I1e−sτ,B=−X0−Y0−Z0+bβ2I1M2−bβ2Λ2μ2−V1M1e−sτ.The inverse Laplace transformation ofLXt, LYt, and LZt will have terms which exponentially increase with time if LXt, LYt, and LZt have poles with positive real parts. Thus, E1 is locally asymptotically stable if and only if all the poles of LXt, LYt, and LZt have negative real parts.By the method of [17] and the Nyquist criterion, the local asymptotic stability of E1 needs to satisfy the following two conditions:(54)ReFiμ0=0,(55)ImFiμ0>0,where(56)Fs=s3+A2s2+A1s+A0+B2s2+B1s+B0e−sτ,where μ0 is the smallest positive root of (54). Thus (54) and (55) can be written as(57)A2μ02−A0=B0−B2μ02cosμ0τ+B1μ0sinμ0τ,−μ03+A1μ0>B0−B2μ02sinμ0τ−B1μ0cosμ0τ.In order to estimate the length of delay to preserve stability, under the premise of ensuring stability, the following conditions need to be satisfied:(58)A2μ2−A0=B0−B2μ2cosμτ+B1μsinμτ,(59)−μ3+A1μ>B0−B2μ2sinμτ−B1μcosμτ.If (58) and (59) are satisfied simultaneously, they are sufficient conditions to guarantee stability. Our aim is to find an upper bound μ+ to μ0 independent of τ and then to estimate τ so that (59) holds true for all values of 0≤μ≤μ+ and in particular at μ=μ0.Sincecosμτ≤1 and sinμτ≤1, from equation (58), we have(60)A2μ2≤B0+B2μ2+B1μ+A0.Let(61)μ+=B1+B12+4A2−B2A0+B02A2−B2,obviously μ+ meets (60) and μ+≥μ0.From equation (59), we obtain(62)μ2<B1cosμτ+B2μ−B0μsinμτ+A1.SinceE1 is locally asymptotically stable for τ=0, the inequality (62) will continue to hold for sufficiently small τ and μ=μ0.On the basis of (58) and (62), we have(63)B1μ−A2B2μ+A2B0μsinμτ+B2μ2+A2B1−B01−cosμτ<B2μ2+A1A2+A2B1−A0−B0.Note the left-hand side of (63) is Φτ,μ and the right-hand side is ρ. By using inequality sinμτ≤μτ and 1−cosμτ=2sin2μτ/2≤μ2τ2/2, we can obtain(64)Φτ,μ≤Φτ,μ≤B1−A2B2μ+A2B0μμτ+B2μ2+A2B1−B0μ2τ22=B1−A2B2μ2+A2B0τ+B2μ2+A2B1−B0μ2τ22.Note the right-hand side of (64) is φτ,μ. Clearly, Φτ,μ≤φτ,μ≤φτ,μ+ when μ∈0,μ+. Thus, if φτ,μ+≤ρ≤K3, we have Φτ,μ+≤ρ≤K3. Let τ∗ be the positive root of φτ,μ+=K3, that is,(65)τ∗=12K1−K2+K22+4K1K3,where(66)K1=B2μ+2+A2B1−B0μ+22,K2=B1−A2B2μ+2+A2B0,K3=B2μ+2+A1A2+A2B1−A0−B0.Summarizing the above discussions, we have the following theorem.Theorem 4.1.
If0<τ<τ∗, then the Nyquist criterion holds true and τ∗ estimates the maximum length of the delay preserving the stability, where τ∗ satisfies (65).


## 5. Direction and Stability of the Hopf Bifurcation

We have obtained the conditions under which the Hopf bifurcation occurs atE1 of the system (7). This section will use the normal form theory and the center manifold theory to give the direction of the Hopf bifurcation and the stability of the bifurcating periodic solutions of system (7). We suppose that system (7) undergoes Hopf bifurcation at E1 for τ=τ˜τ˜=τkj. Let ±iω be a pair of conjugate pure virtual roots at E1 when τ=τ˜.Define(67)x1t=Sτt−S1,x2t=Iτt−I1,x3t=Vτt−V1,τ=τ˜+μ.Thus, system (7) is equivalent to the following functional differential equation in C=C−τ,0,R3.(68)dxdt=Lμxt+fμ,xt,where xt=x1t,x2t,x3tT∈R3. And Lμ:C⟶R3 and f:R×C⟶R3 satisfy(69)Lμϕ=τ˜+μB1ϕ0+τ˜+μB2ϕ−1,(70)fμ,ϕ=τ˜+μf1,f2,f3T,where(71)B1=−bβ1V11+αV1+μ10−bβ1S11+αV12bβ1V11+αV1−σbβ1V1−σbβ1V1+μ1+γbβ1S11+αV12+σbβ1Λ1μ1−S1−I100−μ2,B2=0000000bβ2Λ2μ2−V1−bβ2I1,ϕθ=ϕ1θϕ2θϕ3θT∈C,f1=−bβ11+αV12ϕ10ϕ30+αbβ11+αV13ϕ10ϕ320+αbβ1S11+αV13ϕ320−α2bβ11+αV14ϕ10ϕ330−α2bβ1S11+αV14ϕ330+⋯,f2=−σbβ1ϕ10ϕ30−σbβ1ϕ20ϕ30+bβ11+αV12ϕ10ϕ30−αbβ11+αV13ϕ10ϕ320−αbβ1S11+αV13ϕ320+α2bβ11+αV14ϕ10ϕ330+α2bβ1S11+αV14ϕ330+⋯,f3=−bβ2ϕ2−1ϕ3−1.Applying the Riesz representation theorem, there exists a3×3 matrix-valued function η⋅,μ:−1,0⟶R3×3, such that Lμϕ=∫−10dηθ,μϕθ,ϕ∈C. We choose(72)ηθ,μ=τ˜+μB1δθ+τ˜+μB2δθ+1,where δ is the Dirac delta function, meeting δθ=0θ≠0 and ∫−∞+∞δθdθ=1.We define forϕ∈C−1,0,R3,(73)Aμϕ=dϕθdθ,θ∈−1,0,∫−10dηs,μϕs,θ=0,Rμϕ=0,θ∈−1,0,fμ,ϕ,θ=0.Thus, (68) becomes(74)x˙t=Aμxt+Rμxt,where xt=xt+θ,θ∈−1,0.In order to construct coordinates to describe the integral manifold near the origin, we need to define inner product and the adjoint operatorA∗=A∗0 of A as follows:(75)A∗ψs=−dψsds,s∈0,1,∫−10dηTt,0ψ−t,s=0,ψ,ϕ=ψ¯T0ϕ0−∫θ=−10∫ξ=0θψ¯Tξ−θdηθϕξdξ,where ψ∈C0,1,R3 and ηθ=ηθ,0. Form the discussion in Section 3, we know that ±iωτ˜ are eigenvalues of A0. Thus they are also eigenvalues of A∗. Define qθ=1,q1,q2Teiθωτ˜ and q∗s=D1,q1∗,q2∗Te−isωτ˜ to be the eigenvectors of A0 and A∗ corresponding to the eigenvalues iωτ˜ and −iωτ˜, then(76)A0qθ=iωτ˜qθ,A∗q∗s=−iωτ˜q∗s.We can calculate that(77)q1=I1μ2+bβ2I1e−iωτ˜+ωibβ1V11+αV1+1+αV12ωibβ1S1μ2V1e−iωτ˜,q2=1+αV1V1S1+1+αV12ωibβ1S1,q1∗=bβ1V1+μ1−ωi1+αV1bβ1V1−σbβ1V11+αV1,q2∗=I1σbβ1V1+μ1+γ−ωibβ1V1+μ1−ωi1+αV1μ2V1eiωτ˜bβ1V1−σbβ1V11+αV1.According to (75), we know that(78)q∗s,qθ=D¯1,q¯1∗,q¯2∗1,q1,q2T−∫−10∫ξ=0θD¯1,q¯1∗,q¯2∗e−iξ−θωτ˜dηθ1,q1,q2Teiξωτ˜dξ=D¯1+q¯1∗q1+q¯2∗q2−∫−101,q¯1∗,q¯2∗θeiθωτ˜dηθ1,q1,q2T=D¯1+q¯1∗q1+q¯2∗q2+τ˜q1q¯2∗bβ2Λ2μ2−V1e−iωτ˜−τ˜q2q¯2∗bβ2I1e−iωτ˜=1,thus(79)D=11+q1∗q¯1+q2∗q¯2+τ˜eiωτ˜q¯1q2∗bβ2μ2V1I1−q¯2q2∗bβ2I1.Next using the same notation as in [11] and a computation process similar to that in [18], we compute the center manifold C0 at μ=0.Letxt be the solution of (74).Define(80)zt=q∗,xt,Wt,θ=xtθ−2Reztqθ.On the center manifoldC0, we have(81)Wt,θ=Wzt,z¯t,θ,where(82)Wzt,z¯t,θ=W20θz22+W11θzz¯+W02θz¯22+⋯.For the solutionxt∈C0 with μ=0, we can get(83)z′t=iωτ˜z+q¯∗θ,f0,Wz,z¯,θ+2Rezqθ=iωτ˜z+q¯∗0f0,Wz,z¯,θ+2Rezq0=iωτ˜z+q¯∗0f0z,z¯≜iωτ˜z+gz,z¯,where(84)gz,z¯=q¯∗0f0z,z¯=g20z22+g11zz¯+g02z¯22+g21z2z¯2+⋯.From (80) and (81), we obtain that(85)xtθ=x1tθ,x2tθ,x3tθT=Wt,θ+2Rezqθ=1,q1,q2Teiωτ˜θz+1,q¯1,q¯2Te−iωτ˜θz¯+W20θz22+W11θzz¯+W02θz¯22+⋯.From (70) and (84),(86)gz,z¯=q¯∗0f0z,z¯=q¯∗0f0,xt=τ˜D¯2k1+q¯1∗k2+q¯2∗k3z22+l1+q¯1∗l2+q¯2∗l3zz¯+2m1+q¯1∗m2+q¯2∗m3z¯22+2n1+q¯1∗n2+q¯2∗n3z2z¯2.It follows(87)g20=2τ˜D¯k1+q¯1∗k2+q¯2∗k3,g11=τ˜D¯l1+q¯1∗l2+q¯2∗l3,g02=2τ˜D¯m1+q¯1∗m2+q¯2∗m3,g21=2τ˜D¯n1+q¯1∗n2+q¯2∗n3,where(88)k1=−bβ1q21+αV12+αbβ1S11+αV13q22,k2=−σbβ1q21+q1+bβ1q21+αV12−αbβ1S11+αV13q22,k3=−bβ2q1q2e−2iωτ˜,l1=−bβ11+αV12q¯2+q2+2αbβ1S11+αV13q2q¯2,l2=−σbβ1q¯2+q2+q1q¯2+q¯1q2+bβ11+αV12q¯2+q2−2αbβ1S11+αV13q2q¯2,l3=−bβ2q1q¯2+q¯1q2,m1=−bβ1q¯21+αV12+αbβ1S11+αV13q¯22,m2=−σbβ1q¯21+q¯1+bβ1q¯21+αV12−αbβ1S11+αV13q¯22,m3=−bβ2q¯1q¯2e2iωτ˜,n1=−bβ11+αV12W20102q¯2+W20302+q2W1110+W1130+αbβ11+αV13q22+2q2q¯2+αbβ1S11+αV13W2030q¯2+2q2W1130−3α2bβ1S11+αV14q22q¯2,n2=−σbβ1W20102q¯2+W20302+q2W1110+W1130+W20202q¯2+W20302q¯1+W1120q2+W1130q1+bβ11+αV12W20102q¯2+W20302+q2W1110+W1130−αbβ11+αV13q22+2q2q¯2−αbβ1S11+αV13W2030q¯2+2q2W1130+3α2bβ1S11+αV14q22q¯2,n3=−bβ2W202−12q¯2eiωτ˜+W203−12q¯1eiωτ˜+W112−1q2e−iωτ˜+W113−1q1e−iωτ˜.By using (74) and (80), we can know(89)W′=xt′−z′q−z¯′q¯=AW−2Req¯∗0f0qθ,θ∈−1,0,AW−2Req¯∗0f0qθ+f0,θ=0,≜AW+Hz,z¯,θ,where(90)Hz,z¯,θ=H20θz22+H11θzz¯+H02θz¯22+⋯.From (81), (89), and (90), we can get(91)A−2iωτ˜W20θ=−H20θ,AW11θ=−H11θ.From (89), (91), and the definition of A, we can calculate(92)W20θ=ig20ωτ˜q0eiωτ˜θ+ig¯023ωτ˜q¯0e−iωτ˜θ+Γ1e2iωτ˜θ,(93)W11θ=−ig11ωτ˜q0eiωτ˜θ+ig¯11ωτ˜q¯0e−iωτ˜θ+Γ2,where Γi=Γi1,Γi2,Γi3T∈R3 i=1,2 are the three-dimensional vectors.Now, we determineΓ1 and Γ2. By (91) and the definition of A, we have(94)∫−10dηθW20θ=2iωτ˜W20θ−H20θ,∫−10dηθW11θ=−H11θ,where ηθ=ηθ,0. By (89), when θ=0,(95)Hz,z¯,θ=−q¯∗0f0q0−q∗0f¯0q¯0+f0=−gz,z¯q0−g¯z,z¯q¯0+f0.That leads to(96)H200=−g20q0−g¯02q¯0+2τ˜−bβ11+αV12+αbβ1S11+αV13q22G0−bβ2q1q2e−2iωτ,H110=−g11q0−g¯11q¯0+2τ˜−bβ11+αV12Req2+αbβ1S11+αV13q2q¯2G1−bβ2Req1q¯2,where(97)G0=−σbβ1q2+q1q2+bβ1q21+αV12−αbβ1S11+αV13q22,G1=bβ11+αV12−σbβ1Req2−σbβ1Req1q¯2−αbβ1S11+αV13q2q¯2.Since(98)iωτ˜I−∫−10eiθωτ˜dηθq0=0,−iωτ˜I−∫−10e−iθωτ˜dηθq¯0=0,we have(99)2iωτ˜I−∫−10e2iθωτ˜dηθE1=2τ˜−bβ11+αV12+αbβ1S11+αV13q22G0−bβ2q1q2e−2iωτ.We can calculate(100)Γ11=2Δ1−bβ11+αV12+αbβ1S11+αV13q220bβ1S11+αV12G02iω+σbβ1V1+μ1+γ−bβ1S11+αV12+Q1−bβ2q1q2e−2iωτ˜−bβ2Λ2μ2−V1e−2iωτ˜2iω+μ2+bβ2I1e−2iωτ˜,Γ12=2Δ12iω+bβ1V11+αV1+μ1−bβ11+αV12+αbβ1S11+αV13q22bβ1S11+αV12−bβ1V11+αV1+σbβ1V1G0−bβ1S11+αV12+Q10−bβ2q1q2e−2iωτ˜2iω+μ2+bβ2I1e−2iωτ˜,Γ13=2Δ12iω+bβ1V11+αV1+μ10−bβ11+αV12+αbβ1S11+αV13q22−bβ1V11+αV1+σbβ1V12iω+σbβ1V1+μ1+γG00−bβ2Λ2μ2−V1e−2iωτ˜−bβ2q1q2e−2iωτ˜,where(101)Q1=−σbβ1Λ1μ1−S1−I1,G0=−σbβ1q2+q1q2+bβ1q21+αV12−αbβ1S11+αV13q22,Δ1=2iω+bβ1V11+αV1+μ10bβ1S11+αV12−bβ1V11+αV1+σbβ1V12iω+σbβ1V1+μ1+γ−bβ1S11+αV12−σbβ1Λ1μ1−S1−I10−bβ2Λ2μ2−V1e−2iωτ˜2iω+μ2+bβ2I1e−2iωτ˜.Similarly, we have(102)∫−10dηθΓ2=−2τ˜−bβ11+αV12Req2+αbβ1S11+αV13q22G1−bβ2Req1q¯2.That leads to(103)Γ21=2Δ2−bβ11+αV12Req2+αbβ1S11+αV13q220bβ1S11+αV12G1σbβ1V1+μ1+γ−bβ1S11+αV12+Q1−bβ2Req1q¯2−bβ2Λ2μ2−V1μ2+bβ2I1,Γ22=2Δ2bβ1V11+αV1+μ1−bβ11+αV12Req2+αbβ1S11+αV13q22bβ1S11+αV12−bβ1V11+αV1+σbβ1V1G1−bβ1S11+αV12+Q10−bβ2Req1q¯2μ2+bβ2I1,Γ23=2Δ2bβ1V11+αV1+μ10−bβ11+αV12Req2+αbβ1S11+αV13q22−bβ1V11+αV1+σbβ1V1σbβ1V1+μ1+γG10−bβ2Λ2μ2−V1−bβ2Req1q¯2,where(104)G1=bβ11+αV12−σbβ1Req2−σbβ1Req1q¯2−αbβ1S11+αV13q2q¯2,Δ2=bβ1V11+αV1+μ10bβ1S11+αV12−bβ1V11+αV1+σbβ1V1σbβ1V1+μ1+γ−bβ1S11+αV12−σbβ1Λ1μ1−S1−I10−bβ2Λ2μ2−V1μ2+bβ2I1.We can calculateW200 and W110 from (92) and (93). The above analysis shows that gij in (84) can be represented by using the parameters in the system (7). Thus, we can get the following equation:(105)c10=i2ωτ˜g20g11−2g112−g0223+g212,μ2=−Rec10Reλ′τ˜,β2=2Rec10,T2=−Imc10+μ2Imλ′τ˜ωτ˜,where the signs of μ2 determine the direction of the Hopf bifurcation of the system (7), the signs of β2 determine the stability of the bifurcated periodic solution, and T2 determine the period of the bifurcation periodic solution. We can summarize the following theorem:Theorem 5.1.
For the model (7), the direction and the stability of the periodic solution of Hopf bifurcation is determined by using (105) when τ=τ˜.(i)
Ifμ2>0, then the Hopf bifurcation is supercritical and the bifurcation periodic solution exists for τ>τ˜ in a τ˜ neighborhood; if μ2<0, the Hopf bifurcation is subcritical and the bifurcation periodic solution exists for τ<τ˜ in a τ˜ neighborhood.(ii)
Ifβ2<0, then the bifurcation periodic solution is stable; if β2>0, the bifurcation periodic solution is unstable.(iii)
IfT2>0, then the period of the bifurcation periodic solution increases; if T2<0, the period decreases.


## 6. Numerical Simulations

In this section, we use numerical simulations to illustrate our result about the existence of Hopf bifurcation.The following parameters are selected:Λ1=8, b=0.29, β1=0.0033, α=0.06, μ1=0.0029, σ=0.48, γ=0.56, Λ2=9, β2=0.0059, and μ2=0.03. We can calculate that R0=80.2459 and E1=444.9472,450.8920,288.7707.Equation (37) has a positive root x=0.5779. We have τ∗=2.1386. From Theorem 3.4, we know that the endemic equilibrium E1 is locally asymptotically stable for 0<τ<τ∗ (Figure 1), and system (7) undergoes a Hopf bifurcation at E1 when τ=τ∗. At this time, we can calculate c10=−2.3788−0.8427i, μ2=9.8850>0,β2=−4.7576<0, and T2=1.7152>0. According to Theorem 5.1, system (7) can produce a stable periodic solution, the Hopf bifurcation is supercritical at τ∗, and the period of the bifurcation periodic solution increases (Figure 2).Figure 1
E
1 is stable when τ=2.12.
(a)
(b)Figure 2
E
1 undergoes Hopf bifurcation when τ=2.33.
(a)
(b)


## 7. Conclusions

In this paper, we discuss the dynamics of the vector-borne disease model with delay-saturated infection rate and reinfection. By calculation, we have the basic reproductive numberR0. Through R0, we determined the existence of disease-free equilibrium E0 and the endemic equilibrium E1. According to the characteristic equation of the equilibrium points and using the Routh–Hurwitz criterion, we obtained that if R0<1 the disease-free equilibrium will be stable, and the endemic equilibrium is locally asymptotically stable if R0>1 and in the absence of time delay. Furthermore, by the fluctuation lemma and the limit theory, we analyzed the global stability of the disease-free equilibrium. We find that the time delay does not affect the stability of the boundary equilibrium but can change the stability of E1 and lead to the occurrence of Hopf bifurcation. Then by using the Nyquist criterion, we get the maximum length of delay to preserve stability. Next, we found that the conditions for determining the direction and stability of bifurcating periodic solutions. Finally, the correctness of the main conclusion is verified by numerical simulation.



---

*Source: 1352698-2019-06-09.xml*