library(ramp.falciparum)
library(deSolve)
Superinfection & Queuing
An Introduction to Queuing Theory for Malaria Epidemiology
A basic question in malaria epidemiology is how long it would take to clear an infection. Initial investigations suggested the answer would be counted in months (see the vignette on Infection Duration). In some early studies of exposure, people were getting more than one infective bite, per person, per day [1]. These sorts of observations gave rise questions about the importance of superinfection, or getting infected while already infected.1 Ross’s model had assumed that infected people would clear parasites at a constant rate, and that they could then become infected again. Given the high rates of exposure, that assumption could be wrong.
If people were getting exposed at a higher rate than they were getting infected, how long would it take before they cleared all their infections? If all the parasites in an infecting bite are called a brood, how many broods would be infecting a person? Over the decades, studies of the distribution of the number of broods per person have used the term multiplicity of infection (MoI), a topic we take up in a deep dive vignette (see Multiplicity of Infection).
Walton’s Answer
On the question of the MoI, GA Walton gave this answer in 1947 [3]:
Assume that mosquito-bites are distributed at random. Then if \(m\) be the number of infective bites per person per year, the number of person bitten \(0, 1, 2, 3,\) etc. times per year will be in the proportion of the Poisson series: \(e^{-m}(1, m^2/2!, ...)\)
Walton’s answer, at this stage, is framed only in terms of the number of times a person would be infected in a year, but then he goes a bit deeper. Suppose \(q\) is the proportion of those infections that would be detected, then the fraction of people who would test positive is:
\[1-e^{-mq}\] Notably, Walton’s answer involving detection implicitly includes clearance. The formula assumes (correctly) that a Poisson distribution with mean \(m\) that has been compounded with a binomial with probability of detection \(q\) is another Poisson with mean \(mq.\) The fraction who would test negative is the zero term, \(e^{-mq};\) everyone else would test positive.
It would, perhaps, be useful to formulate the underlying mathematics in a dynamical equation.
Macdonald’s Dynamics
NOTE: In the remainder of this article, we let \(m\) denote the MoI, and we will let \(h\) denote the daily FoI. The notation comes from Ross, who called the FoI the “happenings” rate.
We would like to write down a new formula for the dynamics of malaria with superinfection. In general, we let \(h\) denote the force of infection, and \(r\) the clearance rate, the clearance rate is denoted \(R(h,r),\) and we can write: \[ \frac{dx}{dt} = h (1-x) - R(h,r) x \]
In 1950, George Macdonald published a dynamical model of malaria superinfection [4], where the clearance rate depended on the rate of exposure. In Macdonald’s model, Macdonald describes parasites clearing independently – the same assumptions as Walton – but his mathematical formulas match a different model.
\[
R(h,r) =
\left\{
\begin{array}{rl}
h-r & \mbox{if } h < r \\
0 & \mbox{if } h \geq r \\
\end{array}
\right.
\] Paul Fine’s essay on the problematic formulation of Macdonald’s model for superinfection is highly recommended reading [5].
It’s worth mentioning briefly that Macdonald tended to avoid analysis where the superinfection assumption would have been important. For example, when Macdonald first published his formula for the basic reproductive number, \(R_0,\) the analysis focuses on the dynamics when prevalence is low and superinfection can be ignored [6]. Similarly, when Macdonald led a paper describing the timelines for malaria elimination following the interruption of transmission, he assumed \(h=0\) [7]. In both cases, assumption about clearance ended up being identical to the assumption that Ross made in his simple model (an SIS compartmental model, see The Ross-Macdonald Model).
The Garki Approximation
In 1974, Kaus Dietz and collaborators presented slightly different solution to the problem of clearance with superinfection in a mathematical model developed for the Garki Project [8]. In the Garki model, the mean MoI is assumed to be \(h/r.\) The explanation is:
Inoculations “arrive” according to a Poisson process with rate \(h.\) Each inoculation has an exponentially distributed duration with mean \(r^{-1}.\) Then, in inquilibrium, the number of inoculations present at any time is a Poisson variable with mean \(h/r\).
In this model, clearance can only occur if a person is infected with an MoI of one and that infection clears. If the distribution is Poisson with mean \(m = h/r,\) the clearance rate expressed as a function of \(m\) is the probability of having an MoI equal to one, conditioned on being infected at all. Letting \(\zeta\) denote the random variable.: \[
R(m) = r \frac{\mbox{Pr}(\zeta=1 | m)}{\mbox{Pr}(\zeta>0 | m)} = r \frac{m e^{-m}}{1-e^{-m}} = r \frac{m}{e^{-m}-1}
\] Substituting \(m=h/r\)
\[
R(h,r) = \frac{h}{e^{h/r}-1}.
\] The Garki model solution traces is based on a general solution to the problem published in 1957, by Norman Bailey [9]. Bailey had solved used queuing theory, a framework that provides a general way of modeling both superinfection and the distribution. Dietz’s function was based on a steady state solution to a model from queuing theory, called \(M/M/\infty\). In effect, the Garki model assumes that the MoI instantaneously tracks the steady state of a system with the current FoI.
We will deal with queuing theory below, but before we do that, we present a correct solution to Macdonald’s problem.
Exact Dynamics
Knowing that the MoI, \(m,\) is the mean of a Poisson, and the prevalence of infection is the complement of the zero term from a Poisson, we get that: \[x = 1 - e^{-m}\] Using \(m,\) the change in prevalence can be described by the equation:
\[ \frac{dx}{dt} = h(1-x)-R(m, r)x \]
We can solve for \(m\), \[m = -\ln(1-x),\] and substitute, such that: \[ R(x, r) = - r \ln (1-x) \frac{(1-x) }{x}\] and noting that \[\lim_{x \rightarrow 0} R(x,r) = r,\] a simple formula for the change in prevalence would be: \[ \frac{dx}{dt} = h(1-x) + r \ln (1-x) \frac{(1-x) }{x} \] This equation is, perhaps, only of academic interest. It involves a strong assumption and a clever trick.
While the equation is correct, it is hiding a lot of complex math. Knowing \(x,\) and knowing the distribution is Poisson, we can always recover \(m.\) Without understanding that math, we’re stuck with this model and the assumption. Since we don’t want to get boxed in, we will want to understand a bit more about how all this works, as well as when the equations would not be appropriate. To do that, we will need to learn a bit about queuing theory. In the following, we formulate some queuing models. We do not know if queuing theory will ever be useful for a malaria analyst, but these ideas are woven through malaria theory, and it is never a bad thing to have an extra tool in the box.
In the next vignette, we present Hybrid Models, we cover some of this math, including an elegant solution to the problem of superinfection developed by Ingmar Nåsell [10].
Queueing Theory
The model Macdonald described is called \(M/M/\infty\).
It was correctly formulated and analyzed by Norman Bailey [9]. In this model, each new infection increases the MoI by one; the transition rate to a higher MoI is the force of infection (FoI), denoted \(h\). Each parasite can clear at some rate, \(r\). The following diagram shows how the fraction of the population with a given MoI, denoted \(\zeta_i=i\), changes over time.
\[\begin{equation} \begin{array}{ccccccccc} \zeta_0 & {h\atop \longrightarrow} \atop {\longleftarrow \atop r} & \zeta_1 & {h\atop \longrightarrow} \atop {\longleftarrow \atop {2r}} & \zeta_2 & {h \atop \longrightarrow} \atop {\longleftarrow \atop {3r}} & \zeta_3 & {h \atop \longrightarrow} \atop {\longleftarrow \atop {4r}}& \ldots \end{array} \end{equation}\]
Changes in the fraction uninfected (MoI=0) are described by:
\[\frac{d \zeta_0}{dt} = r \zeta _1 - h \zeta_0.\]
For those who are infected (with MoI\(=i\)), the dynamics are:
\[\frac{d \zeta_i}{dt} = h \zeta_{i-1} + r(i+1) \zeta_{i+1} - (r+h) \zeta_i\] The software in ramp.falciparum
includes functions to set up and solve the queuing model \(M/M/\infty.\)
Numerical Solutions
The following uses a software package called
ramp.falciparum
.
First, we need a model for the force of infection as a function of age and time:
plot(aa, FoI(aa, foiP3), type = "l",
xlab = "a - age (in days)", ylab = expression(FoI[tau](a)))
While \(M/M/\infty\) is an infinite system of differential equations, the function solveMMinfinity
solves a finite system of differential equations. The maximum MoI is chosen to be large enough that it doesn’t affect the results.
<- solveMMinfty(5/365, foiP3, Amax=1095) MMinf
with(MMinf, plot(time, m, type = "l", ylab = expression(m[tau](a)), xlab = "a - cohort age (in days)"))
CAB
The following is a chronological, annotated bibliography (CAB) for superinfection and queuing theory.
1914 – McKendrick’s paper…
1947 – On the control of malaria in Freetown, Sierra Leone; I. Plasmodium falciparum and Anopheles gambiae in relation to malaria occurring in infants. [3]
1950 – The analysis of infection rates in diseases in which superinfection occurs, by George Macdonald [4]
1957 – The Mathematical Theory of Epidemics, by Norman Bailey [9] introduced queuing models as a general way of solving the superinfection problem.
1975 – Superinfection - a problem in formulating a problem, by Paul Fine [5].