Sensitivity Analysis of Dengue Model with Saturated Incidence Rate

M. M. Ojo¹, B. Gbadamosi², Adebimpe Olukayode³, Ogundokun R. Oluwaseun²

¹Department of Ecology and Evolutionary Biology, University of Kansas, Lawrence, USA.
²Department of Computer Science Landmark University, Omu-Aran, Nigeria.
³Department of Physical Sciences, Landmark University, Omu-Aran, Nigeria.

DOI: 10.4236/oalib.1104413                

Abstract

Dengue is a flavivirus, transmitted to human through the bites of infected Aedes aegypti and A. albopictus mosquitoes. In this paper, we analyze a new system of ordinary differential equations which incorporates saturated incidence function, vector biting rate and control measures at both the aquatic and adult stages of the vector (mosquito). The stability of the system is analysed for the dengue-free equilibrium via the threshold parameter (reproduction number) which was obtained using the Next generation matrix techniques. Routh Hurwitz criterion along together with Descartes’ rule of signs change established the local asymptotically stability of the model whenever R₀＜1and unstable otherwise. Furthermore, the sensitivity analysis was carried out and the numerical simulation reveals that increasing the proportion of human antibody and putting into place a control strategy that minimize the vector biting rate are enough to reduce the infection of the disease in the population to its barest minimum.

Keywords

Aedes aegypti, A. albopictus, Dengue Fever, Reproduction Number, Control Measures, Sensitivity Analysis, Aquatic Stage, Adult Stage

1. Introduction

Over the years, mathematical models and computer simulations have been known to be useful experimental tools which are used in building and examing theories, evaluating quantitative speculations, giving answers to particular questions and determining sensitivities to changes in parameter values. Understanding the epidemiology of emerging and re-emerging infectious diseases in a population produces a healthy environment for living. Mathematical models are used in likening, designing, implementing, evaluating and optimizing several detection, prevention and control plans.

Dengue fever is one of the infectious diseases that have continued to be a subject of major concern to the public health. It is known to be a mosquito-borne viral infection which is endemic in more than a hundred countries in the world [1] [2] [3] [4] , usually in a tropic and sub-tropical regions of the world [5] . In recent years, dengue transmission preponderates in urban and semi-urban areas [1] [5] where a figure of 50 to 390 million people worldwide in a year are infected which leads to half a million hospitalizations [6] [7] [8] with an approximate of 25,000 deaths [4] [6] [9] .

The dengue disease has been well known clinically for over 2 centuries, but the etiology of the disease remains unknown until year 1944 [10] [11] . It was first recognized in the Philippines in 1953 and Thailand in 1955 [10] [11] [12] . The threat of the outbreak now exists in Europe which its first local transmission was reported in France and Croatia in 2010, while cases have occurred in Florida (USA) and Yunnan (province of China) in 2013 [1] [5] .

Dengue Hemorrhagic fever being an infectious tropical disease is caused by an infective agent called dengue virus, of the family flaviviridae which has four distinguished serotypes denoted by I, II, III and IV [12] . The virus is transmitted to humans by bites of Aedes mosquitoes [Aedes aegypti and A. albopictus are the primary transistors]. The infection remains in mosquito till death [8] .

Dengue infection causes a range of illness in humans, from clinically in apparent, to severe and fatal hemorrhagic disease [11] [12] . The incubation period; which is the time between infection and appearance of the symptoms in the body is from 3 - 14 days, but often times it ranges from 4 - 7days [4] [6] [13] , and is generally observed clearly in older children and adults [11] . Dengue fever is characterized by sudden onset of fever, frontal headache, nausea, vomiting and some other symptoms.

The use of mathematics in explaining the epidemiology of dengue fever has been extensively studied by many researchers over years. Notable among these studies are [2] [6] [10] [13] [14] [15] [16] . In this study, since dengue fever is spread between two-interacting populations (human-vector), we design and analyses a mathematical compartmental model that considers the human population and the vector population (mosquito). We extended the earlier model [6] by incorporating a “Standard force of infection” with the proportion of an antibodies produced by human in response to the incidence of infection caused by mosquito and vice-versa. Also, an extension of the work is to consider some control effects or precautionary measures of the vector in the absence of vaccination. These measures includes: Larvicides for the Aquatic stage of the vector which prevents the vector from breeding, Naled and EPA-registered insects’ repellants to prevent getting bitten against the adult stage of the vectors.

2. Formulation of the Model

The formulation of dengue model requires the interaction between two-interacting populations (human-vector). The total human population at continuous-time t denoted by $N_h\left(t\right)$ is subdivided into six compartments namely: susceptible humans ( $S_h$ ), exposed humans ( $E_h$ ), infectious humans ( $I_h$ ), migrated population ( $M_h$ ), treatment class ( $T_h$ ), recovered humans ( $R_h$ ). Hence, the total human population $N_h\left(t\right)$ is given by

$N_h\left(t\right)=\left(S_h\right)+\left(E_h\right)+\left(I_h\right)+\left(M_h\right)+\left(T_h\right)+\left(R_h\right)$ (1)

Similarly, the total vector population at continuous-time t denoted by $N_v\left(t\right)$ is subdivided into four compartments namely: aquatic class ( $A_v$ ), susceptible mosquitoes ( $S_v$ ), exposed mosquitoes ( $E_v$ ), infectious mosquitoes ( $I_v$ ). Hence, the total vector population $N_v\left(t\right)$ is given by

$N_v\left(t\right)=\left(A_v\right)+\left(S_v\right)+\left(E_v\right)+\left(I_v\right)$ (2)

The dynamics of the dengue considered here is formulated and studied under the following assumption:

1) the model assumes a homogeneous mixing of the human and vector (mosquito) populations, so that each mosquito bite has equal chance of transmitting the virus to susceptible in the population (or acquiring infection from an infected human);

2) considering saturated incidence rate (Non-linear incidence) which incorporate the production of antibodies in response to parasites causing Dengue in both human and vector population $\left({\upsilon }_h,{\upsilon }_v\right)$ respectively.

3) the model consider the vector-aquatic class so as to investigate on the effect of the control strategies such as Larvicides at the aquatic stage;

4) that the infectious mosquitoes remain infectious until death;

5) there is loss of immunity for the recovered human population;

6) incorporating the controlling rate parameters which will monitor the effects of control strategies at the aquatic stage ( $A_v$ ) and adult stages ( $S_v,E_v,I_v$ ).

In summary, following the assumptions above the transmission dynamics of dengue in a population is given by the following ten compartmental system of non-linear differential equation below:

$\begin{array}{l}{\stackrel{\dot{}}{S}}_h\left(t\right)={\pi }_h-\frac{b{\beta }_{hv}{S}_h\left(t\right)I_v\left(t\right)}{1+{\upsilon }_h{I}_v\left(t\right)}-{\mu }_h{S}_h\left(t\right)+\omega R_h\left(t\right)\\ {\stackrel{\dot{}}{E}}_h\left(t\right)=\frac{b{\beta }_{hv}{S}_h\left(t\right)I_v\left(t\right)}{1+{\upsilon }_h{I}_v\left(t\right)}+{\mu }_1{M}_h\left(t\right)-\left({\mu }_h+{\sigma }_h\right)E_h\left(t\right)\\ {\stackrel{\dot{}}{I}}_h\left(t\right)={\sigma }_h{E}_h\left(t\right)+{\mu }_2{M}_h\left(t\right)-\left({\mu }_h+{\tau }_h+{\delta }_h\right)I_h\left(t\right)\\ {\stackrel{\dot{}}{M}}_h\left(t\right)={\pi }_{mh}-\left({\mu }_1+{\mu }_2+{\mu }_h\right)M_h\left(t\right)\\ {\stackrel{\dot{}}{T}}_h\left(t\right)={\tau }_h{I}_h\left(t\right)-\left({\mu }_h+{\gamma }_1\right)T_h\left(t\right)\\ {\stackrel{\dot{}}{R}}_h\left(t\right)={\gamma }_1{T}_h\left(t\right)-{\mu }_h{R}_h\left(t\right)-\omega R_h\left(t\right)\end{array}\}$

$\begin{array}{l}{\stackrel{\dot{}}{A}}_v\left(t\right)={\pi }_v-\left({\gamma }_m+{\mu }_v+C_a\right)A_v\left(t\right)\\ {\stackrel{\dot{}}{S}}_v\left(t\right)={\gamma }_m{A}_v\left(t\right)-\frac{b{\beta }_{vh}{S}_v\left(t\right)I_h\left(t\right)}{1+{\upsilon }_v{I}_h\left(t\right)}-\left({\mu }_v+C_m\right)S_v\left(t\right)\\ {\stackrel{\dot{}}{E}}_v\left(t\right)=\frac{b{\beta }_{vh}{S}_v\left(t\right)I_h\left(t\right)}{1+{\upsilon }_v{I}_h\left(t\right)}-\left({\theta }_c+{\sigma }_v+{\mu }_v+C_m\right)E_v\left(t\right)\\ {\stackrel{\dot{}}{I}}_v\left(t\right)=\left({\theta }_c+{\sigma }_v\right)E_v\left(t\right)-\left({\delta }_v+{\mu }_v+C_m\right)I_v\left(t\right)\end{array}\}$ (3)

where a dot is representing differentiation with respect to time.

Figure 1 shows the schematic illustration of the dengue model.

Figure 1. The Schematic illustration of dengue model.

Table 1 shows the description of the parameters of the model.

2.1. Basic Properties of the Model

It is important to explore the basic dynamical feature of the model. For the model (3) formulated above to be epidemiologically meaningful, it is very important to prove that all the states variables non-negative for all time (t). In other words, the solution of the model (3) with positive initial values of data will remain positive at all time $t\ge 0$ .

Positivity and Boundedness of Solutions

Since model (3) describe interaction between human and vector population, it is important to state that all the parameters and variables involved are non-negative with respect to time. The dengue model (3) will be consider in the biologically-feasible region $\Im ={\Im }_h\times {\Im }_v\subset {\Re }_{+}^6\times {\Re }_{+}^4$ with

${\Im }_h=\left\{S_h,E_h,I_h,M_h,T_h,R_h\in {\Re }_{+}^6:N_h\le \frac{{\pi }_h}{{\mu }_h}\right\}$ (4)

and

Table 1. Description of the parameters of dengue model (3).

${\Im }_v=\left\{A_v,S_v,E_v,I_v\in {\Re }_{+}^4:N_v\le \frac{{\pi }_v}{{\mu }_v}\right\}$ (5)

It can be shown that the set $\Im$ is a positively invariant set and global attractor of this system. This implies, any phase trajectory initiated anywhere in the non-negative region ${\Re }_{+}^{10}$ of the phase space eventually enters the feasible region $\Im$ and remains in $\Im$ thereafter.

Lemma 2.1 The region
$\Im =\left\{S_h,E_h,I_h,M_h,T_h,R_h,A_v,S_v,E_v,I_v\in {\Re }_{+}^{10}:N_h\le \frac{{\pi }_h}{{\mu }_h},N_h\le \frac{{\pi }_h}{{\mu }_h}\right\}$ is positively-invariant for the model (3).

Proof: The rate of change of the human total population is given by

$\frac{\text{d}{N}_h}{\text{d}t}={\pi }_h-{\mu }_h{N}_h-{\delta }_h{I}_h+{\pi }_{mh}$ (6)

and,

$\frac{\text{d}{N}_v}{\text{d}t}={\pi }_v-{\mu }_v{N}_v-C_a{A}_v-C_m\left(S_v+E_v+I_v\right)-{\delta }_v{I}_v$ (7)

where $N_h=S_h+E_h+I_h+M_h+T_h+R_h$ and $N_v=A_v+S_v+E_v+I_v$ .

Since $\frac{\text{d}{N}_h}{\text{d}t}\le {\pi }_h-{\mu }_h{N}_h+{\pi }_{mh}$ and $\frac{\text{d}{N}_v}{\text{d}t}\le {\pi }_v-{\mu }_v{N}_v$ for special case ${\delta }_h={\delta }_v=C_a=C_m=0$ , it follows that whenever $N_h\left(t\right)>\frac{{\pi }_h+{\pi }_{mh}}{{\mu }_h}$ and $N_v\left(t\right)>\frac{{\pi }_v}{{\mu }_v}$ , then $\frac{\text{d}{N}_h}{\text{d}t}<0$ and $\frac{\text{d}{N}_v}{\text{d}t}<0$ respectively.

Thus, since it follows from the right hand side of Equations ((6) and (7)) that $\frac{\text{d}{N}_h}{\text{d}t}$ is bounded by ${\pi }_h-{\mu }_h+{\pi }_{mh}$ and $\frac{\text{d}{N}_v}{\text{d}t}$ is bounded by ${\pi }_v-{\mu }_v$ , the standard comparison theorem [17] [18] can be used to show that:

$N_h\left(t\right)\le N_h\left(0\right){\exp }^{-{\mu }_{h}t}+\frac{{\pi }_h+{\pi }_{mh}}{{\mu }_h}\left[1-{\exp }^{-{\mu }_{h}t}\right]$ if $N_h\left(0\right)\le \frac{{\pi }_h+{\pi }_{mh}}{{\mu }_h}$ (8)

and

$N_v\left(t\right)\le N_v\left(0\right){\exp }^{-{\mu }_{v}t}+\frac{{\pi }_v}{{\mu }_v}\left[1-{\exp }^{-{\mu }_{v}t}\right]$ if $N_v\left(0\right)\le \frac{{\pi }_v}{{\mu }_v}$ (9)

Thus, $\Im$ is positively invariant under the flow described by (3) so that no solution path leaves through any boundary of $\Im$ . Hence, in the region $\Im$ , the model (3) is recognized to be mathematically and epidemiologically well-posed. Thus, it is sufficient to consider the dynamics of the model in the domain $\Im$ .

2.2. Asymptotic Stability of the Disease Free Equilibrium

The dengue disease-free equilibrium is a point at which the population is free from dengue fever. The disease-free equilibrium of the model (3) exists and is given by

$\begin{array}{l}{\epsilon }_0=\left\{S_h^{*},E_h^{*},I_h^{*},M_h^{*},T_h^{*},R_h^{*},A_v^{*},S_v^{*},E_v^{*},I_v^{*}\right\}\\ =\left\{\frac{{\pi }_h}{{\mu }_h},0,0,0,0,0,\frac{{\pi }_v}{\left({\gamma }_m+{\mu }_v+C_a\right)},\frac{{\gamma }_m{\pi }_v}{\left({\gamma }_m+{\mu }_v+C_a\right)\left({\mu }_v+C_m\right)},0,0\right\}\end{array}$ (10)

2.3. Basic Reproduction Number

A very important concern about an infectious disease is its ability to invade a population [15] [17] . The threshold condition known as the basic reproduction number (usually written as $R_0$ ) is used in determining whether the disease will persist in the population or dies out as time increases; if $R_0<1$ , then the disease free equilibrium (DFE) will be locally asymptotically stable and the disease cannot invade the population while when $R_0>1$ , then the DFE is unstable and invasion is possible which could leads to an endemic equilibrium state [13] [19] .

The linear stability of ${\epsilon }_0$ is studied using the next generation operator technique [20] [21] [22] on the system (3). Using the notations in [20] [22] , it follows that the matrices F and V, for the new infection terms (transmission) and the remaining transfer terms (transition) are respectively given by:

$F=\left(\begin{array}{cccccc}0& 0& 0& 0& 0& b{S}_h^{*}{\beta }_{hv}\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\\ 0& b{S}_v^{*}{\beta }_{vh}& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0\end{array}\right)$

and

$V=\left(\begin{array}{cccccc}{\kappa }_1& 0& -{\mu }_1& 0& 0& b{S}_h^{*}{\beta }_{hv}\\ -{\sigma }_h& {\kappa }_2& -{\mu }_2& 0& 0& 0\\ 0& 0& {\kappa }_3& 0& 0& 0\\ 0& -{\tau }_h& 0& {\kappa }_4& 0& 0\\ 0& 0& 0& 0& {\kappa }_8& 0\\ 0& 0& 0& 0& -{\kappa }_9& {\kappa }_{10}\end{array}\right)$ (11)

where $\begin{array}{l}{\kappa }_1=\left({\mu }_h+{\sigma }_h\right),{\kappa }_2=\left({\mu }_h+{\tau }_h+{\delta }_h\right),{\kappa }_3=\left({\mu }_1+{\mu }_2+{\mu }_h\right),{\kappa }_4=\left({\mu }_h+{\gamma }_1\right),\\ {\kappa }_8=\left({\theta }_c+{\sigma }_v+{\mu }_v+C_m\right),{\kappa }_9=\left({\theta }_c+{\sigma }_v\right)\text{and}{\kappa }_{10}=\left({\delta }_v+{\mu }_v+C_m\right).\end{array}$

It follows that the basic reproduction number of the model (3) denoted by $R_0$ , is given by $R_0=\rho \left(F{V}^{-1}\right)$ where $\rho$ is the spectral radius (maximum eigenvalues) [17] . Hence,

$R_0=\sqrt{\frac{b^2{S}_h^{*}{\beta }_{hv}{S}_v^{*}{\beta }_{vh}{\sigma }_h{\kappa }_9}{{\kappa }_1{\kappa }_2{\kappa }_8{\kappa }_{10}}}$

So that;

$R_0=\sqrt{\frac{b^2{S}_h^{*}{\beta }_{hv}{S}_v^{*}{\beta }_{vh}{\sigma }_h\left({\sigma }_v+{\theta }_c\right)}{\left({\mu }_h+{\sigma }_h\right)\left({\mu }_h+{\tau }_h+{\delta }_h\right)\left({\theta }_c+{\sigma }_v+{\mu }_v+C_m\right)\left({\delta }_v+{\mu }_v+C_m\right)}}$

Hence

$R_0=\sqrt{R_h{R}_v}$ (12)

where

$R_h=\frac{b{S}_h^{*}{\beta }_{hv}{\sigma }_h}{\left({\mu }_h+{\sigma }_h\right)\left({\mu }_h+{\tau }_h+{\delta }_h\right)}$ (13)

and,

$R_v=\frac{b{S}_v^{*}{\beta }_{vh}\left({\sigma }_v+{\theta }_c\right)}{\left({\theta }_c+{\sigma }_v+{\mu }_v+C_m\right)\left({\delta }_v+{\mu }_v+C_m\right)}$ (14)

In (13), $R_h=\frac{b{S}_h^{*}{\beta }_{hv}{\sigma }_h}{\left({\mu }_h+{\sigma }_h\right)\left({\mu }_h+{\tau }_h+{\delta }_h\right)}$ describes the number of humans that just one infectious vector infects over its expected infection period in a completely susceptible humans population. Also, $\frac{{\sigma }_h}{\left({\mu }_h+{\sigma }_h\right)}$ is the probability that a human will survive the exposed state to become infectious, while $\frac{1}{\left({\mu }_h+{\tau }_h+{\delta }_h\right)}$ is the average duration of the infectious period of human.

In (14), $R_v=\frac{b{S}_v^{*}{\beta }_{vh}\left({\sigma }_v+{\theta }_c\right)}{\left({\theta }_c+{\sigma }_v+{\mu }_v+C_m\right)\left({\delta }_v+{\mu }_v+C_m\right)}$ signifies the number of vectors infected by one infectious human during the period of infectiousness in a completely susceptible vector population. Also, $\frac{\left({\sigma }_v+{\theta }_c\right)}{\left({\theta }_c+{\sigma }_v+{\mu }_v+C_m\right)}$ is the probability that a vector will survive the exposed state to become infectious while $\frac{1}{\left({\delta }_v+{\mu }_v+C_m\right)}$ is the average duration of the infectious period of vector.

2.4. Stability of the Disease-Free Equilibrium Point

We will use the basic reproduction obtained for the model (3) to analyze the stability of the equilibrium point in the following result.

Theorem 2.1. The disease-free state ${\epsilon }_0$ of the dengue model considered, is locally asymptotically stable if $R_0<1$ and unstable if $R_0>1$ .

Proof: The Jacobian matrix of the system (3) evaluated at the disease-free equilibrium point ${\epsilon }_0$ , is obtained as

$J\left({\epsilon }_0\right)=\left(\begin{array}{cccccccccc}{J}_{11}& 0& 0& 0& 0& J_{16}& 0& 0& 0& J_{1,10}\\ 0& J_{22}& 0& J_{24}& 0& 0& 0& 0& 0& J_{2,10}\\ 0& J_{32}& J_{33}& J_{34}& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& J_{44}& 0& 0& 0& 0& 0& 0\\ 0& 0& J_{53}& 0& J_{55}& 0& 0& 0& 0& 0\\ 0& 0& 0& 0& J_{65}& J_{66}& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& J_{77}& 0& 0& 0\\ 0& 0& J_{83}& 0& 0& 0& J_{87}& J_{88}& 0& 0\\ 0& 0& J_{93}& 0& 0& 0& 0& 0& J_{99}& 0\\ 0& 0& 0& 0& 0& 0& 0& 0& J_{10,9}& J_{10,10}\end{array}\right)$ (15)

where $\begin{array}{l}{J}_{11}=-{\mu }_h,J_{22}=-{\kappa }_1,J_{32}={\sigma }_h,J_{33}=-{\kappa }_2,J_{53}={\tau }_h,J_{2,10}=b{\beta }_{hv}{S}_h^{*},\\ J_{24}={\mu }_1,J_{34}={\mu }_2,J_{44}=-{\kappa }_3,J_{83}=-b{\beta }_{vh}{S}_v^{*},J_{93}=b{\beta }_{vh}{S}_v^{*},\\ J_{55}=-{\kappa }_4,J_{64}={\gamma }_1,J_{17}=\omega ,J_{66}=-{\kappa }_5,J_{77}=-{\kappa }_6,J_{87}={\gamma }_m,\\ J_{88}=-{\kappa }_7,J_{99}=-{\kappa }_8,J_{10,9}={\kappa }_9,J_{10,10}=-{\kappa }_{10},J_{1,10}=-b{\beta }_{hv}{S}_h^{*}\end{array}$

From (15), it is sufficient for us to show that all the eigenvalues of $J\left({\epsilon }_0\right)$ are negative. The first and eight columns contains only the diagonal terms which form the two negative eigenvalues, $-{\mu }_h$ and $-\left({\mu }_h+C_m\right)$ , so that the other eight eigenvalues can be obtained from the sub-matrix $J_1\left({\epsilon }_0\right)$ , formed by excluding the first and eight rows and columns of $J\left({\epsilon }_0\right)$ . Hence, $J_1\left({\epsilon }_0\right)$ is written as

$J_1\left({\epsilon }_0\right)=\left(\begin{array}{cccccccc}{J}_{22}& 0& J_{24}& 0& 0& 0& 0& J_{2,10}\\ J_{32}& J_{33}& J_{34}& 0& 0& 0& 0& 0\\ 0& 0& J_{44}& 0& 0& 0& 0& 0\\ 0& J_{53}& 0& J_{55}& 0& 0& 0& 0\\ 0& 0& 0& J_{65}& J_{66}& 0& 0& 0\\ 0& 0& 0& 0& 0& J_{77}& 0& 0\\ 0& J_{93}& 0& 0& 0& 0& J_{99}& 0\\ 0& 0& 0& 0& 0& 0& J_{10,9}& J_{10,10}\end{array}\right)$ (16)

In the same way, the fifth and sixth column of $J_1\left({\epsilon }_0\right)$ contains only the diagonal term which forms a negative eigenvalues $-\left({\mu }_h+\omega \right)$ and $-\left({\gamma }_m+{\mu }_v+C_a\right)$ . The remaining six eigenvalues can be obtained from the sub-matrix $J_2\left({\epsilon }_0\right)$ written as

$J_2\left({\epsilon }_0\right)=\left(\begin{array}{cccccc}{J}_{22}& 0& J_{44}& 0& 0& J_{2,10}\\ J_{32}& J_{33}& J_{34}& 0& 0& 0\\ 0& 0& J_{44}& 0& 0& 0\\ 0& J_{53}& 0& J_{55}& 0& 0\\ 0& J_{93}& 0& 0& J_{99}& 0\\ 0& 0& 0& 0& J_{10,9}& J_{10,10}\end{array}\right)$ (17)

Using the same approach, the fourth column and third column of $J_2\left({\epsilon }_0\right)$ contains only the diagonal term which forms a negative eigenvalue, $-\left({\mu }_h+{\gamma }_1\right)$ and $-\left({\mu }_1+{\mu }_2+{\mu }_h\right)$ . The remaining four eigenvalues can now be obtained by the characteristics equation of the sub-matrix $J_2\left({\epsilon }_0\right)$ written as

$J_3\left({\epsilon }_0\right)=\left(\begin{array}{cccc}-\left({\mu }_h+{\sigma }_h\right)& 0& 0& b{\beta }_{hv}{S}_h^{*}\\ {\sigma }_h& -\left({\mu }_h+{\tau }_h+{\delta }_h\right)& 0& 0\\ 0& b{\beta }_{vh}{S}_v^{*}& -\left({\theta }_c+{\sigma }_v+{\mu }_v+C_m\right)& 0\\ 0& 0& \left({\theta }_c+{\sigma }_v\right)& -\left({\mu }_v+{\delta }_v+C_m\right)\end{array}\right)$ (18)

Hence, the eigenvalues of the matrix $J_3\left({\epsilon }_0\right)$ are the roots of the characteristics equation;

$\begin{array}{l}\left(\lambda +{\mu }_h+{\sigma }_h\right)\left(\lambda +{\mu }_h+{\tau }_h+{\delta }_h\right)\left(\lambda +{\theta }_c+{\sigma }_v+{\mu }_v+C_m\right)\left(\lambda +{\mu }_v+{\delta }_v+C_m\right)\\ -b^2{\sigma }_v{\beta }_{hv}{S}_h^{*}{\beta }_{vh}{S}_v^{*}\left({\theta }_c+{\sigma }_v\right)=0\end{array}$

Let ${\kappa }_1=\left({\mu }_h+{\sigma }_h\right)$ , ${\kappa }_2=\left({\mu }_h+{\tau }_h+{\delta }_h\right)$ , ${\kappa }_8=\left({\theta }_c+{\sigma }_v+{\mu }_v+C_m\right)$ , ${\kappa }_{10}=\left({\mu }_v+{\delta }_v+C_m\right)$ and ${\kappa }_9=\left({\theta }_c+{\sigma }_v\right)$ then the equation above becomes:

$B_4{\lambda }^4+B_3{\lambda }^3+B_2{\lambda }^2+B_1\lambda +B_0=0$ (19)

where

$\begin{array}{l}{B}_4=1\\ B_3={\kappa }_1+{\kappa }_2+{\kappa }_8+{\kappa }_{10}\\ B_2=\left({\kappa }_1+{\kappa }_2\right)\left({\kappa }_8+{\kappa }_{10}\right)+{\kappa }_1{\kappa }_2+{\kappa }_8{\kappa }_{10}\\ B_1=\left({\kappa }_1+{\kappa }_2\right)\left({\kappa }_8{\kappa }_{10}\right)+\left({\kappa }_8+{\kappa }_{10}\right)\left({\kappa }_1{\kappa }_2\right)\\ B_0={\kappa }_1{\kappa }_2{\kappa }_8{\kappa }_{10}-b^2{\sigma }_v{\beta }_{hv}{S}_h^{*}{\beta }_{vh}{S}_v^{*}{\kappa }_9\end{array}\}$ (20)

Further perturbation on $B_0$ in terms of reproduction number, $R_0$ yields

$B_0={\kappa }_1{\kappa }_2{\kappa }_8{\kappa }_{10}\left(1-R_0^2\right)$ (21)

We employ the Routh-Hurwitz criterion, [23] [24] [25] , which states that all roots of the polynomial (19) have negative real parts if and only if the coefficients of $B_i$ are positive and matrices $H_i>0$ , for $i=0,1,2,3,4$ from (20), it is obvious that $B_1>0$ , $B_2>0$ , $B_3>0$ , $B_4>0$ , since all $B_i$ ’s are positive. Moreover, if $R_0<0$ , it then follows from (21) that $B_0>0$ . Also, the Hurwitz matrices for the polynomial (19) are found to be positive. That is, $H_1=B_3>0$,

$H_2=\left|\begin{array}{cc}{B}_3& B_4\\ B_1& B_2\end{array}\right|>0$ , $H_3=\left|\begin{array}{ccc}{B}_3& B_4& B_0\\ B_1& B_2& B_3\\ 0& B_0& B_1\end{array}\right|>0$ and $H_4=\left|\begin{array}{cccc}{B}_3& B_4& 0& 0\\ B_1& B_2& B_3& B_4\\ 0& B_0& B_1& B_2\\ 0& 0& 0& B_0\end{array}\right|>0$ .

Hence, all the eigenvalues of the Jacobian matrix $J\left({\epsilon }_0\right)$have negative real parts whenever $R_0<1$, and the disease-free equilibrium point is said to be locally asymptotically stable. However, if $R_0>1$ we deduce that $B_0<0$ and by Descartes’ rule of signs [23] [25] , there exist exactly one sign change in the sequence $B_4,B_3,B_2,B_1,B_0$ of coefficients of the polynomial (19). So, there is one eigenvalue with non-negative real part and hence the disease-free equilibrium point is said to be Unstable which proclaims an existence of an endemic state of equilibria.

3. Sensitivity Analysis

The necessity of conducting an investigation which tells how sensitive the threshold quantity basic reproduction number is with respect to its parameters cannot be over-emphasized. This will help us to know the parameters having the most significant impact on the outcome of the numerical simulations of the model. Sensitivity analysis informs us the importance of each parameter to the disease transmission, and this will help the public health authorities to place priority on a well posed intervention strategy for preventing and controlling the spread of the disease in the population.

Following [19] [26] the normalized forward sensitivity index also called elasticity was employed. The normalized forward sensitivity index of the reproduction number $R_0$ with respect to the parameter k is given by:

$Z_k^{R_0}=\frac{\partial R_0}{\partial k}\times \frac{k}{R_0}$ (22)

Using the elasticity formula presented in (22) and the baseline parameters in Table 2 below to obtain the numerical values for the elasticity’s. The elasticity index was estimated with respect to each parameter, and the results are displayed in Table 3.

Table 2. Parameter values for the dengue model (3).

Table 3. Sensitivity indexes of the dengue model’s parameters with respect to R₀.

Table 2 shows the baseline value of parameters used for the numerical simulation.

Table 3 present the sensitivity index of the dengue model with respect to $R_0$ .

The result of the Sensitivity Index presented above in Table 3 shows that the vector biting rate has the highest positive sensitivity index with value (1), which indicates that increase (decrease) by ?% in the biting rate b will be followed by an immediate increase (decrease) in reproduction number ( $R_0$ ) by ?%. Similarly, increase in the transmission rate from host to vector, transmission rate from vector to human and recruitment rate of humans and vector $\left({\beta }_{hv},{\beta }_{vh},{\pi }_h,{\pi }_v\right)$ respectively will increase the reproduction number. The immediate conclusion is that, at the disease free equilibrium the most effective control strategy is the vector control (Aquatic and Adult stage).

The vector death rate has the highest negative sensitivity index with the value (−0.93) which indicates that any increase by ?% in ${\mu }_v$ will decrease the reproduction number ( $R_0$ ) by ?% and vice-versa. This Furtherly buttress our claim that vector control is the most control strategy to maintain a disease free environment. Also, the sensitivity index of the treatment rate of dengue infectious individuals ( ${\tau }_h$ ) is −0.76, this calls for an increase in treatment rate in other to reduce the reproduction number.

In summary, the dengue sensitivity analysis identifies the most important parameters driving the transmission mechanism of the disease. The results suggest that a control strategy that reduces the vector biting rate, transmission rate from host to vector, transmission rate from vector to human, recruitment rate of humans and vector $\left(b,{\beta }_{hv},{\beta }_{vh},{\pi }_h,{\pi }_v\right)$ respectively and control strategy that increases the treatment rate and the death rate of the vector, will effectively curtail the spread of Dengue Virus in the population.

4. Numerical Simulations and Discussion

In this section, we study numerically the behaviour of the Dengue model (3), using the Rungi-Kutta method with the values presented in Table 2 (unless otherwise stated). The numerical simulations are conducted using computational machine (Maple 18), and the results are presented below.

The behaviour of Susceptible human as the antibody ${\nu }_h$ increases in proportion over a period of time is presented in Figure 2(a) above. It was observed that, the Susceptible human population drastically dropped due to the infection by infectious vectors and increase in the proportion of the human antibodies reduces the decrease in the susceptible population. Also, Figure 2(b) & Figure 2(c) presents the order of magnitude of the exposed human populations and infectious human populations. It was observed that the exposed human populations and infectious human populations decrease as the proportion of human antibody increases.

The behaviour of Susceptible human as the biting rate b decreases over a period of time was presented in Figure 3(a). It was observed that, the Susceptible

(a) (b)(c)

Figure 2. Simulations of Dengue model (3) showing the varying effect of the proportion of antibody ${\nu }_h$ on the human population $\left(S_h,E_h,I_h\right)$, as a function of time when $R_0<1$. Parameter values used are presented in Table 2 with $C_a=C_m=0$ .

(a) (b)

Figure 3. Simulations of Dengue model (3) showing the effect of vector biting rate (b) on the human population $\left(S_h,I_h\right)$ , as a function of time when $R_0<1$ . Parameter values used are presented in Table 2 with $C_a=C_m=0$ .

human population increases due to the decrease in the vector biting rate. Similarly Figure 3(b) presents the effect of vector biting on the infectious human population. As the vector biting rate decreases, the infectious human populations’ also decreases. Thus, minimizing the vector biting rate will reduce the infection rate of the diseases in the population.

The behaviour of the vector population as the Aquatic stage control $C_a$ increases over a period of time was presented in Figure 4(a). It was observed that, the vector population in the Aquatic stage dropped due to increase in the Aquatic stage control. Similarly in Figures 4(b)-(d), it was observed that the order of magnitude of the susceptible vector populations, exposed vector populations and infectious vector populations decrease as the Aquatic stage control increases. Thus, early control of the vectors in their Aquatic stage will reduce the population of the vector to its barest minimum, and this will reduce the probability of

(a) (b)(c) (d)

Figure 4. Simulations of Dengue model (3) showing the effect of Aquatic stage control measure $C_a$ on the vector population $\left(A_v,S_v,E_v,I_h\right)$ , as a function of time when $R_0<1$ . Parameter values used are presented in Table 2 with $C_a\ne 0$ and $C_m=0$ .

humans getting infected by vector infections due to their biting rate.

The behaviour of the vector population $\left(E_v,I_v\right)$ as the Adult stage control $C_m$ increases over a period of time was presented above in Figure 5(a). It was observed that, the exposed vector population reduces due to increase in the Adult stage control. Similarly Figure 5(b) shows that the order of magnitude of the infectious vector populations’ decrease as the Adult stage control increases. Thus, increasing the control of vectors in their Adult stage will reduce the population of the vector to its barest minimum.

The behaviour of the infectious human population ( $I_h$ ) was presented in Figure 6. The vector biting rate decreases with increase in the proportion of

(a) (b)

Figure 5. Simulations of Dengue model (3) showing the effect of Adult stage control measure $C_m$ on the vector population $\left(S_v,E_v,I_h\right)$ , as a function of time when $R_0<1$ . Parameter values used are presented in Table 2 with $C_a=0$ and $C_m\ne 0$ .

Figure 6. Simulation of Dengue model (3) showing the behaviour of the infected human population ( $I_h$ ) with precautionary measures. Parameter values used are presented in Table 2 with $b=0.12,{\nu }_h=0.00$ and $b=0.0012,{\nu }_h=1.00$ .

human antibody over a period of time. It was also observed that, drastic decrease in the vector biting rate and increase in the proportion of human antibody will reduce the population of the infectious human as time increases. Thus, reducing the vector biting rate and increasing the proportion of human antibody are enough to reduce the infection of the disease in the population to its barest minimum.

5. Conclusion

A new system of ordinary differential equations which incorporates saturated incidence function, vector biting rate and control measures at the aquatic and adult stages of the vector (mosquito) was presented and analysed. The disease-free equilibrium represented by ( ${\epsilon }_0$ ) was shown to be locally asymptotically stable whenever the reproduction number ( $R_0$ ) is less than unity. The sensitivity analysis shows that the dominant parameters are the vector biting rate, transmission rate from host to vector, transmission rate from vector to human, recruitment rate of humans and vector, treatment rate and the death rate of the vector $\left(b,{\beta }_{hv},{\beta }_{vh},{\pi }_h,{\pi }_v,{\tau }_h,{\mu }_v\right)$ respectively. Furtherly, the numerical simulation shows that increasing the proportion of human antibody and putting into place a control strategy that minimize the vector biting rate are enough to reduce the infection of the disease in the population to its barest minimum.

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	Bowman, C., Gumel, A.B., van den Driessche, P., Wu, J. and Zhu, H. (2005) Mathematical Model for Assessing Control Strategies against West Nile Virus. Bulletin of 682 Mathematical Biology, 67, 1107-1133. https://doi.org/10.1016/j.bulm.2005.01.002
[2]	Garba, S.M., Gumel, A.B. and Abu Bakar, M.R. (2008) Backward Bifurcations in Dengue Transmission Dynamics. Mathematical Biosciencees, 201, 11-25. https://doi.org/10.1016/j.mbs.2008.05.002
[3]	Ranjit, S. and Kissoon, N. (2011) Dengue Hemorrhagic Fever and Shock Syndromes. Pediatric Critical Care Medicine, 12, Article ID: 90100. https://doi.org/10.1097/PCC.0b013e3181e911a7
[4]	World Health Organization (2017) Dengue. WWW.WHO.int/denguecontrol/faq/en/index6.html
[5]	Iurii, B. (2015) A survey of Mathematical Model of Dengue Fever. Electronic Thesis & Dissertations, 1236. http://digitalcommons.georgiasouthern.edu/etd/1236
[6]	Hossain, Md.S., Nayeem, J. and Podder, Dr.C. (2015) Effects of Migratory Population and Control Strategies on the Transmission Dynamics of Dengue Fever. Journal of Applied Mathematics & Bioinformatics, 5, 43-80.
[7]	Whitehorn, J. and Farrar, J. (2010) Dengue. British Medical Bulletin, 95, Article ID: 161173. https://doi.org/10.1093/bmb/ldq019
[8]	Wilder-Smith, A. and Schwartz, E. (2005) Dengue in Travelers. The New England Journal of Medicine, 353, 924-932. (Effect of Migratory Population and Control Strategies.)
[9]	Varatharaj, A. (2010) Encephalitis in the Clinical Spectrum of Dengue Infection. Neurology India, 58, Article ID: 585591. https://doi.org/10.4103/0028-3886.68655
[10]	Chowell, G., Diaz-Duenas, P., Miller, J.C., Alcazar-Velazco, A., Hyman, J.M., Fenimore, P.W. and Castillo Chavez C. (2007) Estimation of the Reproduction Number of Dengue Fever from Spatial Epidemic Data. Mathematical Biosciences, 208, 571-589. https://doi.org/10.1016/j.mbs.2006.11.011
[11]	Lakshmikantham, V., Leela, S. and Martynyuk, A.A. (1989) Stability Analysis of Nonlinear Systems. Marcel Dekker, Inc., New York and Basel.
[12]	Jelinek, T. (2000) Dengue Fever in International Travelers. Clinical Infectious Diseases, 31, 144-147. https://doi.org/10.1086/313889
[13]	Esteva, L. and Vargas, C. (1999) A Model for Dengue Disease with Variable Human Population. Journal of Mathematical Biology, 38, 220-240. https://doi.org/10.1007/s002850050147
[14]	Coutinho, F.A.B., Burattini, M.N., Lopez, L.F. and Massad, E. (2006) Threshold Conditions for a Non-Autonomous Epidemic System Describing the Population Dynamics of Dengue. Bulletin of Mathematical Biology, 68, 2263-2282. https://doi.org/10.1007/s11538-006-9108-6
[15]	Esteva, L. and Vargas, C. (1998) Analysis of a Dengue Disease Transmission Model. Mathematical Biosciences, 150, 131-151. https://doi.org/10.1016/S0025-5564(98)10003-2
[16]	Esteva, L. and Vargas, C. (2000) Influence of Vertical and Mechanical Transmission on the Dynamics of Dengue Disease. Mathematical Biosciences, 167, 51-64. https://doi.org/10.1016/S0025-5564(00)00024-9
[17]	Akinpelu, F.O. and Ojo, M.M. (2016) A Mathematical Model for the Dynamic Spread of Infection Caused by Poverty and Prostitution in Nigeria. International Journal of Mathematics and Physical Sciences Research, 4, 33-47.
[18]	Esteva, L., Gumel, A. and Vargas, C. (2009) Qualitative Study of Transmission Dynamics of Drug-Resistant Malaria. Mathematical and Computer Modelling, 50, 611-630. https://doi.org/10.1016/j.mcm.2009.02.012
[19]	Akinpelu, F.O. and Ojo, M.M. (2017) Sensitivity Analysis of Ebola Model. Asian Research Journal of Mathematics, 2, 1-10. https://doi.org/10.9734/ARJOM/2017/30642
[20]	Iboi, E. and Okuonghae, D. (2016) Population Dynamics of a Mathematical Model for Syphilis. Applied Mathematical Modelling, 40, 3573-3590. https://doi.org/10.1016/j.apm.2015.09.090
[21]	Hethcote, H.W. and Thieme, H.R. (1985) Stability of the Endemic Equilibrium in Epidemic Models with Sub Populations. Mathematical Biosciences, 75, 205-227. https://doi.org/10.1016/0025-5564(85)90038-0
[22]	Derouich, M. and Boutayeb, A. (2006) Dengue Fever: Mathematical Modelling and Computer Simulation. Applied Mathematics and Computation, 177, 528-544. https://doi.org/10.1016/j.amc.2005.11.031
[23]	Olaniyi, S. and Obabiyi, O.S. (2013) Mathematical Model for Malaria Transmission Dynamics in Human and Mosquito Populations with Nonlinear Forces of Infection. International Journal of Pure and Applied Mathematics, 88, 125-156.
[24]	Takahashi, L.T., Maidana, N.A., Ferreira Jr, W.C., Pulino, P. and Yang, H.M. (2005) Mathematical Models for the Aedes aegypti Dispersal Dynamics: Travelling Waves by Wing and Wind. Bulletin of Mathematical Biology, 67, 509-528. https://doi.org/10.1016/j.bulm.2004.08.005
[25]	Vanden Driessche, P. and Watmough, J. (2002) Repro-duction Numbers and Sub-Threshold Endemic Equilibria for Compartmental Models of Disease Transmission. Mathematical Biosciences, 180, 29-48. https://doi.org/10.1016/S0025-5564(02)00108-6
[26]	Agusto, F.B. and ELmojtaba, I.M. Optimal Control and Cost-Effective Analysis of Malaria/Visceral Leishmaniasis Co-Infection. PLoS ONE, 12, e0171102. https://doi.org/10.1371/journal.pone.0171102