Abstract
Malaria remains by far the world's most important tropical disease, killing more people than any other communicable disease. A number of preventive and control measures have been put in place and most importantly drug treatment. The emergence of drug resistance against the most common and affordable antimalarials is widespread and poses a key obstacle to malaria control. A mathematical model that incorporates evolution of drug resistance and treatment as a preventive strategy is formulated and analyzed. The qualitative analysis of the model is given in terms of the effective reproduction number, . The existence and stability of the disease-free and endemic equilibria of the model are studied. We establish the threshold parameters below which the burden due to malaria can be brought under control. Numerical simulations are done to determine the role played by key parameters in the model. The public health implications of the results are twofold; firstly every effort should be taken to minimize the evolution of drug resistance due to treatment failure and secondly high levels of treatment and development of immunity are essential in reducing the malaria burden.
1. istant parasites, and temporary immune humans and Introduction
Malaria is one of the most devastating infectious diseases in the world, infecting approximately 300–500 million people annually worldwide, resulting in over a million deaths [1, 2]. The high risk groups include young children, pregnant women, and nonimmune travelers to malarious regions. It is reported that one child dies every 30 seconds from this disease. Malaria is caused by Plasmodium parasites. There are 4 main types: Plasmodium vivax, Plasmodium ovale, Plasmodium malaria, and Plasmodium falciparum. It is the Plasmodium falciparum that causes most of the severe diseases and deaths and is most prevalent in Africa [3].
Malaria parasites are transferred by female Anopheles mosquitoes while they bite humans for a blood meal for the development of their eggs. During the blood meal, the mosquito injects sporozoites into the blood stream. In few minutes, the sporozoites enter the liver cells where each sporozoite develops into a tissue schizont that contains 10,000 to 30,000 merozoites. After 1-2 weeks the schizont ruptures and releases the merozoites into the blood stream which then invade the red blood cells. The clinical symptoms of malaria are due to the rupture of the red blood cells and release of the parasites’ waste and cells’ debris into the blood stream. After several asexual cycles, some of the merozoites invade red blood cells and there develop into either male or female gametocytes. When a female mosquito bites an infected person, it ingests blood containing male and female gametocytes. These may be transmitted to another human host during the mosquito’s next blood meal. The viable female and male gametocytes and the mosquito must survive long enough for sporozoite development and subsequent inoculation into another human host. The main symptoms of malaria include fatigue, chill, headache, abdominal and back pain, diarrhea, sometimes vomiting, and fever. The malaria fever is known to undulate approximately every 2-3 days, depending on the Plasmodium species.
There have been a number of control measures against the transmission of malaria and most importantly antimalarial drugs that have played a mainstream role through the treatment of clinical cases, prophylaxis of the high risk groups (infants, nonimmune travelers, and pregnant women). Early and prompt treatment of suspected cases with adequate drugs continue to be the main control strategy of malaria in the sub-Saharan Africa. However, despite intensive control efforts, global incidence of malaria is increasing especially in the sub-Saharan Africa. This has been attributed to poverty, war, and the collapse of health care systems following abandonment of WHO’s Eradicate Malaria Campaign [4] in the 1960s. Another important factor has been largely due to emergence of drug resistant strains of Plasmodium to cheap and affordable antimalarial drugs such as chloroquine (CQ) and sulfadoxine-pyrimethamine (SP, trade name: Fansidar) that have been the mainstay of malaria treatment in disease endemic regions [5–7]. In Uganda, high malaria awareness and access to cheap and effective antimalarials have been able to ensure that the disease has been reasonably controlled. However, it is reported in [8–10] that resistance to these drugs is now widespread in Uganda and many parts of East Africa and this poses a formidable obstacle to malaria control.
Drug resistance has been defined as the ability of the parasite strain to survive and/or multiply despite the administration and absorption of a drug given in doses equal to or higher than those usually recommended but within the limit of tolerance of the subject [18]. The factors associated with drug resistance include longer half-life, single mutation for resistance, poor compliance, host immunity, and the number of people using the drugs. The characteristics of a drug that make it vulnerable to the development of resistance are a long terminal elimination half-life, a shallow concentration-effect relationship, and mutations that confer marked reduction in susceptibility. In order to slow down the spread of drug resistance, the World Health Organization [19] recommend the use of combination therapy to provide effective treatment against malaria where either the pharmacokinetics of the two or three components must be well matched or one of the drug components that rapidly reduce sufficiently parasite biomass with a partner drug that can remove any residual parasites [20]. This has been because the use of drug combinations with different mechanisms of action prevents further selection of resistance and results in higher cure rates [13]. The use of Artemisinin based combination therapies (ACTs) is an alternative strategy that offers promise of both increased efficacy and retardation in the rate of development of parasite drug resistance [12, 15, 16, 21] by delaying resistance for a much longer time period and is seen as an important tool in an effort to roll back malaria (RBM). The ACTs also have an advantage of inhibiting gametocytogenesis which would curtail the role of the vector in the emergence and spread of drug resistant parasites especially in mass drug administration (MDA) programmes. The effects of drug resistant malaria are higher malaria morbidity and mortality, increased cost of malaria case management, increased burden on the health care facilities, and increased relative prevalence of Plasmodium falciparum. All these affect the design of surveillance programmes [22] and the likely impact of malaria control programmes.
When malaria parasites are exposed to subtherapeutic concentrations of a single drug, resistant strains are more likely to develop. The use of combination therapies significantly reduces the risk of the development of resistance to a particular drug [23]. The coexistence of drug sensitive and resistant malaria parasites within an infected host results in competition for space and red blood cells as the malaria infection of strain progresses within the infected human host. The two scenarios likely to occur are that, in the absence of drug treatment, the resistant parasites are at a competitive disadvantage and are rarely noticed in the absence of drug treatment. They only emerge after treatment. According to Koella and Antia [24], sensitive malaria parasites limit the proliferation of the drug resistant parasites and when an infected individual is treated the drug kills off the sensitive strains, leaving the resistant strains to proliferate.
It is important to understand the complex dynamics of the host-vector interactions that enhance the evolution of drug resistance as this would help in designing appropriate control measures of the malaria transmission. The use of mathematical models is a valuable tool in the study of the dynamics of diseases. They help us to understand and predict epidemiological patterns and processes of the disease. The models for the transmission dynamics of resistant and sensitive strains of the parasite can be used to evaluate and compare the factors that influence the spread of antimalarial drug resistance and enable us to figure out strategic decisions that are of significant importance in controlling the spread of drug resistance.
Mathematical models have been used before in the study of the transmission and spread of malaria (see [25–28]). Their use in the study of the evolution of drug resistance of the malaria parasites has been considered in [12, 22, 24, 29]. However, most of these models have been formulated as a simple mass action and assume constant human host and mosquito vector populations in which births and deaths are equal. They also ignore disease induced mortality which seems to be unrealistic especially in the sub-Saharan Africa where most deaths are due to malaria related cases.
In the present paper we extend Tumwiine et al. [30] model that is based on the susceptible-infective-immune (SIRS) in the human population and susceptible-infective (SI) in the mosquito population. The infected human population consists of individuals with drug sensitive and resistant strains, respectively. Infection with the drug sensitive strains will give rise to the drug resistant strains in the event of treatment failure.
This paper is organized as follows. In Section 2, the variables and parameters are defined and the model is formulated. In Section 3, we establish the existence and uniqueness of solutions for the model equations and then carry out the stability analysis. In Section 4, numerical simulation of the model results is presented. The discussion of the results and conclusion is done in Section 5.
2. The Model
In this section, we formulate a mathematical model for the transmission and spread of malaria. We use a compartmental model in which individuals move between susceptible, infected, and immune classes in the human population and between susceptible and infected classes in the mosquito population, respectively. This is an extension of Tumwiine et al. [30] model in which we now divide the infected human population into individuals infected with drug sensitive and resistant malaria strains. Our model ignores superinfections and the latent periods of the disease.
The total human and female mosquito populations are denoted by and , respectively. The human population is divided into the epidemiological classes: susceptible, infected, and temporary immune denoted by , , and , respectively. The infected human population consists of individuals with drug sensitive malaria strains, and, drug resistant malaria strains, . The subscripts and are added to the model variables and parameters in order to specify the sensitive and resistant strains, respectively. Thus, . The mosquito vector population, , has the epidemiological classes denoted by and for the susceptible and infected classes, respectively. Thus, .
The recruitment into the human population is assumed to be at a constant rate and corresponds to births and/or immigration. The interaction of humans and female mosquitoes is modelled by standard incidence [30], with the terms denoting the rate at which susceptible humans get infected by infected mosquitoes and denoting the rate at which the susceptible mosquitoes are infected by infected human hosts .
We let a fraction of humans be infected with drug sensitive strains and the remaining fraction of individuals are infected with drug resistant strains. The drug administered kills off the sensitive malaria strains at a constant rate and these individuals return to the susceptible class. The evolution of drug resistance due to treatment failure is incorporated into the model as a measure of the rate at which humans infected with drug sensitive strains progress to the drug resistant class due to treatment failure. Humans infected with sensitive and resistant strains acquire temporary immunity at constant rates and , respectively. The drug resistant strains are cleared at a faster rate than the sensitive infections, and thus resistant parasites transmit for a relatively short period. Since immunity against malaria infection is temporary, is lost at a constant rate , and the immune humans become susceptible again. All humans are subjected to a nondisease related per capita natural death rate and an additional disease related death rate constant in the sensitive and resistant infective classes.
There are constant per capita birth and death rate and , respectively, in the mosquito population. The susceptible female mosquitoes become infected after biting infected humans at a constant rate , is the fraction of infective bites inflicted by infectious mosquitoes and, of these bites lead to infection with as the probability that a mosquito becomes infectious. We use (susceptible-infectious-recovered-susceptible) model for the human population and (susceptible-infectious) model for the mosquito vector population. The mosquito vector population has no immune class (mosquito vector never recovers from malaria). This is because of their short life cycle.
The dynamics of malaria is given by the following coupled system of nonlinear ordinary differential equations: Since and , this implies that The behavior of this model is examined by means of an equilibrium and stability analysis. This can easily be done by rewriting the equations in terms of the fractions of individuals in each class. Define , , , and , respectively, as the fraction of susceptible humans, humans infected with sensitive parasites, humans infected with resistant parasites, and temporary immune humans and and as the fraction of susceptible and infected mosquitoes, respectively. Given that the total number of bites made by the mosquitoes must equal the number of bites received by the humans, is a constant; see [31]. denotes the number of female mosquitoes per human host; that is, . Then, differentiating the fractions with respect to time gives the following system of differential equations: It is noted that the total human and mosquito population sizes are variable and, in the absence of the disease, the total human population size approaches the carrying capacity .
Using the equations and to eliminate the variables and , respectively, then system (3) is reduced to the 5-dimensional system of differential equations given by Here, we note that the region of biological interest is positively invariant with respect to system (4), where denotes the nonnegative cone of including its lower dimensional faces. We denote the boundary and the interior of by and , respectively.
3. Analysis of the Model
The equilibria are obtained by putting the right hand sides of system (4) equal to zero and then substituting for into the rest of the four expressions of system (4) to give
3.1. Existence of Equilibrium Points
We need to verify that the disease-free equilibrium (DFE) point exists when the proportions of the diseased groups and for the endemic equilibrium point (EEP) we have , , and that are nonzero.
It can be easily established that system (4) has two possible equilibria in the nonnegative cone : the disease-free equilibrium point, , and a unique endemic equilibrium point, , with the coordinates expressed implicitly in system (6).
The boundary equilibrium for the drug sensitive strains only is given by , with the coordinates The boundary equilibrium for the drug resistant strains only is given by with the coordinates For the coexistence equilibrium, , the coordinates are given by system (6). The stability of this equilibrium point is carried out in Section 3.4.
3.2. The Effective Reproduction Number
The effective reproduction number, , is a key parameter that determines the behaviour of the model in the presence of treatment. The effective reproduction number term is used to distinguish it from the basic reproduction number , when there is no treatment. In order to analyze the stability of system (4) we need to first obtain the threshold condition for the establishment of the disease. Then the effective reproduction number is computed using the next generation operator method that is developed by van den Driessche and Watmough [32]. A reproduction number obtained this way determines the local stability of the disease-free equilibrium point with local asymptotic stability for and instability for . We first rearrange the system beginning with the infected classes as follows: And then we distinguish the new infections from all other class transitions in the population. The infected classes are , , and in the human host and mosquito vector, respectively. is the vector of rates of the appearance of new infections in each compartment; , where is the vector of rates of transfer of individuals into the particular compartment; and is the vector of rates of transfer of individuals out of the particular compartment. Thus, from our model, we haveIt is clear that we have three compartments for the infected classes (two from the human population and one from the mosquito vector). The Jacobian matrices and are obtained by linearizing the reduced form about the disease-free equilibrium to obtain from the partial derivatives of and with respect to the infected classes the following matrices: Note that is nonnegative, is a nonsingular -matrix whose inverse, , is nonnegative, and the next generation matrix (see [32, 33]) is nonnegative. The dominant eigenvalue corresponding to the spectral radius of the matrix is the effective reproduction number given by , where The effective reproduction number, , is defined as the number of secondary infections derived from a single primary infection in a population of susceptibles [32–35]. Its biological meaning is readily interpreted from sum of the terms denoted by and in which In the definition is the contribution of sensitive drug strains to the effective reproduction number for model and is the contribution to the effective reproduction number for model by the resistant drug strains. Our major aim is to study the transmission dynamics and spread of malaria within a population in which we have both sensitive and resistant strains. We can hopefully interpret and as the number of secondary human infections caused by one infected with the sensitive and resistant strains, respectively.
It should be noted that the parameters , , , and , from the expression for the effective reproduction number, play important roles in the dynamics of the disease. We can find the basic reproduction number, , in the absence of any intervention. This value of can be obtained from the value of the effective reproduction number when there is no treatment; that is, , , , and . Thus, in the absence of treatment, we have the basic reproduction number, , where It is clear that from the two expressions for the reproduction numbers for all . This implies that treatment in the presence of drug resistant population has a negative impact on the reduction of the spread of malaria. Numerical simulations for the threshold parameters , , , and are shown Section 4.
Using Theorem 2 of van den Driessche and Watmough, [32] we have the following Theorem.
Theorem 1. The disease-free equilibrium point, , is locally asymptotically stable (LAS) if and unstable if .
3.3. Global Stability of the Disease-Free Equilibrium
The disease-free equilibrium has global asymptotic stability (GAS) for and , which is established by proving the following theorem.
Theorem 2. If and , then the disease-free equilibrium is globally asymptotically stable.
Proof. For the Lyapunov function , we have the Lyapunov orbital derivative given by where and . It is clear that in when and . The equality holds when . It can be seen by inspection of system (4) that is the only positively invariant subset of the set where , and, by LaSalle-Lyapunov theorem [36], it follows that all trajectories starting in approach as . Thus, the theorem is proved.
We can also illustrate the result numerically, by considering multiple initial conditions. We note that all the trajectories tend to the disease-free equilibrium.
3.4. Stability of the Endemic Equilibrium
When the disease-free equilibrium point is unstable, then there exists an endemic equilibrium point, . It can be noted that the total human population size, , at the endemic equilibrium point is less than the disease-free carrying capacity, . The stability of this point is established by evaluating the Jacobian matrix at this point and is given bywhere The characteristic equation of the Jacobian matrix (17) at the endemic equilibrium point, , is a fifth-degree polynomial given by where , are the coefficients. It can be shown that all the coefficients are positive. The necessary and sufficient conditions for the local asymptotic stability of endemic equilibrium point are that the Hurwitz determinants, , are all positive for the Routh-Hurwitz criteria. For a fifth-degree polynomial [37], these criteria are given by from which we can conclude whether the endemic equilibrium point is locally asymptotically stable or unstable.
While we did not show the global stability of the endemic equilibrium, numerical results show that, irrespective of the initial conditions, the trajectories tend to an endemic steady state. We only show here the change in the infected human population with similar results being possible for the remaining subpopulations. Numerical simulations of the model showing the time series plots of the infected individuals are shown in Figure 2. It is observed that the infected populations settle to an endemic equilibrium point.
4. Numerical Simulations and Analysis
To investigate the role played by some key epidemiological parameters, we use a fourth-order Runge-Kutta numerical scheme in Matlab for the numerical simulations. The parameters used, their estimated values, and appropriate sources are given in Table 1. The illustrations of some of the numerical results are given in Figures 1 and 2.
We begin by numerically comparing the model reproduction numbers. For the parameter values in Table 1, as drug resistance increases, the change in the reproduction numbers can be tracked by Figure 3. In conformity with the expectation, increased drug resistance leads to an increase in , consequently increasing the overall effective reproduction number, .
The following initial conditions are hypothetically assumed for the proportions , , , , , and in our simulations. We begin by investigating the relationships between various parameters and how they affect the model reproduction number, which ultimately determines the dynamics of malaria. This is done in the form of contour plots given in Figure 4.
(a)
(b)
(c)
(d)
In Figure 4(a), contours of are plotted as a function of the treatment rate of sensitive strain and rate of acquiring immunity of the individuals infected with the sensitive strain, and , respectively. For the parameter values used in the simulations, the figure shows that the malaria epidemic is decreased if there is an increase in the treatment rate, , and rate of acquiring immunity by those infected with the sensitive strain, . For and , . From Figure 4(b), if the treatment rate exceeds 0.35 (i.e., ), then the epidemic can be controlled irrespective of the value of ; otherwise the epidemic will persist if . It is thus envisaged that treatment of individuals with the sensitive strains should always be high enough to control the epidemic and avoid the evolution of drug resistance due to treatment failure. Figure 4(c) shows that a decrease in the rate at which individuals acquire immunity always leads to the persistence of the epidemic. From the contour plot, the reproduction number is however more sensitive to than . In the public health context, this can be interpreted to mean that increased acquired immunity for those with the resistant strains is likely to have a greater effect in the reduction of malaria burden. The contour plot showing the variation of the model reproduction number as and vary is given in Figure 4(d). An increase in and a decrease in can lead to the persistence of malaria. Increased drug resistance is not good for communities and decreased immunity levels for those with the resistant strains are also disastrous. Thus, every effort should be taken to minimize evolution of drug resistance due to treatment failure while development of immunity for those infected by the resistant strain will be a welcome development.
4.1. Effects of Variation on Treatment on the Infected Population
Figure 5 shows the effects of varying on the infected individuals. The parameter is varied from to with a step size of . From Figure 5, it can be seen that increasing the treatment rate reduces significantly the number of infected humans with sensitive strains. The value of increases in the direction of the arrow depicted on the graph.
4.2. Effects of Evolution of Drug Resistance due to Treatment
Figure 6 shows the effects of increasing on the infected humans with drug resistant strains. The parameter is varied from to with a step size of . Increasing leads to an increase in the number of humans infected with drug resistant strains. The value of also increases in the direction of the arrow shown on the graph.
5. Discussion
In this paper, a model for the transmission and spread of sensitive and resistant strains of malaria parasites within a population is formulated and analyzed. The model incorporates that a fraction of the human population are infected with the sensitive strains and the remaining are infected with the resistant strains. The model also captures the effect of drug resistant strains through the treatment of a portion of drug sensitive individuals who progress to the drug resistant epidemiological class because of treatment failure. The model analysis included the determination of equilibrium points and carrying out their stability analysis in terms of the effective reproduction number, which gives the threshold condition for the establishment of the disease. We obtained the effective reproduction number of disease as the square root of the sum of the squares of the effective reproduction number due to sensitive and resistant strains within the human population. Our primary finding is that the effective reproduction number, which obtained in terms the demographic and epidemiological parameters. It is noted from this threshold parameter that the effectiveness of the treatment in the sensitive parasite plays a crucial role in the spread of the disease.
In addition, numerical simulations were carried out to confirm some of the analytical results that were obtained. The effects of varying key epidemiological parameters and their corresponding effects on the model reproduction number were investigated, with the view of determining the implications to the disease. Increased treatment efforts on individuals with the sensitive strains are crucial in malaria control while making sure that evolution of drug resistance is kept to the minimum.
This work can provide insights into several aspects of the malaria control in the presence of drug resistant parasites. It is important to acknowledge that the dynamics of malaria transmission are complex. Malaria transmission, especially in the sub-Saharan Africa, has been affected by funding and service delivery, political instability, poverty, and drug and insecticide resistance. The evolution of these biological and social-economic systems is complex, making it impossible for mathematical models to provide accurate predictions of the malaria dynamics. Mathematical models can however be used to forecast and project disease outcomes in time, but only in fairly gross terms. Overall, mathematical models can assist in shaping control strategies and can provide comprehensive evaluation of model assumptions that influence decisions.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
Julius Tumwiine is grateful for the support from Mbarara University of Science and Technology. Senelani D. Hove-Musekwa acknowledges with thanks the support from African Institute for Mathematical Sciences (AIMS) and National University of Science and Technology (NUST). Farai Nyabadza would like to thank the Mathematics Department of the University of Stellenbosch for its support in the production of this paper.