Global dynamics of an impulsive vector-borne disease model with time delays

: In this paper, we investigate a time-delayed vector-borne disease model with impulsive culling of the vector. The basic reproduction number R 0 of our model is first introduced by the theory recently established in [1]. Then the threshold dynamics in terms of R 0 are further developed. In particular, we show that if R 0 < 1, then the disease will go extinct; if R 0 > 1, then the disease will persist. The main mathematical approach is based on the uniform persistent theory for discrete-time semiflows on some appropriate Banach space. Finally, we carry out simulations to illustrate the analytic results and test the parametric sensitivity on R 0 .


Introduction
Vector-borne disease involves a vector (an organism) that transmits infectious pathogens from the infected host to the uninfected host.For example, the malaria parasites are transmitted through the bite of infected female mosquitoes, spreading the disease from mosquitoes to humans and back to mosquitoes [2].West Nile virus (WNv) is another vector-borne disease caused by Flavivirus and transmitted primarily by mosquitoes (Culex species) to vertebrate hosts, such as humans, birds and horses [3].Those vector-borne diseases have caused a significant threat to public health as well as wildlife worldwide.
Mathematical models give insights to predict the spread of vector-borne diseases, and to test control strategies.Impulsive (delayed) differential systems are widely used in population biology and epidemiology to study the various factors: birth pulses [4,5], effectiveness of vaccination [6][7][8][9][10], and elimination of vector insects [11,12], see also [13,14] for stochastic effects.To investigate the control strategies against vector-borne diseases, Gourley et al. [11] used two time-delayed models with impulses to evaluate the effectiveness of age-structured culling strategies.Yang et al. [10] studied the effect of impulsive controls in multiple patch models with age-structure.Xu and Xiao [9] considered the impacts of periodic impulsive culling of mosquitoes for control of WNv transmission.Based on [9], the incubation period of mosquitoes was furthered incorporated in [15].
To describe the disease transmission between vector and host as well as the strategy of periodically culling the vector, we use a set of impulsive delayed differential equations to incorporate the extrinsic incubation period (EIP) of vector and host.By the theory recently developed in [1], the basic reproduction ratio R 0 of the model is first established.The global dynamics of the model are further investigated in terms of R 0 .The uniform persistence theory based on dynamical systems is used extensively in population biology and epidemiology models, see for instance [16][17][18][19] and references therein.However, it seems like this dynamical system approach is rarely adapted to time-delayed impulsive systems to study the threshold dynamics.There are two possible issues: first, it is not an easy task to find an appropriate phase space for time-delayed impulsive systems so that the evolution operators are well-defined, which was recently solved by [1] for a large class of time-delayed impulsive systems; second, the solutions of the time-delayed impulsive systems may not be continuous in t, and hence, such systems might not generate continuous periodic semiflows (see [20]), which brings trouble for direct application of the uniform persistence theory in [21] to investigate the sharp permanence of infectious compartments when R 0 > 1.Meanwhile, the analytic method developed in [1,6,7] for the uniform persistence of a time-delayed impulsive system is quite useful for a form like Susceptible-Exposed-Infected-Removed (SEIR), but it would not be easy to apply to the vector-borne model (including multiple infected compartments).We refer to [9,15] for attempts in this direction.Combing the ideas in [1] and [21], we establish the sharp persistence of infectious compartments, by using the uniform persistence theory of discrete-time semiflows on some appropriate phase space.It is worth pointing out that our current approach would be easily applied to the models where the period of model coefficients is the same as that of fixed impulsive moments.
The rest of this paper is organized as follows: In Section 2, we formulate a time-delayed impulsive differential model including several factors: EIP and the strategy of periodically culling the vector, EIP and vertical transmission of host.In Section 3, we first introduce the basic reproduction number R 0 for the model and then study the threshold dynamics in terms of R 0 .In Section 4, we use numerical simulations to test the differences with and without culling, and further explore the impacts of parameters such as time delays and vertical transmission rate on R 0 .A brief discussion then concludes the paper.

Model formulation
Inspired by the WNv infection process in [9,15,18], we consider an impulsive delay differential equation model to describe the periodic culling of the vector with two differnt EIPs.Let S v (t), I v (t), S h (t), E h (t) and I h (t) be the the total numbers of the susceptible adult vector, infected adult vector, susceptible host, exposed host and infected host at time t, respectively.The total number of the adult vector is given by N h (t) = S h (t) + E h (t) + I h (t).Let Λ v and Λ h be the recruitment rate of the vector and (susceptible) host, respectively, µ v and µ h be the natural death rate of vector and host, and d h be the disease-induced death rate of the host.Similar to [7], we suppose that the birth rate of the vector equals its natural death rate.Vertical transmission of the virus in the host population is incorporated by a fraction r vt ∈ (0, 1).For the virus transmission, we assume that the susceptible host becomes infectious given contact with the infected vector, and the susceptible vector can receive the infection by biting a virus-carrying host.By similar arguments to those in [2,17,19], we simply assume that the numbers of newly occurred infectious host and newly occurred infected vector per unit time at time t are given by respectively, where β is the average biting rate of the vector, b and c are the transmission probabilities of the virus from vector to host and from host to vector, respectively.However, the newly infected vector and infected host need to survive the EIP to become infectious.We denote the finite constants τ 1 and τ 2 to represent the length of the EIP in the vector and host, respectively.The probability that the vector and host survive the EIP is e −µ v τ 1 and e −µ h τ 2 , then of those vectors and hosts infected τ 1 and τ 2 unit times ago, only the proportions are infectious at time t.Throughout this paper, we suppose that culling occurs at certain particular times, and denote T to be the period of culling.We also assume that spraying reduces both the susceptible and infected vector, p ∈ [0, +∞) is the culling intensity of those killed, where p = 0 means that there is no culling.The transmission of cross-infection between vector and host is shown in Figure 1.Consequently, our model takes the form: All the parameters in model (2.1) are positive.In view of the biological meaning of τ 2 , we impose the following compatibility condition: To investigate the long time behavior of system (2.1) from the point view of dynamical systems, we first introduce a few notations related to the phase space (see also [1,22]).Given two constants a < b, and ϕ(t + ) = ϕ(t) for all but at most a finite number of points t ∈ [a, b) , and ϕ(t + ) = ϕ(t) for all but at most a finite number of points t ∈ [a, b) , where ϕ(t + ) = lim s→t + ϕ(s) and ϕ(t − ) = lim s→t − ϕ(s).Then for any r > 0, PC([−r, 0], R m ) is a Banach space with the norm ∥ • ∥ r given by ∥ϕ∥ r = sup and ϕ(t + ) = ϕ(t) for t ∈ J}, and PC + J := {ϕ ∈ PC J : ϕ(t) ≥ 0, ∀t ∈ [−τ, 0]}.It then follows that PC J is a Banach space endowed with the norm ∥ • ∥ τ, and PC + J is a closed cone of PC J , which induces a partial ordering on PC J .Clearly, PC J ⊂ PC([−τ, 0], R).For any u = (u 1 , u 2 , ..., u 5 ) ∈ PC([−τ, η), R 5 ) with η > 0, define for any ϵ ∈ 0, Λ h µ h +d h .Now we are ready to state the well-posed result on the solutions of system (2.1).Theorem 2.1.For any ϕ ∈ D ϵ , system (2.1) exists a unique non-negative solution u (t, ϕ) on [0, +∞) with the initial value u 0 = ϕ, u nT ∈ D ϵ for any n ∈ N.Moreover, u t+T (ϕ) = u t (u T (ϕ)) for all t ≥ 0 and Φ := u T admits a global attractor in D ϵ .

Threshold dynamics in terms of R 0
In this section, we first introduce the basic reproduction number R 0 for system (2.1) and then investigate its threshold dynamics in terms of R 0 .

Basic reproduction number
The basic reproduction number R 0 is defined as the expected number of secondary infections produced by a typical infectious individual in a completely susceptible population.We will use the theory developed in [1] to give the definition of R 0 .In order to obtain the disease-free periodic solution of system (2.1), we recall a lemma (see [6, Lemma 1] for detail).Lemma 3.1.Consider the following impulsive differential equation: where a, b > 0, θ ∈ (0, 1).Then system (3.1)admits a unique positive periodic solution ûe (t) : . Letting I v = I h = 0 in system (2.1), we then get the following disease-free system: which is globally asymptotically stable.
It then follows that X T is a Banach space.Linearizing system (2.1) at its disease-free periodic solution E 0 (t) = ( Ŝ v (t), 0, Ŝ h (t), 0, 0), we obtain the following linear system for the infectious compartments: where Set Y(t, s), t ≥ s as the Cauchy matrix [24, Section 1.2] of where each t k denotes the impulsive point on [s, t) and . It is easy to check that F(t) and Y(t, s) satisfy the following properties: (a) for each t ∈ R, F(t) is a positive operator; (b) the matrix −V is cooperative and r(Y(T, 0)) < 1, where r(Y(T, 0) is the spectral radius of Y(T, 0).Consequently, the linear impulse periodic differential system (3.5) could be rewritten as: where P = diag ( 1 1+p , 1).After the above settings, we can apply the theory and method in [1] to define basic the reproduction number for system (2.1).Let v(t) be the initial distribution of infected individuals with period T .For t ≥ s, then F(t−s)v t−s is the distribution of newly infected individuals at time t−s.It is produced by infected individuals introduced during the time interval [t − s − τ, t − s], which means it is the distribution of newly infected at time t − s and still infected at time t.Define the linear operator L on X T : We define the basic reproduction number as where r(L) is the spectral radius of the linear operator L. Denote For any ϕ ∈ E J , let w(t, ϕ) = (w 1 (t, ϕ), w 2 (t, ϕ)) be the unique solution of (3.4) with the initial value w 0 = ϕ, and be the solution map of (3.4), then In view of the definition of PC J , we infer that Q := Q(T ) is the operator from E J to E J (in other words, they have the same number of discontinuous points), and set r(Q) to be the spectral radius of Q on E J .By [1, Theorem 1], we have the following observation: Lemma 3.2.R 0 − 1 has the same sign as r(Q) − 1.
where PC τ i is understood as PC J with τ replaced by τ i , and then (Z, Z + ) is an ordered Banach space.For any ψ ∈ Z + , system (3.4) admits a unique non-negative solution z(t, ψ) with the intial value z 0 = ψ.For all t ≥ 0, we have z t (ψ) = (z 1t (ψ), Z 2t (ψ)), for all θ i ∈ [−τ i , 0], and then Z it (ψ)(θ i ) = z i (t + θ i , ψ), i = 1, 2. For any given t ≥ 0, let Q(t) be the solution map of system (3.4) on Z, that is, We remark that z 1t ∈ PC([−τ 2 , 0], R) and z 2t ∈ PC τ 1 for all t ≥ 0 as z 2 (t, ψ) is continuous for all t ≥ 0 even if ψ is piecewise continuous.Therefore, Q(T ) is a map from Z to Z. Next, we will illustrate that Q(t) is eventually strongly positive on Z + .
By the similar argument to that in Theorem 2.1, we obtain that z i (t) ≥ 0 for any t ≥ 0, i = 1, 2.
In the case that φ 2 > 0. We have Then by the integral form of the second equation of system (3.4),we have Denote r( Q) as the spectral radius of Q := Q(T ).One might perform the same argument as in [16, Lemma 3.8] to obtain r(Q) = r( Q).Moreover, we have the following observation: Then there exists a positive T-periodic function ṽ(t) = (ṽ 1 (t), ṽ2 (t)) such that e µt v(t) is a solution of system (3.4) with the feasible domain either Z + or E + J for any t ≥ 0. Proof.Resembling the arguments in [1, Lemma 4], we infer that for each t > τ, Q(t) and Q(t) are compact on E J and Z, respectively.
In the case that the feasible domain is Z + , fix an integer n 0 > 0 such that n 0 T > 3τ.It follows that Qn 0 = Q(n 0 T ) is compact and strongly positively on Z.By [25, Lemma 3.1], we obtain that there is a strongly positive eigenvector φ * = (φ * 1 , φ * 2 ) such that Q(φ * ) = r( Q)φ * .Let v(t, φ * ) be the solution of system (3.4) with the initial value Since φ * ≫ 0, it is easy to see that v t (φ * ) ≫ 0 for any t ≥ 0. Let v(t) = e −µt v(t, φ * ).Then v(t) ≫ 0 satisfies In the case that the feasible domain is E + J , since e µt v(t) given in the above is a function for t ∈ R. We denote φ = ( φ1 , φ2 ) as: By the uniqueness of solutions, we have that w(t, φ) := e µt v(t) satisfies (3.4) for all t ≥ 0 with the initial value w 0 = φ ∈ E + J .

Threshold dynamics
We are now in a position to prove a threshold-type result of system (2.1) in terms of R 0 .
Proof.By the first and sixth equations of system (2.1), we have Consider the following auxiliary system According to Lemma 3.1 and the comparison theorem on impulsive differential equations, we have S v (t, a) ≤ x(t, a) with the same initial condition S v (0) = x(0) = a > 0, ∀t ≥ 0 and |x(t) − Ŝ v (t)| → 0 as t → ∞ with Ŝ v (t) given in (3.3).Then for any δ ∈ (0, Λ h µ h ), there exists t1 > 0 such that When µ h ≥ d h , in view of system (2.1), we have It then follows that for the above δ, there exists t2 > t1 such that Thus, we have We consider the following impulsive differential equations with parameter δ: (3.9) From Lemma 3.2, we know that R 0 < 1 if and only if r(Q) < 1.Let Q δ be the Poincaré map of system (3.9).Observe that lim δ→0 + r(Q δ ) = r(Q) < 1 (as r(Q δ ) is non-decreasing and upper semi-continuous in δ ≥ 0 but close to 0 (see [26,27] for example), we can fix a small positive number δ such that r(Q δ ) < 1.By Lemma 3.4, there is a positive T -periodic function vδ (t) = (v δ 1 (t), vδ 2 (t)) such that e µ δ t vδ (t) is a positive solution of system (3.9),where µ δ = ln r(Q δ ) T < 0. Choose a positive constant K δ such that: Applying the comparison principle for cooperative impulsive delay differential systems [28, Lemma 2.2], we obtain that Letting t → ∞ in (3.10), we have (I v (t), I h (t)) → (0, 0).It then follows from the fourth equation or (2.4) that E h (t) → 0 as t → ∞.In view of the theories of asymptotically periodic semiflows and internally chain transitive sets [21], we further deduce from the third equation of system (2.1) that Then for any given δ 1 ∈ (0, δ), there exists t3 > t2 + τ such that and hence, For this inequality, we consider the following impulsive differential equation: By Lemma 3.1, system (3.11)admits a globally stable T -periodic solution, denoted as Ŝ δ 1 v (t), which depends continuously on δ 1 .It follows that for the above δ, there exists a t4 > t3 such that This, together with (3.8), yields for any δ 1 ∈ (0, δ) Letting δ 1 → 0 + , we obtain for any small δ > 0, that there holds which implies that lim In the remainder of this section, we investigate the uniform persistence of system (2.1).
In the case that I v (t, ψ) ≡ 0 for each t ≥ 0, it follows from the second equation of system (2.1) that I h (t − τ 1 )S v (t − τ 1 ) = 0 for any t ≥ 0. By the first and sixth equations of (2.1), we have This implies there exists t1 > 0, such that S v (t) > 0 for all t > t1 + τ 1 , and hence, I h (t) = 0 and E h (t) = 0 for all t ≥ 0. Since disease-free system (3.2) admits a globally stable T periodic solution ( Ŝ (•), Ŝ h ), we obtain that ω(ψ) = M.
The above claims indicate that M cannot form a cycle for Φ in D ϵ and W s (M) ∩ X 0 = ∅, where W s (M) is the stable set of M for Φ.Now we define a continuous distance function p : D ϵ → R + by p(ϕ) = min{ϕ 2 (0), ϕ 5 (0)}.In view of Claim 1, we have W s (M) ∩ p −1 (0, ∞) = ∅.Now by [21, Lemma 1.2.1 and Theorem 1.3.2]and Theorem 2.1, it then follows that Φ is uniformly persistent with respect to (X 0 , ∂X 0 , p), that is, there exists Now for any t ≥ n 2 T , letting t = nT + t with n = [ t T ] and t ∈ [0, T ), we obtain from the equation of Similarly, we see from > 0, we see that lim inf t→∞ I i (t) ≥ ρ * , i = v, h.

Numerical simulations
In this section, we perform illustrative numerical simulations to verify theoretical results and explore the influences of key model parameters on the disease transmission.Note that the following simulations were based on the WNv disease transmission process studied in [9,15], where vector and host are mosquito and bird, respectively, and some of parameters were chosen only for the test of parameter sensitivity on the basic reproduction number.
Similar to [15], choose T = 10, p = 0.6 while R 0 = 1.1490.Figure 2(a) shows that the number of infected vectors tends to the periodical oscillation when culling occurs, whereas it tends to a equilibrium without culling.In addition, we see from Figure 2(a) that culling significantly reduces the amount of infections of the vector compared to the case of no culling.Choose T = 10, p = 2 while R 0 = 0.7142.As shown in the Figure 3(b), the number of infected hosts that are not culled will continue to increase, and after the action of culling, the number of infected hosts will decrease, and the disease will eventually die out.Next, we use PRCCs (partial rank correlation coefficients) to obtain the sensitivity analysis of R 0 .We take Λ v , Λ h , µ v , µ h , β, b, c, p, d h and r vt as the input variables and the value of R 0 as the output variable.Figure 4 illustrates that parameters Λ h , µ v , µ h , p and d h are negatively correlated with R 0 and the others are positively correlated.We also see that R 0 is more sensitive to µ h , β, b, c and p.Thus, the corresponding control measures should be taken for these sensitive parameters.For instance, to reduce mosquito-bird contact and further control the disease spread, we could burn repellent plants in the bird habitat or at their water sources.In Figure 5, we investigate the joint effects of β, p on R 0 .It shows that when the biting rate β is large, it requires very strong culling effect.In Figure 6 we also choose different r vt = 0.01, r = 0.05, r = 0.1, while β = 0.05, and get R 0 increasing in r vt (compare (a)-(c)).Both Figures 5 and 6 indicate that when all parameters are equal, R 0 is decreasing in p.All the results are consistent with the conclusion in sensitivity analysis.Finally, we examine the impact of the EIP on R 0 .We let τ 1 vary from 0 to 10 and τ 2 vary from 0 to 5. Figure 7 describes the dependence of R 0 on τ 1 and τ 2 for three different frequencies of culling (i.e. 1 T ).As we see from this figure, R 0 is decreasing with respect to τ 1 as well as τ 2 and it increases as T increases for fixed τ 2 and τ 1 .This suggests that culling during the EIP and prolonging the EIP would be beneficial for disease control, while infrequent culling might be counterproductive.

Discussion and conclusions
In this paper, we have formulated an impulsive vector-borne disease model with time delays to investigate the joint effects of the EIP and impulsive intervention.The basic reproduction number R 0 is first derived by the theory in [1], which serves as a threshold value to determine the extinction and uniform persistence of the disease.Unlike most existing works [1,6,9,15], we utilize the dynamical system approach to show the sharp uniform persistence as R 0 > 1.As is well known, the uniform persistence theory often gives rise to the existence of a positive periodic solution (see [21,Theorem 1.3.10]).As a complement, we could also show the existence of a positive periodic solution in [1,9,15] via our strategy.However, we emphasize that since our feasible domain D ϵ in (2.3) is non-convex, we could not directly utilize it to verify the existence of the positive periodic solution, and we will leave

Figure 2 .
Figure 2.While R 0 > 1, comparison of the long-term behavior of infectious vector (the left plot) and host (the right plot) in different scenarios: culling and without culling.

Figure 3 .
Figure 3.While R 0 < 1, comparison of the long-term behavior of infectious vector (the left plot) and host (the right plot) in different scenarios: culling and without culling.

Figure 5 .
Figure 5.The contour plot of R 0 with respect to p and β with r vt = 0.001.

2 Figure 7 .
Figure 7.The curve of R 0 with respect to τ 1 (the left panel) and τ 2 (the right panel) for different culling intervals.