## Emis.de

Computational and Mathematical Methods in MedicineVol. 9, No. 2, June 2008, 143–163
Modelling immune response and drug therapy in human malaria
C. Chiyakaa*, W. Gariraa and S. Dubeb
aDepartment of Applied Mathematics, National University of Science and Technology, P.O. Box
AC 939, Ascot, Bulawayo, Zimbabwe; bDepartment of Applied Biology/Biochemistry, National
University of Science and Technology, P.O. Box AC 939, Ascot, Bulawayo, Zimbabwe
( Received 14 December 2006; final version received 12 December 2007 )
A new intra-host model of malaria that describes the dynamics of the blood stages ofthe parasite and its interaction with red blood cells and immune effectors is proposed.

Local and global stability of the disease free equilibrium are investigated. Conditionsfor existence and uniqueness of the endemic equilibrium are derived. An intra-hostbasic reproductive number is identified. We deduce that drugs based on inhibitingparasite production are more effective than those based on inhibiting merozoiteinvasion of erythrocytes. We extend the model to incorporate, in addition to immuneresponse, drug therapy, following treatment with antimalarial drugs. Using stabilityanalysis of the model, it is shown that infection can be eradicated within the host if thedrug efficacy level exceeds a certain threshold value. It will persist if the efficacy isbelow this threshold. Numerical simulations are done to verify the analytic results andillustrate possible behaviour of the models.

Keywords: intra-host; malaria; global stability; reproductive number; efficacy
AMS Subject Classification: 34D05; 34D23; 93D20; 92B05
Plasmodium falciparum malaria is a major cause of morbidity and mortality, largelyattributable to asexual parasitaemia [33]. Transmission begins when an infected femaleAnopheles mosquito takes a blood meal. As the blood meal is ingested the mosquitosimultaneously injects saliva containing plasmodial sporozoites into the human host. About30 minutes later, the motile, threadlike sporozoites migrate to the liver where they infect theliver cells [39,31,3]. They develop into schizonts (a developmental structure that containsmerozoites) which rupture and the merozoites are released, then enter the blood stream.

These merozoites infect red blood cells (RBCs) and undergo asexual reproduction, which issimilar to but quicker and less prolific than that in the liver cells. This occurs within theparasitophorous vacuole in the RBC [14]. After about 48 hours for P. falciparum the infectederythrocyte ruptures releasing daughter parasites that quickly invade fresh erythrocytes torenew the cycle [7,22]. This (erythrocytic) cycle maintains infection and directly generatesdisease symptoms [17]. Some of the merozoites that invade RBCs, instead of developingasexually, differentiate into sexual forms called gametocytes [7]. They develop through
*Corresponding author. Email: [email protected]
ISSN 1748-670X print/ISSN 1748-6718 online
*q *2008 Taylor & Francis
C. Chiyaka et al.

morphologically distinct stages, designated I – V, within the host RBC [38]. Mature (stage V)gametocytes circulate in the host's bloodstream, available to feeding Anopheles mosquitoes.

Successive erythrocytic cycles result in an increase in parasitaemia and disease
unless they are brought under control by the host's protective immune responses, drugtherapy or until the host dies. The human immune system is a remarkably sophisticateddefender of the body. It has an array of protective cells that can be mobilized to tackle aninvader, with cells of the innate immune system forming the first line of defence andthose of the adaptive immune system arriving later, but with extremely specificweaponry. During asexual reproduction in the RBC, the intracellular parasite, which isnot free within the erythrocyte cytoplasm but resides in a parasitophorous vacuole [14],then begins to modify both the biochemical and physiological processes of the cell [24].

In the process it digests the haemoglobin which makes up around 95% of the erythrocytecytosol. The rapid increase in the parasite biomass then activates innate immunemechanisms [14]. When the infected erythrocyte bursts, the foreign materials alsoactivate monocytes and macrophages, releasing cytokines that stimulate other cells ofthe immune system and presenting foreign antigens to the rest of the immune system[5,26].

Because RBCs do not express the critical molecule major histocompatibility complex
class 1 (MHC-1) on their surface, the blood stage malaria parasites cannot induce thecluster of differentiation 8 (CD8þ) T cell response and, even though the response isinduced by other infected tissue such as the liver, are probably not vulnerable to itscytotoxic action [40,37]. However, CD4þ T cells are essential for immune protectionagainst asexual blood stages in human malaria. CD4þ T cells respond to malaria antigensby proliferation and/or secretion of cytokines, e.g. interferon-gamma (INF-g) orinterleukin-4 (IL-4) [37]. On reinfection, the malaria primed T cells produce greatlyincreased amounts of INF-g which synergize with malaria glycosyphosphatidylinositol(GPI) to upregulate the production of tumour necrosis factor-alpha (TNF-a). Immunity istherefore associated with the ability to regulate the production of pro-inflammatorycytokines to an intermediate level, which mediates parasite clearance whilesimultaneously avoiding severe pathology. Falling antigen concentration leads to aswitch in the predominant T cell phenotype from Th1 (INF-g producing) to a regulatory Tcell phenotype (IL-10 and transforming growth factor-beta (TGF-b) producing). Thesecytokines mediate anti-inflammatory response. A dynamic equilibrium is required with apro-inflammatory effector mechanism targeting and controlling the parasite, and anti-inflammatory cytokines suppressing immunopathology.

Considerable work has been done on mathematical modelling of Plasmodium
falciparum infection [33,7,22,8,40,15,16,21,30,43,18,29,19]. Some of the intra-hostmodels have been reviewed by Molineaux and Dietz [32]. These models do not incorporatethe effect of the immune response and treatment explicitly. That is, they do not show theeffect of immune effectors on merozoite invasion of erythrocytes and suppression ofparasite production by antibodies. Another drawback in these models is that they do nottake into account the accelerated rate of production of RBCs during a malaria infection andthe loss of uninfected RBCs. From the evidence of the effects of the immune system ondisease progression [43,9,12,13,42,45,4], we formulate an intra-host model of malariainfection with immune response and drug therapy. The rest of the paper is organized asfollows: in the following subsection we start by restating the basic intra-host model ofmalaria infection [16,1], which is the foundation of our model. In Section 2, we extend thebasic intra-host model to incorporate the effects of immune response and further show thepositivity of solutions from this new model. The intra-host basic reproductive number
Computational and Mathematical Methods in Medicine
is derived in Section 3 where local and global stability analysis of the disease freeequilibrium is also considered. We show the form of the endemic equilibrium anddetermine its stability numerically in Section 3. We further extend the model to includetreatment with an antimalarial drug in Section 4. We perform the analysis of the model inSection 4 and show that the disease reproduction number of the system, when the infectionis subject to both immune response and drug pressure (Rg), reduces R0 by some factor.

Implications for control of the infection, basing our arguments on the basic reproductionnumber, are also discussed. Numerical simulations are performed in Section 5. A briefdiscussion of the results concludes the paper.

Intra-host model without immunity
Intra-host models of malaria infection describe the dynamics of the blood stages of theparasite and their interaction with host cells which are RBCs and immune effectors [32].

One of the earliest models used to describe intra-host dynamics of malaria infection is thatof Anderson et al. [1]. The mathematical model is given as
X 2 bXðtÞMðtÞ 2 mX XðtÞ;
dYðtÞ ¼ bXðtÞMðtÞ 2 m
Y Y ðtÞ 2 mMMðtÞ 2 bXðtÞMðtÞ;
where X(t) and Y(t) are the concentrations of uninfected RBCs and infected red blood cells(IRBCs) respectively. M(t) are the merozoites (parasites) that infect RBCs. lX is the sourceof RBCs from the bone marrow, b is rate of infection of RBCs, mX is the death rate ofRBCs, mY is the death rate of IRBCs, r is the average number of merozoites produced pereach bursting IRBC and mM is the death rate of merozoites.

This model does not take into account the effect of the immune system and therefore
describes the worst case scenario. The intra-host basic reproductive number for a simplemalaria model is defined as the number of merozoites (r), released per IRBC £ the initialfraction of successful merozoites [32]. This definition gives the basic reproductive numberfor model (1) as R0 ¼ rbXN=ðmM þ bXNÞ. The reproductive number considers the fate ofa single productively infected cell in an otherwise healthy individual with normal targetcell levels XN ¼ lX=mX [6]. In a healthy individual the target cell population is regulatedaccording to the equation dXðtÞ=dt ¼ lX 2 mXXðtÞ and homeostasis is maintained at somesteady state XN ¼ lX=mX . 0. The disease free equilibrium of model system (1) caneasily be shown to be globally stable. We will show in Appendix A that the endemicequilibrium of system (1) is globally asymptotically stable. We note from Appendix A thatthe endemic equilibrium of (1) is globally stable when blX , minðmX; mY ÞmY . Thiscondition is possible for a model without immune response such as (1), where the deathrate of IRBCs can be very high.

Intra-host model with immune response
We build our model from the model system (1) by including the response of theimmune system. In developing the model, we also make the following assumptions:
(1) The model assumes five interacting populations at any given time (t). These
are RBCs X(t), IRBCs Y(t), merozoites M(t), immune cells B(t) and antibodiesA(t).

C. Chiyaka et al.

(2) Due to the many different types of immune cells elicited by the presence of
malaria parasites, the model assumes them to be lumped together in onepopulation for simplicity.

(3) The RBCs are supplied from the bone marrow at a constant rate. We assume
that the supply rate is accelerated by the presence of IRBCs. However, theparticular mechanisms involved in the acceleration of RBCs during a malariainfection are still poorly understood [46]. They are reduced through infectionby merozoites and natural death at a constant rate. The RBCs are alsodestroyed through phagocytosis of erythrocytes bound to merozoites.

(4) IRBCs die at a constant rate and are also killed by the presence of immune
effectors. They produce free merozoites by bursting.

(5) Production rate of merozoites is reduced by immune cells. These free parasites
suffer a natural death, are eliminated from circulation by immune cells and arealso reduced through infecting RBCs.

(6) Immune cells are assumed to have a per capita rate of production and their
production is stimulated by the presence of IRBCs and merozoites. They arereduced (deactivated) at a constant rate.

(7) Antibodies that block invasion of RBCs proliferate in the presence of
merozoites. Proliferation of antibodies is dependent only on merozoites sincewe assume only antibodies that block infection. Antibodies decay at a constantrate.

(8) All parameters are positive.

These assumptions lead to the following system of differential equations which
describe the interaction between uninfected RBCs, IRBCs, merozoites, immune cells(which include CD4þ T cells, dendritic cells, macrophages, etc.), and antibodies:
¼ lX þ sYðtÞ 2
2 mXXðtÞ 2 vXðtÞMðtÞBðtÞ;
2 mY YðtÞ 2 kY BðtÞYðtÞ;
2 mMMðtÞ 2 kMBðtÞMðtÞ 2
¼ lB þ BðtÞ rY
The dynamics of RBCs, X(t) is described by Equation (2.1). The first term on the right
side of Equation (2.1) represents the rate of supply of RBCs from the bone marrow [31].

The second term represents the recruitment of RBCs at a rate s proportional to IRBCswhere 0 # s , 1. RBC production is accelerated during a malaria infection, but theparticular mechanisms involved are still poorly understood [46]. One analysis reported
Computational and Mathematical Methods in Medicine
an average of 37% increase in RBC production in adult first-time P. falciparum patients[23]. The third term represents infection of RBCs by merozoites. bXðtÞMðtÞ models therate at which free merozoites infect RBCs where b is rate of infection. Antibodies specificto malaria parasites (merozoites) inhibit invasion of erythrocytes by merozoites [42,2].

The term f ðAðtÞ; c0Þ ¼ 1=ð1 þ c0AðtÞÞ represents the role of antibodies in controllingparasitaemia, where c0 is the efficiency of antibodies in reducing erythrocytic invasion.

Antibodies bind to the epitopes on the Plasmodium parasite, thereby disabling the parasite,making its entry into the RBC difficult. As the antibody level increases (A(t) ! 1),f ðAðtÞ; c0Þ ! 0. This means that an increase in the number of antibodies reduces the rate ofinfection of RBCs by merozoites. Conversely, low numbers of antibodies increase the rateof infection of RBCs by merozoites. That is as AðtÞ ! 0, f ðAðtÞ; c0Þ ! 1. The fourth term inEquation (2.1) is the natural death rate of RBCs, since their average life span is 1=mX.

The last term describes the destruction of RBCs through phagocytosis of erythrocytesbound to merozoites [23]. Equation (2.2) models the dynamics of infected (parasitized)RBCs, Y(t). The first term on the right hand side of Equation (2.2) is a gain term for IRBCsfrom the loss term in Equation (2.1). The second term is the rate at which IRBCs are lostthrough death. The last term represents the rate at which IRBCs are killed directly byimmune cells [12,42]. The parameter kY represents immunosensitivity of IRBCs [32].

The rate of change of merozoites is described by Equation (2.3). The first term on the right
hand side of Equation (2.3) is the source of free merozoites which are released when IRBCsburst [7,22]. Their source is the second term in Equation (2.2). IRBCs burst during death,hence the number of merozoites produced depends on the death rate of IRBCs. An average of rmerozoites are released per each bursting IRBC. Rate of parasite production (multiplication)by infected cells is suppressed by immune cells with a factor f ðBðtÞ; c1Þ ¼ 1=ð1 þ c1BðtÞÞ.

Immune cells such as gd þ T cells expand during the early stages of infection by P. falciparum[10]. The expanded or activated T-cell populations express mRNA for TNF and IFN-g [43]and can inhibit parasite growth [43,10]. c1 is the efficiency of immune cells in suppressingparasite production. So as B(t) ! 1, f ðBðtÞ; c1Þ ! 0 meaning as immune cells increase, rateof parasite production is reduced. Conversely, as immune cells decrease, parasite productionby IRBCs increases, thus as B(t) ! 0, f ðBðtÞ; c1Þ ! 1. We note however that for malaria, adynamic equilibrium is required between pro-inflammatory and anti-inflammatory cytokinesto suppress immunopathology. The second term is the loss of merozoites through naturaldeath and the third term is the rate at which merozoites are removed by immune cells. IFN-gand CD4þ T cells activate macrophages to phagocytose intra-erythrocytic parasites and freemerozoites [42]. kM is immunosensitivity of merozoites. The last term is the loss ofmerozoites through infection of RBCs.

Equation (2.4) models the population of the immune cells (e.g. macrophages, natural
killer cells, CD4þ T cells). The first term on the right hand side of Equation (2.4) is thesource term for immune cells (combined term for different immune cells). This is followedby the stimulation term for immune cells in the presence of IRBCs and merozoites. rY andrM represent the immunogenicity of IRBCs and merozoites respectively [32]. k0 is thedensity of IRBCs at which immune cells grow at rate rY/2, that is half their maximum rate[13] in the absence of merozoites and k1 is the density of merozoites [25] at which immunecells grow at rate rM/2 in the absence of IRBCs. mB is the natural death rate of the immunecells. The population of antibodies is described by Equation (2.5). The first term on theright hand side of Equation (2.5) is the stimulation term for antibodies. Antibodies thatinhibit invasion of RBCs by merozoites are secreted by immune cells in the presenceof merozoites. These antibodies against the merozoite surface proteins (MSPs) are a majorcomponent of the invasion inhibitory response in humans [36]. h is the maximum rate
C. Chiyaka et al.

of increase of antibodies, k1 is the density of merozoites at which antibodies reach halftheir maximum value (h/2) and mA is the rate at which the antibodies decay.

Unlike Equation (2.4) where we include a natural source term, lB for immune cells, we
do not include a source term for antibodies because we assume that the only source ofantibodies against merozoites in a non-immune human (a human who has no previousexposure to malaria infection) are immune cells.

The solutions of the model (2) are shown to be always positive, which means that the
model is biologically feasible. Positivity of solutions is shown in Appendix B. It can easilybe shown that all solutions of system (2) initiating in IR5þ are bounded and eventually enterthe attracting set
Analysis of the model
In model system (2), some parameters have biological meaning and can be well estimatedand others are mathematical constructs imposed by our understanding of the immunesystem. To examine such a model for its qualitative behaviour and estimate predictedoutcome we examine the following:
(a) the system's equilibria and/or conditions for equilibria,
(b) the stability of the equilibrium points,
(c) the effect of treatment.

Intra-host basic reproductive number and stability analysis
The intra-host basic reproductive number R0 of the malaria parasite is defined as thenumber of secondary IRBCs produced per primary IRBC in a host at the onset of infection[32]. If R0 , 1, then on average an IRBC produces less than one new IRBC and theinfection cannot grow. If however R0 . 1, then on average each IRBC produces more thanone new IRBC and infection is maintained. Following the method of the next generatorapproach [44], the matrices F and V 21 are given as:
mY þkY ðlB=mBÞ
The product FV 21 is called the next generation matrix [8]. The reproductive number
R0, is the dominant eigenvalue of FV 21. For the model system (2), we get
ðmB þ c1lBÞðmBmY þ kY lBÞðmBmMmX þ mXkMlB þ mBblXÞ
Computational and Mathematical Methods in Medicine
The expression for R0 in (3.1) can be written as
Þ21 £ ðmY þ kYBÞ21;
where the first term in (3.2) is the rate at which a merozoite introduced into a completelysusceptible RBC population (where X ¼
A) infects a RBC. The second
term is the number of merozoites produced by the IRBC when it bursts. The third term isthe life span of each produced merozoite and the last term is the life span of the IRBC.

The model system given by Equation (2) has two steady states, the disease-free state
and the endemic state denoted by
To bring the infection (parasitaemia) under control we seek conditions on the
parameters of the transmission process that will guarantee the existence of a stable disease-free state. These conditions are summarized in Appendix C and Appendix D.

The endemic equilibrium and its stability
When R0 . 1, the condition for the stability of the disease-free equilibrium is violated andbesides the disease-free equilibrium, the system of Equation (2) has an endemicequilibrium point, Ee ¼ ðX *; Y *; M *; B*; A*Þ. The explicit form of the endemicequilibrium is quite cumbersome to be produced here, therefore we shall only show itsexistence numerically for a certain range of parameter values and perform its stabilityanalysis. Numerical simulations in Section 5 also confirm the existence of an endemicequilibrium for the parameters used. We show numerically that when R0 , 1, the diseasefree state is stable and when R0 . 1, the system attains an endemic equilibrium which isstable and the disease-free state becomes unstable. The bifurcation diagram fromnumerical simulations to illustrate this is shown in Figure 1 and it was produced by usingthe fourth order Runge Kutta method coded in the C þ þ programming language. It isobtained by varying r, the average number of merozoites produced per each burstingIRBC, while the other parameters are fixed. This shows that there exists a critical value,
ðmB þ c1lBÞðmBmY þ kYlBÞðmBmMmX þ mXkMlB þ mBblXÞ
above which the infection will persist. This is obtained when R ¼
1 in Equation (3.1).

Therefore the uninfected steady state is stable if r , r c.

Intra-host model with immune response and drug therapy
We extend the model to incorporate effects of the drug. Antimalarial drugs takenprophylactically or during infection (blood schizonticides) concentrate particularly inparasitized erythrocytes. The drug diffuses into parasite lysosomal compartments andbecomes protonated in the acidic environment within, so it cannot pass out through the
C. Chiyaka et al.

Bifurcation diagram showing the variation in Y with R0 as bifurcation parameter.

Stability is shown by bold lines and dashed lines indicate unstable equilibrium. The disease freeequilibrium is stable when R0 , 1 and unstable otherwise. The endemic equilibrium exists and isstable when R0 . 1. The bifurcation diagram was obtained by varying r, the average number ofmerozoites produced per IRBC while keeping all the other parameters fixed at numerical valuesgiven as l ¼
41664, s ¼ 0.009, b ¼ 7, c0
h ¼ 0.8 and v ¼ 0.0025.

membrane. It raises the pH of lysosome, inhibiting the polymerase that converts toxic freehaem to a harmless by-product. It prevents digestion of haemoglobin by parasites,reducing its supply of amino acids and therefore makes the parasite survival anddevelopment difficult. If the drug is administered, then the burst size r becomes (1 2 g)rwhere g is drug efficacy and is assumed to lie between zero, meaning totally ineffective,and one, meaning 100% effectiveness. The drug efficacy in this context is the probabilitywith which chloroquine inhibits parasite growth inside IRBCs. The model system thenbecomes
¼ lX þ sYðtÞ 2
2 mXXðtÞ 2 vXðtÞMðtÞBðtÞ;
2 mY YðtÞ 2 kY BðtÞYðtÞ;
Computational and Mathematical Methods in Medicine
Y Y ðtÞ 2 mMMðtÞ 2 kMBðtÞMðtÞ 2
¼ lB þ BðtÞ rY
The Equations (4.1), (4.2), (4.4) and (4.5) are the same as (2.1), (2.2), (2.4) and (2.5)
respectively. Equation (4.3) is a modification of Equation (2.3) to incorporate the effects ofantimalarial drugs.

The model system (4) has two equilibrium states, the disease free state
E 0 ¼ ðX; Y; M; B; AÞ ¼
and the endemic equilibrium E e ¼ ð
AÞ; and due to the complexity of our model,
the endemic equilibrium is shown numerically so that it exists for the parameter valuesused.

For the disease-free equilibrium to be locally stable all the eigenvalues of the Jacobean
matrix evaluated at the disease-free state should be negative. The eigenvalues correspondto the roots of the characteristic Equation jJ 2 zIj ¼ 0, where J is the Jacobian matrix, z isthe eigenvalue and I is a 5 £ 5 identity matrix. The characteristic equation is:
X 2 zÞð2mB 2 zÞð2mA 2 zÞ ¼ 0:
X ð1 þ c1ðlB=mBÞÞ
Clearly three of the roots, 2mX, 2mB and 2mA, are negative. We use Routh – Hurwitz
condition to establish that all roots are negative if the condition
is satisfied. We then deduce the reproductive number Rg when antimalarial drugs areadministered to be
ð1 2 gÞrmY blXm3
B þ c1lBÞðmBmY þ kY lBÞðmBmMmX þ mX kMlB þ mBblX Þ
Furthermore, it can be seen that at least one of these eigenvalues has a positive real part
if Rg . 1. Thus we have established the following result:
Lemma. The disease free equilibrium E 0 of (4) is locally asymptotically stable if Rg , 1and unstable if Rg . 1.

C. Chiyaka et al.

Implications for control
Knowledge of factors that limit parasite numbers offers hope of better-designed treatmentand intervention strategies [15,30,32], as well as providing information on selective forcesthat have moulded parasite life-history strategies [20]. The intra-host basic reproductivenumber R0 is a key parameter of asexual parasitaemia, crucial to calculations concerningits control by any mechanism, natural or artificial. From the formula of R0 in (3.1), wededuce the following: (a) reduction in the transmission rate b reduces the reproductivenumber. (b) Increasing the rate at which immune cells suppress parasite production ismore effective in the control of the disease than increasing the rate at which antibodiesinhibit invasion of parasites into erythrocytes, since the parameter c1, efficiency ofimmune cells appears in R0 and c0, efficiency of antibodies does not. (c) Since mB is acubed term in the expression for R0, increasing mB greatly increases R0.

Biologically speaking, Rg measures the number of secondary infections generated by a
single parasitized RBC in an environment where antimalarial drugs are used as a controlstrategy. The primary focus of drug therapy is on the possibility of clearing the parasites.

If R
g , 1, then parasites are cleared. Using the critical point Rg
1, we find the critical
drug efficacy g c that is required for parasite clearance to be
R0 is the reproductive number when there is immune response only. If immune
response fails to clear the parasites (R0 . 1), the drug is then administered, with the aim ofreducing R0. Rg is always less than R0 for all values of g in the range 0 , g , 1. Rg ¼ R0only when g ¼ 0. This implies that the drug is completely ineffective. Rg , 1 wheneverg . g c, thus malaria parasites can be eradicated from an infected individual if g . g c.

Numerical simulation
To observe the dynamics of the model over time we integrated the system of equations,using the fourth order Runge – Kutta methods in the C þ þ programming language. Forcomputer runs we set the initial densities of RBCs, IRBCs, merozoites, immune cells andantibodies at 500, 5, 50, 30, and 10 respectively. Time is in days. The fixed parametersused are shown in Table 1. The value of the reproductive number R ¼
Choosing parameter values characteristics of the in vivo situation is difficult. Many of
the parameters in our model, as shown in Table 1, are estimated. Even those that have beenquoted may not be as accurate as we need for quantitative predictions. However, our mainthrust is on the overall effect they have on the basic reproductive number, which is thebasis for an infection to persist or be eradicated.

The graphs in Figure 2 show the behaviour of model system (2). It shows that the
populations of uninfected RBCs, IRBCs and merozoites reach a steady endemic state sinceR0 . 1. Figure 3 shows the behaviour of model system (4) where g ¼ 0.97. The graphsshow that the infection can be eradicated if the efficacy of the drug g . g c. If the immunesystem alone fails to clear the infection then a drug of suitable efficacy can accomplish thetask of bringing the infection to disease-free levels. Figure 4 shows that as the death rate ofimmune cells, mB, increases there is an increase in the number of IRBCs and merozoites,while the uninfected RBCs decrease.

Figure 5 shows the behaviour of (a) RBCs (b) uninfected RBCs and (c) merozoites of
model system (4.1) – (4.5) as the rate at which RBCs are eliminated due to the effect
Computational and Mathematical Methods in Medicine
Table showing numerical values of parameters used in the simulations performed where c
are cells, d are days and m are merozoites.

Supply rate of RBC
Rate of recruitment of RBC
Rate of infection
Efficiency of antibodies
Death rate of RBCs
Death rate of IRBCs
Immunosensitivity of IRBC
Rate at which RBCs are eliminated
Merozoites released per each bursting
Rate at which parasite production
Death rate of merozoites
Supply rate of immune cells
Immunogenicity of IRBCs
Immunogenecity of merozoites
Stimulation constant for immune cells
Stimulation constant for immune cells
due to merozoites
Death rate of immune cells
Immunosensitivity of merozoites
Deterioration rate of antibodies
Maximum rate of increase of antibodies
of involvement v is varied from 0.0 to 0.08 in steps of 0.02. The arrows show the directionof increase of v. The graphs show that as v increases, there is a corresponding decrease inthe populations of RBCs, IRBCs and antibodies. A decrease in RBC population gives acorresponding decrease in the number of IRBCs which makes the number of immune cellsand antibodies stimulated decrease as well.

A model of the immune response to Plasmodium falciparum infection was developed andextended to include treatment with antimalarial drugs. As a first step towards more realisticmodelling of immune regulation, we introduce the terms for the effects of immune effectorson merozoite invasion of RBCs and parasite production within a RBC. In addition ,we allowthe RBCs to be eliminated through the effect of involvement [35]. Our mathematicalanalysis yielded a generalization of the intra-host basic reproductive number from which themost effective control strategies that should aid in assessing interventions (drugs andvaccines) are deduced. Our results show that without any treatment the most effective part ofthe immune response in its mission to clear parasites is its ability to inhibit parasite growth inerythrocytes. It is also more effective to increase the death rate of IRBCs than the death rateof merozoites by immune cells. This is because one IRBC is capable of producing about 16merozoites, so killing it will prevent the survival of 16 merozoites. Most of the availablemalaria drugs act by retarding development of parasites. Drugs that kill IRBCs directly,together with those that retard the development of parasites, should also be made availableto make the treatment more effective. For the extended model, with treatment, we have
C. Chiyaka et al.

Graphs that show the behaviour of the model with immune response only. The dynamics
of (a) uninfected RBCs, (b) IRBCs, (c) merozoites, (d) immune cells and (e) antibodies.

Computational and Mathematical Methods in Medicine
Graphs that show the behaviour of the model with immune response and drug therapy.

The dynamics of (a) uninfected RBCs, (b) IRBCs, (c) merozoites, (d) immune cells and (e)antibodies.

C. Chiyaka et al.

Graphs that show the behaviour of (a) uninfected RBCs and (b) merozoites. They were
obtained by varying the death rate of immune cells mB from 0.5 to 2.5 in steps of 0.5, while keepingthe other parameters constant. The arrows show the direction of increase of mB.

shown that parasites can only be cleared from an infected individual if g . g c. This impliesthat for a very large reproductive number, a drug with an efficacy of g < 1 is required to treatthe infection. This might be one of the reasons why for a person from a non endemic area theinfection is more virulent and, if treatment is delayed, then death might be inevitable.

We have also observed that if there are mechanisms in the body that increase the death
rate of immune cells then there is an increase in the severity of the infection. This shows thatfor a person infected by pathogens that kill some of the immune cells, like the humanimmunodeficiency virus which kills the CD4þ T cells [41], or Mycobacterium tuberculosiswhich kills the macrophages, then that person is infected by P. falciparum parasite, theimmune system is compromised and the person suffers severe malaria if not treated in time.

People who suffer from diseases such as acquired immunodeficiency syndrome ortuberculosis should not delay in seeking medical attention when symptoms of malaria startto appear. Analysis also yielded a well-known result that the effect of involvement decreasesthe number of RBCs, which is believed to be one of the causes of anaemia in malaria.

Our mathematical analysis is useful in explaining some of the observed patterns during
a malaria infection. This gives our model the ability to provide a basic representation ofthe complex web that accompanies an intruding malaria parasite during the erythrocyticstages. For people moving to endemic areas, prophylactic drugs should be taken a weekbefore departure so that an adequate concentration is attained in the blood before exposureto parasites and continued a month after return to kill any parasites incubating in the liveras they emerge.

Drug resistance has been a major cause of resurgence of the disease and will be
considered elsewhere. The complexity of the parasite, the different forms it assumes in thebody and the possibility of the patient being infected by several strains concurrently makemodel predictions more difficult. Sporozoites are a target for vaccine development but areable to replace their coats, varying the antigens they present. Relevant data should be madeavailable to be used in model validation. Thus one role of modelling is to point out wherefurther quantitative measurements can improve our understanding of the malaria diseaseprocess.

Computational and Mathematical Methods in Medicine
Graphs that show the behaviour of (a) uninfected RBCs, (b) IRBCs and (c) antibodies.

The graphs are obtained by varying the rate at which uninfected RBCs are eliminated, v from 0.0 to0.08 in steps of 0.02, while keeping the other parameters constant. The arrows show the direction ofincrease of v.

The authors would like to thank the two anonymous referees for their constructive comments whichhave greatly improved the paper. The authors also acknowledge financial support from EagleInsurance Company, Zimbabwe. C. Chiyaka would like to acknowledge financial support byNational University of Science and Technology through a Staff Development Scholarship. We alsowant to thank Professors Hagai Ginsburg and Klaus Dietz for sending us literature on mathematicalmodelling of within-host dynamics of malaria.

[1] R.M. Anderson, R.M. May, and S. Gupta, Non-linear phenomena in host – parasite interactions,
Parasitology 99 (1989), pp. S59 – S79.

[2] K. Artvanis-tsakonas, J.E. Tongren, and E.M. Riley, The war between the malaria parasite and the immune
system: Immunity, immunoregulation and immunopathology, Clin. Exp. Immunol. 133 (2003), pp. 145 – 152.

C. Chiyaka et al.

[3] N.T.J. Bailey, The biomathematics of malaria, in Malaria. Charles Griffin & Co. Ltd, London, 1982,
pp. 1 – 93.

[4] C.A. Bate, J. Tarvene, and J.H.L. Playfair, Malarial parasites induce TNF production by macrophages,
Immunology 64 (1988), pp. 227 – 231.

[5] I.A. Clark, J.L. Virelizier, E.A. Carswell, and P.R. Wood, Possible importance of macrophage-derived
mediators in acute malaria, Infect. Immun. 32 (1981), pp. 1058 – 1066.

[6] P. De Leenheer and H.L. Smith, Virus dynamics: A global analysis, SIAM J. Appl. Math. 63 (2003),
pp. 1313 – 1327.

[7] H.H. Diebner, M. Eincher, L. Molineaux, W.E. Collins, G.M. Jeffert, and K. Dietz, Modelling the transition
of asexual blood stages of Plasmodium falciparum to gametocytes, J. Theor. Biol. 202 (2000), pp. 113 – 127.

[8] O. Diekmann, J.A.P. Heesterbeek, and J.A.J. Metz, On the definition and computation of the basic
reproduction ratio R0 in models for infectious diseases in heterogeneous population, J. Math. Biol. 28 (1990),pp. 365 – 382.

[9] M.M. Elloso, H.C. van der Heyde, J.A. vande Waa, D.D. Manning, and W.P. Weidanz, Inhibition of
Plasmodium falciparum in vitro by human gamma delta T cells, J. Immunol. 153 (1994), pp. 1187 – 1194.

[10] C.R. Engwerda, and M.F. Good, Interactions between malaria parasites and the host immune system, Curr.

Opin. Immunol. 17 (2005), pp. 381 – 387.

[11] M. Fan, M.Y. Li, and K. Wang, Global stability of an SEIS epidemic model with recruitment and a varying
total population size, Math. Biosci. 170 (2001), pp. 199 – 208.

[12] G. Fritsche, C. Larcher, H. Schennach, and G. Weiss, Regulatory interactions between iron and nitric oxide
metabolism for immune defense against Plasmodium falciparum infection, J. Infect. Dis. 183 (2001),pp. 1388 – 1394.

[13] V.V. Ganusov, C.T. Bergstrom, and R. Antia, Within-host dynamics and the evolution of microparasites in a
heterogeneous host population, Evolution 56(2) (2002), pp. 213 – 223.

[14] M.F. Good, H. Xu, M. Wykes, and C.R. Engwerda, Development and regulation of cell-mediated immune
responses to the blood stages of malaria: Implications for vaccine research, Annu. Rev. Immunol. 23 (2005),pp. 69 – 99.

[15] M.B. Gravenor, and D. Kwiatkowski, An analysis of the temperature effects of fever on the intra-host
population dynamics of Plasmodium falciparum, Parasitology 117 (1998), pp. 97 – 105.

[16] M.B. Gravenor, A.R. McLean, and D. Kwiatkowski, The regulation of malaria parasitaemia: Parameter
estimates for a population model, Parasitology 110 (1995), pp. 115 – 122.

[17] M.B. Gravenor, A.L. Lloyd, P.G. Kremsner, M.A. Missinou, M. English, K. Marsh, and D. Kwiatkowski,
A model for estimating total parasite load in Falciparum malaria patients, J. Theor. Biol. 217 (2002),pp. 137 – 148.

[18] D. Gurarie and F.E. McKenzie, Dynamics of immune response and drug resistance in malaria infection,
Malar. J. 5 (2006), p. 86.

[19] D. Gurarie, P.A. Zimmerman, and C.H. King, Dynamic regulation of single- and mixed-species malaria
infection: Insights to specific and non-specific mechanisms of control, J. Theor. Biol. 240 (2006),pp. 185 – 199.

[20] D.T. Haydon, L. Matthews, R. Timms, and N. Colegrave, Top-down or bottom-up regulation of intra-host
blood stage malaria: Do malaria parasites most resemble the dynamics of prey or predator, Proc. R. Soc.

Lond. 270 (2003), pp. 289 – 298.

[21] C. Hertzel and R.M. Anderson, The within-host dynamics of blood stage malaria: Theoretical and
experimental studies, Parasitology 113 (1996), pp. 25 – 38.

[22] M.B. Hoshen, R. Heinrich, W.D. Stein, and H. Ginsburg, Mathematical modelling of the within-host
dynamics of Plasmodium falciparum, Parasitology 121 (2000), pp. 227 – 235.

[23] G.N. Jakeman, A. Saul, W.L. Hogarth, and W.E. Collins, Anaemia of acute malaria infections in non-
immune patients primarily results from destruction of uninfected erythrocytes, Parasitology 119 (1999),pp. 127 – 133.

[24] K. Kirk, Membrane transport in the malaria-infected erythrocyte, Physiol. Rev. 81 (2001), pp. 495 – 537.

[25] R.J. Kosinski, Antigenic variation in trypanosomes: A computer analysis of variant order, Parasitology
80 (1980), pp. 343 – 357.

[26] D. Kwiatkowski, Malaria toxins and the regulation of parasite density, Parasitol. Today 11 (1995),
pp. 206 – 212.

[27] J.P. LaSalle, The Stability of Dynamical Systems. SIAM, Philadelphia, PA, 1976.

[28] W. Ma, Y. Takeuchi, T. Hara, and E. Beretta, Permanence of an SIR epidemic model with distributed time
delays, Tokohu Math. J. 54 (2002), pp. 581 – 591.

[29] D.P. Mason and F.E. McKenzie, Blood-stage dynamics and clinical implications of mixed Plasmodium
vivax – Plasmodium falciparum infections, Am. J. Trop. Med. Hyg. 61 (1999), pp. 367 – 374.

[30] F.E. McKenzie and W.H. Bossert, The dynamics of Plasmodium falciparum blood stage infection, J. Theor.

Biol. 188 (1997), pp. 127 – 140.

[31] P.G. McQueen and F.E. McKenzie, Age-structured red blood cell susceptibility and the dynamics of malaria
infections, Proc. Natl. Acad. Sci. USA 101 (2004), pp. 9161 – 9166.

[32] L. Molineaux and K. Dietz, Review of intra-host models of malaria, Parassitologia 41 (1999), pp. 221 – 231.

Computational and Mathematical Methods in Medicine
[33] L. Molineaux, H.H. Diebner, M. Eincher, W.E. Collins, G.M. Jeffery, and K. Dietz, Plasmodium falciparum
parasitaemia described by a new mathematical model, Parasitology 122 (2001), pp. 379 – 391.

[34] J.S. Muldowney, Compound matrices and ordinary differential equations, Rocky Mount. J. Math. 20 (1990),
pp. 857 – 872.

[35] A. Murase, T. Sasaki, and T. Kagiwara, Stability analysis of pathogen – immune interaction dynamics,
J. Math. Biol. 51 (2005), pp. 247 – 267.

[36] R.A. O'Donnell, T.F. de Koning-Ward, R.A. Burt, M. Bockarie, J.C. Reeder, A.F. Cowman, and B.S. Crabb,
Antibodies against merozoite surface protein (MSP)-1(19) are a major component of the invasion-inhibitoryresponse in individuals immune to malaria, J. Exp. Med. 193 (2001), pp. 1403 – 1412, (247 – 263).

[37] P. Perlmann and M. Troye-Blomberg, Malaria and the immune system in humans, Malaria Immunol.

80 (2002), pp. 229 – 242.

[38] K.P. Piper, R.E. Hayward, M.J. Cox, and K.P. Day, Malaria transmission and naturally acquired immunity
to PfEMP-1, Infect. Immun. 67 (1999), pp. 6369 – 6374.

[39] I.M. Rouzine and F.E. McKenzie, Link between immune response and parasite synchronization in malaria,
Proc. Natl. Acad. Sci. USA 100 (2003), pp. 3473 – 3478.

[40] A. Saul, Models for the in-host dynamics of malaria re-visited: Errors in some basic models lead to large
overestimates of growth rate, Parasitology 117 (1998), pp. 405 – 407.

[41] T. Shiri, W. Garira, and S.D. Musekwa, A two strain HIV-1 mathematical model to assess the effects of
chemotherapy on disease parameters, J. Math. Biosci. Eng. 2 (2005), pp. 811 – 832.

[42] M.M. Stevenson and E.M. Riley, Innate immunity to malaria, Nat. Rev. Immunol. 4 (2004), pp. 169 – 180.

[43] M. Troye-Blomberg, S. Worku, P. Tangteerawatana, R. Jamshaid, K. Soderstrom, G. Elghazali, L. Moretta,
M. Hammarstrom, and L. Mincheva-Nilsson, Human gamma delta T cells that inhibit the in vitro growth ofthe sexual blood stages of the Plasmodium falciparum parasite express cytolytic and proinflammatorymolecules, Scand. J. Immunol. 50 (1999), pp. 642 – 650.

[44] P. van den Driessche and J. Watmough, Reproduction numbers and subthreshold endemic equilibria for
compartmental models of disease transmission, Math. Biosci. 180 (2002), pp. 29 – 48.

[45] N.J. White, Malaria pathophysiology, in Malaria, Parasite Biology, Pathogenesis, and Protection,
I.W. Sherman, ed., ASM Press, Washington, DC, 1998, pp. 371 – 385.

[46] S.N. Wickramasinghe, and S.H. Abdalla, Blood and bone marrow changes in malaria, Bailliers Best Pract.

Res. Clin. Haematol. 13 (2000), pp. 277 – 299.

Appendix A: Global asymptotic stability of the endemic equilibrium state
Let the endemic equilibrium state of model system (1) be denoted by E e ¼ X e; Y e; M e.

A.1. Suppose that R0 . 1 and blX 2 minðmX; mY ÞmY , 0, then E e is a globally
asymptotically stable steady state for system (1) with respect to initial conditions not on theboundary.

The key to verifying the global stability of E e is to rule out the existence of periodic solutions
[11]. This is achieved by showing that any periodic solution to (1) is orbitally asymptotically stable.

Let pðtÞ ; ðp1ðtÞ; p2ðtÞ; p3ðtÞÞT denote the periodic solution and suppose that its minimal period isv . 0.

Theorem A.2. A sufficient condition for a periodic orbit q ¼ {pðtÞ : 0 # t , v} of (1) to beasymptotically stable with asymptotic phase is that the periodic linear system
is asymptotically stable.

›f ½2=›x is the second compound matrix of system (1). (For detailed discussions of compound
matrices we refer the reader to [34]). We will show that
Vðz1; z2; z3; pðtÞÞ ¼ sup
ðjz2ðtÞj þ jz3ðtÞjÞ
is a Lyapunov function for system (1). The right-hand derivative of VðtÞ, DþV, exists. Directcalculations lead to the following differential inequalities
C. Chiyaka et al.

p1ðtÞp3ðtÞ p2ðtÞ
j 1ðtÞjÞ # 2ðbp3ðtÞ þ mX þ mY Þ z
rmY þbp3ðtÞÞ jz1ðtÞj2
2ðtÞðrmY þbp3ðtÞÞ jz1ðtÞj
M 2 bp1ðtÞ 2 minðmX ;mY Þ
Relations (A2) and (A3) lead to
DþV # sup {g1ðtÞ; g2ðtÞ};
1ðtÞ ¼ 2ðmX þ bp3ðtÞÞ þ
; g2ðtÞ ¼ bp2ðtÞ þ
2 minðmX; mY Þ:
It follows from (A5) that g1ðtÞ # 2mX þ p =p
2ðtÞ, and that g1ðtÞ # g2ðtÞ, then we have
DþVðtÞ # g2ðtÞ:
then it will follow from (A6) that V is a Lyapunov function for system (1) and this will conclude theproof of the theorem. Since pðtÞ is a solution of (1), we see that
bp1ðtÞp3ðtÞdt ¼
ðlX 2 mXp1ðtÞÞdt # lXv:
ðbp2ðtÞ 2 min ðmX; mY ÞÞdt # b
2 minðmX; mY Þ v:
It follows that (A7) holds under the assumption that blX 2 minðmX; mY ÞmY , 0 holds.

Computational and Mathematical Methods in Medicine
Appendix B: Positivity of solutions
Theorem B.1. Let the initial data be Xð0Þ ¼ X0 . 0, Y ð0Þ ¼ Y 0 . 0, Mð0Þ ¼ M0 . 0, Bð0Þ ¼B0 . 0 and Að0Þ ¼ A0 . 0. Then solutions ðXðtÞ; YðtÞ; MðtÞ; BðtÞ; AðtÞÞ of (2) are always positive forany t . 0.

Proof. We will perform the proof following ideas by [28]. It is easy to see that XðtÞ . 0 for all t . 0.

If not, we assume that there exists a first time tu such that XðtuÞ ¼ 0, X ðtuÞ # 0 and XðtÞ . 0,YðtÞ . 0, MðtÞ . 0, BðtÞ . 0, AðtÞ . 0 for 0 , t , tu. It follows from Equation (2.1) that we have
X ðtuÞ ¼ lX þ sYðtuÞ 2
2 mXXðtuÞ 2 vXðtuÞMðtuÞBðtuÞ . 0;
which is a contradiction. Similarly we can also prove by contradiction from Equation (2.4) thatdBðtÞ=dt ¼ lB . 0 ; t . 0 when there exists a first time tv such that BðtvÞ ¼ 0.

Therefore XðtÞ and BðtÞ are always positive.

Assume that there exists some time t1 . 0 such that Yðt1Þ ¼ 0, other variables are positive and
YðtÞ . 0 for t [ ½0; t1Þ. Integrating Equation (2.2) from 0 to t1 we have
Yðt1Þ ¼ Yð0Þ exp 2
ðmY þ kY BðtÞÞdt þ exp 2
ðmY þ kY BðtÞÞdt
ðmY þ kY BðuÞÞdu dt . 0:
which contradicts Yðt1Þ ¼ 0.

Assume there is some time t2 . 0 such that Aðt2Þ ¼ 0 and AðtÞ . 0. Then integrating Equation
(2.5) from 0 to t2, we have
ðt2 hBðtÞMðtÞ
Aðt2Þ ¼ Að0Þe 2mAt2 þ e 2mAt2
which contradicts Aðt2Þ ¼ 0. Similarly for MðtÞ, assume that there is some time t3 . 0 such thatMðt3Þ ¼ 0 and MðtÞ . 0. Then integrating Equation (2.3) from 0 to t3, we see that
Mðt3Þ ¼ Mð0Þ exp 2
mM þ kMBðtÞ þ
mM þ kMBðtÞ þ
mM þ kMBðuÞ þ
which contradicts Mðt3Þ ¼ 0.

This implies that the solution ðXðtÞ; YðtÞ; MðtÞ; BðtÞ; AðtÞÞ is always positive for t $ 0.

Appendix C: Local stability of the disease-free equilibrium
Theorem C.1. If R0 , 1, then the disease free state E0 is a locally asymptotically stable state ofsystem (2); if R0 . 1, then it is unstable.

C. Chiyaka et al.

Proof. To prove local stability of the disease-free equilibrium we need to show that all theeigenvalues of the Jacobian matrix of the system evaluated at the disease-free equilibrium arenegative. We linearize the system of Equations (2) at the disease-free state E0, and find the Jacobianmatrix. By inspection, three of the eigenvalues 2mX; 2mB; 2mA are easily deduced and theremaining two eigenvalues are obtained from the remaining 2 £ 2 submatrix
Expanding the characteristic equation of this submatrix we get the following characteristic
Using the Routh-Hurwitz stability criterion we determine that the two eigenvalues of the
submatrix (C1) which correspond to the roots of the characteristic Equation (C2) are negative if
Condition (C3) reduces to R0 , 1. Therefore if R0 , 1, then all the eigenvalues are negative
which implies that the disease free state is stable if R0 , 1.

Appendix D Global stability of the disease-free equilibrium
Theorem D.1. The disease-free equilibrium
lX=mX; 0; 0; lB=mB; 0 of system (2) is globally
asymptotically stable, whenever R0 , 1.

Proof. Consider the following Lyapunov function:
We note from system of Equations (2) that the minimum values of A and B are 0 and lB=mB
respectively and since 0 # s , 1, the maximum value of X ¼ lX=mX.

Computational and Mathematical Methods in Medicine
Therefore the Lyapunov derivative is given by
Since all the parameters of the model are non-negative, it follows that L , 0 for R0 , 1 with
L ¼ 0 if Y ¼ M ¼ 0. If M ¼ 0 ) A ¼ 0. Hence L is a Lyapunov function on F. Since F is invariant
and attracting, it follows that the largest compact invariant set in ððX; Y; M; B; AÞ [ F : L ¼ 0Þ isthe singleton {E0}. LaSalle's Invariance Principle [27] then implies that E0 is globallyasymptotically stable in F.

Source: https://www.emis.de/journals/HOA/CMMM/Volume9_2/849418.pdf

For our clients Preventive Care MedicationsTraditional Prescription Drug List (PDL)1,2,3,4 $0 Cost-share Medications & ProductsU.S. Preventive Services Task Force A & B Recommendation Medications and SupplementsThese medications and supplements will be covered at $0 cost-share when: Prescribed by a health care professional Age and/or gender appropriate Filled at a network pharmacy Most medications are available over-the-counter (OTC) except for prescription fluoride. All brands are

Brain Research 1010 (2004) 151 – 155 Stimulation of the superior cerebellar peduncle during the development of amygdaloid kindling in rats Carmen Rubioa, Vero´nica Custodioa, Francisco Jua´rezb, Carlos Paza,* a Departamento de Neurofisiologı´a, Instituto Nacional de Neurologı´a y Neurocirugı´a M.V.S., Insurgentes Sur 3877, Mexico 14269 D.F., Mexico b Instituto Nacional de Psiquiatrı´a R.F., Mexico