Draining the Pontine Marshes EEB/MCB 182H:
Spring 2006
Home Disinfection, Pontine Marshes
Modelling Malaria.
(Revised 04/25/06)


 
Table I. Variables and Parameters
Symbol Quantity Units
XInfected HumansHumans
YInfected MosquitosMosquitos
NTotal HumansHumans
MTotal MosquitosMosquitos
AMosquito Bite RateBites /
(Mosquito-Year)
BNew Human InfectionsHumans /
(Humans-Bites)
CNew Mosquito Infections Mosquitos /
(Humans-Bites)
mMosquito Death RateYear-1
rHuman Recovery RateYear-1

Ross-MacDonald Equations.

The simplest models [r11, r16, m52, m57] of malaria, the so-called Ross-MacDonald (R-M) equations, date to the early 20th century. These equations compress the manifestly complex interactions among man, mosquito and plasmodium into a pair of coupled differential equations that specify the time evolution of two variables: the numbers, X(t) and Y(t), of parasite-laden humans and female mosquitos. In particular, R-M assumes that transmission of the plasmodium - from mosquito to man and from man to mosquito - depends jointly on the numbers of susceptible and infected individuals of the appropriate species, which is to say, upon their product.

Let N and M be the total numbers of humans and female mosquitos. If parasite transfer is from mosquito to man, the product in question is Y(t) × [N - X(t)], which is the number of plasmodium-bearing mosquitos, Y(t), times the number of malaria-free (and by assumption susceptible) humans, [N - X(t)]. Conversely, the transfer rate from man to mosquito is presumed proportional to the number of infected humans, X(t), times the number of plasmodium-free mosquitos, [M - Y(t)]. This leads us to

dX/dt = A B Y (N - X) - r X                                                              
(1)
dY/dt = A C X (M - Y) - m Y                                                            

where the parameters, A, B, C, m and r are as defined in Table I. In particlar, A is the per mosquito bite rate.

Equations (1) assume constant total numbers, N and M, of humans and mosquitos. They further presume that infected humans recover (and become immediately susceptible) at rate, r, and that sporozite-bearing mosquitos die off at rate mi.e., both loss processes are assumed to be exponential.1   As to the constants, A, B and C, these relate transmission to biting by female mosquitos — the event by which transfer occurs — and the probabilities of parasite transfer from vector to man and man to vector. Note that the products, A × B and A × C, are respectively the per mosquito production rate of newly infected humans and the per human production rate of newly infected mosquitos.



Exercise 1. Verify that Equations (1) "balance," i.e., that the units on both sides of the equality sign are humans/year in the case of the first equation and mosquitos/year in the case of the second.
 



Replacing Absolute Abundance with Proportion.

Inasmuch as Equations (1) presume constant populations, they are customarily simplified by replacing absolute abundances of men and mosquitos with proportions. In this spirit, we introduce the following new variables and parameters:

  x = X/N y = Y/M q = M/N  
 (2)
  a = A C N b = B/C    

Substituting the new quantities into Equations (1) then yields

dx/dt = a b q y (1 - x) - r x                                                              
(3)
dy/dt = a x (1 - y) - m y                                                                  

which is the form most frequently encountered in the literature - see, for example, Aron and May [am82].2



Exercise 2. Verify that the units of Equations (3) - both equations, both sides - are y-1.
 



Isoclines and Equilibria.

Click to enlarge.
Click to enlarge.

Figure 1. Zero-growth isoclines of Equations (3). Stable (unstable) equilibria are shown as solid (hollow) circles.

One way of studying systems of differential equations is to integrate them numerically. Especially if nonlinearities are involved, this is often the only approach possible - i.e., most of the equations arising in mathematical biology do not admit to closed form solutions. In the present case, however, the equations' simplicity, allows us to take an alternative tack that yields geometric insight. Specifically, we study the zero-growth isoclines,

dx/dt = 0:   y = (r/abq)[x/(1-x)]                    
(4)
dy/dt = 0:   y = ax/(ax + m)                         

that are obtained by setting the right-hand sides of Equations (3) equal to 0 and solving for y as a function of x. We now plot these functions in the x-y plane, as shown in Figure 1. Note the essential points:

  1. Given parameter values, a, b, m, q and r, coordinate pairs (x, y) uniquely specify the vector field, F = [(dx/dt), (dy/dt)], that describes motion in the x-y plane.

  2. Plotting the zero-growth isoclines allows us to make qualitative inferences about the nature of F. This is because

    1. Points on the x-isocline, do not move in the x-direction; points on the y-isocline do not move in the y-direction.

    2. Points at which the isoclines intersect are equilibria.

    3. Points to the right of the x-isocline move left; points to the left of the x-isocline move right.

    4. Points above the y-isocline move down; points below the y-isocline move up.

Negative numbers of humans and mosquitos being biologically implausible, we restrict our attention to the plane's positive quadrant, (x ≥ 0; y ≥ 0). As shown in Figure 1, there are two possibilities:

  1. Top graph. The so-called no-disease state, (x, y) = (0, 0), is the sole, non-negative equilibrium, and it is stable.

  2. Bottom graph. There are two non-negative equilibria. One is the origin, which is now unstable; the other, the so-called endemic state, (x, y) = (x*, y*), is stable.

Which possibility obtains depends on the value of a quantity, R, called the basic reproductive rate of the disease. In the present case,3

R = a2 b q/m r .                                                           (5)

Click to enlarge.
Figure 2. Transcritical bifurcation (TC) in Equations (3).
In greater detail, if R < 1, we have the first possibility, and the disease dies out. Conversely, if R > 1, the disease persists, with the proportions of infected humans and mosquitos tending to their endemic values.

How one goes back and forth between these alternatives can be understood by observing the following:

  1. Figure 1 notwithstanding, Equations (3) always admit to the existence of trivial and non-trivial equilibria, the origin and (x*, y*).

  2. When R < 1, the non-trivial equilibrium corresponds to negative values of x and y, and is therefore not observed if we restrict our attention to non-negative solutions. Conversely, and as noted above, for R > 1, the non-trivial equilibrium corresponds to positive values of x and y.

  3. Exactly at R = 1, (x*, y*) = (0, 0), i.e., the non-trivial equilibrium is coincident with the origin.

It follows that as R is varied through R = 1, the equilibria collide, and there is an exchange of stability as shown in Figure 2. This qualitative change in dynamics is called a transcritical bifurcation.4



Exercises.
 
3. Verify that the non-trivial equilibrium of Equations (3) is given by

x* = (R-1)/[R + (a/m)],                                                                                   
(6)
y* = (a/m)(R-1)/{R[1 + (a/m)]}                                                                      

4. For the parameter values in Figure 1, compute the value of a at which the transcritical bifurcation occurs. Compare with Figure 2.

5. Equations (3) are "partially non-dimensionalized," i.e., only time has units. Complete the non-dimensionalization by replacing t, the units of which is years, with the dimensionless quantity, t = r t and by introducing new, dimensionless parameters, a = a/m and b = bmq/r. Write the resulting equations in the form

(dx/dt) = ...                                                                                                
(7)
e (dy/dt) = ...                                                                                                

where e = r/m. Why might e be small? If it is, what are the consequences?



Latency.

Click to enlarge.
Click to enlarge.
Figure 3. Simulated endemic malaria without (Top) and with (Bottom) the inclusion of latency.
Equations (1) and (3) are, of course, mere caricatures. Among the many biological details they omit is the fact that parasites undergo a period of development within the mosquito's gut before entering her salivary glands. In other words, there is a latent period, t, during which the mosquito is infected but not yet infectious. Significantly, t is on the order of 10 days, which is comparable to the life expectancy of adult mosqiuitos. As a result, there is a non-neglible probablity of an infected mosquito's dying before it can pass the disease back to a human host, a circumstance at considerable variance with the assumptions of Equations (1) and (3).

To incorporate latency into our model, we distinguish a third class of female mosquitos, z, which we will say are exposed, and replace our previous equation for (dy/dt) with two equations.5   We further note that the rate at mosquitos enter the infectious class is proportional to the product, x˜ (1- y˜ - z˜), of infected humans and susceptible mosquitos, t time units in the past, i.e., x˜ = x(t-t), etc. Moreover, this product must be multiplied by e-mt, which is the fraction of newly infected mosquitos that survive to become infectious. Accordingly, we write6

dx/dt = a b y (1 - x) - r x                                                                     
 
dz/dt = a x (1 - y - z) - a x˜ (1 - y˜ - z˜) e-mt - m z                              (8)
 
dy/dt = a x˜ (1 - y˜ - z˜) e-mt - m y                                                        

From these equations, we calculate

R˜ = [a2 b q /(r m)] e-mt                                                                (9)

where the symbol, R˜, is used to indicate the inclusion of latency in the model, and

x* = (R˜-1) / [R˜ + (a/m)]                                                               
 
y* = {[(a/m)(R˜-1)/R˜] / [1 + (a/m)]}e-mt                                    (10)
 
z* = [(1 - e-mt) / e-mt] y*                                                                   

From these exprerssions, we conclude that latency diminishes both the basic reproductive rate, = Re-mt, and the equilibrium density of infectious mosquitos, y*, by a factor of e-mt (Figure 3). As pointed out by Aron and May [am82], the latter change brings the model into better accord with observation - i.e., empirically, the fraction of female mosquitos with parasites in the salivary glands is on the order of a few per cent.


Implications.

The most important implication of the foregoing analyses, i.e., with and without the time delay, is that you don't have to kill all the mosquitos or immunize all of the human population to eradicate malaria. What is is necessary is to disrupt propagation to the point that R < 1. This having been accomplished, the disease will then die out of its own accord, as the system tends to the "no-disease" state at the now-stabilized origin.7   What does this tell us about possible control strategies?

In terms of the parameters listed in Table I,

R˜ = [(A B M) / m] × [(A C N) / r] × e-mt                                       (11)

where the quantities in square brackets can be viewed respectively as the per mosquito production of newly infected humans (the term, [(A B M) / m]) and the per case production of newly parasitized mosquitos (the term, [(A C N) / r]). Obviously, one wants to reduce the quantities, A, B, C and M, while increasing m and r. In addition, we note the following:
 

Table II. Control Strategy Efficacy.
  R / RD
D m-Control
(mt = 1)
A-
Control
Non-A/m
Control
1 1 1 1
2 .18 .25 .50
5 .001 .04 .20
10 4.5×10-6 .01 .10

  1. The mosquito mortality rate, m, enters twice, once in proportion to (1/m), the other time, exponentially as e-mt. As a result, reducing mosquito longevity, as might be accomplished by spraying, provides the surest approach to reducing .

  2. The bite rate, A, also enters twice - in proportion to A2. This makes sense, since completing the plasmodium's life cycle requires two bites. Accordingly, a 50% reduction in the number of bites divides R by a factor of 4, etc. It follows that, even without killing mosquitos, control strategies (screens, pesticides, netting) that keep insects out of peoples' houses, or at least out of their beds, yield favorable results. Together with point # 1 above, this fact explains the historical success [a00] of spraying programs (now largely curtailed) in which houses were dusted with DDT once or twice a year.

  3. The remaining parameters enter once: Vaccines that reduce human susceptibility reduce B; drugs that reduce the infectious period by killing plasmodia within the human body increase the recovery rate, r; interupting the parasite's life cycle within mosquitos, e.g., via the creation of transgenic mosquitos, reduces C; reducing the overall mosquito population by spraying or draining marshes and other breeding sites, reduces M; reducing local human population density (relocation programs) reduces N.

In sum, control strategies that focus on the death rate, m, of adult mosquitos or the bite rate, A, are inherently more efficient than "non-A/m" strategies. This is not to deny that complete interuption the plasmodium life cycle at any point would end the scourge. After all, the basic reproductive rate, R, depends multiplicatively on all the parameters, and multiplying anything by 0 gives you 0. At the same time, non-A/m control strategies must, of necessity, be more effective than m- or A- strategies to achieve the same results. Illustrative calculations are given in Table II. Here, D is the value of that obtains upon multiplying or dividing the parameter in question by a factor of D.

Caveat.

Equations (3), ignoring as they do, enormous amounts of biological complexity, are highly simplified.8   Still, something like the partitioning of the basic reproductive rate into two components, and the multiplication of their product by e-mt, is likely to fall out of any reasonable mathematical description of the disease. In effect,

R˜ = R˜hm × R˜mh × e-mt                                                     (12)

where R˜hm is the per mosquito production of newly infected humans, and R˜mh, the per case production of newly parasitized mosquitos, both of which depend on the bite rate, A. This, we believe, is the biological essential, and it guarantees the distinction between "m-", "A-" and "non-A/m" control strategies noted above.



Exercise 6. Enumerate additional aspects of the biology omitted from Equations (8) and discuss their probable conbsequence to the conclusions of the preceding section. In particular, think about the selective consequences to mosquito and plasmodium life histories of various control programs.
 



References.

[am82] Aron, J. L. and R. M. May. The population dynamics of malaria. Pp. 139-179. In, Anderson, R. 1982. The Population Dynamics of Infectious Diseases: Theory and Applications. Chapman and Hall. London.
 
[a00] Attaran. A., Roberts, D. R., Curtis, C. F. and W. L. Kilama. 2000. Balancing risks on the backs of the poor. Nature Medicine. 6: 729-731.
 
[ka03] Koella, J. C. and R. Antia. 2003. Epidemiological models for the spread of antibiotic resistance. Malaria J. 2: 3-14.
 
[ly76] Lajmanovich, A. and J. A. Yorke. 1976. A deterministic model for gonorrhea in a nonhomogeneous population. Math. Biosci. 28: 221-236.
 
[m52] MacDonald, G. 1952. The analysis of equilibrium in malaria. Trop. Diseases Bull. 49: 813-828.
 
[m57] MacDonald, G. 1957. The Epidemiology and Control of Malaria. Oxford Univ. Press. London.
 
[r11] Ross, R. 1911. The Prevention of Malaria. John Murray. London.
 
[r16] Ross, R. 1916. An application of the theory of probabilities to the study of a priori pathometry. I. Proc. Roy. Soc. London Ser. A. 92: 204-230.



Notes.

1. According to this, the mean recovery time in humans is (1/r), and the mean longevity of mosquitos, (1/m). Since (1/r) > (1/m), we say that the time scale of Y is "fast" relative to that of X.
 
2. For readers attempting an exact comparison, our parameter, q = M/N, is their m, and our m, their m.
 
3. We emphasize that the formulaic dependence of R on the various parameters is model-dependent.
 
4. As observed by Lajmanovich and Yorke [ly76], there is a wide range of epidemiological models manifesting this property, i.e., the existence of a function, R, the passage of which through 1 corresponds to a transcritical bifurcation involving the origin.
 
5. We now have three classes of mosquitos: exposed, of which the propotion is z, infectious, of which the proportion is y and susceptible, of which the proportion is (1 - y - z).
 
6. We leave it as an exercise for the interested reader to derive Equations (8) from expressions corresponding to Equations (1), i.e., for a system of equations in which the variables are numbers as opposed to proportions.
 
7. With R < 1, the host population is said to manifest "herd immunity." We emphasize that herd immunity, along with the transcritical bifurcation whence it arises, is a property of most epidemiological models and is, in fact, the essential conceptual underpinning of disease eradication programs across the board.
 
8. For discussion of more complex models see Aron and May [am82] and Koella and Anita [ka03].