New mathematical model of vertical transmission and cure of vector-borne diseases and its numerical simulation

Abstract

In this research article, a new mathematical model for the transmission dynamics of vector-borne diseases with vertical transmission and cure is developed. The non-negative solutions of the model are shown. To understand the dynamical behavior of the epidemic model, the theory of basic reproduction number is used. As this number increases, the disease invades the population and vice versa. The effect of vertical transmission and cure rate on the basic reproduction number is shown. The disease-free and endemic equilibria of the model are found and both their local and global stabilities are presented. Finally, numerical simulations are carried out graphically to show the dynamical behaviors. These results show that vertical transmission and cure have a valuable effect on the transmission dynamics of the disease.

1 Introduction

Vector-borne diseases are infectious diseases transmitted to humans and animals by blood-feeding arthropods. Some common vector-borne diseases are West Nile virus, dengue fever, Rift Valley fever, malaria, and viral encephalitis caused by pathogens such as bacteria, viruses, and parasites. The arthropods are blood sucking insects and arachnids such as ticks, mosquitoes, biting flies, and lice called vectors [1]. The vectors receive pathogens from an infected host and transmit them to a human host, as humans are the major host, or animals. However, direct transmissions, such as transplantation related transmission, transfusion related transmission, and needle-stick-related transmission, are also possible [2]. In case of some diseases such as AIDS and Hepatitis B, it is possible for the offspring of infected parents to be born infected. This type of transmission is called vertical transmission. Now it is found that vector-borne diseases can also be transmitted vertically [3, 4]. Also new research shows that virus is transmitted from female mosquitos to their eggs at a high rate [5], which causes vertical transmission of the disease.

Vector-borne diseases are prevalent in hot areas, such as tropics and subtropics, and are relatively rare in temperate zones. Vector-borne infectious diseases remain amongst the most important cause of global health illness and are major killers, particularly of children. The World Health Organization reports the numbers of deaths in different regions of the world annually. Nearly 700 million people get mosquito-borne illnesses that cause about one million deaths each year. Worldwide, malaria is the leading cause of premature mortality, particularly in children under the age of five. Nearly half of the world’s population is at risk of malaria, and every year 198 million cases (uncertainty range: 124–283 million) and 584,000 deaths (range: 367,000–755,000) occur according to the World Malaria Report 2014 [6]. According to WHO, an estimated of 3.3 billion people in 97 countries are at risk of malaria. Currently, dengue threatens up to 40% of the world’s population, and there may be 50–100 million infections annually [7]. More than 2.5 billion people over 40% of the world’s population are now at risk of dengue.

From the above discussion it is clear that it is necessary to control such epidemic diseases. Control measures for vector-borne diseases are important because most are zoonoses. For the control measure, it is necessary to understand the dynamical features of diseases and treat the infected hosts. Therefore, deciphering the mechanisms and modeling of such diseases are of great interest. Our paper involves such an epidemic model for the transmission dynamics of vector-borne diseases that incorporates both horizontal and vertical transmission in the vector–host population.

Up to date, many mathematical models have been investigated to understand the mechanism of real world phenomena. Researchers investigate different methods to solve these models both analytically and numerically (e.g., see [8–21]). Several models of infectious diseases have been developed in the literature [22–27]. The model first proposed by Ross [28] and subsequently modified by Macdonald [29] has influenced both the modeling and the application of control strategies to a vector-borne disease. The model presented in [30] studied the analysis of a simple vector–host epidemic model with horizontal transmission. We extend their model by including vertical transmission in both vector and host populations, and treatment class in the host population with different interaction rates.

The structure of this paper is as follows: Section 1 represents the introductory remarks with a brief history. Section 2 is about the derivation of SITR epidemic model and shows the non-negative solutions of the proposed model. In Section 3, we find the disease-free and endemic equilibria and prove their local stability. In Section 4, we use mathematical analysis to establish global stability results for the proposed model. We use Lyapunov function theory to show global stability of both disease-free and endemic equilibria. Parameter estimation and numerical results are discussed in Section 5. Finally, we give conclusion.

2 Model framework

The total population sizes at time t for human hosts and vectors are denoted by $N_{1} (t)$ and $N_{2}(t)$, respectively. The population of size $N_{1} (t)$ is divided into four distinct classes: the susceptible population of size $S (t)$, the infectious population of size $I (t)$, the population under treatment of size $T(t)$, and the recovered population of size $R(t)$. Thus $N_{1} (t) = S (t) + I (t) + T(t)+R (t)$. The vector population $N_{2}(t)$ has the subclasses denoted by $V(t)$ and $W(t)$ for the susceptible and infected classes, respectively. Thus, $N_{2}(t) = V(t) + W(t)$. The mathematical model can be represented by the following nonlinear system of ordinary differential equations:

$$ \begin{gathered} \frac{dS}{dt} =(1- \epsilon_{1}I)b_{1}-\beta_{1}S I - \beta_{2}S W-\mu_{1}S , \\ \frac{dI}{dt} =\epsilon_{1}b_{1}I+\beta_{1}S I +\beta_{2}S W-\alpha I -\eta I-\delta_{1} I- \mu_{1}I, \\ \frac{dT}{dt} =\alpha I -\gamma T-\delta_{1} T - \mu_{1}T , \\ \frac{dR}{dt} =\eta I+\gamma T-\mu_{1}R, \\ \frac{dV}{dt} =(1-\epsilon_{2}W)b_{2}- \beta_{3}VI -\mu_{2}V, \\ \frac{dW}{dt} =\epsilon_{2}b_{2}W+ \beta_{3}VI -\delta_{2}W-\mu_{2}W, \end{gathered} $$

(1)

with the initial conditions

$$ S(0)\geq0,\qquad I(0)\geq0,\qquad T(0)\geq0,\qquad R(0)\geq0,\qquad V(0)\geq0,\qquad W(0)\geq0. $$

(2)

The human host population is recruited at a constant birth rate $b_{1}$ in which a fraction $\epsilon_{1}$ were born infected from their infected parents. $\beta_{1}$ is the rate of direct transmission of the disease, $\beta_{2}$ is the vector mediated transmission rate, $\mu_{1}$ is the natural mortality rate of a human. Infectious humans are treated at a rate α, recover naturally at a rate η, and suffer disease-induced death at a rate $\delta_{1}$. Treated humans recover at a rate γ. It is assumed that recovered individuals acquire lifelong immunity against re-infection. Similarly, $b_{3}$ is the constant recruitment rate of vector population in which the ratio $\epsilon_{2}$ are infected by birth from their infected parents. Susceptible mosquitoes become infected by biting infected human at a rate $\beta_{3}$, $\mu_{2}$ is the natural mortality rate of vector population. Infectious vectors die due to disease at a rate $\delta_{2}$. The complete dynamics of the proposed model is represented by the flow chart in Figure 1.

Figure 1

The flow chart represents the interactions and transfer of a vector-borne disease in both human and vector populations

2.1 Properties of solutions

The proposed model (1) is a system of nonlinear ordinary differential equations with the initial conditions (2). To be epidemiologically and mathematically meaningful, it is important to prove that all the solutions with the given initial conditions will remain non-negative and bounded for all finite time. The model shall be analyzed in a biologically meaningful feasible region governed by a positive invariant set.

Theorem 2.1

There exists a unique and bounded solution of the system of equations (1), in a positively invariant set, that remains for all finite time $t\geq0$.

Proof

The right-hand side of each equation is continuous in the convex domain $E=(t,S(t),I(t), T(t),R(t),V(t),W(t))$ of $(6+1)$-dimensional space $R_{+}^{6+1}$ with continuous partial derivatives. So problem (1) has a unique solution in $R_{+}^{6}$ which exists for a given finite time $t\in[0,\infty)$ and initial conditions (2).

As the total population sizes are $N_{1}=S+I+T+R$ and $N_{2}=V+W$, so from (1) we get

$$ \frac{dN_{1}}{dt}=b_{1}-\mu_{1}N_{1}- \delta_{1}(I+T) \quad \text{and}\quad \frac {dN_{2}}{dt}=b_{2}- \mu_{2}N_{v}-\delta_{2}I_{v}. $$

(3)

Then

$$ \begin{gathered} \frac{dN_{1}}{dt}\leq b_{1}-\mu_{1}N_{1} \quad \text{and}\quad \frac{dN_{2}}{dt}\leq b_{2}-\mu_{2}N_{v}. \\ \quad \Rightarrow\quad N_{1}\leq N_{1}(0)e^{-\mu_{1}(t)}+ \frac{b_{1}}{\mu_{1}} \bigl(1-e^{-\mu_{1}(t)} \bigr)\quad \text{and}\\ \hphantom{\quad \Rightarrow\quad{}}N_{2}\leq N_{2}(0)e^{-\mu_{2}(t)} +\frac{b_{2}}{\mu_{2}} \bigl(1-e^{-\mu_{2}(t)} \bigr), \end{gathered} $$

which shows that

$$ \lim_{t\rightarrow\infty}\sup N_{1}\leq\frac{b_{1}}{\mu_{1}} \quad \text{and}\quad \lim_{t\rightarrow\infty} \sup N_{2}\leq \frac{b_{2}}{\mu_{2}}. $$

(4)

The given initial conditions (2) make sure that $N_{1}(0)\geq0$ and $N_{2}(0)\geq0$. Thus the feasible region for system (1) is

$$\Phi= \biggl\{ (S,I,T,R,V,W)\in R^{6}_{+}, N_{1}\leq \frac{b_{1}}{\mu_{1}}, N_{2}\leq\frac{b_{2}}{ \mu_{2}} \biggr\} . $$

Thus the total populations and each population class remain bounded for all finite time $t\geq0$. □

The above theorem shows that model (1) is well posed epidemiologically and mathematically in a positively invariant set Φ. We shall study the dynamics of this basic model in Φ, so, all the solutions of system (1) start and remain in Φ for all $t\geq0$. All the parameters and state variables for the model should be non-negative for all time because they represent the number of the population sizes of humans and vectors.

3 Equilibrium points

3.1 Disease-free equilibrium

The ability to invade a population is an important concern of an infectious disease. The steady state solutions of an epidemiological model at which the population remains in the absence of disease is called disease-free equilibrium point. In order to find the disease-free equilibrium of the proposed model (1), we set the right-hand side of all equations equal to zero and set $I=T=0$ and $W=0$. Also there is no infected recruitment in the populations, so we put the parameters $\epsilon_{1}=\epsilon_{2}=0$, which implies that $(1-\epsilon _{1})b_{1}=b_{1}$ and $(1-\epsilon_{2})b_{2}=b_{2}$ mean that the total recruited population is only susceptible. By direct calculations, we get the disease-free equilibrium point $E_{1}$ in the feasible region Φ, which is given by

$$E_{1}=(S_{1},I_{1},T_{1},R_{1},V_{1},W_{1})= \biggl(\frac{b_{1}}{\mu_{1}},0,0,0,\frac {b_{2}}{\mu_{2}},0\biggr). $$

The dynamics of model (1) is analyzed by a dimensionless number called basic reproduction number denoted by $R_{0}$, defined as “The expected number of secondary cases produced by a typical infected individual during its entire period of infectiousness in a completely susceptible population” [31]. Mathematically, $R_{0}$ is defined as

$$R_{0}\propto \biggl(\frac{\text{infection}}{\text{contact}} \biggr)\cdot \biggl( \frac{\text{contact}}{\text{time}} \biggr)\cdot \biggl(\frac{\text{time}}{\text{infection}} \biggr). $$

More precisely,

$$R_{0}\propto T\cdot C\cdot D, $$

where T is the transmissibility (i.e., probability of infection given contact between a susceptible individual and an infected one), C is the average rate of contact between susceptible and infected individuals, and D is the duration of infectiousness. This quantity serves as a threshold parameter that predicts whether a disease will spread in a community or will simply die out. It can be calculated by the method of next generation matrix given in [32]. In the vector–host model (1), infected states are I, T, and W and uninfected states are S, R, and V. The matrices $\mathcal{F}$ and $\mathcal{V}$ are the rate of production of new infections and the transition rates between states, respectively, which are given by

$$\mathcal{F}= \left ( \begin{matrix} \epsilon_{1}b_{1}I+\beta_{1}SI+\beta_{2}SW\\ 0\\ 0 \end{matrix} \right ),\qquad \mathcal{V}= \left ( \begin{matrix} (\alpha+\eta+\delta_{1}+\mu_{1})I\\ -\alpha I+(\gamma+\delta_{1}+\mu _{1})T\\ -\epsilon_{2}b_{2}W-\beta_{3}VI+(\delta_{2}+\mu_{2})W \end{matrix} \right ). $$

At the disease-free equilibrium $S=N_{1}=\frac{b_{1}}{\mu_{1}}$, $I=T=0$, $V=N_{2}=\frac{b_{2}}{\mu_{2}}$, and $W=0$. The Jacobian matrices at the disease-free equilibrium of $\mathcal{F}$ and $\mathcal{V}$ are F and V, respectively, where

$$\begin{gathered} F=\left ( \begin{matrix} \epsilon_{1}b_{1}+\beta_{1}N_{1}&0&\beta_{2}N_{1}\\ 0&0&0\\ 0&0&0 \end{matrix} \right ),\\ V=\left ( \begin{matrix} \alpha+\eta+\delta_{1}+\mu_{1}&0&0\\ -\alpha&\gamma +\delta_{1}+\mu_{1}&0\\ -\beta_{3}N_{2}&0&-\epsilon_{2}b_{2}+\delta_{2} +\mu_{2} \end{matrix} \right ). \end{gathered} $$

F and V are the rates for new infections and transitions near the equilibrium. We used MATLAB(R2010A) to find $V^{-1}$ and $FV^{-1}$, which gives the times spent in each state and the total production of new infections over the course of an infection, respectively. The largest eigenvalue of $F{V}^{-1}$ is the basic reproduction number $R_{0}$, given by

$$R_{0}=\frac{\epsilon_{1}b_{1}+\beta_{1}N_{1}}{k}+\frac{\beta_{2}\beta_{3}N_{1}N_{2}}{ mk}, $$

where $k=\alpha+\delta_{1}+\mu_{1}+\eta$ and $m=\delta_{2}+\mu_{2}-\epsilon_{2}b_{2}$. When there is no vertical transmission, $\epsilon_{1}=\epsilon_{2}=0$, then $R_{0}$ is the basic reproductive number for the model with only horizontal transmission. Geometrically it means that the number of new infections comes from both direct and indirect transmission. In the presence of vertical transmission, $\epsilon_{1},\epsilon_{2} > 0$, $R_{0}$ increases as these vertical transmission parameters increase, because vertical transmission directly increases the number of infectious populations. Also we can see the inverse relation of treatment strategies with $R_{0}$ and the direct relation with new infections and total population.

The basic reproduction number $R_{0}$ has a significant effect on the dynamics of infection. As we can see from the first and second equations of model (1),

$$ \frac{dS}{dt}=b_{1}-kR_{0}I- \mu_{1}S,\qquad \frac {dI}{dt}=k(R_{0}-1)I. $$

(5)

When $R_{0}<1$, it means that each infected individual infects less than one other individual averagely by ever kind of transmission, then the change in the number of infected population is negative, so the disease simply dies out. On the other hand, when $R_{0}>1$, it means that each infected individual infects more than one other individual, then the change is positive and invasion is always possible (see the survey paper by Hethcote [33]). For $R_{0}=1$, it means that each infectious individual infects one other individual as a whole, then there is no change in the infected population, so the infection constantly remains in the population. Also the effect of $R_{0}$ on the susceptible population is shown in the first equation of (5). All these facts are shown in Figures 2 and 3.

Figure 2

The first figure shows that $R_{0}$ decreases with increasing cure rate. The second figure shows that $R_{0}$ increases as vertical transmission increases

Figure 3

The figures show the threshold behavior of $R_{0}$ and its critical value $R_{0}=1$

Theorem 3.1

The disease-free equilibrium point $E_{1}$ is locally asymptotically stable if $R_{0}<1$, otherwise unstable.

Proof

This can be proved by linearizing system (1) around $E_{1}$, which gives the following Jacobian matrix:

$$J_{1}= \begin{bmatrix} -\mu_{1} & -\epsilon_{1}b_{1}-\beta_{1}\frac{b_{1}}{\mu_{1}} &0& 0 & 0 & -\beta_{2}\frac{b_{1}}{\mu_{1}}\\ 0 & \epsilon_{1}b_{1}+\beta_{1}\frac{b_{1}}{\mu_{1}}-k & 0 & 0 & 0 & \beta _{2}\frac{b_{1}}{\mu_{1}} \\ 0 & \alpha& -l & 0 & 0 & 0 \\ 0 & \eta& \gamma& -\mu_{1} & 0 & 0 \\ 0 & -\beta_{3}\frac{b_{2}}{\mu_{2}}& 0 & 0 &-\mu_{2}&-\epsilon_{2}b_{2}\\ 0 & \beta_{3}\frac{b_{2}}{\mu_{2}} & 0 & 0 & 0& -m \end{bmatrix} , $$

where $l=\gamma+\delta_{1}+\mu_{1}$.

The characteristic equation of $J_{1}$ is

$$ (x+\mu_{1}) (x+\mu_{1}) (x+\mu_{2}) (x+l) \bigl(c_{0}x^{2}+c_{1}x+c_{2} \bigr)=0, $$

(6)

where

$$\begin{aligned}& c_{0}=\mu_{1}\mu_{2}, \\& c_{1}=k\mu_{1}\mu_{2}+m\mu_{1} \mu_{2}-\beta_{1}b_{1}\mu_{2}-b_{1} \epsilon_{1}\mu _{1}\mu_{2}, \\& c_{2}=km\mu_{1}\mu_{2}(1-R_{0}). \end{aligned}$$

Four eigenvalues $-\mu_{1}$, $-\mu_{1}$, $-\mu_{2}$, and −l out of six have a negative real part. The remaining two eigenvalues are the roots of the equation $c_{0}x^{2}+c_{1}x+c_{2}=0$. For $R_{0}<1$ and $k+m>\beta_{1}N_{1}+b_{1}\epsilon_{1}$, we have $c_{1}>0$ and $c_{1}c_{2}>0$. So, according to the Routh–Hurwitz criteria [34], these two eigenvalues have a negative real part.

Since each eigenvalue of the characteristic equation (6) has a negative real part when $R_{0}<1$, according to the Routh–Hurwitz method [34], system (1) is locally asymptotically stable at the disease-free equilibrium point $E_{2}$ and unstable when $R_{0}>1$. The dynamical behaviors of the model at disease-free equilibrium are shown in Figure 4.

Figure 4

The plots show the dynamical behavior of the model at disease-free equilibrium

 □

3.2 Endemic equilibrium

The constant presence of a disease or an infectious agent within a given geographic area is called endemic. The endemic equilibrium state is the state where the disease cannot be totally eradicated but remains in the population. In order to find positive solutions of system (1), let $E_{2} = (S_{2},I_{2},T_{2}, R_{2}, V_{2},W_{2})$ represent any arbitrary endemic equilibrium. Setting left-hand side equal to zero and solving the equations simultaneously at steady state, we obtain

$$\begin{aligned}& S_{2}=\frac{b_{1}-kI_{2}}{\mu_{1}},\qquad T_{2}=\frac{\alpha I_{2}}{l},\qquad R_{2}=\frac{(l\eta+\gamma\alpha)I_{2}}{\mu_{1}l}, \\& V_{2}= \frac{mW_{2}}{\beta_{3}I_{2}},\qquad W_{2}= \frac{\mu_{2}\beta_{3}N_{2}I_{2}}{\beta _{3}(\delta_{2}+\mu_{2})I_{2}+\mu_{2}m}. \end{aligned}$$

Theorem 3.2

The endemic equilibrium point $E_{2}$ is locally asymptotically stable if $R_{0} > 1$, otherwise unstable.

Proof

To show these results, we linearize system (1) around $E_{2}$, which gives the following Jacobian matrix:

$$J_{2}= \begin{bmatrix} -Q &-T & 0 & 0 & 0 & -\beta_{2}S_{2}\\ Q-\mu_{1} & T-K & 0 & 0 &0& \beta_{2}S_{2}\\ 0 & \alpha& -l & 0 & 0 & 0 \\ 0 & \eta& \gamma& -\mu_{1} & 0 & 0 \\ 0 &-\beta_{3}V_{2} &0&0& -\beta_{3}I_{2}-\mu_{2} & -\epsilon_{2}b_{2}\\ 0 &\beta_{3}V_{2}& 0 & 0 & \beta_{3}I_{2}& -m \end{bmatrix} , $$

where

$$Q=\beta_{1}I_{2}+\beta_{2}W_{2}+ \mu_{1},\qquad T=\epsilon_{1}b_{1}+\beta_{1}S_{2}. $$

Two of the eigenvalues are $-\mu_{1}$ and −l. The remaining eigenvalues are the eigenvalues of the following matrix:

$$J_{2}^{*}= \begin{bmatrix} -Q &-T & 0 & -\beta_{2}S_{2}\\ Q-\mu_{1} & T-K & 0& \beta_{2}S_{2}\\ 0 &-\beta_{3}V_{2} & -\beta_{3}I_{2}-\mu_{2} & -\epsilon_{2}b_{2}\\ 0 &\beta_{3}V_{2}& \beta_{3}I_{2}& -m \end{bmatrix} . $$

We make an elementary row operation for the Jacobian matrix $J_{2}^{*}$ to obtain the following matrix:

$$J_{2}^{*}= \begin{bmatrix} -Q &-T & 0 & -\beta_{2}S_{2}\\ 0 & \frac{\mu_{1}T}{Q}-K & 0& \frac{\mu_{1}}{Q}\beta_{2}S_{2}\\ 0 &0 & -\mu_{2} & -\epsilon_{2}b_{2}-m\\ 0 &0& 0& -M \end{bmatrix} , $$

where

$$M=m+L+\frac{\beta_{3}I_{2}}{\mu_{2}}(m+\epsilon_{2}b_{2}) \quad \text{and}\quad L=m+\frac{\mu_{1}\beta_{2}\beta_{3}S_{2}V_{2}}{\mu_{1}T-KQ}. $$

$J_{2}^{*}$ is a lower triangular matrix and its eigenvalues are the elements of the main diagonal which are given by −Q, $\frac{\mu_{1}T}{Q}-K$, $-\mu_{2}$, and −M. Three of the eigenvalues have a negative real part. The second eigenvalue $\frac{\mu_{1}T}{Q}-K$ has a negative real part if and only if $\frac{\mu_{1}T}{Q}-K<0$. Using the value of Q and T, we can rewrite this equation by rearranging it as follows:

$$ -2\beta_{1}\beta_{3}K(\delta_{2}+ \mu_{2})I_{2}^{2}+\bigl[\beta_{3}( \delta_{2}+\mu _{2})\mu_{1}K(1-R_{0}) \bigr]I_{2}+\mu_{2}m\mu_{1}K(1-R_{0}). $$

(7)

All the coefficients of this equation are negative if $R_{0}>1$. Thus all the eigenvalues have negative real parts, which shows that the endemic equilibrium point $E_{2}$ is locally asymptotically stable iff $R_{0}>1$. □

4 Global stability analysis

In this section, we study the global analysis of the disease-free and endemic equilibria using the direct Lyapunov method which requires the construction of a function with specific properties. In order to do this, we derive the following results.

Theorem 4.1

When $R_{0}<1$, then the disease-free equilibrium $E_{1}$ of system (1) is globally asymptotically stable on Φ.

Proof

To show the global stability of the disease-free equilibrium $E_{1}$, we construct the following Lyapunov function, following the method used in [35]:

$$ U(t)=I+\frac{\beta_{2}b_{1}}{m\mu_{1}}W,\quad \mbox{with time derivative } U'(t)= \dot{I}+\frac{\beta_{2}b_{1}}{m\mu_{1}}\dot{W}. $$

(8)

Then U is $C^{1}$ on the interior of Φ, $E_{1}$ is the global minimum of U on Φ, and $U(t)=0$ at $E_{1}$. Putting the values from model (1), we obtain

$$ \begin{aligned}[b] U'(t)&= \epsilon_{1}b_{1}I+\beta_{1}S I + \beta_{2}S W-\alpha I -\eta I-\delta _{1}I-\mu_{1}I\\ &\quad {}+ \frac{\beta_{2}b_{1}}{m\mu_{1}}(\epsilon_{2}b_{2}W+ \beta_{3}VI - \delta_{2}W-\mu_{2}W), \\ &\leq\epsilon_{1}b_{1}I+\beta_{1}N_{1} I +\beta_{2}N_{1}W-kI\\ &\quad {}+\frac{\beta_{2}b_{1}}{m\mu_{1}}(\beta_{3}N_{2}I-mW), \quad \text{since } S\leq N_{1}, \text{ and } V\leq N_{2} \\ &=(R_{0}-1)I. \end{aligned} $$

(9)

Equation (9) shows that $U'(t)$ is negative if $R_{0}<1$. Also $U'(t)=0$ at $E_{1}$. Substituting $I=T=R=W=0$ in the equations for $S(t)$ and $V(t)$ of model (1) shows that $S(t)\rightarrow \frac{b_{1}}{\mu_{1}}$ and $V(t)\rightarrow\frac{b_{2}}{\mu_{2}}$ as $t\rightarrow\infty$. Similarly, substituting in the equations for $T(t)$ and $R(t)$ shows that $(T(t),R(t))\rightarrow(0,0)$ as $t\rightarrow\infty$. Therefore the largest compact invariant set in $\{(S_{h}, E_{h}, I_{h}, N_{h}, S_{v}, E_{v}, I_{v})\in\Phi:U'(t)= 0\}$ is the singleton disease-free equilibrium point $\{E_{f}\}$. Therefore, from LaSalle’s principle [36], the disease-free equilibrium $E_{f}$ is globally asymptotically stable in Φ. □

Theorem 4.2

For $R_{0}>1$, the endemic equilibrium $E_{2}$ is globally asymptotically stable.

Proof

For the global stability of the endemic equilibria, we construct the following Lyapunov function:

$$ Y(t)=\frac{1}{\beta_{1}S_{2}}(S-S_{2}\log S)+\frac{1}{\beta _{3}V_{2}}(V-V_{2} \log V)+\frac{1}{\beta_{1}S_{2}}I+\frac{1}{\beta_{3}V_{2}}W. $$

(10)

Taking the time derivative of W, we get

$$ \begin{aligned} Y'(t)&= \frac{1}{\beta_{1}S_{2}}(S-S_{2}) \biggl[\frac{b_{1}}{S}- \frac {\epsilon_{1}b_{1}I}{S}-\beta_{1}I-\beta_{2}W-\mu_{1} \biggr]\\ &\quad {}+ \frac{1}{\beta_{3}V_{2}}(V-V_{2}) \biggl[\frac{b_{2}}{V}- \frac{\epsilon _{2}b_{2}W}{V}-\beta_{3}I-\mu_{2} \biggr] \\ &\quad{}+\frac{1}{\beta_{1}S_{2}}[\beta_{1}S I +\beta_{2}S W-K_{1}I], \end{aligned} $$

(11)

where $K_{1}=\alpha+\delta_{1}+\mu_{1}+\eta-\epsilon_{1}b_{1}$. Let us consider

$$ \begin{gathered}\mu_{1}=\frac{b_{1}}{S_{2}}\quad \Rightarrow\quad b_{1}= \mu_{1}S_{2},\qquad \mu_{2}=\frac {b_{2}}{V_{2}}\quad \Rightarrow\quad b_{2}=\mu_{2}V_{2},\\ K_{1}=2 \beta_{1}S_{2},\quad \mbox{and}\quad m=\frac{\beta_{2}\beta_{3}V_{2}}{\beta_{1}}. \end{gathered} $$

(12)

Rearranging equation (11), we get

$$ Y'(t)=-\frac{\mu_{1}}{\beta_{1}} \biggl(\frac{S}{S_{2}}+ \frac {S_{2}}{S}-2 \biggr)-\frac{\mu_{2}}{\beta_{3}} \biggl(\frac{V}{V_{2}}+ \frac{V_{2}}{V}-2 \biggr). $$

(13)

Since

$$ \frac{S}{S_{2}}+\frac{S_{2}}{S}\geq2 \quad \text{and}\quad \frac {V}{V_{2}}+ \frac{V_{2}}{V}\geq2, $$

(14)

because the arithmetic mean is greater than or equal to the geometric mean. Thus $Y'(t)\leq0$ for all $(S,I,T,R,V,W)\in\Phi$ and the equality $(Y'(t)=0)$ holds for $E_{2}$. The proof is completed as in the proof of Theorem (4.1). □

5 Numerical simulation and graphs

We collect data from different sources and use the Runge–Kutta fourth order scheme to solve the model. Some of the parameter values are based on reality, for example, the death rate of humans by nature, corresponding to life expectancy of a 70-year-old human, is $\mu_{1} = 0.000039$ per day, and the death rate of mosquitoes is $\mu _{2} = 0.1$ per day corresponding to mosquito’s average life span of 10 days. Some of the parameter values are chosen from [25, 35]. The human’s and vector’s recruitment rates are $b_{1}=20$ and $b_{2}=100$ per day, respectively. The disease-induced death rates of humans and mosquitoes are $\delta_{1}=0.01$ and $\delta_{2}=0.21$, respectively. $\beta_{1}= 0.00001$ and $\beta_{2} = 0.0012$ are the transmission probabilities of dengue from human to human and vector to human population, respectively, $\beta_{3}=0.001$ is the transmission probability of dengue from human to vector population. Given different values to the treatment parameter $0\leq \alpha\leq1$ to check the treatment effects. The natural recovery rate is $\eta=0.01$, and the recovery rate due to treatment is $\gamma=0.4$. We suppose the values of $\epsilon_{1}$, $\epsilon_{2}$ and the initial population sizes. In rare cases the new offspring of infected parents are infected so take $\epsilon _{1}=0.001$ and the vertical transmission rate for mosquitos is $\epsilon_{2}=0.002$. For initial values, let $S(0)=100$, $I(0)=30$, $T(0)=25$, $R(0)=10$, $V(0)=600$, and $W(0)=100$. After solving we draw the results graphically and show the effect of cure rate and vertical transmission. Figure 5 shows the effect of cure rate on each population class, and Figure 6 shows the effect of vertical transmission. Figures 7 and 8 show the phase portraits of susceptible population versus infected population of human and vector populations, respectively.

Figure 5

The plots show the dynamical behaviors of population sizes with increasing cure rate

Figure 6

The plots show the dynamical behavior of population sizes with increasing vertical transmission rate

Figure 7

The plots show the phase portrait of the susceptible human population versus the infected human population

Figure 8

The plots show the phase portrait of the susceptible vector population versus the infected vector population

6 Conclusion

The spread of different infectious diseases causes very high mortality rates in a population. Vector-borne diseases are infectious diseases transmitted to humans and animals through vectors. These diseases propagate from the infected to the susceptible population in different ways. This paper formulated an epidemic model for the transmission dynamics of vector-borne diseases with both vertical and horizontal transmissions with treatment strategy. The equilibrium points and the basic reproduction of the model are found. The basic reproduction number, which is a threshold quantity, has an important role in the epidemiology of the disease. As this number increases the disease invades the population, and as it decreases the disease simply dies out. Figure 2 shows that $R_{0}$ decreases as treatment strategies increase and increases as vertical transmission increases. Figure 3 shows the threshold behavior of $R_{0}$ and the critical value $R_{0}=1$. As $R_{0}$ increases, the infected population increases with time. For $R_{0}<1$, the number of infected population decreases; for $R_{0}=1$, the infected population remains constant; and for $R_{0}>1$, the number of infected population increases. It is also shown that when $R_{0}<1$ the disease-free equilibrium is locally and globally asymptotically stable; and for $R_{0}>1$, the positive endemic equilibrium is locally and globally asymptotically stable.

Numerical simulations are carried out graphically to show the dynamical behavior of the diseases. Figure 5 shows the effect of cure rate on the transmission dynamics of the disease. As treatment strategy increases, the susceptible population and the recovered human population increase while the infected population decreases. Figure 6 shows the effect of vertical transmission. As vertical transmission increases, the susceptible population decreases and the infected population increases. Finally, Figures 7 and 8 show the phase portraits of the susceptible populations versus the infected populations which move towards the stable points.