|Home | About | Journals | Submit | Contact Us | Français|
Building a mathematical model of population dynamics of pathogens within their host involves considerations of factors similar to those in ecology, as pathogens can prey on cells in the host. But within the multicellular host, attacked cell types are integrated with other cellular systems, which in turn intervene in the infection. For example, immune responses attempt to sense and then eliminate or contain pathogens, and homeostatic mechanisms try to compensate for cell loss. This review focuses on modeling applied to malarias, diseases caused by single-cell eukaryote parasites that infect red blood cells, with special concern given to vivax malaria, a disease often thought to be benign (if sometimes incapacitating) because the parasite only attacks a small proportion of red blood cells, the very youngest ones. However, I will use mathematical modeling to argue that depletion of this pool of red blood cells can be disastrous to the host if growth of the parasite is not vigorously check by host immune responses. Also, modeling can elucidate aspects of new field observations that indicate that vivax malaria is more dangerous than previously thought.
The online version of this article (doi:10.1007/s12551-010-0034-3) contains supplementary material, which is available to authorized users.
Population biology is concerned with the dynamics of a population of organisms and influences on those dynamics. For example, consider a herd of zebras on a plain. Lions and hyenas prey on the zebras. If there is a drought, fewer zebras may survive, thus fewer lions and hyenas survive. Lions and hyenas compete against each other, and parasites infect both zebras and predators. A less nutritious grass may invade the plain, etc. From a biophysical standpoint, interactions among the species facilitate the exchange of energy along with nutrients necessary to use that energy. The plain is a substrate on which natural selection acts on the component species and thus shapes their genomes. Mathematical modeling of an ecosystem lets us investigate “what if” situations: given a proposed interaction or process in the ecosystem, what is the logical outcome of the presence of that interaction or process upon the dynamics of the populations in that system? Will these outcomes be robust against perturbations on the system? To quote John Maynard Smith (Maynard Smith 1968), modeling (1) “provides a rapid way of discovering the kind of effect various features—for examples, predation, parasitism, territorial behavior, a long or short period of immaturity—may have on the behavior of a population” and (2) “suggests what needs to be measured before the behavior of any particular species can be understood.”
A group of pathogenic agents in a multicellular host also forms a population. The pathogen is a predator if it attacks a cell type. But multicellular hosts are usually proactive against an invader, mounting immune responses to kill or limit the pathogen and starting processes to repair the tissue damage. From a modeler’s perspective, perhaps the ecosystem is a morass: we may not know all the components of the host’s response, and of those we do know, we may not know which are most important. Host responses may be self-amplifying, self-regulating, and possibly self-harming. They may differ greatly among individuals of the same species due to genetic background or contingency of environmental exposure. The pathogen might have mechanisms to evade or interfere with host responses. But again, we can use mathematical modeling to answer “what-if” questions, and this in turn can guide researchers in the laboratory or clinic. There are two levels of population dynamics in a host-pathogen system: (1) those of individual cells or virions of the invader as well as host cells and responses, and (2) those of both host population and disease incidence, the subject of epidemiology. In this review, I discuss the application of mathematical modeling to a type of human malaria, with the focus on the within-host dynamics of the pathogen. I will also discuss the significance of the results from both public health and evolutionary standpoints.
Malarias are diseases of red blood cells (RBCs) in vertebrates caused by 100+ species of single-cell eukaryotes of the genus Plasmodium, order Haemosporida of the class Apicomplexa. Birds, reptiles, and mammals are targeted (Garnham 1966; Perkins and Austin 2009). Apicomplexans have a set of specialized organelles, together called the apical complex, that enables them to penetrate the cells of their hosts. All known Apicomplexans are intracellular parasites, and some are human disease agents (Cryptosporidium parvum, Toxoplasma gondii) or serious livestock pathogens (Babesia spp., Theileria spp.). At least five Plasmodium spp. infect humans regularly, causing disease characterized by repeating bouts of severe fever and anemia [see Cox-Singh and Singh (2008) for commentary]. Part of the life cycle of these species is spent in mosquitoes of genus Anopheles. The mosquitoes serve as vectors that spread the disease from human to human by picking up the gametocytes, the sexual forms of the parasite that mate within the mosquito gut. Mosquitoes transmit to humans the haploid descendants of the mated forms, the sporozoites, which then enter liver cells. The form of parasite that infects RBCs is the merozoites, which are released initially from the liver as the descendants of sporozoites but then later by the bursting of infected RBCs. [Baer et al. (2007) is an illustrative study of the initial merozoite release.] After a merozoite enters an RBC, its nucleus undergoes rounds of division by schizogony. The infection is sustained among the host’s RBCs in this manner, but some asexual forms become gametocytes in a process that is not well understood. The disease pathology is due to the asexual forms that cannot be transmitted to mosquitoes, although blood transfusion can spread them.
Two species, Plasmodium falciparum and Plasmodium vivax, cause the bulk of the public health burden of the malarias. [For review, see White and Breman (2008) and Miller et al. (2002).] P. vivax attacks mainly reticulocytes, which are RBCs that have aged 36 h or less after differentiation from erythroblasts and that still retain traces of a nucleus. This parasite infects ≈108 people per year in both temperate and tropical areas. Infection is temporarily incapacitating, but as I discuss below, its true public health burden is only now becoming apparent. P. falciparum can infect RBCs of any age. This species infects ≈5×108 people each year and causes >106 deaths per year, mostly in tropical regions and especially in sub-Saharan Africa. Severe anemia and organ complications, especially in the brain, cause death or long term disability. Figure 1 shows a human RBC parasitized with this species.
In 2003, Ellis McKenzie of the National Institutes of Health, Bethesda, approached me with a conundrum: It was assumed that P. vivax caused less severe disease than P. falciparum because a mere 1.25% of RBCs are reticulocytes (Krogstad 1995). Yet data from neurosyphilis patients who were infected with P. vivax showed that this parasite has a rather large growth rate early in infection, with parasites from a single infected RBC infecting up to 15 other RBCs (McKenzie et al. 2002). (Before penicillin came into use, the induction of antibacterial fever by malaria was the standard of treatment for neurosyphilis.) This was similar to the initial growth rate of parasitemia seen in P. falciparum infections in neurosyphilis patients (Collins and Jeffery 1999) and other patient groups (White 1997). Could P. vivax cause severe anemia by depleting the reticulocyte population, thus choking off replacements for aging RBCs removed by the spleen? If so, what makes P. vivax supposedly less dangerous than P. falciparum? We set about modeling the dynamics of the parasite population in an hypothetical patient to address this question. This review describes the results of our efforts.
Incidentally, there are many fine theoretical studies of the within-host population dynamics of malaria parasites that attempt to address several aspects of malaria infection. For a selection, see Gravenor et al. (1995), Hoshen et al. (2000), Gravenor et al. (2002), Rouzine and McKenzie (2003), Haydon et al. (2003), Maire et al. (2006), Chiyaka et al. (2009), and Thibodeaux (2010).
We first considered the following question: what would happen if parasitemia grew to the maximum level constrained only by availability of RBCs, with no immune response from the host against the parasite? This situation resembles Malthusian growth of the within-host parasite population, although Mathus was primarily concerned with the effects of human population growth on economies. But the basic idea is that a population grows unchecked until limited by availability of resources, and the growth can lead to catastrophe for the population if the host dies due to anemia. (However, if a mosquito has picked up some of the gametocytes, the catastrophe is not for the Plasmodium species itself. The within-host population is sacrificed.) The 1798 essay of Malthus had a great impact on biologists in the 19th century and is considered one of the founding works in population biology (Malthus 2003). In some ways, the models that I describe are akin to the logistics equation introduced by Malthus’s younger contemporary P.F. Verhulst to account for population growth limited by a resource (Hoppensteadt 1982). Before introducing the models, I must introduce two important concepts.
The lifetime τ μ of a free merozoite is ≈10 min (Johnson et al. 1980), much shorter than the development time of the parasite within RBCs (≈2 days for both P. vivax and P. falciparum). The difference between these time scales allows us to derive, by statistical arguments, the average number of descendants an infected RBC has after one development cycle, R, early in the infection before an immune response can respond (or if there is no immune response). Assuming that all vulnerable RBCs are equally likely to be invaded, we introduce a merozoite-RBC binding affinity ζ, so that the rate of infection of RBCs per unit volume of blood, Inf r, is
where μ is the density of merozoites in the blood, and V is the density of vulnerable RBCs (reticulocytes for P. vivax, all RBCs for P. falciparum). Making an additional assumption that early in the infection V changes slowly during one cycle of the within-RBC parasite, we obtain
Here p is the average number of merozoites released from an infected RBC. [The details of the derivation are given in McQueen and McKenzie (2004)]. Direct photographic evidence of bursting infected RBCs suggests that p≈16 (Glushakova et al. 2005). As already mention in the section “A conundrum of malaria biology” above, R for both species at the beginning of infection is ≈15. Thus, the asexual forms of both species are very effective at propagating themselves before host immune reactions can intervene, with a characteristic time to invade an RBC ζV≈24 s. I call the initial value of R the basic reproductive rate, R 0. Given V, τ μ, and p, knowing R 0 determines ζ. As infection progresses, R declines from R 0 because (1) the vulnerable RBC population is depleted and (2) the host mounts an immune response against the parasite.
To develop a reasonable mathematical model, one has to treat the various constituents of the pathogen-host system (the erythocytes, the parasite, etc.) as populations of aging individuals: that is, on a population average, individuals spend D amount of time in a given stage of development with some standard deviation σ. If simple rate equations were used to model the population levels of the constituents, then population level changes would be like that of radioactive nuclei decaying at random. Yet the development time of an individual parasite within an RBC is 48 h for the two species, with maybe a few hours of variation. Reticulocytes last only 36 h, plus or minus a few hours. Individuals in such populations are not like radioactive nuclei that randomly “pop” with probability per unit time of D −1. To account realistically for aging, we used the technique of compartment ordinary differential equations in which a stage of development of an organism (with duration D) is described as a set of N substages with duration DN −1 that decay into one another. Letting P n be the population in substage n (and ignoring interactions between individuals for a moment):
It can be shown that the standard deviation about D for individuals in the development stage is
I stress that D, σ, and the total population ΣP n are the tangible quantities connected to measurement, not N. The N substage populations are just a fictional construct. This technique was developed in the 1960s for population modeling, with Lloyd adapting it to epidemiology and within-host pathogen population dynamics (Lloyd 2001a, b). [See also supplementary material in McQueen and McKenzie (2004).] We wrote a set of modular code libraries in Java (now in C++) that allowed for rapid modification of the numerical model as needed. It was found that a fourth-order Runge-Kutta method with adaptive time stepping (Cash and Karp 1990) worked very well to solve the ordinary differential equation system arising from the application of Eqs. 3 and 4, even though such systems might contain hundreds of dependent variables.
This model was first described in McQueen and McKenzie (2004); it can thought as two developmental chains, one for RBCs and the other for the parasite. The stages of each chain are indicated in Fig. 2. For each such stage, a system of differential equations similar to Eq. 3 describes the population changes; the corresponding values of D and σ are also shown in Fig. 2. The two chains interact through the binding of merozoites to the vulnerable RBCs, as described in Eq. 1. For convenience of the reader, the full system of differential equations that describe this model is given in the Electronic supplementary material (ESM).
The simulation is started with a small number of merozoites to mimic the primary release of merozoites from the liver. The normal basal production rate for human RBCs is about 1,736 μl−1 h−1; we considered models in which this rate could be amplified to compensate for the RBC loss due to infection. We also considered models with dyserythropoiesis, a decline in the rate of production due to factors of the parasite. [For a review of observed effects of malaria upon erythropoiesis, see Wickramasinghe and Abdalla (2000).] For some sets of simulations, we varied the σ for the erythrocytic stage of the parasite to test the effects of changes in this model parameter upon the results. This will be discussed more in the section “Sensitivity to the standard deviation of the life cycletime of the intraerythrocyte parasite.”
The erythrocyte-limited model of parasite dynamics predicts that even for a basic reproduction rate R 0≤2, P. vivax infection would eventually deplete the reticulocyte pool enough to kill the patient: there would not be enough replacements for aged RBCs removed by the spleen (McQueen and McKenzie 2004). Figure 3 shows some typical simulation results. [The results in the figure are slightly different from those presented in McQueen and McKenzie (2004) because (1) σ of the infected RBC stage is 2.4 h instead of 4.8 h as in that paper, and (2) the level of lethal anemia was set lower than in that paper.] One effect apparent in Fig. 3 is that a compensatory response of RBC production would make matters worse for the host with P. vivax infection by lowering the value of R 0 at which death via severe anemia could occur, as well as increasing the overall parasitemia. Conversely, dyserythropoiesis would tend to lower the number of parasites produced in lethal infections. These latter effects would also occur in the resource-limited model of P. falciparum infection, but then there is less effect on the survival of the host (see Fig. S1 in the ESM).
Since R 0≈15 in actual infections, even in patients who survive the disease, human immune constraints on the growth of P. vivax must be very tight indeed! Even though the erythropoietic function affects some of the details of the erythrocyte-limited model infection outcome, infection with P. vivax in this model would be lethal even if R 0 were just above 1 regardless of the erythropoietic function.
The importance of the age-structured nature of the attack of P. vivax becomes apparent when examining Fig. 4, which shows results for simulations using parameters for the pathogen Plasmodium malariae, an attacker of senescent human RBCs. [The duration of the infected RBC stage for this pathogen is 72 h instead of 48 h. See McKenzie et al. (2001) for a review of knowledge about this parasite.] No simulated P. malariae infection ended in host death unless there was dyserythropoiesis, even though P. malariae attacks a slightly larger subset of RBCs in this model. Since no other stage of RBC development depends on the senescent cells, depletion of this population would not cause severe anemia unless a secondary effect disrupted the production of RBCs or killed uninfected RBCs.
If both P. vivax and P. falciparum infect the host in an erythrocyte-limited model, the results depend on which infects first, and their respective values of R 0. Basically, we found that the species with the largest value of R 0 dominates. See McQueen and McKenzie (2006) for more detailed discussion.
In real malaria infections, hosts mount aggressive immune responses, both innate and adaptive. For example, the report of Karunaweera et al. (1992) illustrates the quick production (and then decay) of cytokines associated with innate responses in patients infected by P. vivax. [Some good reviews of immune responses to malaria are Engwerda and Good (2005), Karunaweera et al. (2003), and Plebanski and Hill (2000).] To address how the presence of immune responses might affect the virulence of Plasmodium spp., we added additional components to the model to incorporate two types of immune reactions: (1) a fast acting, fast decaying response to mimic the intense innate reactions seen in malarias, and (2) a slowly activating, slowly decaying response to mimic an antibody response. Since some of the innate immune reaction is believed to be triggered by the bursting of the schizonts, the model innate response is triggered when the level of late-stage-infected RBCs is above a threshold (that varies from simulation to simulation). Then a factor is produced that attacks some subset of the parasite; the rate of removal of parasites, C, would be
Here, Att Inn is the attacker level in the blood for the model innate reaction, ξ Inn is the binding affinity of the attacker to the targeted parasites, and P is the subset of parasite population that is attacked. (For example, it could be only infected RBCs that are trophozoites. In another simulation it might be only schizonts.) In a sense, it is the product ξ Inn Att Inn that is the tangible quantity, as it is a rate of removal. The model antibody response is similar, except that a delay is incorporated before the level of the attacker factor grows; a formalism similar to that of Eq. 3 above enforces the delay. The dynamics of both types of immune responses are modeled with a system of differential equations; for details see McQueen and McKenzie (2008). Also, the immune components of the model are expounded in the ESM.
The targets, trigger levels, and binding affinities of both types of responses were varied so that simulations were performed for thousands of points in the parameter space of the models. (In addition, we consider several values of σ for the infected erythrocyte stage; see the section “Sensitivity to the standard deviation of the life cycletime of the intraerythrocyte parasite” below.) The main conclusion is that given an immune response with the same capacity to detect and clear parasites, P. vivax infection would be just as severe in terms of risk of lethal outcome as P. falciparum. This can be seen in Fig. 5, which shows the results of simulations with basic reproduction rate R 0=15 but with immune responses present. In most of the phase space of the model immune parameters, the immune responses protect the host from death by severe anemia. A compensatory increase in the rate of erythropoiesis to blood cell loss due to infection tends to boost overall parasitemia, especially for P. vivax, although the fraction of phase space in which the host dies is about same as for the case when the erythropoiesis is fixed at the basal rate. On the other hand, dyserythropoiesis tends to reduce on average the burden of parasitemia to which the host can mount an immune response, but increases the fraction of phase space that is lethal to the host, especially for P. vivax. One wonders if the dyserythropoiesis observed in some patients has been shaped by an evolutionary struggle between humans and Plasmodium spp.
Since even models with plausible immune responses predict that P. vivax can cause infections as severe as those of P. falciparum, why does the latter species cause more severe disease on average than the former, especially in people who have no previous exposure to either pathogen? One possible answer is that actual human immune responses to P. vivax are in a part of the immune system phase space in which infections are controlled more readily, either by being detected at a lower threshold than P. falciparum or by more efficient processes of removal of P. vivax-infected cells. Furthermore, immune responses are immensely complex, and perhaps some important component has been left out. For example, we did not consider the effects of antigenic variation in P. falciparum infection, in which the parasite modifies over time some of its proteins exposed to the immune system (Hommel et al. 1983). And RBCs infected with P. falciparum can sequester from the blood by adhering to vascular endothelium (Udomsangpetch et al. 2002), which perhaps gives that species an extra ability to evade immune attack. Although both species produce factors to suppress host responses (Kubata et al. 1998; del Portillo et al. 2004), perhaps P. falciparum can achieve suppression more readily. Yet there is another possibility: maybe severe P. vivax infection is more common than thought.
As we were doing our theoretical work on the population dynamics of P. vivax, empirical evidence mounted that this pathogen can cause severe disease after all. In Korea, vivax malaria was associated with retinal hemorrhaging (Choi et al. 2004) and splenic infarction (Kim et al. 2007). Kochar et al. (2005) reported 11 patients with severe malaria (renal and cerebral complications, severe anemia, fatal outcomes) from northwestern India in which polymerase chain reaction (PCR) analysis of parasite DNA as well as pathology examination showed P. vivax present, but not P. falciparum [see also Kochar et al. (2007)]. A follow-up report from that region (Kochar et al. 2009) showed that out of ≈1,100 patients admitted to a hospital in Bikaner for malaria, 40 of them had severe malaria in which PCR analysis could only find P. vivax and not P. falciparum. If further investigation can confirm this work, then P. vivax is a more severe public health problem than realized. The question is not academic, as the drug regime for P. vivax disease is somewhat different from that for P. falciparum disease. There are other questions, too, that need investigating: (1) Has severe vivax malaria always been present but cryptic, or are the cases seen in India and Korea due to the emergence of new strains? (2) Do patients with severe vivax malaria have some other health problem that exacerbates the disease? Only further research will show.
Empirical observation is what ultimately defines the truth of the matter, but the construction of a plausible mathematical model of malaria infection can help us understand why P. vivax can cause severe or even lethal disease. Of course, not all the pathology is due to anemia, as it is in the models discussed in this review. Some disease symptoms will be due to self-damage by the host immune responses, and perhaps by toxins released by the parasite (Elsheikha and Sheashaa 2007). Nonetheless, the modeling shows that there is no a priori reason for P. vivax to be harmless just because it attacks only a small subset of RBCs, the reticulocytes. The fact that all RBCs must pass through the reticulocyte stage means that P. vivax has the capability to do immense harm, especially if it can interfere with erythropoiesis in the bone marrow. Many of the vivax malaria patients described by Kochar et al. (2009) had severe anemia, a result understandable from our models.
In addition to elucidating how severe anemia can come about, mathematical modeling can help frame other questions concerning basic malaria biology. One interesting phenomenon in human malaria is the synchronization of merozoite release from bursting infected RBCs. Various mechanisms have been proposed for this: (1) different effects of fever during the earlier and late parts of the erythrocytic stage (Kwiatkowski 1995), (2) modulation of the life cycle of the parasite by immune responses (Hoshen et al. 2000; Rouzine and McKenzie 2003), and (3) modulation by the melatonin cycle of the host (Hotta et al. 2000). We found in our simulation that a quickly acting, quickly decaying innate response would induce synchronized release of merozoites as illustrated in Fig. 6; otherwise, bursting of individual parasitized cells becomes out of phase (McQueen and McKenzie 2008). One question we investigated was the sensitivity of simulation results to the standard deviation in the time that parasites spend within RBCs, σ IBC. The value of this quantity by itself does not determine if synchronization occurs, but it does effect the duration of merozoite release from the synchronized bursting of schizonts. In particular, we wanted to know if short duration of release would enhance the parasite levels by allowing an evasion of immune pressures.
Simulation results indicate that a small value of σ IBC (“tight” synchronization) might not necessarily benefit the parasite; in fact, as is apparent in Fig. 7, it might actually suppress overall parasitemia. Furthermore, the fraction of the model parameter space in which the host dies of anemia is smaller as the standard deviation in the within-host duration decreases, especially for P. vivax infection. From a purely dynamical point of view, a tight synchronization of merozoite release might well benefit the host, not the parasite. However, one effect not included is the possible release of a factor that might suppress immune responses: a quick, overwhelming release of such a factor could boost parasitemia. Another issue is that for the propagation of the parasite genome, it is the survival of the gametocytes and their ability to be taken up by mosquitoes that is important, and we do not know what bearing tight synchronization has on that process. Both of these issues merit further empirical and theoretical investigation.
The production of gametocytes is an area of active research, as much remains to be discovered. For example, are gametocytes produced at a constant rate throughout infection or is their production triggered by some stressor on the within-host parasite population (Carter and Miller 1979)? Indeed, Diebner et al. (2000) used mathematical modeling results to argue that the rate of gametocyte production is highly variable during P. falciparum infection and is modulated by host innate immunity. Furthermore, the actual transition steps from asexual forms to gametocytes are still not known, although careful gene expression analysis is giving us a clearer answer (Eksi et al. 2008). Can modeling help with the understanding of that process? The development of gametocytes differs between P. vivax and P. falciparum (McKenzie et al. 2007); what are the implications of that for their population dynamics (and perhaps on their transmissibility to mosquitoes)? Further modeling could well prove fruitful in understanding the genesis and development of gametocytes.
Finally, in the discussions above, I used severe anemia as a proxy for severity of infection in the hypothetical patient. But one could easily use additional criteria to define an infection in which the host is at major risk of dying or suffering long term disability. For example, the level of cytokines in a model immune response could be used as a proxy for possible cerebral damage. One could also consider a build-up of toxins from the parasite as well. “What-if” questions could be used to elucidate the consequences of plausible processes that may damage the host in cases of malaria.
I thank Dr. F. Ellis McKenzie, Forgarty International Center, National Institutes of Health, for insightful discussion. Some of the studies utilized the high-performance computational capabilities of the Biowulf Linux cluster at the National Institutes of Health, Bethesda, MD (http://biowulf.nih.gov). This research was supported by the Intramural Research Program of the NIH, Center of Information Technology.
This article is distributed under the terms of the Creative Commons Attribution Noncommercial License which permits any noncommercial use, distribution, and reproduction in any medium, provided the original author(s) and source are credited.