A mathematical model for effective control and possible eradication of malaria

In this paper, a deterministic mathematical model for the transmission and control of malaria is formulated. The main innovation in the model is that, in addition to the natural death rate of the vector (mosquito), a proportion of the prevention efforts also contributes to a reduction of the mosquito population. The motivation for the model is that in a closed environment, an optimal combination of the percentage of susceptible people needed to implement the preventative strategies ( α ) and the percentage of infected people needed to seek treatment can reduce both the number of infected humans and infected mosquito populations and eventually eliminate the disease from the community. Prevention effort α , was found to be the most sensitive parameter in the reduction of R 0 . Hence, numerical simulations were performed using different values of α to determine an optimal value of α that reduces the incidence rate fastest.It was discovered that an optimal combination that reduces the incidence rate fastest comes from about 40% of adherence to the preventive strategies coupled with about 40% of infected humans seeking clinical treatment as this will reduce the infected human and vector populations considerably.


Introduction
There are many infectious diseases for which there is a cure, yet they remain endemic in many communities around the world. Malaria is one of such infectious diseases. In spite of the fact that there is not yet a vaccine for malaria, there are a number of preventive measures which, when properly administered, can play the role of "vaccination", by reducing considerably, the size of the susceptible class. The key to eradicating any infectious disease lies in effective treatment for the infected population, plus a rapid reduction in the susceptible class below a certain threshold.
The Global Malaria Control Strategy is a concerted effort meant to bring about changes in the way malaria problem is addressed. As a result, this strategy stresses the selective use of preventive measures wherever they can lead to sustainable results. The measures are aimed at halting the deteriorating effects of the malaria situation, minimizing the wasteful use of resources and contributing appropriately to the development of health services, intersectoral cooperation and community participation. Malaria affects the health and wealth of nations and individuals alike. In Africa today, malaria is understood to be both a disease of poverty and a cause of poverty [14,31]. Malaria is said to have significant measurable direct and indirect costs, and has been proven to be a main constraint to economic development [31]. This means the gap in prosperity between countries with malaria and countries without malaria has become wider every single year. It has been stated that where malaria has been eliminated, economic growth has increased substantially [13]. Hence we need to find cost-effectiveness of the intervention strategies.
Mathematical modelling uses ordinary, partial, delay and other forms of differential equations in understanding the dynamics of natural occurrences. [25] both qualitative and quantitative features of delay differential equations were applied to real-life problems. Investigation into the sensitivity of model solutions due to perturbing the parameters appearing in delay differential systems, using variational and direct approaches were done in [24].
Mathematical modelling of infectious diseases and their constructive analytical study together with numerical simulations have greatly contributed to understanding diseases dynamics and their possible control measures. It plays an essential role in informing policy and management decisions during outbreak diseases. Mathematical models play an important role in understanding the dynamics and tracking tumour and immune populations over time. [27] applied delay differential model with control variables that describe the interactions of immune cells, tumour cells, normal cells, and immunochemotherapy treatment with control variables, it was confirmed that the optimal treatment strategies reduce the tumour cells load and increase the effector cells after few days of therapy. [28] presented a delay differential model to describe the interactions between the effector and tumour cells. [4] addressed the critical role tick play in Q fever transmission in livestock, a nonlinear integer order mathematical model was developed to represent the spread of this infectious disease in livestock. Caputo, Caputo-Fabrizio and Atangana-Baleanu derivatives and integrals were used for the analysis.
Over a century, mathematical models have been used to produce an explicit framework for understanding malaria transmission dynamics in human population.
The mathematical modelling of malaria transmission started with the work of [29] and [16] and was continued by [3]. Sir Ronald Ross is seen as the originator of modern mathematical epidemiology. In his pioneering work on malaria, he observed that it is transmitted from humans to mosquitoes and vice versa while working at the Indian Medical Service in 1911. R. Ross formulated in [30] a simple model, which is called the classical "Ross-Macdonald model". This model explained the relationship between the number of mosquitoes and incidence of malaria in humans. Different researchers like [23], [3], [34] and others have formulated different models from the Ross's model. They incorporate different factors such as latent or exposed classes in mosquitoes and humans, age-related differential susceptibility to malaria in human population, acquired immunity and genetic heterogeneity of host and parasite. [7] stated that planning and modelling efforts should focus on preserving the efficacy of our current drugs and insecticides as the minimum requirement for any chance of success in the control and eventual eradication of malaria.
[1] proposed a new deterministic model for the dynamics of malaria and the effectiveness of drugs. The human population was portioned into three distinct compartments of susceptible, infected and the recovered, while the mosquito population was portioned into two distinct Compartments of susceptible and the infected. [10] modelled malaria as a 7-staged state Ordinary Differential Equations (ODE's) to compare intervention strategies for malaria control for two representative areas of high and low transmission focusing on the most effective prevention strategies. [32] formulated a mathematical model of nonautonomous ordinary differential equations describing the dynamics of malaria transmission with age structure for the vector population.It was shown that the transmission of malaria is highly influenced by the dynamics of immature mosquitoes. [19], considered SEIR − SEI model of malaria transmission between humans and mosquitoes. They concluded that Malaria may be controlled by reducing the contact rate between human and mosquito, the use of active malaria drugs, insecticides and mosquito treated nets can help to reduce mosquitoes population and malaria transmission respectively. In this study, we formulate a malaria model similar to that of [32,10,12]. We incorporate prevention strategies (such as sleeping under insecticide-treated bed nets, intermittent prophylactic treatment, the use of mosquito repellent and all other prevention strategies) that move people directly from the susceptible class to the recovered class. We also model the effect of the preventive strategies on the vector population by assuming that a certain proportion, c where 0 ≤ c ≤ 1, of these prevention strategies increases the death rate of the mosquitoes. By this, we mean prevention strategies that kill mosquitoes in the process such as spraying (space and indoor residual spraying), larviciding and the others. Numerical simulations were performed to come up with the best levels of α needed to reduce or eradicate malaria.
The paper is organized as follows: A deterministic mathematical model for the transmission of malaria is formulated in Section 2. Local stability of the equilibrium points is analyzed in Section 3. Parameter estimation and sensitivity analysis is done in Section 4. In Section 5, numerical simulations are performed and a conclusion given in Section 6.

Proposed Model of malaria transmission 2.1 Model Formulation
Many researchers have formulated different mathematical models in an attempt to combat malaria from our communities but none has been able to come up with a model that suggest the optimal combinations (percentage) of preventive strategies and the number of infected humans needed to seek clinical treatment in order to reduce both the infected human and vector populations. In [33], the modelling exercise demonstrates that combining interventions with suitable delivery strategies can protect drugs from selective pressure through indirect effects on the transmission dynamics. [33] used a simple mathematical model to demonstrate the interrelationships between different interventions and the potential for malaria elimination when they are combined. Based on this, we formulate and analyze a malaria model as an SEIRS for the host population and SEI for the vector population similar to that of [10,32,19,12] by combining biting rate and transmission rate asβ h for host andβ v for vector, describing the transmission of malaria.The human population is divided into: susceptible, S h (t), exposed or latent, E h (t), infectious, I h (t), and recovered (immune), R h (t) classes. People are born (at a constant per capita rate) or immigration (at a constant rate) into the susceptible class.
The mosquito population is divided into: susceptible, S v (t), exposed or latent, E v (t), and infectious, I v (t) classes. Female mosquitoes (we do not include male mosquitoes in our model because only female mosquitoes bite human hosts for blood meals) are born or migrate to the susceptible class through birth.
The following assumptions are made in order to formulate the equations of the model: • A certain proportion of the prevention strategies contribute to the death rate of the mosquitoes. bites the human host.
• The infective human population recovers with temporary immunity with clinical treatment.
• Mosquitoes do not die from the infection.
• All new humans and vectors are born susceptible.
The transition diagram for the proposed model is presented in Figure 1.  Table 1 shows the parameters used in the model 1 and their interpretations. The systems of equations (1) is obtained from the transition diagram Figure 1: Let All population compartments are non negative for all t > 0 in the feasible region Ω, . By summing up the first four equations and the three last equations of (1) it can be shown that all the solutions are bounded in Ω, ∀t > 0 such that For positivity of solutions, we apply the theorem that follows. Theorem 2.1 Let the initial data be (1) is positive for all t > 0.

Using separation of variables and integrating both sides gives
using the initial conditions: From the second equation we have, Using separation of variables and integrating both sides gives it can be shown that for all t > 0, the remaining equations of the system (1) are positive. Therefore, the region Ω is positively-invariant ( i.e. solutions of the system (1) generated by nonnegative initial conditions remain inside Ω for all time t ≥ 0 and the model (1) is biologically, epidemiologically meaningful and mathematically well-posed in the domain Ω. Let It must be noted thatN v N h gives an average number of female mosquitoes per one human host (or average vectorial density). This term is important because the same human host may get bitten by various mosquitoes during one day [5]. Similarly, let Then system (1) can be written as the equivalent seven-dimensional non-linear system of ODEs:

Equilibrium points of the model and Basic Reproduction Number (R 0 )
The system (2) has two equilibrium points namely; • The Disease Free equilibrium, and • The Endemic equilibrium,

The Disease Free equilibrium
At the disease free equilibrium all the disease classes are zero. It a scenario which depicts an infectionfree state in the community or society. The Disease Free equilibrium of the system (2) is given by In the next section, we compute the basic reproduction number. The endemic equilibrium will expressed in terms of the basic reproduction number (R 0 ).

The basic reproduction number (R 0 )
The basic reproduction number (R NG 0 ) of the model was computed using the next generation matrix approach discussed in [11]. It is defined as the number of expected secondary cases that is produced by a single infected individual in a completely disease-free or susceptible population over the time of its infectious period. An R 0 < 1 implies the disease cannot invade the population while for an R 0 > 1, it means that each individual in the infectious classes (I h or I v ) produces more than one secondary infected individual. The basic reproduction number here is the spectral radius of the product derived from new infection matrix and the inverse of the transition or transfer matrix. The basic reproduction for (2) is given as .
According to [18], the square of the reproduction number obtained via the next-generation approach [11]: is the basic reproduction number (R 0 ) of the model because only R 0 gives the number of secondary infections that one infective host individual will produce in an entirely susceptible host population during its lifespan as infective. Therefore, the basic reproduction number for (2) is given as .
Now, for a population where there is no prevention effort, α = 0, the basic reproduction number is represented by . where: means the probability that a human will survive the exposed state to become infectious; is the probability that a vector will survive the exposed state to become infectious; is the number of humans that one vector infects during its infectious lifetime, provided all humans are susceptible; is the number of vectors that one human infects during his/her lifetime provided all vectors are susceptible;

Endemic equilibrium
The solution for the endemic equilibrium is obtained in terms of the infected humans which is also expressed in terms of R 0 and after some algebraic manipulation, we have; 3 Local stability analysis of the equilibrium points In this section, a corollary of Gershgorin circle theorem stated in [2] is used to investigate the local stability of the equilibrium points of the model. The stability of an equilibrium point (stationary states) of a mathematical model for an infectious disease helps to determine whether the solutions remain near the equilibrium point or get further away or not. The equilibrium point can be either stable or unstable or a saddle point [15,21]. The Gershgorin's theorem provides sufficient conditions for the eigenvalues to lie in the left half of the complex plane [2,22,20,6]. That is, the local stability can be established without the need to calculate the eigenvalues, instead the basic reproduction number which also gives a condition for an equilibrium point to be stable is used for the analysis. determines the sign of the constant term [17].
Corollary 3.1 (Corollary of Gershgorin Circle Theorem) Let A be an n × n matrix with real entries. If the diagonal elements a ii of A satisfy for i = 1, ..., n, then the eigenvalues of A are negative or have negative real parts [2].
To use the corollary of Gershgorin circle theorem to investigate the local stability of the system, we have to obtain an n × n matrix from the system that is the jacobian matrix. The jacobian matrix for the system (2) is given as where n * = β h I v − (α + µ h ).

Local stability analysis of the disease free equilibrium point
Now, we investigate the local stability of the disease free equilibrium point with α = 0, since the basic reproduction for the model without prevention is always greater than the basic reproduction for the model with prevention, that is, R 0 < R ′ 0 .
Theorem 3.1 The disease free equilibrium (E 0 ) with α = 0 is locally asymptotically stable if R ′ 0 < 1. Proof 3.1 Evaluating the jacobian matrix at the disease equilibrium point (E 0 ) yields From rows 1, 4 and 5, we get three of the eigenvalues λ 1 = −µ h , λ 2 = −(µ h + ρ) and λ 3 = −µ c which are negatives. Eliminating rows 1, 4 and 5 from J 0 , we end up with the sub-matrix The matrix J s satisfies the corollary of Gershgorin's circle theorem, if the following inequalities hold; Combining Eq. (9) to Eq. (12) gives, Therefore, R 0 < 1, ensures local stability of the disease free equilibrium.

Global stability conditions for the Disease Free Equilibrium
We have established that the disease-free equilibrium is locally asymptotically stable whenever R 0 < 1 and unstable when R 0 > 1. Now, we want to establish the global stability of the disease-free equilibrium.
According to [8], for the system (2) with disease-free equilibrium . can be expressed as:ẋ E = G(x, I) G(x, 0) = 0 where x ∈ R 3 denotes (its components) the number of susceptible human, susceptible vector and recovered human(that is, uninfected individuals) and I ∈ R 4 denotes (its components) the number of latent human, infectious human, latent vector and infectious vector (that is, infected individuals including latent and infectious). E 0 = (x * , 0)denotes the disease-free equilibrium of this system. And the conditions (1) and (2) below must be met to guarantee global asymptotic stability.  0) is an M-matrix (the off diagonal elements of A are nonnegative) and D is the region where the model makes biological sense. If System (13) satisfies the above two conditions then the following theorem holds: Theorem 3.2 The equilibrium point E 0 is globally asymptotically stable (g.a.s.) equilibrium of (13) provided that R 0 < 1 (l.a.s.) and that assumptions (l) and (2) are satisfied.
Considering the system (2), It is also clear that, 0). Hence, by the above theorem E 0 is globally asymptotically stable.

Local stability analysis of endemic equilibrium point
where Using the corollary of Gershgorin's circle theorem, the matrix (J e ) will have negative eigenvalues if the following inequalities are satisfied. That is, From Eq. (16) and Eq. (17) we get, Let Then Eq. (22) becomes Adding Eq. (19) and Eq. (20) yields Since h * = −d *

Estimated parameters
A suitable curve for the model is obtained using the data by applying Least Squares Estimation method. The method of least squares is about estimating parameters by minimizing the squared discrepancies between observed data and their expected values. Non-linear least squares is a form of least square analysis used to fit a set of m observations with a model that is non-linear in n unknown parameters.   Table 2 shows the estimated parameters and their sources for model (1). The rates are given per day.

Sensitivity analysis
The assessment of the changes associated with model outcomes as the model parameters are varied is the main goal of sensitivity analysis. This method is very important because it helps us to identify the weak points in the model [26]. To know how best to reduce or eradicate an infectious disease such as malaria, and to minimize it mortality, it is essential to know the relative importance of different factors responsible for its transmission and prevalence. A parameter is said to be sensitive if small change in the parameter value produce a large changes in the solution of the differential equations. The derivative of the solution for the parameter is used to examine a change in the solution of the model with respect to the particular parameter [9]. From Table 3, the sensitivity indices for each parameter involved in the basic reproduction number R 0 of model (2) were computed. The numerical values were obtained after solving the respective partial derivative for the parameters in the R 0 .
The sensitivity index with negative signs indicate that, for an increase in the corresponding parameters there is a decrease in the value of the reproduction number R 0 and vice versa. It is noted that the sensitivity index of R 0 with respect to β v and β h do not depend on the values of the parameters because are always exactly equal to 1. The most sensitive parameter is α which was obtained by ) , the basic number R 0 decreases as α increases which support the claim that α has effect on both human and vector populations.The second sensitive parameter is c, which is the proportion of prevention efforts that increase the death rate of the vector. From the sensitive analysis, the R 0 decreases as c increases. So, as the proportion of prevention efforts increases, the eradication of malaria becomes much easier. The third sensitive parameter is γ. For almost all parameters, the sign of the sensitivity indices of R 0 (i.e., whether R 0 increases or decreases when a parameter increases) correspond with an intuitive expectation.

Numerical simulations
A numerical simulation of the model is conducted to find out the dynamics of the disease in the human and vector populations. The simulations were conducted using MATLAB's ode45. The initial conditions used are proportions of the total humans and mosquito populations. The total human population was N h = 2567369. For instance, for the initial value of I h we have I h = 59000 2567369 , = 0.023, the initial values for the other state variables for the human classes are given as S h = 0.977, E h = 0, R h = 0. Assumed the mosquito population is about 10 times that of the human population, the proportions will be the same. Therefore, S v = 0.977, E v = 0 and I v = 0.023. The main strategy considered in controlling malaria is to find the right combinations of; 1. the percentage of people needed to carry out the preventative strategies (α) in order to eradicate malaria.
2. the percentage of infected people needed to seek proper treatment to help eradicate malaria.

Effect of α on the infected human and vector populations
The effect of α on the infected populations is examined by finding the basic reproduction numbers for each value of α for α = 0.1, 0.2, 0.3, 0.4, 0.6, and 0.8. where α = 0.1 corresponds to 10% of individuals needed to apply the prevention strategies in quest of combating malaria. The same analogy goes for the other values of α. Table 4 shows the R 0 values for different values of α. It can be observed from Table 4 that some R 0 values are less than 1, which imply that if we have those levels of α there is the possibility of reducing or eradicating malaria.  Table 4. That is, the R 0 values corresponding to α = 0.1, 0.2, 0.3 are all greater than 1 whereas for α = 0.8, the curve drops to about 0.002 (which is about 5,135 infected humans) in just 50 days. The R 0 values for α ≥ 0.4 are less than 1. This means that for α ≤ 0.3 they will all converge to the endemic equilibrium because they never touch the time axis but for each α ≥ 0.4, the curves converge to the disease free equilibrium. As α increases, the time taken for the graph to touch the time axis decreases. This means the time taken for the infected humans to reach the disease free equilibrium (for the disease to die out) decreases. In Figure 4, we examine the effect of increasing the value of α on the infected mosquito population,

Effect of α on the infected vector population
It can be observed that, as α increases from 0.1 through to 0.8 the number of infected mosquitoes decreases. The implication is that, the susceptible mosquitoes find it difficult to find an infected human to bite. This buttresses our claim that preventive measures such as spraying (space and indoor residual spraying), larviciding and environmental management contribute to the death of mosquitoes in all the vector class. If the cycle continues, it will decrease their population in the long run and the end result will be a malaria-free society.

Effect of varying cα on the infected vector population
The

Effect of increasing the proportion of infected individuals needed to seek clinical treatment on infected human and vector populations
In this section, we want to look at 10%, 20%, 40%, 60%, 80% of individuals needed to be clinically treated coupled with fixed percentage of individuals needed to apply the prevention strategies in combating malaria. The figures below show the effect of keeping α constant and varying proportions of infected humans needed to seek clinical treatment on the infected human population. The figures show that increasing the proportions of infected humans needed to seek clinical treatment decreases the infected mosquito population. With fixed α values, increasing the proportions of infected humans needed to seek clinical treatment yield a significant drop in the proportion of infected mosquitoes number. This means that there will be less infected mosquitoes and less susceptible humans to propagate the disease.

Conclusion
This paper formulated a deterministic mathematical model for the transmission and control of malaria. We incorporate prevention strategies (such as sleeping under insecticide-treated bed nets, intermittent prophylactic treatment, the use of mosquito repellent and all other prevention strategies) that move people directly from the susceptible class to the recovered class. We also model the effect of the preventive strategies on the vector population by assuming that a certain proportion, c where 0 ≤ c ≤ 1, of these prevention strategies increases the death rate of the mosquitoes.The equilibrium points and basic reproduction number (R 0 ) for the model were found. Stability analysis was done for the disease-free and endemic equilibrium points of the model. From the data on confirmed malaria cases used, the unknown parameters in the model were estimated. The estimated values and known values from literature for the parameters were used to perform numerical simulations using MATLAB's ode45 and sensitivity analysis of the model was also performed.
From the sensitivity analysis, it was shown that the most sensitive parameter in the reduction of R 0 was found to be α. Based on this result, numerical simulations were performed using different percentages of susceptible humans needed to apply the prevention measures (α) to help reduce the incidence rate fastest. From Figure 3, it was observed that as α value increases, the number of infected humans decreases and the basic reproduction number also decreases. From Figure 3, it takes about 300 days for the infected human population to reach the disease free equilibrium when α = 0.4 while it takes about 100 days for the infected human population to reach the disease free equilibrium when α = 0.8. For α ≤ 0.3 all the infected human populations converge to the endemic equilibrium and for α ≥ 0.4 the infected human population converge to the disease free equilibrium. Also, as the number of recovered humans and cα increases, the time taken to reach the disease free equilibrium decreases. This can be seen in Figures 11 to 14. Also, for the infected vector population I v , it can be observed from Figure 4 that as α increases from 0.1 through to 0.8 and also from Figures 15 to 18, as the number of recovered humans and cα increases, the number of infected vector population decreases. The implication is that, the susceptible mosquitoes find it difficult to find an infected human to bite. It can be observed from Figures 7 to 10 that as the value of c increases the infected vector population decreases and the curves touch the time axis (the disease free equilibruim) faster. For instance, for c = 0.1 and α = 0.8, it takes about 100 days before the infected vector population decreases and touches the time axis whiles in Figure 10, for c = 0.25 and α = 0.8, it takes just 50 days for the infected vector population to decrease and touch the time axis. Also, the value of R 0 decreases as the value of c increases, which confirms that cα contribute to the death rate of mosquitoes. This means that preventive measures such as spraying (space and indoor residual spraying), larviciding and environmental management contribute to the death of mosquitoes in all the vector classes. If the cycle continues, it will decrease their population in the long run and the end result will be a malaria-free society. It is recommended that since malaria eradication continues to pose a big challenge to national malaria control programmes in most developing countries, there is the need to strengthen the control strategies at hand as well as looking for some new ones and early clinical treatment campaign should be intensified. Based on the result of the simulations, it is suggested that if by a concerted effort about 40% and above of the people in the closed neighbourhood (that is, the communities) apply due diligence to preventive measures such as space and indoor residual spraying, larviciding and environmental management and about 40% of the infected humans as seen in Figures 11 to 18 are treated well, there will be no infected human to infect the mosquitoes. Also, since the mosquito's life expectancy is just about 30 days, all the infected mosquitoes will die out and this will help humans and mosquitoes to co-exist in the communities without malaria. This implies that if the prevention efforts are implemented by all contiguous communities in the region, acting in concert, then a marked reduction would be seen in the infected vector population leading to a complete eradication of malaria in the region, by implication, through out Ghana and the world as well.