# A Mathematical Model of Malaria Transmission with Structured Vector Population and Seasonality

Authors: Bakary Traoré; Boureima Sangaré; Sado Traoré
Journal: Journal of Applied Mathematics
 (2017)

---

## Abstract

In this paper, we formulate a mathematical model of nonautonomous ordinary differential equations describing the dynamics of malaria transmission with age structure for the vector population. The biting rate of mosquitoes is considered as a positive periodic function which depends on climatic factors. The basic reproduction ratio of the model is obtained and we show that it is the threshold parameter between the extinction and the persistence of the disease. Thus, by applying the theorem of comparison and the theory of uniform persistence, we prove that if the basic reproduction ratio is less than1, then the disease-free equilibrium is globally asymptotically stable and if it is greater than 1, then there exists at least one positive periodic solution. Finally, numerical simulations are carried out to illustrate our analytical results.


---

## Body

## 1. Introduction

Malaria is an infectious disease caused by plasmodium parasite which is transmitted to humans through the bites of infectious female mosquitoes. According to the estimations of World Health Organization (WHO) in 2015, 3.2 billion persons were at risk of infection and 2.4 million new cases were detected with 438,000 cases of deaths. However sub-Saharan Africa remains the most vulnerable region with high rate of deaths due to malaria.To reduce the impact of malaria in the world, many scientific efforts were done including mathematical models construction. The first model of malaria transmission was developed by Ross [1]. According to Ross, if the mosquito population can be reduced to below a certain threshold, then malaria can be eradicated. Later, Macdonald did some modifications to the model and included superinfection. He showed that reducing the number of mosquitoes has little effect on the epidemiology of malaria in areas of intense transmission [2]. Nowadays, several mathematical models have been developed in order to reduce the malaria death rate in the world [3, 4]. In spite of the efforts made, it is still difficult to predict future malaria intensity, particularly in view of climate change.It must be noticed that transmission and distribution of vector-borne diseases are greatly influenced by environmental and climatic factors. Seasonality and circadian rhythm of mosquito population, as well as other ecological and behavioural features, are strongly influenced by climatic factors such as temperature, rainfall, humidity, wind, and duration of daylight [5]. Moreover, in most mathematical models, the mosquito life cycle is generally ignored because eggs, larvae, and pupae are not involved in the transmission cycle. That is a useful simplification of the system but unfortunately the results of these models do not predict malaria intensity in most endemic regions. Thus, it is necessary to consider the life cycle of mosquitoes and the seasonality effect, which are very important aspects of the dynamics of malaria transmission.Recently, Moulay et al. [6] have formulated a mathematical model describing the mosquito population dynamics which takes into account autoregulation phenomena of eggs and larvae stages. They have defined a threshold and proved that the growth of the mosquito population is governed by that threshold. Considering the climatic factors and the mosquitoes life cycle, we formulate a mathematical model describing the dynamics of malaria transmission. We analyze the impact of the model describing the mosquito population dynamics on the model of malaria transmission. Besides, by using the comparison theorem and the theory of uniform persistence, we, respectively, study the global stability of the nontrivial disease-free equilibrium [7–10] and the existence of positive periodic solutions.This paper is organized as follows. In Section2, we formulate the mathematical model of our problem. Section 3 provides the mathematical analysis of the model. Computational simulations are performed in Section 4 in order to illustrate our mathematical results. In the last section, Section 5, we conclude and give some remarks and future works.


## 2. Model Formulation

Motivated by the compartmental models in [6, 11], we derive an age-structured malaria model with seasonality to account for the cross infection between mosquitoes and humans. The human population is divided into four epidemiological categories representing the state variables:the susceptible class Sh,exposed class Eh,infectious class Ih, andrecovered class Rh (immune and asymptomatic, but slightly infectious humans). In the life cycle of anopheles, there are mainly two major stages: mature stage and aquatic stage. Therefore, we divide the mosquitoes population into these stages: immature and mature. The immature stage is divided in two compartments:eggs class E, larvae and pupae class L. In the mature stage, we have three compartments:the susceptible class Sm,exposed class Em, andinfectious class Im. At any time, the total number of humans and mature mosquitoes is given, respectively, by(1)Nht=Sht+Eht+Iht+Rht,(2)At=Smt+Emt+Imt.It is assumed throughout this paper that (H1)
all vector population measures refer to densities of female mosquitoes,(H2)
the mosquitoes bite only humans,(H3)
there is no vertical transmission of malaria,(H4)
all the new recruits are susceptibles.


### 2.1. Interactions between Humans and Mosquitoes

When an infectious mosquito bites a susceptible human, the parasite enters the body of the human with a probabilitycmh and the human moves into the exposed class Eh. Some time after, he leaves from class Eh to class Ih with rate α. Infectious humans migrate into the class Rh after acquisition of their immunity with rate rh. The immunized lose their immunity with rate γ if they do not have continuous exposure to infection. Humans leave the total population through natural death rate dh and malaria death rate dP.Similarly, when a susceptible mosquito bites an infectious human, it enters the classEm with a probability chm. Some time after, it leaves from class Em to infective class Im with rate νm where it remains for life. Mature mosquitoes leave the population through natural mortality dm.Using the standard incidence as in the model of Ngwa and Shu [4], we define, respectively, the infection incidence from mosquitoes to humans, kh(t), and from humans to mosquitoes, km(t):(3)kht=cmhβtImtNht,(4)kmt=chmβtIhtNht+c-hmβtRhtNht.Furthermore, using the above assumptions, we obtain the transfer diagram (Figure 1) of the model.Figure 1
The dashed arrows indicate the direction of the infection and the solid arrows represent the transition from one class to another.


### 2.2. The Mathematical Model

Using the above assumptions and by making a balance of the movements in each class, we obtain the following system:(5)dEdtt=b1-EtKEAt-s+dEt,dLdtt=s1-LtKLEt-sL+dLLt,dShdtt=Λ+γRht-dh+khtSht,dEhdtt=khtSht-dh+αEht,dIhdtt=αEht-dh+dp+rhIht,dRhdtt=rhIht-dh+γRht,dSmdtt=sLLt-dm+kmtSmt,dEmdtt=kmtSmt-νm+dmEmt,dImdtt=νmEmt-dmImt.The growth of the whole human population and mature vector is, respectively, described by the following equations:(6)dNhdtt=Λ-dhNht-dpIht,dAdtt=sLLt-dmAt.Using (2), we get Sm(t)=A(t)-Em(t)-Im(t) and then the model can be rewritten as follows:(7)dEdtt=b1-EtKEAt-s+dEt,dLdtt=s1-LtKLEt-sL+dLLt,dAdtt=sLLt-dmAt,dShdtt=Λ+γRht-dh+khtSht,dEhdtt=khtSht-dh+αEht,dIhdtt=αEht-dh+dp+rhIht,dRhdtt=rhIht-dh+γRht,dEmdtt=kmtAt-kmtImt-νm+dm+kmtEmt,dImdtt=νmEmt-dmImt.Mathematically model (7) can be written as follows:(8)X˙t=Ft,Xt,where X(t)=(E(t),L(t),A(t),Sh(t),Eh(t),Ih(t),Rh(t),Em(t),Im(t))T. The function F:R+×R9→R9 is C∞(R9) and defined by(9)Ft,Xt=b1-EtKEAt-s+dEts1-LtKLEt-sL+dLLtsLLt-dmAtΛ+γRht-dh+khtShtkhtSht-dh+αEhtαEht-dh+dp+rhIhtrhIht-dh+γRhtkmtAt-kmtImt-νm+dm+kmtEmtνmEmt-dmImt.Let us consider F=(F1,F2)T and X(t)=(X1(t),X2(t))T with X1(t)=(E(t),L(t),A(t))T and X2(t)=(Sh(t),Eh(t),Ih(t),Rh(t),Em(t),Im(t))T. Then system (8) can be rewritten as follows:(10)X˙1t=F1X1t,X2t,(11)X˙2t=F2t,X1t,X2t,with the functions F1 and F2 defined as follows: (12)F2t,X1t,X2t=Λ+γRht-dh+khtShtkhtSht-dh+αEhtαEht-dh+dp+rhIhtrhIht-dh+γRhtkmtAt-kmtImt-νm+dm+kmtEmtνmEmt-dmImt,F1X1t,X2t=b1-EtKEAt-s+dEts1-LtKLEt-sL+dLLtsLLt-dmAt.System (10) describes the maturation cycle of mosquitoes and system (11) describes the dynamics of malaria transmission. System (10) is biologically well defined in (13)Δ≔E,L,A∈R+3E≤KE,L≤KL,A≤sLdmKLand system (11) is biologically well defined in (14)Ω=Sh,Eh,Ih,Rh,Em,Im∈R+6∣Sh+Eh+Ih+Rh≤Λdh,Em+Im≤sLKLdm;then model (7) is biologically well defined in Γ≔Δ×Ω.


## 2.1. Interactions between Humans and Mosquitoes

When an infectious mosquito bites a susceptible human, the parasite enters the body of the human with a probabilitycmh and the human moves into the exposed class Eh. Some time after, he leaves from class Eh to class Ih with rate α. Infectious humans migrate into the class Rh after acquisition of their immunity with rate rh. The immunized lose their immunity with rate γ if they do not have continuous exposure to infection. Humans leave the total population through natural death rate dh and malaria death rate dP.Similarly, when a susceptible mosquito bites an infectious human, it enters the classEm with a probability chm. Some time after, it leaves from class Em to infective class Im with rate νm where it remains for life. Mature mosquitoes leave the population through natural mortality dm.Using the standard incidence as in the model of Ngwa and Shu [4], we define, respectively, the infection incidence from mosquitoes to humans, kh(t), and from humans to mosquitoes, km(t):(3)kht=cmhβtImtNht,(4)kmt=chmβtIhtNht+c-hmβtRhtNht.Furthermore, using the above assumptions, we obtain the transfer diagram (Figure 1) of the model.Figure 1
The dashed arrows indicate the direction of the infection and the solid arrows represent the transition from one class to another.


## 2.2. The Mathematical Model

Using the above assumptions and by making a balance of the movements in each class, we obtain the following system:(5)dEdtt=b1-EtKEAt-s+dEt,dLdtt=s1-LtKLEt-sL+dLLt,dShdtt=Λ+γRht-dh+khtSht,dEhdtt=khtSht-dh+αEht,dIhdtt=αEht-dh+dp+rhIht,dRhdtt=rhIht-dh+γRht,dSmdtt=sLLt-dm+kmtSmt,dEmdtt=kmtSmt-νm+dmEmt,dImdtt=νmEmt-dmImt.The growth of the whole human population and mature vector is, respectively, described by the following equations:(6)dNhdtt=Λ-dhNht-dpIht,dAdtt=sLLt-dmAt.Using (2), we get Sm(t)=A(t)-Em(t)-Im(t) and then the model can be rewritten as follows:(7)dEdtt=b1-EtKEAt-s+dEt,dLdtt=s1-LtKLEt-sL+dLLt,dAdtt=sLLt-dmAt,dShdtt=Λ+γRht-dh+khtSht,dEhdtt=khtSht-dh+αEht,dIhdtt=αEht-dh+dp+rhIht,dRhdtt=rhIht-dh+γRht,dEmdtt=kmtAt-kmtImt-νm+dm+kmtEmt,dImdtt=νmEmt-dmImt.Mathematically model (7) can be written as follows:(8)X˙t=Ft,Xt,where X(t)=(E(t),L(t),A(t),Sh(t),Eh(t),Ih(t),Rh(t),Em(t),Im(t))T. The function F:R+×R9→R9 is C∞(R9) and defined by(9)Ft,Xt=b1-EtKEAt-s+dEts1-LtKLEt-sL+dLLtsLLt-dmAtΛ+γRht-dh+khtShtkhtSht-dh+αEhtαEht-dh+dp+rhIhtrhIht-dh+γRhtkmtAt-kmtImt-νm+dm+kmtEmtνmEmt-dmImt.Let us consider F=(F1,F2)T and X(t)=(X1(t),X2(t))T with X1(t)=(E(t),L(t),A(t))T and X2(t)=(Sh(t),Eh(t),Ih(t),Rh(t),Em(t),Im(t))T. Then system (8) can be rewritten as follows:(10)X˙1t=F1X1t,X2t,(11)X˙2t=F2t,X1t,X2t,with the functions F1 and F2 defined as follows: (12)F2t,X1t,X2t=Λ+γRht-dh+khtShtkhtSht-dh+αEhtαEht-dh+dp+rhIhtrhIht-dh+γRhtkmtAt-kmtImt-νm+dm+kmtEmtνmEmt-dmImt,F1X1t,X2t=b1-EtKEAt-s+dEts1-LtKLEt-sL+dLLtsLLt-dmAt.System (10) describes the maturation cycle of mosquitoes and system (11) describes the dynamics of malaria transmission. System (10) is biologically well defined in (13)Δ≔E,L,A∈R+3E≤KE,L≤KL,A≤sLdmKLand system (11) is biologically well defined in (14)Ω=Sh,Eh,Ih,Rh,Em,Im∈R+6∣Sh+Eh+Ih+Rh≤Λdh,Em+Im≤sLKLdm;then model (7) is biologically well defined in Γ≔Δ×Ω.


## 3. Mathematical Analysis

### 3.1. Positivity and Boundedness of Solutions

Lemma 1 (see [6]).
The setΔ is a positive invariant region under the flow induced by (10).We assume that(H5)
β(t) is a ω-periodic positive function with ω=12 months,(H6)
all the parameters of the model are positive except the disease-induced death rate,dp, which is assumed to be nonnegative.Theorem 2.
For any initial conditionϕ∈R+9, system (8) has a unique solution. Further, the compact Γ is a positively invariant set, which attracts all positive orbits in R+9.Proof.
For allϕ∈R+9, the function F is locally Lipschitzian in X(t). It then follows through Cauchy-Lipschitz theorem that system (8) has a unique local solution.
Furthermore, according to (6), we have (15)dNhdtt=Λ-dhNht-dpIht≤Λ-dhNht,dAdtt=sLLt-dmAt≤sLKL-dmAt.It then follows that if Nh(t)>Λ/dh and A(t)>sLKL/dm, then dNh/dt(t)<0 and dA/dt(t)<0.
Let us consider the following differential equations:(16)dNhdtt=Λ-dhNht,dAdtt=sLKL-dmAt with general solutions:(17)Nht=Λdh+Nh0-Λdhe-dht,At=sLKLdm+A0-sLKLdme-dmt.
By applying the standard comparison theorem, we obtain, for allt≥0, Nh(t)≤Λ/dh and A(t)≤sLKL/dm if Nh(0)≤Λ/dh and A(0)≤sLKL/dm. Thus, the set Ω is positively invariant with respect to system (11). Therefore, from Lemma 1, the set Δ is positively invariant with respect to system (10). Then, we conclude that the compact set Γ=Δ×Ω is positively invariant. Thus, all the solutions of system (8) are nonnegative and bounded.


### 3.2. Disease-Free Equilibriums

Let us consider the following threshold parameter:r=(b/s+d)(s/sL+dL)(sL/dm). Then we have the following result.Proposition 3 (see [6]).
System (10) always has the mosquito-free equilibrium P0=(0,0,0). (i)
Ifr≤1, then system (10) has no other equilibrium.(ii)
Ifr>1, there is a unique endemic equilibrium (18)P1=E∗,L∗,A∗=1-1rKEγE,KLγL,sLdmKLγL,where (19)γE=1+s+ddmKEbsLKL,γL=1+sL+dLKLsKE.Lemma 4.
Model (7) has (i)
trivial disease-free equilibriumE0=(0,0,0,Sh∗,0,0,0,0,0) if r≤1,(ii)
nontrivial disease-free equilibriumE1=(E∗,L∗,A∗,Sh∗,0,0,0,0,0) if r>1, where Sh∗=Λ/dh, A∗=Sm∗=sLL∗/dm, and E∗, L∗, and A∗ are given above.Proof.
By solving the systemF2(t,X1(t),X2(t))=0 at the disease-free equilibrium, Eh(t)=Ih(t)=Rh(t)=Em(t)=Im(t)=0,∀t≥0, we get the equilibrium point E1+=(Sh∗,0,0,0,0,0) for system (11), with Sh∗=Λ/dh. Moreover, thanks to Proposition 3, system (10) has a unique mosquito-free equilibrium (0,0,0) if r≤1 and a unique endemic equilibrium (E∗,L∗,A∗) if r>1. Thus, we conclude that system (7) has a trivial disease-free equilibrium E0=(0,0,0,Sh∗,0,0,0,0,0) if r≤1 and a nontrivial disease-free equilibrium E1=(E∗,L∗,A∗,Sh∗,0,0,0,0,0) if r>1.Remark 5.
We will only consider the equilibrium stateE1 because it is more biologically realistic. So, in the rest of the paper, we assume that r>1.


### 3.3. Threshold Dynamics

Linearizing system (8) at the equilibrium state E1, we obtain the following system (here we write down only the equations for the “diseased” classes): (20)dEhdtt=cmhβtImt-dh+αEht,dIhdtt=αEht-dh+dp+rhIht,dRhdtt=rhIht-dh+γRht,dEmdtt=chmβtA∗Sh∗Iht+c¯hmβtA∗Sh∗Rht-νm+dmEmt,dImdtt=νmEmt-dmImt.This system can be rewritten as(21)dZtdt=Ft-VtZt,where Z(t)=(Eh(t),Ih(t),Rh(t),Em(t),Im(t))T and F(t) and V(t) are 5×5 matrix defined as follows: (22)Ft=0000cmhβt00000000000chmβtA∗Sh∗c¯hmβtA∗Sh∗0000000,Vt=dh+α0000-αdh+dp+rh0000-rhdh+γ00000νm+dm0000-νmdm.Let us assume thatY(t,s), t≥s, is the matrix solution of the linear ω-periodic system(23)dydt=-Vty.That is, for each s∈R, the 5×5 matrix Y(t,s) satisfies the equation (24)ddtYt,s=-VtYt,s,∀t≥s,Ys,s=I,where I is the 5×5 identity matrix. Thus, the monodromy matrix Φ-V(t) of (23) is equal to Y(t,0), ∀t≥0.LetCω be the ordered Banach space of all ω-periodic functions from R to R5 which is equipped with the maximum norm · and the positive cone Cω+≔{ϕ∈Cω:ϕ(t)≥0,∀t∈R}. Then, we can define a linear operator L:Cω→Cω by(25)Lϕt=∫0∞Yt,t-aFt-aϕt-ada,∀t∈R,ϕ∈Cω.It then follows from [12] that L is the next infection operator, and the basic reproduction ratio is R0=ρ(L), the spectral radius of L.In order to calculateR0, we consider the following linear ω-periodic system:(26)dwtdt=1λFt-Vtwt,∀t∈R+,λ∈0,∞.Let W(t,s,λ), t≥s, s∈R, be the evolution operator of system (26) on R5. Clearly W(t,0,1)=ΦF-V(t),∀t≥0. The following result will be used in our numerical calculation of the basic reproduction ratio.Lemma 6 (see [12]).
(i) Ifρ(W(ω,0,λ))=1 has a positive solution λ0, then λ0 is an eigenvalue of L, and hence R0>0.
(ii) IfR0>0, then λ=R0 is the unique solution of ρ(W(ω,0,λ))=1.
(iii)R0=0 if and only if ρ(W(ω,0,λ))<1, for all λ>0.


### 3.4. Stability of Equilibrium StateE1

In this section, we will study the asymptotic behaviour of the nontrivial equilibriumE1; thus we have the following result, which will be used in the proofs of our main results.Lemma 7 (see [12]).
The following statements are valid:(i)
R0=1 if and only if ρ(ΦF-V(ω))=1.(ii)
R0<1 if and only if ρ(ΦF-V(ω))<1.(iii)
R0>1 if and only if ρ(ΦF-V(ω))>1.Lemma 8 (see [6]).
Ifr>1, then P1 is globally asymptotically stable in int(Δ), with respect to system (10).Theorem 9.
The nontrivial equilibriumE1 is locally asymptotically stable if R0<1 and unstable if R0>1.Proof.
LetA(t) be the Jacobian matrix of (8) evaluated at E1. Then we have (27)At=A11A12A21A22t,where (28)A12=0000b1-E∗KEb1-E∗KE0000000000-dm-dm,A21=000000000000000000,A11=-s+d-bA∗KE0s+dE∗A∗sL+dLL∗E∗-sL+dL-sE∗KL00sL-dm,A22t=-dhCt0^Ft-Vt with (29)Ct=00γ0-cmhβt,0^=0,0,0,0,0T.E1 is locally asymptotically stable if ρ(ΦA(ω))<1. The matrix A11 is a constant matrix and its characteristic equation is given by π(z)=z3+a1z2+a2z+a3, where (30)a1=1-1rsKEγEKL+bsLKLdmγLKE+sL+dL+s+d+dm,a2=bsLKLdmγLKE1-1rsKEγEKL1-1r+sL+dL+dm+dmsKEγEKL1-1r+dmsL+dL,a3=dm1-1rbssLdmγLγE1-1r+s+dsKEγEKL+sL+dLbsLKLdmγLKE. If r>1, then a1, a2, a3 and a1a2-a3 are clearly positive. So, thanks to Routh-Hurwitz criterion, all eigenvalues of A11 have negative real part. It then follows that ρ(ΦA11(ω))<1. Thus, the stability of E1 depends on ΦA22(ω).
Thus, ifρ(ΦF-V(ω))<1, then ρ(ΦA22(ω))<1 and then E1 is stable. If ρ(ΦF-V(ω))>1 then E1 is unstable. So, thanks to Lemma 7,  E1 is locally asymptotically stable if R0<1 and unstable if R0>1.Lemma 10 (see [13]).
Letθ=1/ωln⁡ρ(ΦA(·)(ω)); then there exists a positive ω-periodic function v(t) such that eθtv(t) is a solution of x˙(t)=A(t)x(t).Theorem 11.
IfR0<1 and dp=0, then E1 is globally asymptotically stable.Proof.
Ifdp=0, we can rewrite (6) as follows: (31)dNhdtt=Λ-dhNht,dAdtt=sLLt-dmAt.Thus, there exists a period ω′ such that ∀t≥ω′, Nh(t)≥Nh∗-ϵ and A(t)≤A∗+ϵ, ∀ϵ>0.
At disease-free equilibrium, we haveNh∗=Sh∗ and Sm∗=A∗. So, A(t)/Nh(t)≤A∗+ϵ/Sh∗-ϵ. It then follows from system (11) that(32a)dEhdtt≤cmhβtImt-dh+αEht,(32b)dIhdtt=αEht-dh+rhIht,(32c)dRhdtt=rhIht-dh+γRht,(32d)dEmdtt≤chmβtA∗+ϵSh∗-ϵIht+c-hmβtA∗+ϵSh∗-ϵRht-νm+dmEmt,(32e)dImdtt=νmEmt-dmImt.Let us consider the following auxiliary system: (33)dE¯hdtt=cmhβtI¯mt-dh+αE¯ht,dI¯hdtt=αE¯ht-dh+rhI¯ht,dR¯hdtt=rhI¯ht-dh+γR¯ht,dE¯mdtt=chmβtA∗+ϵSh∗-ϵI¯ht+c¯hmβtA∗+ϵSh∗-ϵR¯ht-νm+dmE¯mt,dI¯mdtt=νmE¯mt-dmI¯mt,which can be rewritten as follows:(34)dh¯dtt=Mϵth¯t;h¯t=E¯ht,I¯ht,R¯ht,E¯mt,I¯mtTwith(35)Mϵt=-dh+α000cmhβtα-dh+rh0000rh-dh+γ000chmβtA∗+ϵSh∗-ϵc¯hmβtA∗+ϵSh∗-ϵ-νm+dm0000νm-dm.From Lemma 7, if R0<1, then ρ(ΦF-V(ω))<1. Clearly, limϵ→0+ΦMϵ(ω)=ΦF-V(ω) and, by continuity of the spectral radius, we have limϵ→0+ρ(ΦMϵ(ω))=ρ(ΦF-V(ω))<1. Thus, there exists ϵ1>0 such that ρ(ΦMϵ(ω))<1, ∀ϵ∈[0,ϵ1[.
From Lemma10, there exists a positive ω-periodic function v(t) such that h¯(t)=eθtv(t) is a solution of (34). Since ρ(ΦMϵ(ω))<1, θ<0. The ω-periodic function v(t) is bounded and it then follows that limt→∞h¯(t)=0. Applying comparison theorem on system (32a)–(32e), we get limt→∞(Eh(t),Ih(t),Rh(t),Em(t),Im(t))=(0,0,0,0,0). Using the theory of asymptotically periodic semiflow [[14], Theorem 3.2.1], we have limt→∞Sh(t)=Sh∗, limt→∞A(t)=A∗=Sm∗. From Lemma 8, if r>1 then P1 is globally asymptotically stable, so limt→∞E(t)=E∗ and limt→∞L(t)=L∗. Hence, the equilibrium E1 is globally attractive.


### 3.5. Existence of Positive Periodic Solutions

System (8) is constructed by coupling two subsystems. The term coupling these two systems is given by the function sLL(t). The coupling takes place only in one direction because the dynamics of system (11) depend on the dynamics of system (10). The asymptotic behaviour of system (10) is given by Lemma 8. Now we are going to study the existence of positive periodic solutions of system (11): (36)dShdtt=Λ+γRht-dh+khtSht,dEhdtt=khtSht-dh+αEht,dIhdtt=αEht-dh+dp+rhIht,dRhdtt=rhIht-dh+γRht,dEmdtt=kmtAt-kmtImt-νm+dm+kmtEmt,dImdtt=νmEmt-dmImt.Model (11) is well defined in Ω and if r>1 it has a disease-free equilibrium E1+=Sh∗,0,0,0,0,0 with Sh∗=Λ/dh.Let us consider the following sets:(37)X≔R+6,X0≔Sh,Eh,Ih,Rh,Em,Im∈X∣Eh>0,Ih>0,Rh>0,Em>0,Im>0,∂X0≔X∖X0.Let u(t,ψ) be the unique solution of (11) with initial conditions ψ, Φ(t) the periodic semiflow generated by periodic system (11), and P:X→X the Poincaré map associated with system (11); namely, (38)Pψ=Φωψ=uω,ψ,∀ψ∈X,Pmψ=Φmωψ=umω,ψ,∀m≥0.Proposition 12.
The setsX0 and ∂X0 are positively invariant under the flow induced by (11).Proof.
Note that ifX0 is positively invariant, then ∂X0 is positively invariant. Thus we only need to prove that X0 is positively invariant.
For any initial conditionψ∈X0, solving the equations of system (11) we derive that (39)Sht=exp⁡-∫0tkhs+dhdsSh0+∫0tΛ+Ihs+γRhsexp⁡∫0skhc+dhdcds≥exp⁡-∫0tkhs+dhds∫0tΛ+Ihs+γRhsexp⁡∫0skhc+dhdcds>0,∀t>0,Eht=e-dh+αtEh0+∫0tkhsShsedh+αsds≥e-dh+αt∫0tkhsShsedh+αsds>0,∀t>0,Iht=e-dh+dp+rhtIh0+∫0tαEhsedh+dp+rhsds≥e-dh+dp+rht∫0tαEhsedh+dp+rhsds>0,∀t>0,Rht=e-dh+γtRh0+∫0trhIhsedh+γsds≥e-dh+γt∫0trhIhsedh+γsds>0,∀t>0,Emt=e∫0t-kms+dm+νmdsEm0+∫0tkmsAs-Imse∫0skmc+dm+νmdcds≥e∫0t-kms+dm+νmds∫0tkmsAs-Imse∫0skmc+dm+νmdcds>0,∀t>0,Imt=e-dmtIm0+∫0tνmEmsedms≥e-dmt∫0tνmEmsedms>0,∀t>0. Thus, X0 is positively invariant. So, ∂X0 is also positively invariant.Note that, from Theorem2, Ω is a compact set which attracts all positive orbits in X, which implies that the discrete-time system P:X→X is point dissipative. Moreover, ∀n0≥1, Pn0 is compact; it then follows from Theorem 2.9 in [15] that P admits a global attractor in X.Lemma 13.
IfR0>1, there exists η>0 such that when ψ-E1+≤η, ∀ψ∈X0, we have limsupm→∞Pm(ψ)-E1+≥η.Proof.
Suppose by contradiction thatlimsupm→∞Pm(ψ)-E1+<η for some ψ∈X0. Then, there exists an integer n≥1 such that, for all m≥n, Pm(ψ)-M<η. By the continuity of the solution u(t,ψ), we have u(t,Pm(ψ))-ut,E1+≤σ for all t≥0 and σ>0. For all t≥0, let t=mω+t1, where t1∈[0,ω] and m=[t/ω]. [t/ω] is the greatest integer less than or equal to t/ω. If ψ-E1+≤η, then by the continuity of the solution u(t,ψ) we have (40)ut,ψ-ut,E1+=ut1+mω,ψ-ut1+mω,M=Φt1+mωψ-Φt1+mωE1+=Φt1Φmωψ-Φt1ΦmωE1+=Φt1Pmψ-Φt1PmE1+=Φt1Pmψ-Φt1E1+≤σ.It then follows that Sh∗-σ≤Sh(t)≤Sh∗+σ and A∗-σ≤A(t)≤A∗+σ. So, there exists σ∗>0 such that Sh(t)/Nh(t)≥1-σ∗ and A(t)/Nh(t)≥A∗/Nh∗-σ∗.
From (11) we have(41)dEhdtt≥cmhβt1-σ∗Imt-dh+αEht,dIhdtt=αEht-dp+dh+rhIht,dRhdtt=rhIht-dh+γRht,dEmdtt≥βtA∗Nh∗-σ∗chmIht+c-hmRht-νm+dmEmt,dImdtt=νmEmt-dmImt.Let us consider the following auxiliary linear system:(42)dh^dtt=Mσ∗th^t;h^t=E^ht,I^ht,R^ht,E^mt,I^mtTwith(43)Mσ∗t=-dh+α0001-σ∗cmhβtα-dh+rh0000rh-dh+γ000chmβtA∗Sh∗-σ∗c¯hmβtA∗Sh∗-σ∗-νm+dm0000νm-dm.By applying the same method as above, if R0>1 then ρΦMσ∗ω>1. In this case θ is positive, and then h^(t)→∞ as t→∞. Moreover, since X0 is positively invariant, then there exists an integer q≥n and a real number κ>0 such that (44)Ehqω,Ihqω,Rhqω,Emqω,Imqω≥κh^0. Applying the theorem of comparison principle, we get (45)Ehqω+t,Ihqω+t,Rhqω+t,Emqω+t,Imqω+t≥κh^t,∀t≥0. It then follows that limt→∞Eht,Iht,Rht,Emt,Imt=∞, which contradicts the fact that solutions are bounded.Theorem 14.
IfR0>1, then system (7) has at least one positive periodic solution.Proof.
We first prove thatP is uniformly persistent with respect to (X0,∂X0).
We define the following sets:(46)M∂=ψ∈∂X0∣Pmψ∈∂X0,foranym≥0,D=Sh,0,0,0,0,0∈X∣Sh≥0.Let us prove that M∂=D.
It is easy to remark thatD⊂M∂. We only need to prove that M∂⊂D.
Letψ∈∂X0∖D. If (i)
Ih(0)>0, Im(0)>0, and Eh(0)=Em(0)=Rh(0)=0, then we have Sh(t)>0, Ih(t)>0, Im(t)>0, Em(t)>0, Eh(t)>0, Rh(t)>0, ∀t>0,(ii)
Ih(0)=Im(0)=0 and Eh(0)>0, Em(0)>0, Rh(0)>0, then we have Sh(t)>0, Ih(t)>0, Im(t)>0, Em(t)>0, Eh(t)>0, Rh(t)>0, ∀t>0.
For any cases, it follows that(Sh(t),Eh(t),Ih(t),Rh(t),Em(t),Im(t))∉∂X0 for t>0 sufficiently small, which contradicts the fact that ∂X0 is positively invariant. Hence, M∂⊂D. Thus, it then follows that M∂=D.
The equalityM∂=D implies that E1+ is a fixed point of P and acyclic in M∂; every solution in M∂ approaches to E1+. Moreover, Lemma 13 implies that E1+ is an isolated invariant set in X and Ws(E1+)∩X0=∅. By the acyclicity theorem on uniform persistence for maps, Theorem 1.3.1 and Remark 1.3.1 in [14], it follows that P is uniformly persistent with respect to X0. Thus, Theorem 3.1.1 in [14] implies that the periodic semiflow Φ(t):X→X is also uniformly persistent with respect to X0. Thanks to Theorem 1.3.6 in [14], model (11) has at least one ω-periodic solution u~(t,ψ∗) with ψ∗∈X0 and t≥0. Now, we show that u~(t,ψ∗) is positive.
Suppose thatψ∗=0; then, for all t>0, we obtain u~i(t,ψ∗)>0, for i=1,2,3,4,5,6. By using the periodicity of the solution, we have Sh∗(0)=Sh∗(nω)=0, Eh∗(0)=Eh∗(nω)=0, Ih∗(0)=Ih∗(nω)=0, Rh∗(0)=Rh∗(nω)=0, Em∗(0)=Em∗(nω)=0, Im∗(0)=Im∗(nω)=0, ∀n≥1, which contradicts the fact that u~i(t,ψ∗)>0 for i=1,2,3,4,5,6. So, the periodic solution is positive.


## 3.1. Positivity and Boundedness of Solutions

Lemma 1 (see [6]).
The setΔ is a positive invariant region under the flow induced by (10).We assume that(H5)
β(t) is a ω-periodic positive function with ω=12 months,(H6)
all the parameters of the model are positive except the disease-induced death rate,dp, which is assumed to be nonnegative.Theorem 2.
For any initial conditionϕ∈R+9, system (8) has a unique solution. Further, the compact Γ is a positively invariant set, which attracts all positive orbits in R+9.Proof.
For allϕ∈R+9, the function F is locally Lipschitzian in X(t). It then follows through Cauchy-Lipschitz theorem that system (8) has a unique local solution.
Furthermore, according to (6), we have (15)dNhdtt=Λ-dhNht-dpIht≤Λ-dhNht,dAdtt=sLLt-dmAt≤sLKL-dmAt.It then follows that if Nh(t)>Λ/dh and A(t)>sLKL/dm, then dNh/dt(t)<0 and dA/dt(t)<0.
Let us consider the following differential equations:(16)dNhdtt=Λ-dhNht,dAdtt=sLKL-dmAt with general solutions:(17)Nht=Λdh+Nh0-Λdhe-dht,At=sLKLdm+A0-sLKLdme-dmt.
By applying the standard comparison theorem, we obtain, for allt≥0, Nh(t)≤Λ/dh and A(t)≤sLKL/dm if Nh(0)≤Λ/dh and A(0)≤sLKL/dm. Thus, the set Ω is positively invariant with respect to system (11). Therefore, from Lemma 1, the set Δ is positively invariant with respect to system (10). Then, we conclude that the compact set Γ=Δ×Ω is positively invariant. Thus, all the solutions of system (8) are nonnegative and bounded.


## 3.2. Disease-Free Equilibriums

Let us consider the following threshold parameter:r=(b/s+d)(s/sL+dL)(sL/dm). Then we have the following result.Proposition 3 (see [6]).
System (10) always has the mosquito-free equilibrium P0=(0,0,0). (i)
Ifr≤1, then system (10) has no other equilibrium.(ii)
Ifr>1, there is a unique endemic equilibrium (18)P1=E∗,L∗,A∗=1-1rKEγE,KLγL,sLdmKLγL,where (19)γE=1+s+ddmKEbsLKL,γL=1+sL+dLKLsKE.Lemma 4.
Model (7) has (i)
trivial disease-free equilibriumE0=(0,0,0,Sh∗,0,0,0,0,0) if r≤1,(ii)
nontrivial disease-free equilibriumE1=(E∗,L∗,A∗,Sh∗,0,0,0,0,0) if r>1, where Sh∗=Λ/dh, A∗=Sm∗=sLL∗/dm, and E∗, L∗, and A∗ are given above.Proof.
By solving the systemF2(t,X1(t),X2(t))=0 at the disease-free equilibrium, Eh(t)=Ih(t)=Rh(t)=Em(t)=Im(t)=0,∀t≥0, we get the equilibrium point E1+=(Sh∗,0,0,0,0,0) for system (11), with Sh∗=Λ/dh. Moreover, thanks to Proposition 3, system (10) has a unique mosquito-free equilibrium (0,0,0) if r≤1 and a unique endemic equilibrium (E∗,L∗,A∗) if r>1. Thus, we conclude that system (7) has a trivial disease-free equilibrium E0=(0,0,0,Sh∗,0,0,0,0,0) if r≤1 and a nontrivial disease-free equilibrium E1=(E∗,L∗,A∗,Sh∗,0,0,0,0,0) if r>1.Remark 5.
We will only consider the equilibrium stateE1 because it is more biologically realistic. So, in the rest of the paper, we assume that r>1.


## 3.3. Threshold Dynamics

Linearizing system (8) at the equilibrium state E1, we obtain the following system (here we write down only the equations for the “diseased” classes): (20)dEhdtt=cmhβtImt-dh+αEht,dIhdtt=αEht-dh+dp+rhIht,dRhdtt=rhIht-dh+γRht,dEmdtt=chmβtA∗Sh∗Iht+c¯hmβtA∗Sh∗Rht-νm+dmEmt,dImdtt=νmEmt-dmImt.This system can be rewritten as(21)dZtdt=Ft-VtZt,where Z(t)=(Eh(t),Ih(t),Rh(t),Em(t),Im(t))T and F(t) and V(t) are 5×5 matrix defined as follows: (22)Ft=0000cmhβt00000000000chmβtA∗Sh∗c¯hmβtA∗Sh∗0000000,Vt=dh+α0000-αdh+dp+rh0000-rhdh+γ00000νm+dm0000-νmdm.Let us assume thatY(t,s), t≥s, is the matrix solution of the linear ω-periodic system(23)dydt=-Vty.That is, for each s∈R, the 5×5 matrix Y(t,s) satisfies the equation (24)ddtYt,s=-VtYt,s,∀t≥s,Ys,s=I,where I is the 5×5 identity matrix. Thus, the monodromy matrix Φ-V(t) of (23) is equal to Y(t,0), ∀t≥0.LetCω be the ordered Banach space of all ω-periodic functions from R to R5 which is equipped with the maximum norm · and the positive cone Cω+≔{ϕ∈Cω:ϕ(t)≥0,∀t∈R}. Then, we can define a linear operator L:Cω→Cω by(25)Lϕt=∫0∞Yt,t-aFt-aϕt-ada,∀t∈R,ϕ∈Cω.It then follows from [12] that L is the next infection operator, and the basic reproduction ratio is R0=ρ(L), the spectral radius of L.In order to calculateR0, we consider the following linear ω-periodic system:(26)dwtdt=1λFt-Vtwt,∀t∈R+,λ∈0,∞.Let W(t,s,λ), t≥s, s∈R, be the evolution operator of system (26) on R5. Clearly W(t,0,1)=ΦF-V(t),∀t≥0. The following result will be used in our numerical calculation of the basic reproduction ratio.Lemma 6 (see [12]).
(i) Ifρ(W(ω,0,λ))=1 has a positive solution λ0, then λ0 is an eigenvalue of L, and hence R0>0.
(ii) IfR0>0, then λ=R0 is the unique solution of ρ(W(ω,0,λ))=1.
(iii)R0=0 if and only if ρ(W(ω,0,λ))<1, for all λ>0.


## 3.4. Stability of Equilibrium StateE1

In this section, we will study the asymptotic behaviour of the nontrivial equilibriumE1; thus we have the following result, which will be used in the proofs of our main results.Lemma 7 (see [12]).
The following statements are valid:(i)
R0=1 if and only if ρ(ΦF-V(ω))=1.(ii)
R0<1 if and only if ρ(ΦF-V(ω))<1.(iii)
R0>1 if and only if ρ(ΦF-V(ω))>1.Lemma 8 (see [6]).
Ifr>1, then P1 is globally asymptotically stable in int(Δ), with respect to system (10).Theorem 9.
The nontrivial equilibriumE1 is locally asymptotically stable if R0<1 and unstable if R0>1.Proof.
LetA(t) be the Jacobian matrix of (8) evaluated at E1. Then we have (27)At=A11A12A21A22t,where (28)A12=0000b1-E∗KEb1-E∗KE0000000000-dm-dm,A21=000000000000000000,A11=-s+d-bA∗KE0s+dE∗A∗sL+dLL∗E∗-sL+dL-sE∗KL00sL-dm,A22t=-dhCt0^Ft-Vt with (29)Ct=00γ0-cmhβt,0^=0,0,0,0,0T.E1 is locally asymptotically stable if ρ(ΦA(ω))<1. The matrix A11 is a constant matrix and its characteristic equation is given by π(z)=z3+a1z2+a2z+a3, where (30)a1=1-1rsKEγEKL+bsLKLdmγLKE+sL+dL+s+d+dm,a2=bsLKLdmγLKE1-1rsKEγEKL1-1r+sL+dL+dm+dmsKEγEKL1-1r+dmsL+dL,a3=dm1-1rbssLdmγLγE1-1r+s+dsKEγEKL+sL+dLbsLKLdmγLKE. If r>1, then a1, a2, a3 and a1a2-a3 are clearly positive. So, thanks to Routh-Hurwitz criterion, all eigenvalues of A11 have negative real part. It then follows that ρ(ΦA11(ω))<1. Thus, the stability of E1 depends on ΦA22(ω).
Thus, ifρ(ΦF-V(ω))<1, then ρ(ΦA22(ω))<1 and then E1 is stable. If ρ(ΦF-V(ω))>1 then E1 is unstable. So, thanks to Lemma 7,  E1 is locally asymptotically stable if R0<1 and unstable if R0>1.Lemma 10 (see [13]).
Letθ=1/ωln⁡ρ(ΦA(·)(ω)); then there exists a positive ω-periodic function v(t) such that eθtv(t) is a solution of x˙(t)=A(t)x(t).Theorem 11.
IfR0<1 and dp=0, then E1 is globally asymptotically stable.Proof.
Ifdp=0, we can rewrite (6) as follows: (31)dNhdtt=Λ-dhNht,dAdtt=sLLt-dmAt.Thus, there exists a period ω′ such that ∀t≥ω′, Nh(t)≥Nh∗-ϵ and A(t)≤A∗+ϵ, ∀ϵ>0.
At disease-free equilibrium, we haveNh∗=Sh∗ and Sm∗=A∗. So, A(t)/Nh(t)≤A∗+ϵ/Sh∗-ϵ. It then follows from system (11) that(32a)dEhdtt≤cmhβtImt-dh+αEht,(32b)dIhdtt=αEht-dh+rhIht,(32c)dRhdtt=rhIht-dh+γRht,(32d)dEmdtt≤chmβtA∗+ϵSh∗-ϵIht+c-hmβtA∗+ϵSh∗-ϵRht-νm+dmEmt,(32e)dImdtt=νmEmt-dmImt.Let us consider the following auxiliary system: (33)dE¯hdtt=cmhβtI¯mt-dh+αE¯ht,dI¯hdtt=αE¯ht-dh+rhI¯ht,dR¯hdtt=rhI¯ht-dh+γR¯ht,dE¯mdtt=chmβtA∗+ϵSh∗-ϵI¯ht+c¯hmβtA∗+ϵSh∗-ϵR¯ht-νm+dmE¯mt,dI¯mdtt=νmE¯mt-dmI¯mt,which can be rewritten as follows:(34)dh¯dtt=Mϵth¯t;h¯t=E¯ht,I¯ht,R¯ht,E¯mt,I¯mtTwith(35)Mϵt=-dh+α000cmhβtα-dh+rh0000rh-dh+γ000chmβtA∗+ϵSh∗-ϵc¯hmβtA∗+ϵSh∗-ϵ-νm+dm0000νm-dm.From Lemma 7, if R0<1, then ρ(ΦF-V(ω))<1. Clearly, limϵ→0+ΦMϵ(ω)=ΦF-V(ω) and, by continuity of the spectral radius, we have limϵ→0+ρ(ΦMϵ(ω))=ρ(ΦF-V(ω))<1. Thus, there exists ϵ1>0 such that ρ(ΦMϵ(ω))<1, ∀ϵ∈[0,ϵ1[.
From Lemma10, there exists a positive ω-periodic function v(t) such that h¯(t)=eθtv(t) is a solution of (34). Since ρ(ΦMϵ(ω))<1, θ<0. The ω-periodic function v(t) is bounded and it then follows that limt→∞h¯(t)=0. Applying comparison theorem on system (32a)–(32e), we get limt→∞(Eh(t),Ih(t),Rh(t),Em(t),Im(t))=(0,0,0,0,0). Using the theory of asymptotically periodic semiflow [[14], Theorem 3.2.1], we have limt→∞Sh(t)=Sh∗, limt→∞A(t)=A∗=Sm∗. From Lemma 8, if r>1 then P1 is globally asymptotically stable, so limt→∞E(t)=E∗ and limt→∞L(t)=L∗. Hence, the equilibrium E1 is globally attractive.


## 3.5. Existence of Positive Periodic Solutions

System (8) is constructed by coupling two subsystems. The term coupling these two systems is given by the function sLL(t). The coupling takes place only in one direction because the dynamics of system (11) depend on the dynamics of system (10). The asymptotic behaviour of system (10) is given by Lemma 8. Now we are going to study the existence of positive periodic solutions of system (11): (36)dShdtt=Λ+γRht-dh+khtSht,dEhdtt=khtSht-dh+αEht,dIhdtt=αEht-dh+dp+rhIht,dRhdtt=rhIht-dh+γRht,dEmdtt=kmtAt-kmtImt-νm+dm+kmtEmt,dImdtt=νmEmt-dmImt.Model (11) is well defined in Ω and if r>1 it has a disease-free equilibrium E1+=Sh∗,0,0,0,0,0 with Sh∗=Λ/dh.Let us consider the following sets:(37)X≔R+6,X0≔Sh,Eh,Ih,Rh,Em,Im∈X∣Eh>0,Ih>0,Rh>0,Em>0,Im>0,∂X0≔X∖X0.Let u(t,ψ) be the unique solution of (11) with initial conditions ψ, Φ(t) the periodic semiflow generated by periodic system (11), and P:X→X the Poincaré map associated with system (11); namely, (38)Pψ=Φωψ=uω,ψ,∀ψ∈X,Pmψ=Φmωψ=umω,ψ,∀m≥0.Proposition 12.
The setsX0 and ∂X0 are positively invariant under the flow induced by (11).Proof.
Note that ifX0 is positively invariant, then ∂X0 is positively invariant. Thus we only need to prove that X0 is positively invariant.
For any initial conditionψ∈X0, solving the equations of system (11) we derive that (39)Sht=exp⁡-∫0tkhs+dhdsSh0+∫0tΛ+Ihs+γRhsexp⁡∫0skhc+dhdcds≥exp⁡-∫0tkhs+dhds∫0tΛ+Ihs+γRhsexp⁡∫0skhc+dhdcds>0,∀t>0,Eht=e-dh+αtEh0+∫0tkhsShsedh+αsds≥e-dh+αt∫0tkhsShsedh+αsds>0,∀t>0,Iht=e-dh+dp+rhtIh0+∫0tαEhsedh+dp+rhsds≥e-dh+dp+rht∫0tαEhsedh+dp+rhsds>0,∀t>0,Rht=e-dh+γtRh0+∫0trhIhsedh+γsds≥e-dh+γt∫0trhIhsedh+γsds>0,∀t>0,Emt=e∫0t-kms+dm+νmdsEm0+∫0tkmsAs-Imse∫0skmc+dm+νmdcds≥e∫0t-kms+dm+νmds∫0tkmsAs-Imse∫0skmc+dm+νmdcds>0,∀t>0,Imt=e-dmtIm0+∫0tνmEmsedms≥e-dmt∫0tνmEmsedms>0,∀t>0. Thus, X0 is positively invariant. So, ∂X0 is also positively invariant.Note that, from Theorem2, Ω is a compact set which attracts all positive orbits in X, which implies that the discrete-time system P:X→X is point dissipative. Moreover, ∀n0≥1, Pn0 is compact; it then follows from Theorem 2.9 in [15] that P admits a global attractor in X.Lemma 13.
IfR0>1, there exists η>0 such that when ψ-E1+≤η, ∀ψ∈X0, we have limsupm→∞Pm(ψ)-E1+≥η.Proof.
Suppose by contradiction thatlimsupm→∞Pm(ψ)-E1+<η for some ψ∈X0. Then, there exists an integer n≥1 such that, for all m≥n, Pm(ψ)-M<η. By the continuity of the solution u(t,ψ), we have u(t,Pm(ψ))-ut,E1+≤σ for all t≥0 and σ>0. For all t≥0, let t=mω+t1, where t1∈[0,ω] and m=[t/ω]. [t/ω] is the greatest integer less than or equal to t/ω. If ψ-E1+≤η, then by the continuity of the solution u(t,ψ) we have (40)ut,ψ-ut,E1+=ut1+mω,ψ-ut1+mω,M=Φt1+mωψ-Φt1+mωE1+=Φt1Φmωψ-Φt1ΦmωE1+=Φt1Pmψ-Φt1PmE1+=Φt1Pmψ-Φt1E1+≤σ.It then follows that Sh∗-σ≤Sh(t)≤Sh∗+σ and A∗-σ≤A(t)≤A∗+σ. So, there exists σ∗>0 such that Sh(t)/Nh(t)≥1-σ∗ and A(t)/Nh(t)≥A∗/Nh∗-σ∗.
From (11) we have(41)dEhdtt≥cmhβt1-σ∗Imt-dh+αEht,dIhdtt=αEht-dp+dh+rhIht,dRhdtt=rhIht-dh+γRht,dEmdtt≥βtA∗Nh∗-σ∗chmIht+c-hmRht-νm+dmEmt,dImdtt=νmEmt-dmImt.Let us consider the following auxiliary linear system:(42)dh^dtt=Mσ∗th^t;h^t=E^ht,I^ht,R^ht,E^mt,I^mtTwith(43)Mσ∗t=-dh+α0001-σ∗cmhβtα-dh+rh0000rh-dh+γ000chmβtA∗Sh∗-σ∗c¯hmβtA∗Sh∗-σ∗-νm+dm0000νm-dm.By applying the same method as above, if R0>1 then ρΦMσ∗ω>1. In this case θ is positive, and then h^(t)→∞ as t→∞. Moreover, since X0 is positively invariant, then there exists an integer q≥n and a real number κ>0 such that (44)Ehqω,Ihqω,Rhqω,Emqω,Imqω≥κh^0. Applying the theorem of comparison principle, we get (45)Ehqω+t,Ihqω+t,Rhqω+t,Emqω+t,Imqω+t≥κh^t,∀t≥0. It then follows that limt→∞Eht,Iht,Rht,Emt,Imt=∞, which contradicts the fact that solutions are bounded.Theorem 14.
IfR0>1, then system (7) has at least one positive periodic solution.Proof.
We first prove thatP is uniformly persistent with respect to (X0,∂X0).
We define the following sets:(46)M∂=ψ∈∂X0∣Pmψ∈∂X0,foranym≥0,D=Sh,0,0,0,0,0∈X∣Sh≥0.Let us prove that M∂=D.
It is easy to remark thatD⊂M∂. We only need to prove that M∂⊂D.
Letψ∈∂X0∖D. If (i)
Ih(0)>0, Im(0)>0, and Eh(0)=Em(0)=Rh(0)=0, then we have Sh(t)>0, Ih(t)>0, Im(t)>0, Em(t)>0, Eh(t)>0, Rh(t)>0, ∀t>0,(ii)
Ih(0)=Im(0)=0 and Eh(0)>0, Em(0)>0, Rh(0)>0, then we have Sh(t)>0, Ih(t)>0, Im(t)>0, Em(t)>0, Eh(t)>0, Rh(t)>0, ∀t>0.
For any cases, it follows that(Sh(t),Eh(t),Ih(t),Rh(t),Em(t),Im(t))∉∂X0 for t>0 sufficiently small, which contradicts the fact that ∂X0 is positively invariant. Hence, M∂⊂D. Thus, it then follows that M∂=D.
The equalityM∂=D implies that E1+ is a fixed point of P and acyclic in M∂; every solution in M∂ approaches to E1+. Moreover, Lemma 13 implies that E1+ is an isolated invariant set in X and Ws(E1+)∩X0=∅. By the acyclicity theorem on uniform persistence for maps, Theorem 1.3.1 and Remark 1.3.1 in [14], it follows that P is uniformly persistent with respect to X0. Thus, Theorem 3.1.1 in [14] implies that the periodic semiflow Φ(t):X→X is also uniformly persistent with respect to X0. Thanks to Theorem 1.3.6 in [14], model (11) has at least one ω-periodic solution u~(t,ψ∗) with ψ∗∈X0 and t≥0. Now, we show that u~(t,ψ∗) is positive.
Suppose thatψ∗=0; then, for all t>0, we obtain u~i(t,ψ∗)>0, for i=1,2,3,4,5,6. By using the periodicity of the solution, we have Sh∗(0)=Sh∗(nω)=0, Eh∗(0)=Eh∗(nω)=0, Ih∗(0)=Ih∗(nω)=0, Rh∗(0)=Rh∗(nω)=0, Em∗(0)=Em∗(nω)=0, Im∗(0)=Im∗(nω)=0, ∀n≥1, which contradicts the fact that u~i(t,ψ∗)>0 for i=1,2,3,4,5,6. So, the periodic solution is positive.


## 4. Numerical Simulation

In this section, we will present a series of numerical simulations of model (11) in order to support our theoretical results, to predict the trend of the disease, and to explore some control measures.


### 4.1. Initial Conditions and Estimation ofβ(t)

To validate our results, we choose the following initial conditions:E(0)=2400, L(0)=1200, Sh(0)=1500, Eh(0)=50, Ih(0)=200, Rh(0)=50, Sm(0)=3000, Em(0)=100, Im(0)=500, and A(0)=3600. Our numerical simulation will be performed using the MATLAB technical computing software with the fourth-order Runge–Kutta method [16].Using the method developed in [11], we express the biting rate as follows: (47)βt=α0-1.83692cos⁡0.523599t-0.175817cos⁡1.0472t-0.166233cos⁡1.5708t-0.16485cos⁡2.0944t-0.17681cos⁡2.61799t-1.37079sin⁡0.523599t+0.296267sin⁡1.0472t+0.2134sin⁡1.5708t-0.295228sin⁡2.0944t-0.201712sin⁡2.61799t,with α0≥3.


### 4.2. The Model Parameters and Their Dimensions

Numerical values of parameters are given in Table1.Table 1
Values for constant parameters for the malaria model.
Parameter
Description
Value
Reference
Dimension
Λ
Constant recruitment rate for humans
400
Estimated
Humans/month
dh
Human death rate
0.019
Estimated
/month
α
Transmission rate of humans fromEh to Ih
3.04
[17]
/month
dp
Disease-induced death rate for humans
0.0028
[11]
/month
rh
Recovery rate of humans
0.0159
[11]
/month
γ
Per capita rate of loss of immunity for humans
0.0167
[11]
/month
s
L
Transfer rate fromL to adult
15
[6]
/month
d
m
Death rate for adult vectors
3.4038
[11]
/month
ν
m
Transmission rate of mosquitoes fromEm to Im
2.523
[11]
/month
c
m
h
Probability of transmission of infection fromIm to Sh
0.022
[17]
Dimensionless
c
h
m
Probability of transmission of infection fromIh to Sm
0.48
[17]
Dimensionless
c
-
h
m
Probability of transmission of infection fromRh to Sm
0.048
[17]
Dimensionless
K
E
Available breeder sites occupied by eggs
30000
Estimated
Space
K
L
Available breeder sites occupied by larvae
18000
Estimated
Space
s
Transfer rate fromE to L
15
[6]
/month
b
Eggs laying rate
180
[6]
/month
d
Death rate of eggs
6
[6]
/month
dL
Larvae death rate
6
[6]
/month


### 4.3. Numerical Results

Using the above initial conditions, we now simulate model (11) in order to illustrate our mathematical results.By takingα0=7, dp=0.0028, cmh=0.022, chm=0.48, c¯hm=0.048, b=180, s=15, d=6, dL=7.5, sL=15, dm=3.4038 and considering the above initial conditions, we get r=25.1819, R0=1.3310>1 and  Figures 2, 3, and 4.Figure 2
Distribution of infected humans.
(a)
(b)Figure 3
Distribution of infected mosquitoes.
(a)
(b)Figure 4
Distribution of susceptible humans and mosquitoes.
(a)
(b)Figure2 describes the evolution of infected (exposed and infectious) humans. Figure 3 describes the evolution of infected (exposed and infectious) mosquitoes and Figure 4 describe the evolution of susceptible humans and mosquitoes. Figures 2 and 3 show that malaria remains persistent in the two populations. Besides, we observe that system (11) has one positive periodic solution. So, these numerical results illustrate the result of our Theorem 14.In order to understand the model behaviour around the disease-free equilibrium, we consider the same above initial conditions and the following values:α0=4, dp=0, cmh=0.022, chm=0.24, c¯hm=0.024, b=180, s=15, d=6, dL=7.5, sL=15, dm=6. Then we get r=14.2857 and R0=0.2602<1. Figures 5 and 6 illustrate that the disease dies out in both populations. Thus, the numerical results are the same as what we got in Theorem 11.Figure 5
Distribution of susceptible humans and mosquitoes.Figure 6
Distribution of infected humans and mosquitoes.


### 4.4. Parameters of Control of Malaria

Now, we assume that people became more conscious about the malaria disease and they use some efficient methods to reduce the proliferation of mosquitoes. That reduction can perhaps consist in fighting against the development of eggs, larvae, and pupa, firstly, by using chemical application methods (larvicide) or by introducing larvivore fish, and secondly, by using ecological methods (cleaning up the environment) to reduce the breeding sites of eggs and larvae. Letμ1,μ2∈[0,1[, respectively, be the efficiency of both intervention measures. So, we will use r~=(1-μ1)r, K~E=(1-μ2)KE, and K~L=(1-μ2)KL in order to evaluate their impact on the dynamics of malaria transmission.Thus, by considering the above initial conditions and by takingα0=7, dp=0.0028, cmh=0.022, chm=0.48, c¯hm=0.048, dm=3.4038, we obtain the following results.(i) Numerical Results for
μ
1
≃
89
%. For this value, we get r~=2.8204 and R0=0.6414. Moreover, according to Figure 7, we notice that the distribution of infected humans and mosquitoes has highly reduced and the malaria is progressively dying out in the populations.Figure 7
Distribution of infected humans and mosquitoes forb=80, s=10, d=15, sL=6, and dL=14.(ii) Numerical Results for
μ
2
=
80
%. Using μ2=0.8, we get K~E=6000, K~L=3600, and R0=0.5953. Further, Figure 8 clearly shows that the disease is quickly disappearing from the populations.Figure 8
Distribution of infected humans and mosquitoes forb=180, s=15, d=6, sL=15, and dL=7.5.Remark 15.
We must notice that the two parameters are important in the malaria transmission because a little perturbation of those parameters influences the dynamics of malaria transmission. So they can be used to fight against the persistence of the disease. The controlμ1 is efficient but its action is very slow in finite time, but the control μ2 is the best because it is more optimal and its action is very quick. Thus cleaning up the environment can be a very good mean of controlling malaria in the populations.


## 4.1. Initial Conditions and Estimation ofβ(t)

To validate our results, we choose the following initial conditions:E(0)=2400, L(0)=1200, Sh(0)=1500, Eh(0)=50, Ih(0)=200, Rh(0)=50, Sm(0)=3000, Em(0)=100, Im(0)=500, and A(0)=3600. Our numerical simulation will be performed using the MATLAB technical computing software with the fourth-order Runge–Kutta method [16].Using the method developed in [11], we express the biting rate as follows: (47)βt=α0-1.83692cos⁡0.523599t-0.175817cos⁡1.0472t-0.166233cos⁡1.5708t-0.16485cos⁡2.0944t-0.17681cos⁡2.61799t-1.37079sin⁡0.523599t+0.296267sin⁡1.0472t+0.2134sin⁡1.5708t-0.295228sin⁡2.0944t-0.201712sin⁡2.61799t,with α0≥3.


## 4.2. The Model Parameters and Their Dimensions

Numerical values of parameters are given in Table1.Table 1
Values for constant parameters for the malaria model.
Parameter
Description
Value
Reference
Dimension
Λ
Constant recruitment rate for humans
400
Estimated
Humans/month
dh
Human death rate
0.019
Estimated
/month
α
Transmission rate of humans fromEh to Ih
3.04
[17]
/month
dp
Disease-induced death rate for humans
0.0028
[11]
/month
rh
Recovery rate of humans
0.0159
[11]
/month
γ
Per capita rate of loss of immunity for humans
0.0167
[11]
/month
s
L
Transfer rate fromL to adult
15
[6]
/month
d
m
Death rate for adult vectors
3.4038
[11]
/month
ν
m
Transmission rate of mosquitoes fromEm to Im
2.523
[11]
/month
c
m
h
Probability of transmission of infection fromIm to Sh
0.022
[17]
Dimensionless
c
h
m
Probability of transmission of infection fromIh to Sm
0.48
[17]
Dimensionless
c
-
h
m
Probability of transmission of infection fromRh to Sm
0.048
[17]
Dimensionless
K
E
Available breeder sites occupied by eggs
30000
Estimated
Space
K
L
Available breeder sites occupied by larvae
18000
Estimated
Space
s
Transfer rate fromE to L
15
[6]
/month
b
Eggs laying rate
180
[6]
/month
d
Death rate of eggs
6
[6]
/month
dL
Larvae death rate
6
[6]
/month


## 4.3. Numerical Results

Using the above initial conditions, we now simulate model (11) in order to illustrate our mathematical results.By takingα0=7, dp=0.0028, cmh=0.022, chm=0.48, c¯hm=0.048, b=180, s=15, d=6, dL=7.5, sL=15, dm=3.4038 and considering the above initial conditions, we get r=25.1819, R0=1.3310>1 and  Figures 2, 3, and 4.Figure 2
Distribution of infected humans.
(a)
(b)Figure 3
Distribution of infected mosquitoes.
(a)
(b)Figure 4
Distribution of susceptible humans and mosquitoes.
(a)
(b)Figure2 describes the evolution of infected (exposed and infectious) humans. Figure 3 describes the evolution of infected (exposed and infectious) mosquitoes and Figure 4 describe the evolution of susceptible humans and mosquitoes. Figures 2 and 3 show that malaria remains persistent in the two populations. Besides, we observe that system (11) has one positive periodic solution. So, these numerical results illustrate the result of our Theorem 14.In order to understand the model behaviour around the disease-free equilibrium, we consider the same above initial conditions and the following values:α0=4, dp=0, cmh=0.022, chm=0.24, c¯hm=0.024, b=180, s=15, d=6, dL=7.5, sL=15, dm=6. Then we get r=14.2857 and R0=0.2602<1. Figures 5 and 6 illustrate that the disease dies out in both populations. Thus, the numerical results are the same as what we got in Theorem 11.Figure 5
Distribution of susceptible humans and mosquitoes.Figure 6
Distribution of infected humans and mosquitoes.


## 4.4. Parameters of Control of Malaria

Now, we assume that people became more conscious about the malaria disease and they use some efficient methods to reduce the proliferation of mosquitoes. That reduction can perhaps consist in fighting against the development of eggs, larvae, and pupa, firstly, by using chemical application methods (larvicide) or by introducing larvivore fish, and secondly, by using ecological methods (cleaning up the environment) to reduce the breeding sites of eggs and larvae. Letμ1,μ2∈[0,1[, respectively, be the efficiency of both intervention measures. So, we will use r~=(1-μ1)r, K~E=(1-μ2)KE, and K~L=(1-μ2)KL in order to evaluate their impact on the dynamics of malaria transmission.Thus, by considering the above initial conditions and by takingα0=7, dp=0.0028, cmh=0.022, chm=0.48, c¯hm=0.048, dm=3.4038, we obtain the following results.(i) Numerical Results for
μ
1
≃
89
%. For this value, we get r~=2.8204 and R0=0.6414. Moreover, according to Figure 7, we notice that the distribution of infected humans and mosquitoes has highly reduced and the malaria is progressively dying out in the populations.Figure 7
Distribution of infected humans and mosquitoes forb=80, s=10, d=15, sL=6, and dL=14.(ii) Numerical Results for
μ
2
=
80
%. Using μ2=0.8, we get K~E=6000, K~L=3600, and R0=0.5953. Further, Figure 8 clearly shows that the disease is quickly disappearing from the populations.Figure 8
Distribution of infected humans and mosquitoes forb=180, s=15, d=6, sL=15, and dL=7.5.Remark 15.
We must notice that the two parameters are important in the malaria transmission because a little perturbation of those parameters influences the dynamics of malaria transmission. So they can be used to fight against the persistence of the disease. The controlμ1 is efficient but its action is very slow in finite time, but the control μ2 is the best because it is more optimal and its action is very quick. Thus cleaning up the environment can be a very good mean of controlling malaria in the populations.


## 5. Conclusion

In this paper, we have presented a seasonal determinist model of malaria transmission. From the theoretical point of view, we have shown that the basic reproduction ratio,R0, is the distinguishing threshold parameter of the extinction or the persistence of the disease: if R0 is less than 1 malaria disappears in the human and mosquito populations and if it is greater than 1 malaria persists.It also emerges from our study that the transmission of malaria is highly influenced by the dynamics of immature mosquitoes and depends on the regulatory threshold parameter of the mosquito population,r. Thus, the severity of malaria increases with this parameter. So, the life cycle of the anopheles is a very important aspect that must be taken into account in malaria modeling.Moreover, we have shown that malaria transmission can be controlled by fighting against the proliferation of the mosquitoes, namely, by reducing the value ofr or by reducing the value of available breeder sites, KE and KL. We have proved that the reduction of the available breeder sites is a very efficient and more ecological method in fighting against malaria transmission. It then follows that environmental sanitation can be a very good means to control malaria in the endemic regions.However, it must be noticed that our model is limited due to the following reasons: (i) we have not considered the effect of climate change on the life cycle of mosquitoes. (ii) The larva and pupa class were not distinguished.In the future, one can develop a more realistic model by incorporating the above important factors and by considering the general force of infection. In addition, we can also take into account the degree of vulnerability of human populations in the model.



---

*Source: 6754097-2017-06-04.xml*