Modelling the transmission dynamics of two-strain Dengue in the presence awareness and vector control

Journal of Theoretical Biology 443 (2018) 82–91

Contents lists available at ScienceDirect

Journal of Theoretical Biology journal homepage: www.elsevier.com/locate/jtbi

Modelling the transmission dynamics of two-strain Dengue in the presence awareness and vector controlR Ting-Ting Zheng, Lin-Fei Nie∗ College of Mathematics and Systems Science, Xinjiang University, Urumqi 830046, PR China

a r t i c l e

i n f o

Article history: Received 26 July 2017 Revised 30 November 2017 Accepted 17 January 2018

MSC: 34A37 34D23 92D30 Keywords: Dengue virus Two-strain model Stability Sensitivity analysis Optimal control

a b s t r a c t In this paper, a mathematical model describing the transmission of two-strain Dengue virus between mosquitoes and humans, incorporating vector control and awareness of susceptible humans, is proposed. By using the next generation matrix method, we obtain the threshold values to identify the existence and stability of three equilibria states, that is, a disease-free state, a state where only one serotype is present and another state where both serotypes coexist. Further, explicit conditions determining the persistence of this disease are also obtained. In addition, we investigate the sensitivity analysis of threshold conditions and the optimal control strategy for this disease. Theoretical results and numerical simulations suggest that the measures of enhancing awareness of the infected and susceptible human self-protection should be taken and the mosquito control measure is necessary in order to prevent the transmission of Dengue virus from mosquitoes to humans.

1. Introduction Dengue fever is a vector-borne viral disease transmitted to humans through the bite of an infective female mosquito (for example, Aedes aegypti and Aedes albopictus, which are known as the principal vector of Dengue (Esteva and Vargas, 1999, 2003)). In recent decades, Dengue fever has been made an international public health concern because of its high morbidity and mortality, which happens in most tropical, subtropical and temperate countries. The symptoms of this disease are characterized by high fever, frontal headache, pain behind the eyes, joint pains, nausea, vomiting and other symptoms (Derouich and Boutayeb, 2006). Dengue virus is generally observed clearly in children and older adults, and people usually recover from Dengue fever in three to seven days (Derouich and Boutayeb, 2006). One recent estimate indicates 390 million Dengue infections per year (95% credible interval 284–528 million), of which 96 million (67–136 million) manifest clinically (with any severity of disease) (Bhatt et al., 2013).

R This research is partially supported by the National Natural Science Foundation of China (Grant Nos. 11461067 and 11771373), the Natural Science Foundation of Xinjiang (Grant Nos. 2016D01C046 and 2016D03022). ∗ Corresponding author. E-mail addresses: [email protected], [email protected] (L.-F. Nie).

https://doi.org/10.1016/j.jtbi.2018.01.017 0022-5193/© 2018 Elsevier Ltd. All rights reserved.

© 2018 Elsevier Ltd. All rights reserved.

With more than one-third of the world’s population living in areas at risk for infection of Dengue virus, how to prevent and control the spread of this disease has been one of the hot topic from many points of view, including medical scientists and mathematics (Centers for Disease Control and Prevention). So far, many mathematical models have been proposed to study the dynamic behaviors of Dengue transmission. For example, Esteva and Vargas (1999) proposed an SIR model for the transmission of Dengue fever with variable human population size, and found three threshold parameters which govern the existence of the endemic proportion equilibrium, the increase of the human population size, and the behaviours of the total number of human infective. And in Amaku et al. (2013), considered the impact of vector-control strategies on the human prevalence of Dengue virus. Garba et al. (2008) proposed a deterministic model for the transmission dynamics of Dengue fever, which allows transmission by exposed humans and mosquitoes, and obtained the existence and local asymptotical stability of the disease-free equilibrium for the basic reproduction number is less than unity. Besides, authors also discussed the phenomenon of backward bifurcation. Other examples also can be found in Blayneh et al. (2009); Cai and Li (2010); Cai et al. (2017); Esteva and Vargas (1998, 20 0 0); Rodrigues et al. (2010) and the research in this area is still going on.

T.-T. Zheng, L.-F. Nie / Journal of Theoretical Biology 443 (2018) 82–91

It is well known that there are five distinct serotypes of Dengue virus (DEN1, DEN2, DEN3, DEN4 and DEN5) according clinical data collected from the past many years (https://en.wikipedia.org/wiki/Dengue.fever). So, an individual reside in an endemic area can be there are five Dengue viruses infectious during his lifetime, one with each serotype. Epidemiological studies (Stech and Williams, 2008) support the recovery people re-infected with a different serotype face an increased risk of developing Dengue hemorrhagic fever and Dengue shock syndrome. A person recovery from infection by one of the five serotypes will provide lifelong homologous immunity, but this immunity confers only temporary cross-immunity againsts subsequent infection by the four other serotypes. Once infected, a mosquito remains infected for life because of its short life-span, transmitting the virus to susceptible individuals during probing and feeding. Therefore, taking into account Dengue fever pathogen diversity and the transmission mechanism, several investigators have studied the multistrain Dengue fever transmission models. Examples can be found in Feng and Velasco-Hernández (1997) where Feng et al. proposed an SIR vector transmitted disease with two pathogen strains to simulate the spread of Dengue virus in mosquitoes and between people, to study both the epidemiological trends of the disease and conditions that permit coexistence in competing serotypes. Sriprom et al. (2007) introduced a mathematical model describing the sequential transmission dynamics of Dengue virus infection in the presence of two serotypes of Dengue virus, and obtained the conditions governing the stabilities of disease-free equilibrium, boundary equilibria and interior equilibrium. More research can be found in Esteva and Vargas (2003); Kooi et al. (2013); Mishra and Gakkhar (2014); Nuraini et al. (2007) and the references therein. As is known to all, there is no specific treatment for Dengue, although there are a number of vaccines under development, including the tetravalent Dengvaxia which has been approved for use in several countries and recommended for introduction in ares with high endemicity (World Health Organization, 2016). Until the efficacy of this vaccine is properly established, however, Dengue fever control strategies are based on taking appropriate preventive measures. The main measures are mosquitoes reduction mechanisms and personal protection against exposure to mosquitoes. Mosquitoes reduction mechanisms entail the elimination of mosquitoes breeding sites (such as cleaning culverts, roadside ditches etc., disposing of solid waste properly and removing artificial man-made habitats) and adulticiding (killing of adult mosquitoes by spraying insecticides) (Mishra and Gakkhar, 2014; Yang and Ferreira, 2008). On the other hand, personal protection is based on preventing mosquitoes from biting human (by using window screens, long-sleeved clothes, insecticide treated materials, coils and vaporizers). Motivated by the above discussion, this paper is aimed at describes more effectively control strategies to control and eliminate Dengue virus. For all this, we propose a non-linear multistrain model to describe the dynamics of primary and secondary infection of Dengue fever, where two control strategies, susceptible human awareness and vector control are introduced to prevent the spread of Dengue virus. The paper is structured as follows. In Section 2, a Dengue transmission model with two serotypes viruses and control strategies is introduced. In Section 3, we calculate the basic reproduction number and analyze the stability of boundary equilibria for this model. The existence of the interior equilibrium and the persistence of this disease are given in Section 4. The sensitivity of threshold conditions, and the existence and uniqueness of the optimal control strategies are derived in Section 5 and Section 6, respectively. In Section 7, by carrying out some numerical simulations, we try to confirm the main theoretical results, to illustrate some key factors to prevent and control

83

the spread of Dengue fever. Some concluding remarks are present in the last section. 2. Model formulation In this section, we present a mathematical model to examine the transmission dynamics of Dengue virus between mosquitoes and humans. This multi-strain model considers the presence of two serotypes namely, serotype-1 and serotype-2 where 1 and 2 can be DEN1, DEN2, DEN3, DEN4 or DEN5. The populations involved in the transmission are human and mosquitoes, and N(t) and M(t) denote the total human and mosquitoes population size at time t, respectively. Consider the transmission mechanism of Dengue fever, we divide the population into eight compartments: susceptible humans denoted S(t); primary infective humans with serotype-i denoted Ii (t); recovered from serotype-i (and immune serotype-i, susceptible to serotype-j) humans, Ri (t); secondary infective humans with serotype-j, Yj (t); recovered from secondary infection of either of the serotypes, and immune to both serotypes, R(t); susceptible mosquitoes denoted U(t); infected mosquitoes with serotype-i, Vi (t), where i, j = 1, 2 and i = j. Based on the transmission rules of Dengue virus, a two-strains Dengue model with awareness and vector control is given by the following nonlinear differential equations

⎧ dS(t )  = ω − S(t ) 2i=1 αiVi (t ) − μS(t ) − mS(t ), ⎪ dt ⎪ ⎪ dIi (t ) ⎪ ⎪ = αi S(t )Vi (t ) − γi Ii (t ) − μIi (t ), ⎪ dt ⎪ ⎪ dRi (t ) ⎪ ⎪ ⎨ dt = γi Ii (t ) − σ j α jV j (t )Ri (t ) − μRi (t ), i = j, dYi (t ) = σi αiVi (t )R j (t ) − (γi + d + μ )Yi (t ), i = j, dt ⎪ 2 ⎪ dR(t ) ⎪ = ⎪ i=1 γiYi (t ) − μR (t ) + mS (t ), dt ⎪ ⎪ 2 ⎪ dU (t ) ⎪ = ω ⎪ 1 −U (t ) i=1 βi (Ii (t ) +Yi (t )) −cU (t ) − μ1U (t ), ⎪ ⎩ dVdi t(t ) = βi (Ii (t ) + Yi (t ))U (t ) − cVi (t ) − μiVi (t ). dt

(1)

Here, the susceptibility index to serotype-i σ i is a positive real number that may mimic either cross-immunity (0 < σ i < 1) or increased susceptibility (σ i > 1) by Antibody Dependent Enhancement (ADE), where the pre-existing antibodies to previous Dengue infection cannot neutralize but rather enhance the new infection. There is no definite argument about the value of parameters σ i in medicine and mathematics. In terms of σ i > 1, there are many works have been done by some scholars, see Feng and VelascoHernández (1997); Kooi et al. (2013, 2014) and references therein. Throughout this paper, we just investigate the case of 0 < σ i < 1. The meanings and possible values of other parameters of model (1) are given in Table 1. From the above model, we obtain N (t ) = S(t ) + I1 (t ) + R1 (t ) + Y2 (t ) + I2 (t ) + R2 (t ) + Y1 (t ) + R(t ), M (t ) = U (t ) + V1 (t ) + V2 (t ), and

dN (t ) = dt dM (t ) = dt

ω − μN (t ) − d (Y1 (t ) + Y2 (t )) ≤ ω − μN (t ), ω1 − (c + μ1 )M (t ).

It follows that M (t ) → ω1 /(c + μ1 ), N(t) ≤ ω/μ as t → ∞. Therefore, we only in the following areas

 = {(S(t ), I1 (t ), I2 (t ), R1 (t ), R2 (t ), Y1 (t ), Y2 (t ), R(t ), U (t ), V1 (t ), V2 (t ) ) ∈ R11 + | 0 < S (t ) + I1 (t ) + R1 (t ) + Y2 (t ) + I2 (t ) + R2 (t ) + Y1 (t ) + R(t ) ≤ ω/μ, 0 < U (t ) + V1 (t ) + V2 (t ) ≤ ω1 /(c + μ1 )} discuss the dynamic behavior of model (1), where R11 + := {(x1 , x2 , · · · , x11 )|xi ≥ 0, i = 1, 2, · · · , 11}.

84

T.-T. Zheng, L.-F. Nie / Journal of Theoretical Biology 443 (2018) 82–91 Table 1 Definition and possible values of the basic parameter for model (1). Param.

Description

Value

Source

ω ω1 1/μ 1/μ1

Recruitment rate of host Recruitment rate of vector Average lifespan of host Average lifespan of vector Disease induced death rate in host population Degree of awareness in host population Degree of effective vector control Rate of primary transmission of serotype-1, 2 infection to susceptible host Recovery rate of serotype-1, 2 Susceptibility index to serotype-1, 2 Transmission rate of infection from human to mosquito by serotype-1

Variable Variable 70years 14days Variable [0, 1) [0, 1) (1.3035e−5 , 2.0 0 09e−4 ) (0.1521, 0.4440) [0, 5] (2.8239e−5 , 3.6568e−4 )

– – Feng and Velasco-Hernández (1997) Feng and Velasco-Hernández (1997) – – – Pandey et al. (2013) Pandey et al. (2013) Feng and Velasco-Hernández (1997) Pandey et al. (2013)

d m c

α1 , α2 γ 1, γ 2 σ 1, σ 2 β1, β2

3. The analysis of boundary equilibria We attempt to calculate the threshold conditions that allow the invasion and persistence of different serotypes in humans. Analysis of model (1) reveals the existence of boundary equilibria (where no serotype or only one serotype is present), and interior equilibrium (where serotype-1, serotype-2 both exist). In this section, we analyze the stability properties of the former. Through calculation, it is easy to find that model (1) always exists the disease-free equilibrium is E0 (S¯∗ , 0, 0, 0, 0, 0, 0, R¯ ∗ , U¯ ∗ , 0, 0 ), where S¯∗ = ω/(μ + m ), R¯ ∗ = mω/μ(μ + m ) and U¯ ∗ = ω1 /(c + μ1 ). The basic reproduction numbers which determine whether Dengue virus can be invade a human population are threshold values, R0i , i = 1, 2, the average number of secondary infections produced by a single infective of serotype-i in a completely susceptible population. Now, we using the next generation matrix method calculate the basic reproduction number for model (1). Letting

⎛ 0 ⎜ 0 ⎜ 0 ⎜ F =⎜ 0 ⎜β ω ⎝ 1 1

μ1 + c

0

0 0 0 0 0

0 0 0 0

β1 ω 1 μ1 + c

β2 ω 1 μ1 + c

0

0 0 0 0 0

ωα1 μ+ m

0 0 0 0

β2 ω 1 μ1 + c

0

0

⎞

0 0 0

⎟ ⎟ ⎟ ⎟ ⎠

⎛ γ1 + μ 0 0 0 0 0 ⎜ 0 γ2 + μ ⎜ 0 0 γ1 + d + μ 0 V=⎜ ⎜ 0 0 0 γ2 + d + μ ⎝ 0 0

0 0

0 0

ωα2 μ+ m ⎟

⎜ ⎜ ⎜ −1 FV = ⎜ ⎜ ⎝

0 0 0 0

β1 ω 1 (μ1 +c )(μ+γ1 )

0

0 0 0 0 0 β2 ω 1 (μ1 +c )(μ+γ2 )

⎞

0 0 0 0 ⎟ 0 0 ⎟ ⎟, 0 0 ⎟ ⎠ μ1 + c 0 0 μ1 + c

β1 ω 1 (μ1 +c )(μ+d+γ1 )

0

The stability properties of the disease-free equilibrium E0 (S¯∗ , 0, 0, 0, 0, 0, 0, R¯ ∗ , U¯ ∗ 0, 0 ) are given in the next theorem. Theorem 1. For model (1), the disease-free equilibrium E0 is locally asymptotically stable if R0 < 1 and unstable if R0 > 1.

where,



R01 ,





R02 ,

ωα1 (μ+m )(μ1 +c )

0 0 0 0 0

0

⎞

0 0 0

⎟ ⎟ ⎟. ⎟ ⎠

ωα2 (μ+m )(μ1 +c ) ⎟

0 0 0 0

β2 ω 1 (μ1 +c )(μ+d+γ2 )

The basic reproduction number R0 is defined as the spectral radius (dominant eigenvalue) of matrix FV −1 , denoted by ρ (FV −1 ). That is,

R0 = max

(c + μ1 )[ωβ1 + (γ1 + μ )(c + μ1 )] , β1 [ω1 α1 + (μ + m )(c + μ1 )] (R01 − 1 )(μ + m )(c + μ1 )2 γ1 ∗ I1∗ = , R∗1 = I , β1 [ω1 α1 + (μ + m )(c + μ1 )] μ 1 (μ + γ1 )[α1 ω1 + (μ + m )(c + μ1 )] U∗ = , α1 [β1 ω + (c + μ1 )(μ + γ1 )] (R01 − 1 )(μ + m )(μ + γ1 )(c + μ1 ) m V1∗ = , R∗ = S ∗ . α1 [ωβ1 + (γ1 + μ )(c + μ1 )] μ S∗ =

(c + μ1 )[ωβ2 + (γ2 + μ )(c + μ1 )] , β2 [ω1 α2 + (μ + m )(c + μ1 )] (R02 − 1 )(μ + m )(c + μ1 )2 γ2 ∗ I2∗ = , R∗2 = I , β2 [ω1 α2 + (μ + m )(c + μ1 )] μ 2 (μ + γ2 )[α2 ω1 + (μ + m )(c + μ1 )] U ∗∗ = , α2 [β2 ω + (c + μ1 )(μ + γ2 )] (R02 − 1 )(μ + m )(μ + γ2 )(c + μ1 ) m V2∗ = , R∗∗ = S∗∗ . α2 [ωβ2 + (γ2 + μ )(c + μ1 )] μ

0

0 0 0 0

It is easy to solve that the serotype-2 free equilibrium of model (1) is E1 (S∗ , I1∗ , 0, R∗1 , 0, 0, 0, R∗ , U ∗V1∗ , 0 ), where

S∗∗ =

then,

⎛

R02

The state E1 will exist only for R01 > 1. Similarly, the serotype-1 free equilibrium of model (1) is E2 (S∗∗ , 0, I2∗ , 0, R∗2 , 0, 0, R∗∗ , U ∗∗ , 0, V2∗ ) and which exists for R02 > 1, where

and

0 0

α1 β1 ωω1 , (μ + m )(μ + γ1 )(μ1 + c )2 α2 β2 ωω1 = . (μ + m )(μ + γ2 )(μ1 + c )2

R01 =

0

0

Proof. The stability of the disease-free equilibrium E0 of model (1) is determine by the eigenvalues of the corresponding Jacobian matrix are given by, −μ(multipl icity3 ), −(μ + m ), −(γ1 + d + μ ), −(γ2 + d + μ ), λ1 , λ2 , λ3 and λ4 , where

λ1,2 =

−a1 ±



a21 − 4b1

2

,

λ3,4 =

−a2 ±



a22 − 4b2

2

T.-T. Zheng, L.-F. Nie / Journal of Theoretical Biology 443 (2018) 82–91

and

ai = (γi + μ ) + (c + μ1 ) > 0,

bi = (γi + μ )(c + μ1 )(1 − R0i ),

i = 1, 2. It is seen that the eigenvalues λ1 , λ2 , λ3 and λ4 have negative real parts when R0 < 1, and at least one of them have positive real part when R0 > 1. We deduce that E0 is locally asymptotically stable when R0 < 1, and unstable when R0 > 1.  Now, we investigate the global stability if the disease-free equilibrium E0 for special case m = 0. In this case, the disease-free equilibrium is (ω/μ, 0, 0, 0, 0, 0, 0, 0, ω1 /(c + μ1 ), 0, 0 ). To do so, we construct a suitable Liapunov functional, and then obtain the following theorem. Theorem 2. For model (1), if R0 < 1, then the disease-free equilibrium is globally asymptotically stable when m = 0 and σi ≤ 1(i = 1, 2 ).

⎛

−α

⎜ ⎜

(disease-free state to an endemic state). We obtain that if R0 ≤ 1 the disease will eventually disappear, and the disease will be able to invade the human population for R0 > 1. By analyzing the eigenvalues at E1 and E2 , we conclude the following theorems on the local asymptotical stability of these boundary equilibria. Theorem 3. Consider the case when R01 > 1. If

R02 <

R01 , −1 )γ1 σ2 (μ+m )(c+μ1 ) 1 + (Rμ01[ωβ 1 + (γ1 +μ )(c+μ1 )]

− (μ + m ) α1V1∗ 0 0 m 0

J11 = ⎜ ⎜

⎝

⎛−(γ + μ ) 2 γ2 ⎜ J22 = ⎜ β2U ∗ ⎝ 0 0

0 −(γ1 + μ )

γ1 β1U ∗

0 −β1U ∗

0 0 −μ 0 0 0

−α1 S α1 S∗ 0 − ( c + μ1 ) 0 0 ∗

0 0 0 0 −μ 0

Proof. For model (1), consider the boundary equilibrium E1 , the corresponding Jacobian matrix is



J ( E1 ) =

J11 0

0 −σ1 α1V1∗ − μ 0 σ1 α1V1∗ 0

α2 S

0 − ( c + μ1 ) 0 σ2 α2 R∗1

0 0 0 −(γ1 + d + μ ) 0

Proof. We consider a Lyapunov function L(t) := L(I1 (t), I2 (t), Y1 (t), Y2 (t), V1 (t), V2 (t)) as follows.

L(t ) = I1 (t ) + I2 (t ) + Y1 (t ) + Y2 (t ) +

α2 ω + V (t ), μ ( c + μ1 ) 2

α1 ω V (t ) μ ( c + μ1 ) 1

⎞

0 0 β2U ∗ 0 −(γ2 + d + μ )

⎞ ⎟ ⎟. ⎠

The eigenvalues of J(E1 ) are given by the eigenvalues of J11 and J22 . The eigenvalues of J11 are −μ(multipl icity2 ), −(c + μ1 )and the roots of the polynomial

λ3 + c1 λ2 + d1 λ + e1 = 0,

(3)

where

whose orbital derivative is given by dL(t ) dI1 (t ) dI2 (t ) dY1 (t ) dY2 (t ) = + + + + dt dt dt dt dt



J12 . J22

0 0 ⎟ ⎟ 0 ⎟ ∗ ⎟ β1 I1 ⎠ 0 ∗ −β1 I1 − (c + μ1 ) and

∗

(2)

then E1 is locally asymptotically stable. If (2) is reversed, then E1 is unstable.

where ∗ 1V1

85

α1 ω dV1 (t ) × μ ( c + μ1 ) dt

dV2 (t ) α2 ω + × μ ( c + μ1 ) dt   ω ≤ α1 S(t ) + R2 (t ) − V1 (t ) + α2 S(t ) + R1 (t ) μ ω − V (t ) − (γ1 + μ )(I1 (t ) + Y1 (t )) μ 2 α1 β1 ωω1 − (γ2 + μ )(I2 (t ) + Y2 (t )) + (I (t ) + Y1 (t )) μ ( c + μ1 ) 1 α2 β2 ωω1 + (I (t ) + Y2 (t )) μ ( c + μ1 ) 2   ω ω = α1 S(t ) + R2 (t ) − V1 (t ) + α2 S(t ) + R1 (t ) − V (t ) μ μ 2 + (R01 − 1 ) × (γ1 + μ )(I1 (t ) + Y1 (t )) + (R02 − 1 )(γ1 + μ )(I2 (t ) + Y2 (t )).

It is observed that the subset of  where dL(t)/dt < 0 when R01 < 1 and R01 < 1. So, we obtain that for R0 < 1, the diseasefree equilibrium is globally asymptotically stable when σ i ≤ 1 and m = 0.  The stability of the equilibrium points E1 and E2 , determine whether the disease changes from one state to another

c1 = [α1V1∗ + (μ + m )] + [β1 I1∗ + (c + μ1 )] + (γ1 + μ ), d1 = [α1V1∗ + (μ + m )][β1 I1∗ + (c + μ1 )] + [α1V1∗ + (μ + m )](γ1 + μ ) + [β1 I1∗ + (c + μ1 )](γ1 + μ ) − (c + μ1 )(μ + γ1 ) and

e1 = [α1V1∗ + (μ + m )][β1 I1∗ + (c + μ1 )](γ1 + μ ) − (c + μ1 ) × (μ + γ1 )(μ + m ). Observe that the coefficients of Eq. (3) is always positive. Also it can be verified that

c1 d1 > [α1V1∗ + (μ + m )][β1 I1∗ + (c + μ1 )](γ1 + μ ) > e1 . Therefore, the coefficients of Eq. (3) satisfy Routh–Hurwitz conditions, that is, the three roots of Eq. (3) have negative real parts. Eigenvalues of J22 are given by −(γ1 + d + μ ), −(σ1 α1V1∗ + μ ) and by the roots of the polynomial

λ3 + c2 λ2 + d2 λ + e2 = 0, where

c2 = (γ2 + μ ) + (c + μ1 ) + (γ2 + d + μ ) > 0, d2 = (γ2 + μ )(c + μ1 ) + (c + μ1 )(γ2 + d + μ ) + (γ2 + μ ) × (γ2 + d + μ ) − σ2 α2 β2 R∗1U ∗ − α2 β2 S∗U ∗

(4)

86

T.-T. Zheng, L.-F. Nie / Journal of Theoretical Biology 443 (2018) 82–91

ω1 xi , (c + μ1 )(1 + x1 + x2 ) ω1 ∗ = U (x1 + x2 + 1 )(c + μ1 ) ∗ = V i

and

e2 = (γ2 + μ )(c + μ1 )(γ2 + d + μ ) − σ2 α2 β2 R∗1U ∗ (γ2 + μ ) −

α2 β2 S∗U ∗ (γ2 + d + μ ).

where a = (c + μ1 )(μ + m ), b = μ(c + μ1 ), f = ω1 (α1 x1 + α2 x2 ) + ∗ a(1 + x1 + x2 ) and gi = σi αi ωxi + b(1 + x1 + x2 ). Substituting  Ii∗ , Y i ∗ ∗   in xi = βi (Ii + Yi )/(c + μ1 ), we obtain two equations for x1 and x2 as follow,

Hence

c2 d2 > (γ2 + μ )[(c+ μ1 )(γ2 +d+ μ ) − σ2 α

β αβ

αβ (γ

∗ ∗ ∗ ∗ 2 2 R1 U − 2 2 S U ] ∗ ∗ 2 2 2 R1 U 2 +

= (γ2 + μ )(c + μ1 )(γ2 + d + μ ) − σ −

μ)

α2 β2 S∗U ∗ (γ2 + μ ) > e2

[(ω1 α1 + a )x1 + (ω1 α2 + a )x2 ][bx1 + (σ2 α2 ω1 + b)x2 ]

and

e2 ≥ (γ2 + d+ μ )[(γ2 + μ )(c+ μ1 ) − σ2 α2 β2 R∗1U ∗ − α2 β2 S∗U ∗ ]. (5) Obviously, the coefficient c2 is always positive. Substituting S∗ , R∗1 , U∗ in (5), we can obtain that e2 is greater than zero if and only if

R02 <

R01

. −1 )γ1 σ2 (μ+m )(c+μ1 ) 1 + (Rμ01[ωβ 1 + (γ1 +μ )(c+μ1 )]

Further, d2 is greater than zero. Therefore, the coefficients of Eq. (4) satisfy Routh–Hurwitz conditions, that is, the eigenvalues of matrix J22 have negative real parts under the above condition. In conclusion, E1 is locally asymptotically stable if R01 > 1 and inequality (2) is satisfied. If the inequality is reversed, E1 is unstable.  Theorem 4. For R02 > 1, E2 is locally asymptotically stable if the inequality

R02

R01 <

(6)

−1 )γ2 σ1 (μ+m )(c+μ1 ) 1 + (Rμ02[ωβ 2 + (γ2 +μ )(c+μ1 )]

The proof of Theorem 4 is similar to Theorem 3, here omitted Remark 1. The inequality (2) or (6) is sufficient to determine whether it is possible for a serotype to invade a population. 4. The existence of interior equilibrium and persistence In this section, we will discuss the existence of the interior equilibrium and persistence of this disease. In fact, if two boundary equilibria E1 and E2 are unstable, then it is possible for the two serotypes to coexist. That is, model (1) exists a interior equilibrium ∗ , Y ∗ ,  ∗ , Y ∗ , R ∗ , U ∗ , V ∗ , V ∗ ) at which both serotypes reE3 ( S∗ ,  I1∗ , R I∗ , R 1 2 2 2 1 1 2 main endemic. Supposing that

R02 >

R02 1+

(R02 −1 )γ2 σ1 (μ+m )(c+μ1 ) μ[ωβ2 +(γ2 +μ )(c+μ1 )]

1+

(R01 −1 )γ1 σ2 (μ+m )(c+μ1 ) μ[ωβ1 +(γ1 +μ )(c+μ1 )]

f

,

∗ = aγ j σi αi ω1 R0 j (c + μ1 )x1 x2 , Y i β j f gi (γi + d + μ )

(7)

a R 0 i ( c + μ1 ) x i ∗  Ii = , βi f



aβ2 γ1 σ2 α2 ω1 R01 x β1 (γ2 + d + μ ) 1

+ ab(1 − R02 ) = 0, [(ω1 α1 + a )x1 + (ω1 α2 + a )x2 ][(σ1 α1 ω1 + b)x1 + bx2 ] + [b(ω1 α1 + a ) + a(σ1 α1 ω1 + b)(1 − R01 )]x1



+

b(ω1 α2 + a ) + ab(1 − R01 ) −



aβ1 γ2 σ1 α1 ω1 R02 x β2 (γ1 + d + μ ) 2

+ ab(1 − R01 ) = 0. It can be seen that these equations define two hyperbolas which have a unique intersection in the region



β1 ω β2 ω 1 = ( x 1 , x 2 ) | 0 < x 1 < , 0 < x2 < μ ( c + μ1 ) μ ( c + μ1 )



Remark 2. It is not easy to analyze the stability of the interior equilibrium E3 , since one must show that all eigenvalues of the 11 × 11 Jacobian matrix at E3 have negative real parts. This task can carry out by using mathematical software and the method of Reference (Chung and Lui, 2016), however, this progress very complex and no practical significance. Therefore, the stability of the interior equilibrium is not discussed here. Next, we show that model (1) is persistent when condition (7) holds. In the other words, we prove that all solution  x(t ) = (S(t ), I1 (t ), R1 (t ), Y1 (t ), I2 (t ), R2 (t ), Y2 (t ), R(t ), U (t ), V1 (t ), V2 (t )) starting in int() has the character that lim inft→∞  x(t ) is at a positive distance from the boundary of region . We denote

aγ1 γ2 ω1 σ1 α1 R02 (c + μ1 )x1 x2 + + β2 μ f g2 (γ1 + d + μ )

are invariant sets in ∂  and the vector field on ∂  − (N1 ∪ N2 ) points inward . The only ω-limit sets are E0 , E1 , and E2 in N1 ∪ N2 . From the above discuss, we know that when R0i > 1, E0 is unstable with stable manifold the subspace (x1 , 0, 0, 0, 0, 0, 0, x8 , x9 , 0, 0) and Ei is locally asymptotically stable in the region Ni , i=1,2. Proposition 1. If (7) holds, Ei (i = 1, 2 ) cannot be the ω-limit of any orbit in int(). Proof. For the boundary equilibrium E2 , consider the Lyapunov function L2 (t ) = o2 I2 (t ) + p2Y2 (t ) + q2V2 (t ), where

2 ∗ = aγi R0i (c + μ1 ) (1 + x1 + x2 )xi , R i βi f g j

∗ = aγ1 γ2 ω1 σ2 α2 R01 (c + μ1 )x1 x2 R β1 μ f g1 (γ2 + d + μ )

b(ω1 α1 + a ) + ab(1 − R02 ) −

N2 =  ∩ { x ∈ R11 + |x2 = x3 = x4 = x7 = x10 = 0}

R01

ω (c + μ1 )(1 + x1 + x2 )



+

N1 =  ∩ { x ∈ R11 + |x4 = x5 = x6 = x7 = x11 = 0},

,

∗ )/(c + μ ), i = 1, 2. Substituted into model and let xi = βi ( Ii∗ + Y 1 i (1), we have

 S∗ =

+ [b(ω1 α1 + a ) + a(σ2 α2 ω1 + b)(1 − R02 )]x2

if and only if condition (7) is satisfied.

is satisfied. If (6) is reversed, E2 is unstable.

R01 >

(i, j = 1, 2, i = j ),

1 q2 = 2

ω (c + μ1 )(1 + x1 + x2 ) , μf

o2 =





σ2 α2 R∗1 1 α2 S∗∗ + + , β2U ∗∗ (c + μ1 )(μ + γ2 ) (c + μ1 )(γ2 + d + μ )

1

μ + γ2

,

p2 =

1

γ2 + d + μ

whose orbital derivative is given by

,

T.-T. Zheng, L.-F. Nie / Journal of Theoretical Biology 443 (2018) 82–91

dL2 (t ) = dt



α2 S(t ) α2 (S(t ) − S∗∗ ) σ2 α2 (R1 (t ) − R∗1 ) + + 2(μ + γ2 ) 2(μ + d + γ2 ) 2(μ + γ2 )  σ2 α2 R1 (t ) c + μ1 + − V2 (t ) 2(γ2 + d + μ ) 2β2U ∗∗  U (t ) α2 β2U (t ) + + S∗∗ 2U ∗∗ 2(c + μ1 )(μ + γ2 )  σ2 α2 β2U (t ) ∗ + R − 1 (I2 (t ) + Y2 (t )). 2(c + μ1 )(μ + γ2 + d ) 1

Therefore, we can choose small enough positive constant ε and a neighborhood U1 of E2 , such that for any  x(t ) = (S(t ), I1 (t ), R1 (t ), Y2 (t ), I2 (t ), R2 (t ), Y1 (t ), R(t ), V1 (t ), V2 (t )) ∈ U1 ∩ int (), one have

α2 (S(t ) − S∗∗ ) σ2 α2 (R1 (t ) − R∗1 ) + > −ε , 2(μ + γ2 ) 2(μ + d + γ2 )

U (t ) 1 > − ε. 2U ∗∗ 2

However, from (7), it follows that there exists a neighborhood U2 of E2 , such that for  x(t ) ∈ U2 ∩ int (),

c + μ1 α2 S(t ) σ2 α2 R1 (t ) α2 β2U (t ) + − > ε, S∗∗ 2(μ+γ2 ) 2(γ2 + d + μ ) 2β2U ∗∗ 2(c + μ1 )(μ + γ2 ) +

1 σ2 α2 β2U (t ) R∗ − 1 > − + ε . 2(c + μ1 )(μ + γ2 + d ) 1 2

Therefore, if x ∈ U1 ∩ U2 ∩ int(), it follows that dL2 (t)/dt > 0 unless o2 = p2 = q2 = 0. On the other hand, the level sets of L2 (t) are subspaces o2 I2 + p2Y2 + q2V2 = k that go away from N2 as k increases. Since L2 (t) is increasing along orbits starting in U1 ∩ U2 ∩ int(). We conclude that they go away from E2 . We choose a similar Lyapunov function L1 (t ) = o1 I1 (t ) + p1Y1 (t ) + q1V1 (t ), where

q1 =

1 2





σ1 α1 R∗2 1 α1 S∗ + + , ∗ β1U (c + μ1 )(μ + γ1 ) (c + μ1 )(γ1 + d + μ )

1 o2 = , μ + γ1

1 p2 = , γ1 + d + μ

also can prove that the boundary equilibrium E1 cannot be the ωlimit of any orbit in int() if condition (7) holds. This proves the proposition.  Remark 3. The condition (7) is sufficient condition for the coexistence of serotype-1 and serotype-2. 5. Sensitivity analysis We just according to the logic to build model (1), there is only one get results. In fact, there is a lot of uncertainty, the basic reproductive number will change with the change of the parameter values. Sensitivity indices enable us to measure the relative change in a state variable when a model parameter changes. Therefore, in this section, we present a sensitivity analysis for R01 and R02 . First, we introduce in the definition of sensitivity analysis. Definition 1 (Nakul et al.(2008)). The normalized forward sensitivity index of a variable, u, that depends differentiable on a parameter, p, is defined as,

γ pu :=

∂u p × . ∂p u

Table 2 denotes sensitivity indices of model parameters to R0 , and the values of parameters for model (1) are fixed as: c = 0.35, m = 0.08, μ = 0.0 0 0 04 and μ1 = 0.07, respectively. Table 2 shows the results of the sensitivity analysis with respect to R01 and R02 . The results represent the relative amount of variation (expressed in perceptual variation) in the variable if we vary the parameters by 1%. Accordingly, a reduction of 1% in

87

Table 2 Sensitivity indices of R01 and R02 to parameters values for model (1). Variable

Parameter

Sensitivity index

Variable

Parameter

R01

β1 α1 ω1 μ1

1 1 1 −0.39995 −0.99502 −1.99977

R02

β2 α2 ω1 μ1

m c

Sensitivity index

m c

1 1 1 −0.39995 −0.99502 −1.99977

the mosquito control degree c increases R01 and R02 by 1.99977%; a reduction of 1% in the personal protection awareness degree m increases R01 and R02 by 0.99502% and a reduction of 1% in the mosquito death rate μ1 increases R01 and R02 by 0.39995%. In contrast, a reduction of 1% in the mosquito recruitment rate ω1 and in the transmission rates α i , β i decreases R01 and R02 by 1%. Our results indicate that the most effective strategy to reduce R01 and R02 increases the degree of effective vector control c (killing of adult mosquitoes by spraying). The second effective strategy is reducing the mosquito recruitment rate ω1 , which strategy can be obtained by increasing the mosquito control degree, eliminate mosquito breeding sites, or kill larvae before they become adults. In addition, decrease transmission rates α i and β i , through increases the degree of infected individuals personal protection. The next strategy suggested that increases the degree of susceptible humans self-protection m, by avoiding locations where mosquitoes are biting and using barrier methods such as window screens and long-sleeved clothing. Based on above discussion and analysis, we can conclude that Dengue fever control strategies are primarily relied on mosquitoes reduction strategies and personal protection in the absence of an effective vaccine. 6. Optimal control problem Optimal control techniques are most the of major use in developing the optimal strategies to against the spread of Dengue virus. In this section, we present an optimal problem for model (1) to find a suitable compromise between minimize number of the susceptible individuals, the infected individuals, the total mosquitoes and the cost of the campaign. The degrees of awareness in the host population and effective vector control are seen as a control variable to reduce or even eradicate the disease. That is, let m and c be the control variables in model (1) and 0 ≤ m(t), c(t) ≤ 1. The objective function is given by

min J (m, c ) =



0

+

tend

[ξ1 S(t ) + ξ2 I1 (t ) + ξ3 I2 (t ) + ξ4 M (t )

 ξ5 m2 (t ) + ξ6 c2 (t ) dt.

(8)

In equation (8), ξ 1 , ξ 2 , ξ 3 and ξ 4 denote weight constants of S(t), I1 (t), I2 (t) and M(t), respectively; ξ 5 , ξ 6 are weight constants for the degree of personal protective awareness and vector control, respectively and ξ 5 m2 (t), ξ 6 c2 (t) denote the costs associated with personal awareness and vector control, respectively. These costs come from many sources, for instance, cost is related to human originates from public health education to human populations, or cost is related to mosquitoes control originates from the insecticide application. Consider the control strategy set as follows:

= {(m, c )|m, c are Lebesgue measurable on [0, tend ], 0 ≤ m, c < 1}. The integrand of the objective function given by (8), ξ1 S(t ) + ξ2 I1 (t ) + ξ3 I2 (t ) + ξ4 M (t ) + ξ5 m2 (t ) + ξ6 c2 (t ) is convex in which is also convex and closed by definition. Since

88

T.-T. Zheng, L.-F. Nie / Journal of Theoretical Biology 443 (2018) 82–91

Fig. 1. The stability of disease-free equilibrium of model (1) with R01 ≈ 0.1538 < 1 and R02 ≈ 0.1428 < 1: (a) the blue curves and red curves represent infected human I1 (t ) + I2 (t ) + Y1 (t ) + Y2 (t ) and recovered individual R(t), respectively; (b) the three-dimensional diagram of infected serotype-1 human I1 (t ) + Y1 (t ), infected serotype-2 human I2 (t ) + Y2 (t ) and infected mosquitoes V1 (t ) + V2 (t ). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 2. The stability of boundary equilibrium E1 of model (1) with R01 ≈ 1.1147 > 1 and R02 ≈ 0.0073 < 1: (a) the red curves and blue curves represent susceptible human and primary infected human with serotype-1, respectively; (b) (I2 (t ), Y1 (t ) + Y2 (t ), V2 (t )). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

model (1) is linear in the control variables and bounded by a linear system in the state variables, then the conditions for the existence of optional control are satisfied (more details can be found in Fleming and Rishel (1975)). Consider the objective function J with (m, c) ∈ subject to model (1). There exists (m∗ , c∗ ) ∈ such that J (m∗ , c∗ ) = {min(J (m, c ))|(m, c ) ∈ }. Therefore, Pontryagins maximum principle (Pontryagin and Boltyanskii, 1987) can be used for the necessary conditions for this optimal control problem. We begin by forming the Hamiltonian function

H=

ξ1 S(t ) + ξ2 I1 (t ) + ξ3 I2 (t ) + ξ4 M (t ) + ξ5 m2 (t ) + ξ6 c2 (t ) +

11 

λi fi ,

i=1

where fi represent the right-hand side of model (1) of i-th state variable. In the following theorem, we derive the necessary conditions for the optimal control problem. Theorem 5. The problems (1) and (8) admit a unique optimal solution (Sˇ(t ), Iˇ1 (t ), Rˇ1 (t ), Yˇ2 (t ), Iˇ2 (t ), Rˇ2 (t ), Yˇ1 (t ), Rˇ(t ), Uˇ (t ), Vˇ1 (t ), Vˇ2 (t )) associated with an optimal control J∗ (m∗ , c∗ ) on [0, tend ], with a fixed final time tend . Moreover, there are adjoint function λ(t),

i = 1, 2, · · · , 11, satisfying

⎧ dλ1 (t ) ⎪ dt ⎪ ⎪ ⎪ ⎪ ⎪ d λ 2 (t ) ⎪ ⎪ dt ⎪ ⎪ d λ (t ) ⎪ ⎪ d3t ⎪ ⎪ d λ (t ) ⎪ ⎪ d4t ⎪ ⎪ d λ 5 (t ) ⎪ ⎪ dt ⎪ dλ6 (t ) ⎪ ⎨ dt

= −ξ1 + (λ1 − λ2 )α1Vˇ1 (t ) + (λ1 − λ5 )α2Vˇ2 (t ) + (λ1 − λ8 )m∗ + λ1 μ, = −ξ2 + (λ2 − λ3 )γ1 + (λ9 − λ10 )β1Uˇ (t ) + λ2 μ, = (λ3 − λ4 )σ2 α2Vˇ2 (t ) + λ3 μ, = (λ4 − λ8 )γ2 + (λ9 − λ11 )β2Uˇ (t ) + λ4 (d + μ ), = −ξ3 + (λ5 − λ6 )γ2 + (λ9 − λ11 )β2Uˇ (t ) + λ5 μ, = (λ6 − λ7 )σ1 α1Vˇ1 (t ) + λ6 μ, dλ7 (t ) = (λ7 − λ8 )γ1 + (λ9 − λ10 )β1Uˇ (t ) + λ7 (d + μ ), dt ⎪ dλ8 (t ) ⎪ = λ8 μ, ⎪ ⎪ dt ⎪ dλ9 (t ) ⎪ = −ξ4 + (λ9 − λ10 )β1 (Iˇ1 (t ) + Y1 ˇ(t )1 (t )) ⎪ dt ⎪ ⎪ ⎪ + (λ9 − λ11 )β2 (Iˇ2 (t ) + Yˇ2 (t )) + λ9 (c∗ + μ1 ), ⎪ ⎪ ⎪ d λ ( t ) 10 ⎪ = −ξ4 + (λ1 − λ2 )α1 Sˇ(t ) + (λ6 − λ7 )σ1 α1 Rˇ2 (t ) ⎪ dt ⎪ ⎪ ⎪ + λ10 (c∗ + μ1 ), ⎪ ⎪ dλ11 (t ) ⎪ ⎪ ⎩ dt = −ξ4 + (λ∗ 1 − λ5 )α2 Sˇ(t ) + (λ3 − λ4 )σ2 α2 Rˇ1 (t ) + λ11 (c + μ1 ).

(9)

the terminal conditions are λi (tend ) = 0 f or i = 1, 2, · · · , 11. Moreover, optimal control strategies m∗ and c∗ are given by

T.-T. Zheng, L.-F. Nie / Journal of Theoretical Biology 443 (2018) 82–91

89

Fig. 3. The existence of interior equilibrium and persistence of model (1) with R01 ≈ 1.9146 > 1 and R02 ≈ 1.7782 > 1. (a)-(c): the quantities of firstly infected human I1 (t), I2 (t), the quantities of secondly infected human Y1 (t), Y2 (t) and infected mosquitoes V1 (t), V2 (t), where the red curves and blue curves represent serotype-1 infection and serotype-2 infection in model (1), respectively; (d) the stability of interior equilibrium. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)





m∗ = min 1, max 0,





c∗ = min 1, max 0,

(λ1 − λ8 )Sˇ(t ) 2 ξ5



7. Numerical simulation and discussion

,

λ9Uˇ (t ) + λ10Vˇ1 (t ) + λ11Vˇ2 (t ) 2 ξ6

 .

Proof. Adjoint equations and transversality conditions can be obtained using Pontryagin’s Maximum Principle such that

dλ1 (t ) ∂H =− , dt ∂S dλ2 (t ) ∂H =− , dt ∂ I1

λ1 (tend ) = 0, λ2 (tend ) = 0,

···

···

dλ11 (t ) ∂H =− , dt ∂ V2

λ11 (tend ) = 0.

The optimal control strategies m∗ , c∗ can be solved from the optimality conditions

∂H = 0, ∂m

∂H =0 ∂c

This completes the proof.



To illustrate the theoretical results and the feasibility of control strategies, in this section, some numerical simulations are carried out using the Runge–Kutta method in the software MATLAB routines with different parameters values. The possible values of parameters for model (1) are listed in Table 1, and we fixed some basic model parameters as follows: ω = 8, ω1 = 40 0 0 0, μ = 4.0 × 10−5 and d = 0.0 0 0 01. First, according to Table 1, we choose model parameters α1 = 1.405 × 10−5 , α2 = 1.409 × 10−5 , β1 = 8.8 × 10−5 , β2 = 8.15 × 10−5 , σ1 = 0.8, σ2 = 0.6, γ1 = 1/7, γ2 = 1/7, μ1 = 0.2, m = 0.2 and c = 0.1. It is easy to calculate that the threshold values R01 ≈ 0.1538 < 1 and R02 ≈ 0.1428 < 1, then model (1) has a unique globally asymptotically stable disease-free equilibrium (S¯∗ , 0, 0, 0, 0, 0, 0, R¯ ∗ , U¯ ∗ 0, 0 ) by Theorem 1, which is verified here by Figs. 1(a) and 1(b). Next, we choose α1 = 8.475 × 10−5 , α2 = 1.025 × 10−5 , β1 = 8.115 × 10−5 , β2 = 1.025 × 10−5 , σ1 = 0.8, σ2 = 0.6, γ1 = 1/12, γ2 = 1/5, μ1 = 1/20, m = 0.05 and c = 0.635, to discuss the existence and stability of boundary equilibrium of model (1). It’s easy to calculate that parameter values satisfy all conditions of Theorem 3, that is, R01 ≈ 1.1242 > 1 and R02 < R01 /{1 + (R01 − 1 )γ1 σ2 (μ + m )(c + μ1 )/μ[ωβ1 + (γ1 + μ )(c + μ1 )]}. In this case, serotype-1 has a higher reproductive number and it

90

T.-T. Zheng, L.-F. Nie / Journal of Theoretical Biology 443 (2018) 82–91

Fig. 4. The quantities of infected and recovered human of model (1) with different values of m: (a) infected human I1 (t ) + I2 (t ) + Y1 (t ) + Y2 (t ); (b) recovered human R(t).

Fig. 5. The quantities of infected mosquitoes and human in model (1). (a)-(b): infected mosquitoes V1 (t ) + V2 (t ) and human I1 (t ) + I2 (t ) + Y1 (t ) + Y2 (t ) with different values of c, respectively.

needs a lower level of susceptible to start an epidemic, one caused by this serotype will be more probable than one epidemic caused by serotype-2, then serotype-1 will displace serotype-2 after a finite amount of time. Therefore, model (1) has a locally stable boundary equilibrium E1 . The plots in Figs. 2(a) and 2(b) show this theoretical result. In addition, the existence and stability of E2 is similar to E1 hence we omit it here. Further, we choose parameters α1 = 1.405 × 10−5 , α2 = 1.409 × 10−5 , β1 = 8.8 × 10−5 , β2 = 8.15 × 10−5 , σ1 = 0.8, σ2 = 0.6, γ1 = 1/7, γ2 = 1/7, μ1 = 0.07, m = 0.05 and c = 0.1. Numerical calculation follows that R01 ≈ 1.9146 > 1 and R02 ≈ 1.7782 > 1, that is, conditions (2) and (6) are invalid. This is indicating that both of the boundary equilibria E1 and E2 are unstable. As we can see, in Figs. 3(a)-3(c), various types of infected human and mosquitoes completely tend to some positive constants rather than zero, and solution curves from different initial values all tend to a point in the first quadrant, which is shown Fig. 3(d). Therefore, numerical simulations imply that model (1) exists a stable interior equilibrium. Finally, we focus our attention on what the significant effects of control parameters m and c to against the spread of Dengue fever. We fixed parameter values of model (1) as Fig. 3, change the control parameter m to be 0, 0.2, 0.4, 0.6 and 0.8, respectively. The

plots in Fig. 4(a) show that the density of the total infected individuals decreases with increasing personal protective awareness m. At the same time, the number of recovered human increases as the degree of personal protection for m increases, and which has a significant increase in the early phases. This is shown Fig. 4(b). Numerical simulations show that we can control the quantities of infected individuals by adjusting the control parameter m. On the other hand, we fixed the personal protective awareness rate m = 0.05 and changed the force of vector control c to be 0, 0.05, 0.1, 0.2, 0.4 and 0.8, respectively. Numerical simulations show that the quantities of infected mosquitoes and the total infected human was reduced remarkably when the force of vector control c increased. Which are shown in Figs. 5(a) and 5(b), respectively. These numerical simulations illustrate that the value of the control parameters m and c play very important roles in controlling the process of infectious disease. 8. Conclusion The dynamics of two-strains Dengue fever model with vector control and awareness of susceptible humans control strategies are systematically analyzed in this paper. The model, which incorporates essential elements of Dengue fever transmission and enables the assessment of various anti-Dengue preventive strate-

T.-T. Zheng, L.-F. Nie / Journal of Theoretical Biology 443 (2018) 82–91

gies. The threshold values of this model are obtained, R01 , R02 , which govern whether the disease dies out or not. The theoretical results show that the disease-free equilibrium is locally asymptotically stable if R0 = max{R01 , R02 } < 1, and the disease-free equilibrium is globally asymptotically stable when R0 < 1, m = 0 and σi ≤ 1, i = 1, 2. However, our numerical simulation shows that for the disease-free equilibrium is also globally asymptotically stable when R0 < 1, σi ≤ 1, i = 1, 2 and m = 0 (see Fig. 1). But, we can’t still obtain the global stability of the disease-free equilibrium under the condition R0 < 1, σi ≤ 1, i = 1, 2 and m > 0 by theoretical calculation though we try to choose various Lyapunov functions. Therefore, we put forward an interesting open question: If R0 < 1 and σi ≤ 1, i = 1, 2, then the disease-free equilibrium is globally asymptotically stable. The boundary equilibrium E1 is locally asymptotically stable if R01 > 1 and condition (2) is satisfied (see Fig. 2). And we have given sufficient conditions (7) for the coexistence of two different serotypes in the same population. From the expressions of the threshold values, it is easy to observe that there are many factors influence its values, for example, the rates of transmission and the rate of recruitment vector, and so on. Therefore, we put forward a sensitivity analysis of the impact of six parameters on the threshold values. The analysis shows that the more effective control measures in Dengue-endemic areas seem to increase the degree of personal protection and vector control. For individuals personal protection against mosquito bites represents the first line of defence for Dengue fever prevention. To do so, the aggressive control measures, such as, use clothing protection, insect repellents and eliminate mosquito breeding site by improved drainage and prevention of standing water, help in reducing the infection rate and increasing the removal rate of the human population. However, numerical simulation shows that the personal protective awareness of susceptible human and vector control are not obvious for reducing the infected human. Such as, the total infected human has a small degree decrease as the awareness of personal protection increase, from 0 to 0.2 to 0.4 to 0.6 and ultimately to 0.8 (see Fig. 4(a)). It seems that take the above two control strategies to prevent the transmission of Dengue fever in human and mosquitoes only is not sufficient to quickly reduce infected human. Based on the transmission mechanisms of Dengue virus, reducing the source of mosquitoes access to Dengue virus is very critical to prevent the transmission of Dengue fever. Therefore, if we want to quickly and effectively eliminate the infectious disease in Dengue outbreak areas, much taken a three-pronged approach by controlling vector, raising the personal protective awareness of susceptible human and encouraging the infected individuals to take selfprotective measures.

91

References Amaku, M., Coutinho, F.A.B., Raimundo, S.M., Lopez, L.F., Burattini, M.N., Massad, E., 2013. A comparative analysis of the relative efficacy of vector-control strategies against dengue fever. Bull. Math. Biol. 76, 697–717. Bhatt, S., Gething, P.W., Brady, O.J., Messina, J.P., Farlow, A.W., 2013. The global distribution and burden of dengue. Nature 496, 504–507. Blayneh, K., Cao, Y., Kwon, H.D., 2009. Optimal control of vector-borne diseases: treatment and prevention. Discrete Contin. Dyn. Syst. Ser. B. 11 (3), 587–611. Cai, L.M., Li, X.Z., 2010. Global analysis of a vector-host epidemic model with nonlinear incidences. Appl. Math. Comput. 217 (7), 3531–3541. Cai, L.M., Li, X.Z., Fang, B., Ruan, S.G., 2017. Global properties of vector-host disease models with time delays. J. Math. Biol. 74 (6), 1397–1423. Centers for Disease Control and Prevention, https://www.cdc.gov/dengue/index. html. Chung, K.W., Lui, R., 2016. Dynamics of two-strain influenza model with cross-immunity and no quarantine class. J. Math. Biol. 73 (6), 1–23. Derouich, M., Boutayeb, A., 2006. Dengue fever: mathematical modelling and computer simulation. Appl. Math. Comput. 177 (2), 528–544. Esteva, L., Vargas, C., 1998. Analysis of a dengue disease transmission model. Math. Biosci. 150 (2), 131–151. Esteva, L., Vargas, C., 1999. A model for dengue disease with variable human population. J. Math. Biol. 38 (3), 220–240. Esteva, L., Vargas, C., 20 0 0. Influence of vertical and mechanical transmission on the dynamics of dengue disease. Math. Biosci. 167 (1), 51–64. Esteva, L., Vargas, C., 2003. Coexistence of different serotypes of dengue virus. J. Math. Biol. 46 (1), 31–47. Feng, Z., Velasco-Hernández, J.X., 1997. Competitive exclusion in a vector-host model for the dengue fever. J. Math. Biol. 35 (5), 523–544. Fleming, W.H., Rishel, R.W., 1975. Deterministic and Stochastic Optimal Control. Springer, New York. Garba, S.M., Gumel, A.B., Bakar, M.R.A., 2008. Backward bifurcations in dengue transmission dynamics. Math. Biosci. 215, 11–25. Kooi, B.W., Aguiar, M., Stollenwerk, N., 2013. Bifurcation analysis of a family of multi-strain epidemiology models. J. Comput. Appl. Math. 252, 148–158. Kooi, B.W., Aguiar, M., Stollenwerk, N., 2014. Analysis of an asymmetric two-strain dengue model. Math. Biosci. 248 (1), 128–139. Mishra, A., Gakkhar, S., 2014. The effects of awareness and vector control on two strains dengue dynamics. Appl. Math. Comput. 246, 159–167. Nakul, C., James, M.H., Jim, M.C., 2008. Determining important parameters in the spread of malaria through the sensitivity analysis of a mathematical model. Bull. Math. Biol. 70, 1272–1296. Nuraini, N., Soewono, E., Sidarto, K.A., 2007. Mathematical model of dengue disease transmission with severe DHF compartment. Bull. Malays. Math. Sci. Soc. 30 (2), 143–157. Pandey, A., Mubayi, A., Medlock, J., 2013. Comparing vector-host and SIR models for dengue transmission. J. Math. Biol. 246 (2), 252–259. Pontryagin, L.S., Boltyanskii, V.G., 1987. The mathematical theory of optimal process. American Control Conference 27 (3), 36–40. Rodrigues, H.S., Monteiro, M.T.T., Torres, D.F.M., 2010. Dynamics of dengue epidemics using optimal control. Math. Comput. Modelling 52 (9–10), 1667–1673. Sriprom, M., Barbazan, P., Tang, I.M., 2007. Destabilizing effect of the host immune status on the sequential transmission dynamic of the dengue virus infection. Math. Comput. Modelling 45 (9–10), 1053–1066. Stech, H., Williams, M., 2008. Alternate hypothesis on the pathogenesis of dengue hemorrhagic fever (DHF)/dengue shock syndrome (DSS) in dengue virus infection. Exp. Biol. Med. 233 (4), 401–408. Wikipedia, https://en.wikipedia.org/wiki/Dengue.fever. World Health Organization, 2016. Factsheets no. 117. Dengue and severe dengue. Yang, H.M., Ferreira, C.P., 2008. Assessing the effects of vector control on dengue transmission. Appl. Math. Comput. 198 (1), 401–413.