September 20, 2017, Wednesday
University of Colorado at Boulder Search A to Z Campus Map University of Colorado at Boulder CU 
Search Links

MBW:Epizootics of the Plague in Priarie Dog Populations

From MathBio

Jump to: navigation, search
Black tailed prairie dog.jpg

Basic Summary

(a) Math Used:

  • Deterministic Model
  • Stochastic Model
  • Sensitivity Analysis

(b) Disease Modeling, particularly plagues in prairie dog communities

(c) Studied prairie dog systems


In a recent wiki page Epizootics of the Plague in Prarie Dog Populations , we summarized an article describing the transmission of the plague (Y. Pestis) within local prairie dog populations. The fairly recent article, by Colleen Webb, et al [1] demonstrated that the commonly perceived threat to the spread of this bacterium (a certain reaction the dogs' fleas have to the disease), in fact, has very little affect on the rapid extinction due to the plague. The following wiki page is our interpretation of the article along with some analyses we conducted in response to what we learned in reading the newly conducted research.

Article Summary

In the study by Webb, after collecting data in the field as well as from previous literature, numerical and sensitivity analysis is performed on a population of prairie dogs in the Pawnee National Grasslands (PNG) in northern Colorado. The authors analytically and numerically succeed in demonstrating that the plague is spread through a town of these small mammals not by the previously believed "blocked fleas", but by other factors. The authors discover there may be several other factors driving the epizootics that are much more likely to cause extinction of a local population.


"Yersinia pestis" is the bacterium known to cause three different plagues: pneumonic, septicemic, and bubonic, to which both animals and humans are susceptible. These plagues are known to be responsible for highly contagious and deadly epidemics, including the Black Death, which destroyed approximately a third of the European population between the years 1347 and 1353. [2].

It has been long known that the black tailed prairie dog (Cynomys ludovicianus) is one of the rodent populations most susceptible to the plague (Yersinia pestis). Because of its high susceptibility, if an infected host is introduced into a "town" of these mammals, frequently, an epizootic follows in the population, resulting in a 98% mortality rate (based on the stochastic models utilized in this article). Another advantage to researching this particular species is the fact that it is vulnerable to all forms of transmission of the bacteria.

For many years, researchers believed that the mass deaths in prairie dog towns were mainly caused by the bites of "blocked" fleas; or fleas which, after contracting the disease, form a blockage in their midguts, resulting in starvation that leads to an aggressive feeding behavior and frequent regurgitation (in attempt to rid themselves of the blockage). The authors of this article disprove this common misconception, demonstrating with a thorough sensitivity analysis of the observed population, that these "blocked" fleas are not nearly as important in driving the epizootic as some other factors. There are two other ways of transmitting the disease: through airborne methods and contact with a "short-term infectious reservoir." This reservoir is a combination of deceased prairie dogs carcasses, waste from infected prairie dogs and possibly other smaller mammals which live in the vicinity of the prairie dog colony. Another alternative is that there exist certain species of insects that do not form a "blockage", but in fact carry the bacterium on their mouth parts. This article shows that the spread of Y. pestis results from one or a combination of these two alternative transmission routes, and their analysis helps explain how and why this occurs.

History of Mathematical Model

In the study done by Webb and her colleagues, Susceptible-Infected-Removed (SIR) models were used to model the spread of the Y. Pestis bacterium. (More information on the SIR model can be found under Population Biology of Infectious Diseases.) For general SIR depending on the individual topic being researched, parameter values will vary and the dynamics of the system can vary as well. The simplest model contains a susceptible class which moves to the infected class at a certain rate based on the density of the susceptibles and the infected, an infected class that can move to either the removed or the susceptible class, and a removed class which can possibly become part of the susceptible class. This model is very basic and while informative, brushes over many dynamics that a system can also exhibit.

The model under consideration has the same basic principles, but has been extended to more clearly explain the behavior of the system. Instead of only one species being considered, another species is introduced that helps spread the disease. Also, the prairie dogs themselves are broken up further into classes: susceptibles, the exposed class which grows in accordance to the amount of fleas that are "questing" for food, an infected class and a reservoir class. Ignoring the dynamics of the fleas, only an exposed class has been added. Also, the removed class consists of the deceased prairie dogs and other factors which contribute to the spread of the disease as well as the infected class. Thus, in this case, the "removed" class can also infect the susceptible and the exposed classes.

Mathematical Model

Deterministic Model

Parameter Values

There are multiple parameters for the model used in this article. They were derived to best represent actual parameters of prairie dog populations from previous literature on the same topic as well as from the data obtained in the field. Note that time units are in days.

Parameter Value Description
r 0.0866 Intrinsic rate of increase (host)
K 200 Carrying capacity (host)
\mu 0.0002 Natural mortality rate (host)
\beta _{F} 0.09 Blocked vector transmission rate
\beta _{C} 0.073 Airborne transmission rate
\beta _{R} 0.073 Transmission rate from reservoir*
B 20 Number of burrows host enters**
\sigma 0.21 1 / Exposed period (host)
\alpha _{F} 0.5 Blocked vector mortality rate
\alpha _{C} 0.5 Contact mortality rate**
\lambda 0.006 Reservoir decay rates
\delta 0.05 Rate of leaving hosts
a 0.004 Searching efficiency of questing fleas
\mu _{F} 0.07 Natural mortality rate (vector)
r_{F} 1.5 Conversion efficiency (vector)
\gamma 0.28 Transmission rate: hosts to vector
\tau 0.009 1 / Exposed period (vector)
s 0.33 Disease induced mortality rate (vector)

*Assumed to be the same as airborne transmission

**Estimated from field data

Research by the authors in the Pawnee National Grasslands in Northern Colorado and literature from other authors led to the estimates of most parameters in the model. All of the parameters that do not rely on the bacteria were obtained through this research as well as other statistics on prairie dog and flea dynamics. For the disease parameters, laboratory experiments were used with the bacterium Y. Pestis and O. Hirsuta, which portrays similar dynamics with Y. Pestis. For some of the transmission parameters, estimates were also taken from observations and data obtained from similar species of mammals like the California vole (Microtus Californicus).

Differential Equations

The differential equations used are:

The host model:

{\frac  {dS}{dt}}=rS(1-N/K)-{\frac  {\beta _{{C}}S(I_{{C}}+I_{{F}})}{N}}-\beta _{{F}}F_{{IQ}}S(1-e^{{-aN/\beta }})-\beta _{{R}}S{\frac  {M}{B}}-\mu S,

{\frac  {dE_{{F}}}{dt}}=\beta _{{F}}F_{{IQ}}S(1-e^{{-aN/B}})-E_{{F}}(\sigma +\mu ),

{\frac  {dE_{{C}}}{dt}}=\beta _{{R}}S{\frac  {M}{B}}+{\frac  {\beta _{{C}}S(I_{{C}}+I_{{F}})}{N}}-E_{{C}}(\sigma +\mu ),

{\frac  {dI_{{F}}}{dt}}=\sigma E_{{F}}-I_{{F}}\alpha _{{F}},

{\frac  {dI_{{C}}}{dt}}=\sigma E_{{C}}-I_{{C}}\alpha _{{C}},

{\frac  {dM}{dt}}=\alpha _{{C}}I_{{C}}+\alpha _{{F}}I_{{F}}-\lambda M.

The vector submodel:

{\frac  {dF_{{SQ}}}{dt}}=\delta F_{{SH}}+r_{{F}}F_{{O}}({\frac  {N}{1+N+F_{0}}})-F_{{SQ}}[\mu _{{F}}+(1-e^{{-aN/B}})],

{\frac  {dF_{{SH}}}{dt}}=F_{{SQ}}(1-e^{{-aN/B}})-F_{{SH}}[\mu _{{F}}+\delta +\gamma ({\frac  {I_{{C}}+I_{{F}}}{N}})],

{\frac  {F_{{EQ}}}{dt}}=\delta F_{{EH}}-F_{{EQ}}[(1-e^{{-aN/B}})+\tau +\mu _{{F}}],

{\frac  {F_{{EH}}}{dt}}=F_{{SH}}\gamma ({\frac  {I_{{C}}+I_{{F}}}{N}})+F_{{EQ}}(1-e^{{-aN/B}})-F_{{EH}}(\tau +\mu _{{F}}+\delta ),

{\frac  {dF_{{IQ}}}{dt}}=\delta F_{{IH}}+\tau F_{{EQ}}-F_{{IQ}}[s+\mu _{{F}}+(1-e^{{-aN/B}})],

{\frac  {dF_{{IH}}}{dt}}=F_{{IQ}}(1-e^{{-aN/B}})+\tau F_{{EH}}-F_{{IH}}(s+\mu _{{F}}+\delta ).

In this model, the prairie dog population are subdivided into six classes: susceptibles, S; exposed, E_{{F}}, and infectious, I_{{F}} by contamination by the blocked fleas; and those exposed, E_{{C}}, and infectious I_{{C}} through direct contact with the reservoir, M. From this, the total number of prairie dogs is described by N=S+E_{{F}}+E_{{C}}+I_{{F}}+I_{{C}}.

Similarly, the fleas are divided into six classes: susceptible and questing,F_{{SQ}}, susceptible and on the host,F_{{SH}}, exposed and questing,F_{{EQ}}, exposed and on the host, F_{{EH}}, infectious and questing,F_{{IQ}}, and infectious and on the host,F_{{IH}}.

These ODES use probability- and density-dependent contact between different groups, depending on the characteristics of each (i.e. questing fleas and on-host fleas exhibit different types of growth and death). It also seems that prairie dog colonies are usually structured in groups, which affect how diseases are spread throughout the population. Depending on the proximity of different family groups, the transmission rates between prairie dogs can vary. Because of the more random-characteristics of flea transmission, dynamics resulting from fleas are modeled primarily using frequency-dependent methods. To model the transmission caused by the short-term reservoir, density-dependent methods are used. As more prairie dogs die, the population becomes mixed as the structure of the colonies breaks down. Thus, density-dependent methods are more appropriate.

In the prairie dog classes, mortality rates vary due to biological differences. The prairie dog mortality rate, \mu , only affects the susceptible class and the exposed class due to the average length of time from exposure to death from infection is around 2 days. The number of holes that any prairie dog will enter, B is used to estimate area of the prairie dog colony. The amount of time that a prairie dog will remain in the exposed class before it moves to the infected class is given by \sigma ^{{-1}}. Thus, \sigma indicates the rate at which the exposed class moves to the infected class. Also, the infected class contributes to the short-term reservoir proportional to its density by the parameter \lambda .

The death rate of fleas is given by the parameter \mu _{{F}} and the rate at which the fleas die due to blockage is given by the parameter s. The transition from on-host fleas to questing fleas is given by the parameter \delta and the function (1-e^{{iaN/B}}). Thus, the transition from questing to on-host classes is based on the number of prairie dogs, N, the proximity of other burrows to the current one, B, and the efficiency of the fleas in finding a new host, a. Due to the short infection time and the abundant amount of blood needed to reproduce, reproduction of fleas is restricted completely to the on-host class. Since \tau ^{{-1}} is the time is takes for the fleas' proventriculus to get blocked, the rate \tau indicates the rate at which exposed fleas become able to infect the host due to regurgitation after the blockage.

Sensitivity Analysis

Sensitivity analysis was used to find which parameter values had the greatest impact on changing the extinction probability using the following equation:

\Sigma \approx {\frac  {ln(V(P)-ln(V(P_{{0}}))}{ln(P)-ln(P_{{0}}))}}.

This average sensitivity was calculated based on over 1000 simulations. The extinction probability is given by the function V. P_{{0}} is the default parameter whose sensitivity is being analyzed and P is another arbitrary parameter value.

Stochastic Model

The stochastic model is fairly similar to the deterministic model. The only difference is that if the prairie dogs die out (i.e. N = 0), the rate at which infected and exposed classes of fleas grow due to contact with the prairie dogs will become zero. The assumptions of the stochastic model are as follows: events occur only one at a time, all events occur independently of any other event and that the probability of an event occurring per unit time is held constant. For an in-depth study of stochastic epidemic behavior, see MBW:Stochastic Epidemic Modeling.


The authors of the article set up two different types of models; deterministic (figure a), and stochastic (figure b).

Dynamics of the model for the default parameter values starting with one infected host

Deterministic Model

Using numerical analysis, it was found that for the model of ordinary differential equations (ODES), there exist three different equilibria. The first is a stable equilibrium for a population with no existence of the Y. pestis bacterium, the second, an equilibrium where both susceptible and infected coexist, and the third, where all species become extinct.

In both of the models, the prairie dogs quickly become exposed, then infected. The population dies out in a matter of weeks. After the prairie dog population drops to extremely low levels, the fleas similarly begin to die out (or leave the colony) because they cannot reproduce without the necessary food supply. However, even after both populations die out, the short-term reservoir continues to persist at a large density for quite some time.

The only solution that can be solved without numerical techniques is the equilibrium were the plague no longer exists. However, the parameters can only lie in a very small region and most of default parameters did not lie anywhere it. The equilibrium where the plague persists but both of the species still survive had similar characteristics. In most of the parameter space, the plague causes extinction of both species.

In accordance with the sensitivity analysis, for the fleas to cause extinction, the transmission rate from blocked fleas would have to be two to five orders of magnitude higher than the default parameter.

Stochastic Model

Because in most populations, effects are generally density and frequency dependent, they also created a stochastic model to include these varying terms. They found that in both models, when infected individuals are introduced at the start, exctinction of all individuals occurred over a short period of time, and this time was only shortened further when the density of infected individuals increased at t=0. It is in this model the authors were able to predict that if only one infected host is introduced into the local population, the probability of extinction becomes 98%, which would occur within an average of 52 days. This is in agreement with the data taken from the Pawnee National Grasslands, which reported that the prairie dogs from infected sites dropped below detectable levels within 6-8 weeks.

Through the stochastic model, it was calculated that the average time that the reservoir lasted was approximately 2.73 years. According to data, the usual amount of time that a prairie dog colony is recolonized is about 2.59 years.

Sensitivity Analysis

Traditional sensitivity: extinction probability and time

After calculating the sensitivities of all of the parameters, it was found that none of the parameters had much of an impact on the extinction probability. However, it turns out that the blocked fleas do not show as much of an impact on extinction as was previously thought. The parameter that had the highest impact was the reservoir decay rate (\Sigma _{{\lambda }}=0.17). Also, extinction time was relatively sensitive to other parameters dealing with the reservoir: reservoir transmission rate(\Sigma _{{\beta _{{R}}}}=0.27) and the reservoir decay rate (\Sigma _{{\lambda }}=0.10). In most of the other parameters, the sensitivities in respect to extinction probability and extinction time were below 0.04.

Numerical Analysis

Recreating The Model

Using SimBiology (a program that uses MATLAB), we were able to numerically exhibit the dynamics of the ODE. After testing a range of parameters, we verifies the authors' conclusions that over most of the (realistic) parameter space, the outcome is that the entire colony dies with a matter of weeks. Of all of the simulations tested, the most interesting happened when the reservoir parameters were set to zero, thus taking it out of the equation entirely. With the parameters set to the default values in the paper and the reservoir parameters set to zero, the prairie dogs, fleas and plague continued to exist at a non-trivial equilibrium. With the reservoir rates at zero, the range at which the plague completely dies out is much larger. However, looking at data taken from the PNG and other prairie dog sites, this is not a realistic model. In almost 100% of cases studied, the colony is completely taken over by plague, Therefore, in accord with the authors' results, it seems as if the reservoir is what primarily drives the plague in prairie dog colonies.

Model without susceptible fleas Model without influence of reservoir

Sensitivity Analysis

We also performed sensitivity analysis on all of the parameters to see if the reservoir was greatly affected by changes in the values. However, it seems as if the reservoir is impervious to any of the parameters. The sensitivity is approximately zero for the entire set. This is not to say that changing the rates cannot effect the system. By dropping the rate at which the infected class contribute to the reservoir and increasing the rate at which the reservoir decays, the time it took for the population of prairie dogs to reach zero more than doubled, but the result was the same.

A more beneficial approach to trying to change the dynamics of the system was to test for the sensitivities of the susceptible prairie dogs in relation to the parameters. This was much more fruitful. The susceptible prairie dogs seem to be more sensitive to the exposed period, \sigma , the transmission rate from the reservoir, \beta _{R} and the number of burrows that a prairie dog enters, \mathrm{B} . Below are two sensitivity equations for {\frac  {dS}{d\beta _{R}}} and {\frac  {dS}{d\sigma }}.

dS/dBr dS/dSigma

Generalized Sensitivity Functions

From the authors' analysis of the system and our study of the numerics, it seems as if the reservoir drives the prairie dogs to extinction rather than the blocked fleas. However, the current parameters are just estimates of the actual system. In different colonies and with different species of mammals and insects, the time it takes for extinction to happen can vary dramatically. Since the susceptible class is highly sensitive to changes in certain parameters, correct estimates are vital to predicting the behavior of the system.

The usual sensitivity equations are derived from taking the partial derivative of the equations in respect to the parameter concerned. Sensitivity equations are useful to seeing at what time each group is sensitive to changes in parameters. However, they are not as useful at finding when to take field estimates to approximate the parameters because they do not take into account how the change in the parameter will affect the parameter estimates. To do that, another mathematical technique is needed: the Generalized Sensitivity Functions [3]

To generate the Generalized Sensitivity Functions (GSF), a nonlinear regression is performed. The usual method in fitting a model to the data at hand involves minimizing the the weighted residual sum of squares: \sum _{{i=1}}^{{M}}{\frac  {[y(t_{i})-f(t_{i},\theta )]^{2}}{\sigma ^{2}(t_{i})}} where \theta is the set of parameters, f(t_{i},\theta ) is the proposed model andy(t_{i}) is the actual data taken at time t_{i}. Manipulating this equation and assuming that all parameter estimates are unbiased leads to the GSF:

gs(t_{k})=\sum _{{i=1}}^{{k}}[([\sum _{{j=1}}^{M}\nabla _{{\theta }}f(x_{i})\nabla _{{\theta }}f(x_{i})^{T}]^{{-1}})\times \nabla _{{\theta }}f(x_{i}))\cdot \nabla _{{\theta }}f(x_{i})].

In this equation, \nabla _{{\theta }}f(x_{i}) is a matrix of the partial derivatives in respect to set of parameters, \nabla _{{\theta }}f(x_{i})\nabla _{{\theta }}f(x_{i})^{T} is the Fisher Information Matrix and \cdot denotes element-wise multiplication.

The GSF are equations that start at zero and end at one at the end of the time considered. Since these functions are dealing with information, this makes sense. At the first time step, there is no information and at the end, all of the information needed is available. These equations tell us the time where the most information is collected for a certain parameter. When the GSF has the highest rate in increase is the time when an experimentalist should take the data.

Prairie Dog GSF

In this article, the reservoir was a hypothetical idea. Since it seems that the reservoir actually contributes much more to the dynamics of the system than the fleas, the next step of the process is to try and find what the short-term reservoir actually is and what values the parameters take with that case.

Since the reservoir was minimally sensitive to any parameter, the GSF was not very useful in finding when to take data. Next to the reservoir, the susceptible class was of highest importance. Thus, we produced the GSF of the susceptible class in respect to many parameters concerned with the reservoir. Here are graphs for two of the parameters, \lambda and \beta _{R}.

GSdSdLambda.jpg GSdSdBr.jpg

For \lambda and \beta _{R}, the graphs do not change much in the first few weeks of the epidemic. This means that all data for the reservoir should be taken after the force of the flea infection has died down a bit. At the beginning of the epidemic, There is no reservoir, so the fleas are the driving force of the infection. However, as the days pass, the reservoir builds up and takes over the work and drives the prairie dogs to extinction. Therefore, it makes sense that data for the reservoir should be taken near the end of the epidemic. This will be where the information is most useful.


Both the deterministic and the stochastic models indicated that extinction probability and extinction time are relatively more sensitive to parameters dealing with the short-term reservoir than any of the parameters concerning the blocked fleas. Extinction times were more sensitive to changes in the parameters, but the end result is that both species will eventually die out with the default values for the parameters estimated from the data.

This is very useful information in respect to the spread of this bacteria. The reservoir transmission rate has the greatest impact on the extinction probability. This can be used to determine methods that may stop the plague from spreading or killing out the colony. Removal of the reservoir would be the best way to stop the extinction. Since the bubonic plague is a strain of the Y. Pestis bacteria and the dynamics of other strains of bacteria that infect humans are quite similar to the bacteria discussed, smarter methods in the treatment of such diseases could be implemented. Since the spread of the disease by fleas does not easily change the outcome, making sure that the short-term reservoir is removed or treated is a more effective way of slowing the outbreak.

Our numerical studies back up what the authors where trying to prove with their work. The reservoir seems to have a higher impact on the extinction of the prairie dog colony than any other factor in the spread of the disease. Without the reservoir in the equation, the range of parameters that lead to a non-trivial equilibrium was much larger than when the reservoir was introduced.

Many of the parameters that the authors used were from other species similar to the fleas and the prairie dogs. To expand on their work, we studied the sensitivity equations of the system and the Generalized Sensitivity Functions to determine when would be the most opportune time to take data from a colony to better estimate the parameters and get a clearer picture to what is exactly going on in the system. According to our analysis, the best time to take data on the reservoir would be near the end of the epidemic. This is when the reservoir is most actively driving the infection. Once the fleas spread the plague, the reservoir kicks in and drives the prairie dogs to extinction.

Now that it seems as if something else is driving the plague, the next step will be to find what the reservoir exactly is. There are many possibilities including dead prairie dogs, bacteria in the soil and even other animals who are less susceptible to the plague. Once this new reservoir is found, this model gives a good start to modeling the dynamics of the plague in the colony. With our analysis, it also gives them an idea as to when to take data to best fit the system to the actual prairie dog colony.

Recent Citations

R. J. Eisen, A. P. Wilder, S. W. Bearden, J. A. Monteneiri, and K. L. Gage. Early-Phase Transmission of Yersinia pestis by Unblocked Xenopsylla cheopis (Siphonaptera: Pulicidae) Is as Efficient as Transmission by Blocked Fleas. EcoHealth, 2008

M. Buzby, D. Neckels, M. F Antolin, D. Estep. Analysis of the sensitivity properties of a model of vector-borne bubonic plague J R Soc Interface 2008 5 (26) 1099-1107

R. J. Eisen, S. W. Bearden, A. P. Wilder, J. A. Montenieri, M. F. Antolin, and K. L. Gage. Early-phase transmission of Yersinia pestis by unblocked fleas as a mechanism explaining rapidly spreading plague epizootics Proc. Natl. Acad. Sci. USA 2006 103 (42) 15380-15385

Discussion of Recent Citations

The purpose of article by R. J. Eisen et. al is to discuss the transmission of "Yersinia pestis" bacterium. In particular, the C.T. Webb et. al article is cited by R. J. Eisen et. al as evidence that flea-born transmission could not have occurred from flea-borne transmission alone. The Webb article showed that through sensitivity analysis, their model was not sensitive to flea-borne transmission because the blocked fleas die before they become infectious. This result is cited multiple times throughout the Eisen article.


  1. Webb, Colleen T., et al, 2006, Classic flea-borne transmission does not drive plague epizootics in prairie dogs, [], PNAS vol. 103:6236-6241
  2. Yersinia pestis [1], Wikipedia
  3. Thomaseth, Karl et al, 1999, Generalized Sensitivity Functions in Physiological Systems, Annals of Biomedical Engineering vol. 27:6087-616