# The Failure ofR0

Authors: Jing Li; Daniel Blakeley; Robert J. Smith?
Journal: Computational and Mathematical Methods in Medicine
 (2011)

---

## Abstract

The basic reproductive ratio,R0, is one of the fundamental concepts in mathematical
biology. It is a threshold parameter, intended to quantify the spread
of disease by estimating the average number of secondary infections in a wholly
susceptible population, giving an indication of the invasion strength of an epidemic: if
R0<1, the disease dies out, whereas if R0>1, the disease persists. R0
has
been widely used as a measure of disease strength to estimate the effectiveness
of control measures and to form the backbone of disease-management policy.
However, in almost every aspect that matters, R0
is flawed. Diseases can persist
with R0<1, while diseases with R0>1
can die out. We show that the same
model of malaria gives many different values of R0, depending on the method
used, with the sole common property that they have a threshold at 1. We also
survey estimated values of R0
for a variety of diseases, and examine some of the alternatives that
have been proposed. If R0
is to be used, it must be accompanied by caveats
about the method of calculation, underlying model assumptions and evidence
that it is actually a threshold. Otherwise, the concept is meaningless.


---

## Body

## 1. Introduction

The basic reproductive ratio—also known as the basic reproductive number, the basic reproduction number, the control reproduction number, orR0—is one of the foremost concepts in epidemiology [1–3]. R0 is the most widely used epidemiological measurement of the transmission potential in a given population [4]. It is a measure of initial disease spread, such that if R0>1, then the disease can invade an otherwise susceptible population and hence persist, whereas if R0<1, the disease cannot successfully invade and will die out. The concept is defined as the number of secondary infections produced by a single infectious individual in an otherwise susceptible population  [5].Despite its place at the forefront of mathematical epidemiology, the concept ofR0 is deeply flawed. Defining R0 proves to be significantly more difficult than it appears. Few epidemics are ever observed at the moment an infected individual enters a susceptible population, so calculating the value of R0 for a specific disease relies on secondary methods. There are many methods to calculate R0 from mathematical models, few of which agree with each other and few of which produce the average number of secondary infections. Methods to calculate R0 from theoretical models include the survival function, the next-generation method, the eigenvalues of the Jacobian matrix, the existence of the endemic equilibrium, and the constant term of the characteristic polynomial. R0 can also be estimated from epidemiological data via the number of susceptibles at endemic equilibrium, the average age at infection, the final size equation and calculation from the intrinsic growth rate. For an overview, see Heffernan et al. [2].Furthermore, there are many diseases that can persist withR0<1, while diseases with R0>1 can die out, reducing the utility of the concept as a threshold. R0 is also used as a measure of eradication for a disease that is endemic, but issues such as backward bifurcations, stochastic effects, and networks of spatial spread mean that an invasion threshold does not necessarily coincide with a persistence threshold. This results in a reduction of the usefulness of R0. For example, it is possible that a disease can persist in a population when already present but would not be strong enough to invade.  Finally, the threshold value that is usually calculated is rarely the average number of secondary infections, diluting the usefulness of this concept even further.In this paper, we outline the problems withR0 and examine a number of alternatives that have been proposed. We include a worked example of malaria to demonstrate the many different results that the various methods give for the same model. Finally, we survey some of the recent uses of R0 in the literature. The number of articles that use R0 likely numbers in the tens of thousands, so an exhaustive review is not feasible. We have restricted ourselves to articles published since 2005 and which include interesting or novel explorations of R0.


## 2. Methods for CalculatingR0

In this section, we identify some of the more popular methods (although by no means all) used to calculateR0. We also describe the limitations that each method presents and demonstrate one of the core problems with R0.  Specifically, we address a key problem with R0: how do biologists make sense of it from mathematical models? (See, for example, the puzzled discussion in van den Bosch et al. [6].)Although the “R” in R0 is derived from “reproductive”, based on the original formulation of the concept as the average number of secondary infections, many thresholds have been denoted by “R0”, even when they are not related to the average number of secondary infections. Thus, in keeping with the notation, we will use the notation R0,X to denote an R0-like surrogate associated with a particular method, symbolised by X.


### 2.1. The Survival Function

The survival function is given by(1)R0,S=∫0∞(the average number of susceptibleindividuals that an infected individualnewly infects per unit timewhen infected for total timea)×(the probability that anewly infected individualremains infectiousfor at least timea)daThe survival function has the advantage that it always produces the average number of secondary individuals infected by a single infected individual, in the same class. Thus, in Figure1, where one human infects two mosquitoes, who each infect three humans, the survival function produces R0=6. This is the number of humans infected by a single infected human via mosquitoes; or, equivalently, the number of mosquitoes infected by a single mosquito via humans. The survival function is a generalised method of calculating the basic reproductive ratio that is not restricted to ODEs.Figure 1
A example of a two-stage basic reproductive ratio. A single human infectsRH=2 mosquitoes, each of whom subsequently infects RM=3 humans. Thus, a single human results in R0=6 infected humans. The next-generation method instead estimates R0,N=6 infected individuals in the subsequent generation, regardless of whether that generation is human or mosquito; R0,N is the geometric mean of RH=2 and RM=3. (Reproduced with permission from Smith? [7]).However, determining the individual probabilities can be cumbersome, especially if multiple states are involved. For a vector-borne infection such as malaria, with two infection states (human and mosquito), calculating the first probability involves determining the probability that a human infected at time 0 exists at timet, the probability that a human infected for total time t infects a mosquito and the probability that an infected mosquito lives to be age a-t, where 0≤t≤a [2]. For diseases with more states, such as Guinea Worm disease, where there is a waterborne parasite, which can attach itself to copepods, which in turn are ingested by humans and which subsequently grow into an internal nematode, the calculations of the survival probabilities become unwieldy.Thus, although this method always produces the correctR0, in practice, it is difficult to use. This is especially true for models with sufficient complexity, which are often those encountered most frequently.


### 2.2. The Jacobian

The Jacobian matrix is used to linearise a nonlinear system of differential equations. Around the disease-free equilibrium, the linear system will have the same stability properties as the nonlinear system if it is hyperbolic; that is, if no eigenvalues have zero real part. In particular, if all eigenvalues have negative real part, then the equilibrium is stable, whereas if there is an eigenvalue with positive real part, the equilibrium is unstable.It follows that a threshold isλmax=0, where λmax is the largest eigenvalue of the Jacobian matrix (or the largest real part if the eigenvalues are complex). However, for a system of n differential equations, this requires solving an nth order polynomial, which may be impossible. Furthermore, rearranging the condition λmax=0 to produce a threshold R0,J=1 is not a unique process  and does not always produce the average number of secondary infections.


### 2.3. Constant Term of the Characteristic Polynomial

Whenλmax=0, the constant term of the characteristic polynomial will be zero. However, the reverse is not true, as the polynomial could have both zero and positive roots. If the characteristic polynomial is (2)λn+an-1λn-1+⋯+a1λ+a0=0,
then a0=0 is a threshold if aj≥0 for all j. More generally, the Routh-Hurwitz condition allows the coefficients to take on other signs, under certain restrictions, but the nonconstant coefficients all being positive is a sufficient condition. Another sufficient condition is that aj≥0 under the constraint a0=0 (so that the largest eigenvalue at a0=0 is 0).This method is significantly easier to use than finding the largest eigenvalue, although verifying thata0=0 necessarily corresponds to the largest eigenvalue can become complicated for some models. However, similar to the above, rearranging a0=0 to produce R0,C=1 is not a unique process  and does not always produce the average number of secondary infections.


### 2.4. The Next-Generation Method

The next-generation method, developed by Diekmann et al. [8] and Diekmann and Heesterbeek [9], and popularised by van den Driessche and Watmough [10], is a generalisation of the Jacobian method. It is significantly easier to use than Jacobian-based methods, since it only requires the infection states (such as the exposed class, the infected class and the asymptomatically infected class) and ignores all other states (such as susceptible and recovered individuals). This keeps the size of the matrices relatively manageable.However, the next-generation method suffers from a lack of uniqueness. In order to determine the matricesF and V  (where F accounts for the “new” infections and V accounts for the transfer between infected compartments),  biological insight must be used in order to decide which terms count as “new” infections and which terms are transfer terms. While this may seem intuitive for most models, it is easy to construct a counterexample (3)I′=βSI-μI=βSI+5I-5I-μI.
Here, the term +5I might represent, for example, new infections arising from vertical transmission, whereas the term -5I might represent a disease-specific death rate. Although this construction is clearly arbitrary, it demonstrates that identifying “new” infections is not a unique process and relies on the modeller's judgement.We would then have(4)F=βS+5,V=5+μ,
and thus (5)R0,5=βS+55+μ.
This has the same threshold property as R0,N=βS/μ, but is clearly not the average number of secondary infections. In essence, the next-generation method is a mathematical generalisation of putting the negative values on one side and dividing (this is what V-1 is) so that the eigenvalue threshold at zero is transformed into an R0-like threshold at one. However, as we have seen, this does not produce a unique result.van den Driessche and Watmough [10] note that other decompositions of F and V can be chosen, which lead to different values for R0,N. They claim that only one choice of F is epidemiologically correct. However, this is not true, as the above example shows. Furthermore, it means that the definition of R0 relies upon the judgement of the modeller as to what “epidemiologically correct” means.Furthermore, the next-generation method does not produce the number of humans infected by a single human if there is an intermediate host, but rather the geometric mean of the number of infections per generation.For example, consider a mosquito-borne disease where humans infect two mosquitoes, while mosquitoes infect three humans, as shown in Figure1. For convenience, label these RH=2 and RM=3. Then the number of humans infected from a primary human (via mosquitoes) is R0=2×3=6. (This is also the value calculated by the survival function.)However, the next-generation method would calculateR0,N=6, which is a weighted average (2<6<3) of the number of infectives each individual produces in the next infection event. While mathematically sound, it is questionable whether this is biologically meaningful.This example could be extended. Consider a three-stage disease, such as tularemia, where ticks may transmit between humans and livestock, but humans may also be infected by eating livestock. Suppose that a single human directly infects two ticks (soRH=2), each tick infects four animals (so RT=4) and each animal infects three humans (so RA=3). Then, a single human has resulted in (6)R0=2×4×3=24
infected humans. However, R0,N=243≈2.88, since there are three infection stages. As the number of infection stages increase, R0,N becomes a progressively higher surd.The next-generation method is likely the most frequently used method to calculateR0. It has been used extensively to calculate R0-like values from host-vector models (see Wonham et al. [11], Gaff et al. [12] or Gubbins et al. [13]). However, it does not produce the number of newly infected individuals in the same infection class and does not always produce the average number of secondary infections.


### 2.5. The Graph-Theoretic Method

In de Camino-Beck et al. [14], a graph-theoretic method for calculating R0 is given. Starting from the definition of R0=ρ(FV-1), they derived a series of rules for reducing the digraph associated with Fλ-1-V to a digraph with zero weight, from which λ=R0. The rules are as follows.Rule.
To reduce the loop-aii<0 to −1 at node i, every arc entering i has weight divided by aii.Rule.
For a trivial nodei on a path j→i→k, the two arcs are replaced by j→k with weight equal to the product of weights on arc j→i and i→k. Weights on multiple arcs j→k are added.The graph-theoretic method has the advantage that it avoids calculation ofV-1, which may be complicated for large systems. However, it always produces the same threshold value as the next-generation method.They considered the simple vector-host model(7)dIdt=βsSW-(b+γ)I,dWdt=βmMI-cW,dSdt=b-bS+γI-βsSW,dMdt=c-cM-βmMI,
where I, W, S, and M are proportions of infective hosts, infective vectors, susceptible hosts, and susceptible vectors, respectively, b is the host birth and death rate, γ is the host recovery rate, c is the vector birth and death rate, and βs and βm are disease transmission coefficients. The disease-free equilibrium is (0,0,1,1)T. They found matrices (8)F=(0βsβm0),V=(b+γ00c).
From this, they produced a digraph of Fλ-1-V and concluded that (9)R0=βmβsc(b+γ).However, this method still contains the same issues as the next-generation method. TheR0 value calculated is not the number of humans infected by a single human, but rather the (less biologically meaningful) geometric mean of the number of humans infected by the vector and the number of vectors infected by a human. Furthermore, the creation of F and V is not a unique process, as we have seen above,  and does not always produce the average number of secondary infections.For example, consider the mathematically equivalent model(10)dIdt=(βs+9)SW-9SW-(b+γ)I,dWdt=βmMI-cW,dSdt=b-bS+γI-βsSW,dMdt=c-cM-βmMI,
with matrices (11)F=(0βs+9βm0),V=(b+γ90c).The graph reduction then proceeds as in Figure2. It follows (either by the graph-reduction method or using the conventional next-generation method) that (12)R0=12[9(βs+9)c(b+γ)+81(βs+9)2c2(b+γ)2+4βm(βs+9)c(b+γ)].
Thus, the graph-theoretic method inherits the existing problems from the next-generation method.Figure 2
The graph-theoretical method of de Camino-Beck et al. [14] applied to a vector-host model may not produce a unique R0.


### 2.6. Existence of the Endemic Equilibrium

R0 can also be calculated from the endemic equilibrium, in the case where there is a bifurcation at R0=1 such that the endemic equilibrium does not exist for R0<1. The existence of the endemic equilibrium is, thus, a threshold that can be rearranged to produce an R0-like surrogate.However, for many models, calculating the endemic equilibrium can be quite difficult. Furthermore, in the case of a backward bifurcation, the endemic equilibrium still exists forR0<1 (in fact, two endemic equilibria exist in this case, one of which is stable and the other unstable). It follows that the endemic equilibrium is not a useful general method for calculating R0  and does not always produce the average number of secondary infections.van den Bosch et al. [6] use the endemic equilibrium to determine a threshold given by Î=(1/α)(μq-μ+αK) so that Î>0 if μq-μ+αK>0 (where Î represents infected plants at the endemic equilibrium, α is the transmissibility, μ is the death rate, q is the fraction of infectious seeds, and K is the total plant population density). They show that two different rearrangements of this inequality give (13)αμ(1-q)K>1orαμK+q>1,
and note that either would suffice as an R0 value, but were unable to resolve the question of which was appropriate.


### 2.7. Summary ofR0 Methods

In summary, there are many methods available for calculatingR0, but few of them agree with each other and almost none reliably calculate the average number of secondary infections in a wholly susceptible population. The only method which does is the survival function, but this method is too cumbersome to be used except for the simplest of models.


## 2.1. The Survival Function

The survival function is given by(1)R0,S=∫0∞(the average number of susceptibleindividuals that an infected individualnewly infects per unit timewhen infected for total timea)×(the probability that anewly infected individualremains infectiousfor at least timea)daThe survival function has the advantage that it always produces the average number of secondary individuals infected by a single infected individual, in the same class. Thus, in Figure1, where one human infects two mosquitoes, who each infect three humans, the survival function produces R0=6. This is the number of humans infected by a single infected human via mosquitoes; or, equivalently, the number of mosquitoes infected by a single mosquito via humans. The survival function is a generalised method of calculating the basic reproductive ratio that is not restricted to ODEs.Figure 1
A example of a two-stage basic reproductive ratio. A single human infectsRH=2 mosquitoes, each of whom subsequently infects RM=3 humans. Thus, a single human results in R0=6 infected humans. The next-generation method instead estimates R0,N=6 infected individuals in the subsequent generation, regardless of whether that generation is human or mosquito; R0,N is the geometric mean of RH=2 and RM=3. (Reproduced with permission from Smith? [7]).However, determining the individual probabilities can be cumbersome, especially if multiple states are involved. For a vector-borne infection such as malaria, with two infection states (human and mosquito), calculating the first probability involves determining the probability that a human infected at time 0 exists at timet, the probability that a human infected for total time t infects a mosquito and the probability that an infected mosquito lives to be age a-t, where 0≤t≤a [2]. For diseases with more states, such as Guinea Worm disease, where there is a waterborne parasite, which can attach itself to copepods, which in turn are ingested by humans and which subsequently grow into an internal nematode, the calculations of the survival probabilities become unwieldy.Thus, although this method always produces the correctR0, in practice, it is difficult to use. This is especially true for models with sufficient complexity, which are often those encountered most frequently.


## 2.2. The Jacobian

The Jacobian matrix is used to linearise a nonlinear system of differential equations. Around the disease-free equilibrium, the linear system will have the same stability properties as the nonlinear system if it is hyperbolic; that is, if no eigenvalues have zero real part. In particular, if all eigenvalues have negative real part, then the equilibrium is stable, whereas if there is an eigenvalue with positive real part, the equilibrium is unstable.It follows that a threshold isλmax=0, where λmax is the largest eigenvalue of the Jacobian matrix (or the largest real part if the eigenvalues are complex). However, for a system of n differential equations, this requires solving an nth order polynomial, which may be impossible. Furthermore, rearranging the condition λmax=0 to produce a threshold R0,J=1 is not a unique process  and does not always produce the average number of secondary infections.


## 2.3. Constant Term of the Characteristic Polynomial

Whenλmax=0, the constant term of the characteristic polynomial will be zero. However, the reverse is not true, as the polynomial could have both zero and positive roots. If the characteristic polynomial is (2)λn+an-1λn-1+⋯+a1λ+a0=0,
then a0=0 is a threshold if aj≥0 for all j. More generally, the Routh-Hurwitz condition allows the coefficients to take on other signs, under certain restrictions, but the nonconstant coefficients all being positive is a sufficient condition. Another sufficient condition is that aj≥0 under the constraint a0=0 (so that the largest eigenvalue at a0=0 is 0).This method is significantly easier to use than finding the largest eigenvalue, although verifying thata0=0 necessarily corresponds to the largest eigenvalue can become complicated for some models. However, similar to the above, rearranging a0=0 to produce R0,C=1 is not a unique process  and does not always produce the average number of secondary infections.


## 2.4. The Next-Generation Method

The next-generation method, developed by Diekmann et al. [8] and Diekmann and Heesterbeek [9], and popularised by van den Driessche and Watmough [10], is a generalisation of the Jacobian method. It is significantly easier to use than Jacobian-based methods, since it only requires the infection states (such as the exposed class, the infected class and the asymptomatically infected class) and ignores all other states (such as susceptible and recovered individuals). This keeps the size of the matrices relatively manageable.However, the next-generation method suffers from a lack of uniqueness. In order to determine the matricesF and V  (where F accounts for the “new” infections and V accounts for the transfer between infected compartments),  biological insight must be used in order to decide which terms count as “new” infections and which terms are transfer terms. While this may seem intuitive for most models, it is easy to construct a counterexample (3)I′=βSI-μI=βSI+5I-5I-μI.
Here, the term +5I might represent, for example, new infections arising from vertical transmission, whereas the term -5I might represent a disease-specific death rate. Although this construction is clearly arbitrary, it demonstrates that identifying “new” infections is not a unique process and relies on the modeller's judgement.We would then have(4)F=βS+5,V=5+μ,
and thus (5)R0,5=βS+55+μ.
This has the same threshold property as R0,N=βS/μ, but is clearly not the average number of secondary infections. In essence, the next-generation method is a mathematical generalisation of putting the negative values on one side and dividing (this is what V-1 is) so that the eigenvalue threshold at zero is transformed into an R0-like threshold at one. However, as we have seen, this does not produce a unique result.van den Driessche and Watmough [10] note that other decompositions of F and V can be chosen, which lead to different values for R0,N. They claim that only one choice of F is epidemiologically correct. However, this is not true, as the above example shows. Furthermore, it means that the definition of R0 relies upon the judgement of the modeller as to what “epidemiologically correct” means.Furthermore, the next-generation method does not produce the number of humans infected by a single human if there is an intermediate host, but rather the geometric mean of the number of infections per generation.For example, consider a mosquito-borne disease where humans infect two mosquitoes, while mosquitoes infect three humans, as shown in Figure1. For convenience, label these RH=2 and RM=3. Then the number of humans infected from a primary human (via mosquitoes) is R0=2×3=6. (This is also the value calculated by the survival function.)However, the next-generation method would calculateR0,N=6, which is a weighted average (2<6<3) of the number of infectives each individual produces in the next infection event. While mathematically sound, it is questionable whether this is biologically meaningful.This example could be extended. Consider a three-stage disease, such as tularemia, where ticks may transmit between humans and livestock, but humans may also be infected by eating livestock. Suppose that a single human directly infects two ticks (soRH=2), each tick infects four animals (so RT=4) and each animal infects three humans (so RA=3). Then, a single human has resulted in (6)R0=2×4×3=24
infected humans. However, R0,N=243≈2.88, since there are three infection stages. As the number of infection stages increase, R0,N becomes a progressively higher surd.The next-generation method is likely the most frequently used method to calculateR0. It has been used extensively to calculate R0-like values from host-vector models (see Wonham et al. [11], Gaff et al. [12] or Gubbins et al. [13]). However, it does not produce the number of newly infected individuals in the same infection class and does not always produce the average number of secondary infections.


## 2.5. The Graph-Theoretic Method

In de Camino-Beck et al. [14], a graph-theoretic method for calculating R0 is given. Starting from the definition of R0=ρ(FV-1), they derived a series of rules for reducing the digraph associated with Fλ-1-V to a digraph with zero weight, from which λ=R0. The rules are as follows.Rule.
To reduce the loop-aii<0 to −1 at node i, every arc entering i has weight divided by aii.Rule.
For a trivial nodei on a path j→i→k, the two arcs are replaced by j→k with weight equal to the product of weights on arc j→i and i→k. Weights on multiple arcs j→k are added.The graph-theoretic method has the advantage that it avoids calculation ofV-1, which may be complicated for large systems. However, it always produces the same threshold value as the next-generation method.They considered the simple vector-host model(7)dIdt=βsSW-(b+γ)I,dWdt=βmMI-cW,dSdt=b-bS+γI-βsSW,dMdt=c-cM-βmMI,
where I, W, S, and M are proportions of infective hosts, infective vectors, susceptible hosts, and susceptible vectors, respectively, b is the host birth and death rate, γ is the host recovery rate, c is the vector birth and death rate, and βs and βm are disease transmission coefficients. The disease-free equilibrium is (0,0,1,1)T. They found matrices (8)F=(0βsβm0),V=(b+γ00c).
From this, they produced a digraph of Fλ-1-V and concluded that (9)R0=βmβsc(b+γ).However, this method still contains the same issues as the next-generation method. TheR0 value calculated is not the number of humans infected by a single human, but rather the (less biologically meaningful) geometric mean of the number of humans infected by the vector and the number of vectors infected by a human. Furthermore, the creation of F and V is not a unique process, as we have seen above,  and does not always produce the average number of secondary infections.For example, consider the mathematically equivalent model(10)dIdt=(βs+9)SW-9SW-(b+γ)I,dWdt=βmMI-cW,dSdt=b-bS+γI-βsSW,dMdt=c-cM-βmMI,
with matrices (11)F=(0βs+9βm0),V=(b+γ90c).The graph reduction then proceeds as in Figure2. It follows (either by the graph-reduction method or using the conventional next-generation method) that (12)R0=12[9(βs+9)c(b+γ)+81(βs+9)2c2(b+γ)2+4βm(βs+9)c(b+γ)].
Thus, the graph-theoretic method inherits the existing problems from the next-generation method.Figure 2
The graph-theoretical method of de Camino-Beck et al. [14] applied to a vector-host model may not produce a unique R0.


## 2.6. Existence of the Endemic Equilibrium

R0 can also be calculated from the endemic equilibrium, in the case where there is a bifurcation at R0=1 such that the endemic equilibrium does not exist for R0<1. The existence of the endemic equilibrium is, thus, a threshold that can be rearranged to produce an R0-like surrogate.However, for many models, calculating the endemic equilibrium can be quite difficult. Furthermore, in the case of a backward bifurcation, the endemic equilibrium still exists forR0<1 (in fact, two endemic equilibria exist in this case, one of which is stable and the other unstable). It follows that the endemic equilibrium is not a useful general method for calculating R0  and does not always produce the average number of secondary infections.van den Bosch et al. [6] use the endemic equilibrium to determine a threshold given by Î=(1/α)(μq-μ+αK) so that Î>0 if μq-μ+αK>0 (where Î represents infected plants at the endemic equilibrium, α is the transmissibility, μ is the death rate, q is the fraction of infectious seeds, and K is the total plant population density). They show that two different rearrangements of this inequality give (13)αμ(1-q)K>1orαμK+q>1,
and note that either would suffice as an R0 value, but were unable to resolve the question of which was appropriate.


## 2.7. Summary ofR0 Methods

In summary, there are many methods available for calculatingR0, but few of them agree with each other and almost none reliably calculate the average number of secondary infections in a wholly susceptible population. The only method which does is the survival function, but this method is too cumbersome to be used except for the simplest of models.


## 3. A Worked Example

In this section, we take a sample model and apply the various methods for calculatingR0 to it. We show that each method can produce a different result, all of which have the property that they have an outbreak threshold at R0=1 but otherwise bear little relation to one another.Consider the following model for malaria:(14)HS′=ΛH-βMHMIHS-μHHS,HI′=βMHMIHS-(μH+σ)HI,MS′=ΛM-βHMMSHI-μMMS,MI′=βHMMSHI-μMMI.
Humans may be susceptible or infected (HS and HI, resp.), while mosquitoes may be susceptible or infected (MS and MI, resp.). The birth rates are ΛH for humans and ΛM for mosquitoes. The background death rates are μH and μM for humans and mosquitoes, respectively, while the disease-specific death rate for humans is σ. The transmission rate from mosquitoes to humans is βMH, while the transmission rate from humans to mosquitoes is βHM. See Figure 3.Figure 3
A standard model of malaria. Humans can be susceptible or infected, with birth rateΛH, background death rate μH and disease-specific death rate σ. Mosquitoes can be susceptible or infected, with birth rate ΛM and background death rate μM.The Jacobian matrix at the disease-free equilibrium is(15)J=[-μH00-βMHH¯S0-μH-σ0βMHH¯S0-βHMM¯S-μM00βHMM¯S0-μM],
where H¯S and M¯S are the equilibrium values of susceptible humans and mosquitoes, respectively. Thus, the nontrivial part of the characteristic polynomial satisfies (16)λ2+(μH+σ+μM)λ+μM(μH+σ)-βHMβMHH¯SM¯S=0.Using the Jacobian method, the largest eigenvalue is(17)λmax⁡=12[(μM-σ-μH)2+4βHMβMHH¯SM¯S-μH-σ-μM+(μM-σ-μH)2+4βHMβMHH¯SM¯S].
Rearranging λmax=0, we have(18)R0,J=(μM-σ-μH)2+4βHMβMHH¯SM¯S(μM+σ+μM)2.Conversely, since the nonconstant coefficients of the characteristic polynomial are positive, all eigenvalues will be negative if the constant term of the characteristic polynomial is positive, whereas there will be an eigenvalue with positive real part if the constant term is negative. Rearranging, we have(19)R0,C=βMHβHMH¯SM¯S(μH+σ)μM.This is not the only rearrangement possible, so, like the nonuniqueness of the next-generation method, it is possible to rearrange by adding and subtracting arbitrary constants. Thus, for example,(20)R0,9=βMHβHMH¯SM¯S+9(μH+σ)μM+9
is also a measure of disease spread with the property that the disease persists if R0,9>1 and is eradicated if R0,9<1.Other rearrangements are possible, so long as the threshold at 0 is transformed into a threshold at 1. Thus,(21)R0,e=exp⁡(βMHβHMH¯SM¯S-(μH+σ)μM)
has the same threshold property. A generalised version of this formulation was used in Smith? et al. [15] called T0.The endemic equilibrium satisfies(22)M̂I=βHMM̂SĤIμM,M̂S=ΛMβHMĤI+μM,ĤS=(μM+σ)μM(βHMĤI+μM)βMHβHMΛM,ĤI=ΛHΛMβMHβHM-μH(μH+σ)μM2ΛMβMHβHM(μH+σ)+βHMμHμM(μH+σ).
It follows that there will be a biologically meaningful endemic equilibrium if ĤI>0, or (23)R0,E=βMHβHMH¯SM¯S(μH+σ)μM>1,
since M¯S=ΛM/μM and H¯S=ΛH/μH. As above, this is not the only rearrangement possible.The next-generation matrices are(24)F=[0βMHH¯SβHMM¯S0],V=[μH+σ00μM],
since the next-generation method only considered the two infected classes. Then, using the next-generation method,(25)R0,N=βMHβHMH¯SM¯S(μH+σ)μM.Thus, for the same model,R0,N, R0,J, R0,C, R0,9, and R0,e are all distinct. Although R0,E=R0,C in this case, this does not hold in general. Note that R0,C and R0,E produce the number of infected humans per infected human (or, equivalently, the number of infected mosquitoes per infected mosquito), but this is also not necessarily true in general. Figure 4 illustrates how these R0-like expressions vary with σ, when all other parameters are held constant.Figure 4
The changing face ofR0. All expressions for R0 were generated from the same model, using standard techniques. As the death rate varies, all have the property that if R0>1 then an outbreak will occur, whereas if R0<1, then the disease will be eradicated. However, aside from the common threshold at R0=1, none of these expressions are equal to each other. The curve labelled “Jacobian” illustrates R0,J, as given by (18). The curve labelled “Adding and subtracting 9” illustrates the value R0,9 as given by (20). The curve labelled “Next  generation” illustrates R0,N, as given by (25). The curve labelled “Constant term of the characteristic polynomial” illustrates R0,C, as given by (19). The curve labelled “Exponential” illustrates R0,e, as given by (21). In this case, the values are H¯S=45, M¯S=30, βMH=0.05, βHM=0.03, μH=0.15 and μM=0.6.Anecdotally, we have found that most biologists believe that there is one and only oneR0 value for a given disease and, in order to verify an R0 is correct, all that is required is to check that it has a threshold at 1. We have demonstrated here that R0 is not a unique concept. Indeed, it is straightforward to construct simple variations that satisfy the threshold criterion at 1. Furthermore, since multiple methods calculate multiple thresholds, it follows that the vast majority of them cannot be the average number of secondary infections. However, the nonuniqueness is not the only problem associated with the concept, as we will demonstrate.


## 4. Backward Bifurcations

Backward bifurcations occur when multiple stable equilibria coexist forR0<1. This presents a serious complication when a disease is already endemic, since lowering the basic reproduction number below 1 may no longer be a viable control measure; hence, different prevention and control measures may have to be considered. In particular, a backward bifurcation makes the system more complicated, since the behaviour now depends on the initial conditions.A backward bifurcation atR0=1 may result in persistence of the disease when R0<1. In this case, the disease will always persist for R0>1. However, there is a point Ra<1 such that the endemic equilibrium exists for Ra<R0<1 and a third, unstable, equilibrium also exists between the two stable equilibria. Hence, an endemic disease is only eradicated if 0<R0<Ra. For Ra<R0<1, the outcome depends on initial conditions. If the disease is still in its early stages—that is, if the initial conditions are sufficiently small—then the system will approach the disease-free equilibrium and the disease will be eradicated. However, if the initial conditions are large, then the system will approach the endemic equilibrium and the disease will persist. Thus, a backward bifurcation prevents the system from switching to the disease-free equilibrium as soon as R0<1 is reached. See Figure 5(a).The effects of backward bifurcations. Solid curves indicate stable equilibria, while dashed curves indicate unstable equilibria. (a) A backward bifurcation atR0=1 may result in persistence of the disease when R0<1. There is a point Ra<1 such that the endemic equilibrium exists for Ra<R0<1 and a third, unstable, equilibrium also exists. Hence, the disease-free equilibrium is only globally stable if 0<R0<Ra. (b) Backward bifurcations at other points may also affect the outcome. Although the disease persists for all R0>1 and is eradicated when R0<1 (due to the transcritical bifurcation at R0=1), there is a region Rm<R0<Rn, where three equilibria coexist. In this region, the outcome depends on the initial conditions.
(a)(b)Several mechanisms have been shown to lead to backward bifurcation in epidemic models, but backward bifurcations in compartmental models have only recently attracted serious research attention.Gómez-Acevedo and Li  [16] investigated a mathematical model for human T-cell lymphotropic virus type I (HTLV-I) infection of CD4+ T cells that incorporates both horizontal transmission through cell-to-cell contact and vertical transmission through mitotic division of infected T cells. They assumed that a fraction σ of the infected cells survive the immune system attack after the error-prone viral replication. Under the biologically sound assumptions that the fraction σ should be very low and the rate of the mitotic division should be high, their model has a bifurcation that predicts persistent infection for an extended range of the basic reproduction R0>R0(σ0), where R0(σ0)<1. This model undergoes a backward bifurcation as σ increases: multiple stable equilibria exist for an open set of parameter values, when the basic reproduction number is below one.Safan et al. [17] studied an epidemiological model under the assumption that the susceptibility after a primary infection is r times the susceptibility before a primary infection. They present a method for determining the control effort required to eliminate an infection from a host population when subcritical persistence may occur. This effort can be interpreted as a reproduction number but is not necessarily the basic reproduction number.Forr>1+μ/α, this model exhibits backward bifurcations, where μ is the death rate and α is the recovery rate. For such models, the authors presented a method for determining the minimum effort required to eradicate the infection from the endemic steady state if one concentrates on control measures affecting the transmission rate constant.Garba et al. [18] presented a deterministic model for the transmission dynamics of a single strain of dengue by realistically adopting a standard incidence formulation and allowing dengue transmission by exposed humans and vectors. The model was extended to include an imperfect vaccine for dengue. A backward bifurcation was observed in both models. This makes R0<1 no longer sufficient for effectively controlling dengue in a community. However, this phenomenon can be removed by replacing the standard incidence function in the model with a mass-action formulation.Reluga et al. [19] proposed a series of epidemic models for waning immunity that can be applied in many different settings. With biologically realistic hypotheses, they found that immunity alone never creates a backward bifurcation. However, this does not rule out the possibility of multiple stable equilibria, which can be shown by a counterexample of the forward bifurcation at R0=1. See Figure 5(b).


## 5.R0 in Spatial Contexts

WhenR0 is considered in a spatial context, many of its properties fail to hold. In particular, diseases with R0>1 can fail to persist, depending on the nature of the spatial transmission. Since many diseases are spatially dependent, this further limits the utility of R0. As R0 increases beyond 1, the probability of disease invading the initially infected host group increases, but additional criteria are important to determining the probability of the spreading of the disease to other groups.


### 5.1. Networks

Green et al. [20] related deterministic mean-field models to network models, taking into account the manner in which the contact rate and infectiousness change over time. For a node with exactly m connections, the expected number of secondary infections at infection age u is given by (26)R(m,u)=m(1-e-βu/k),
where β is the transmissibility and k is the mean number of connections per node. This model applies when the rate of infectious contact is independent of k.If there is a constant rate of generation of new cases, then the expected number of secondary infections in an infectious period of lengthu is given implicitly by (27)R(m,u)=∫0u(m-R)CIS(t)ψ(t)dt,
where CIS(t) denotes the contact per susceptible neighbour and ψ(t) denotes the infectiousness at time t after infection.Kao [21] showed that novel pathogens may evolve towards a lower R0, even if this results in pathogen extinction. This is because the presence of exploitable heterogeneities, such as high variance in the number of potentially infectious contacts, increases R0; thus, pathogens that can exploit heterogeneities in the contact structure have an advantage over those that do not. The exploitation of heterogeneities results in a more rapid depletion of the potentially susceptible neighbourhood for an infected host. While the low R0 strategy is never evolutionarily stable, invading strains with higher R0 will also converge to the low R0 strategy if not sufficiently different from the resident strain. This is in contrast to the conventional belief that the emergence of novel pathogens is driven by maximisation of R0.In a randomly mixed epidemiological network,R0 can be approximated by(28)R0=〈linlout〉〈lout〉,
where lin and lout are, respectively, the number of inward and outward “truly infectious” links per node; the angled brackets represent the expected value of the relevant quantity [22].Meyers [23] showed that in a contact network framework, (29)R0=T(〈k2〉-〈k〉〈k〉),
where T is the mean probability of transmission between individuals and 〈k〉 and 〈k2〉 are the mean degree and mean square degree of the network. Here, R0 depends explicitly on the structure of the network (i.e., on 〈k〉 and 〈k2〉). A single pathogen may, therefore, have very different transmission dynamics depending on the population through which it spreads. If two networks have the same mean degree, k, then the one with the larger variance in degree, 〈k2〉-〈k〉2, will be more vulnerable to the spread of disease.In compartment models, infected hosts are assumed to have potentially disease-causing contacts with random individuals from the population according to a Poisson process that yields an average contact rate ofβ per unit time. The mass-action assumption of compartmental models is tantamount to assuming that the underlying contact patterns form a random graph with a Poisson degree distribution. Estimates of R0 that assume a mass-action model may, therefore, be invalid for populations with non-Poisson contact patterns and, in particular, will underestimate the actual growth rate of the disease in highly heterogeneous networks.


### 5.2. Individual-Level Models

Rahmandad and Sterman [24] noted that if R0>1 in an agent-based model then, due to the stochastic nature of interactions, it is possible that no epidemic occurs or that it ends early if, by chance, the few initially contagious individuals recover before generating new cases.Schimit and Monteiro [25] showed that in an individual-level model, R0 cannot be uniquely determined from some features of transient behaviour of the infective group. The value of R0 can be unambiguously determined from the asymptotical stable stationary concentrations, but this relies on waiting for the system to reach its permanent regime, which is not feasible in practice. The same value of R0 can be associated to networks with distinct values of clustering coefficients and average shortest path length. This result can affect the evaluation of the effectiveness concerning different strategies employed for controlling a disease. Because distinct values of topological properties can produce the same value of R0 in a model considering the spatial structure of the contact network, it is difficult to evaluate the effective contribution of each control measure. This is because the correspondences among R0 and the topological properties of the contact network are not one-to-one.


### 5.3. Metapopulation Models

Cross et al. [26] showed that when R0 is based on data collected from simulated epidemics mimicking epidemiological contact-tracing data, R0 can be substantially greater than one and yet not cause a pandemic. In populations with social or spatial structure, a chronic disease is more likely to invade than an acute disease with the same R0, because it persists longer within each group and allows for more host movement between groups.Under the settings where the rate of host population turnover was negligible relative to the rate of disease processes of infection and recovery, they showed thatR0>1 was insufficient for disease invasion when the product of the average group size and the expected number of between-group movements made by each individual while infectious was less than 1.Smith? et al. [15] examined a metapopulation model with travel between two regions, with reproductive ratios R0,1(0) and R0,2(0) for each region in the absence of travel, and R¯0,1 and R¯0,2 when only susceptibles travel, but infectives do not. They showed that if both R0,1(0)<1 and R0,2(0)<1, then there are conditions on the travel of susceptible such that R¯0,1>1 and R¯0,2<1. Thus, a disease which would otherwise be eradicated in both regions could be sustained in one of the regions if there were sufficient travel of the susceptibles (not infectives). Furthermore, if R0,1(0)<1 and R0,2(0)>1, then there are conditions on the travel of susceptibles such that R¯0,1>1 and R¯0,2>1. Thus, if one region sustains the disease on its own, while the other does not, then sufficient travel of susceptibles (not infectives) could sustain the disease in both regions.


### 5.4. Partial Differential Equation Models

Althaus et al. [27] examined an age-dependent partial differential equation model of in-host HIV infection. They showed that (30)R0=1∫0∞e-rag(a)da,
where the denominator is the Laplace transform of the generation time distribution g(a) and r is the growth rate. They found that estimates for R0 were generally smaller than those derived from the standard model when the generation time was taken into account.


## 5.1. Networks

Green et al. [20] related deterministic mean-field models to network models, taking into account the manner in which the contact rate and infectiousness change over time. For a node with exactly m connections, the expected number of secondary infections at infection age u is given by (26)R(m,u)=m(1-e-βu/k),
where β is the transmissibility and k is the mean number of connections per node. This model applies when the rate of infectious contact is independent of k.If there is a constant rate of generation of new cases, then the expected number of secondary infections in an infectious period of lengthu is given implicitly by (27)R(m,u)=∫0u(m-R)CIS(t)ψ(t)dt,
where CIS(t) denotes the contact per susceptible neighbour and ψ(t) denotes the infectiousness at time t after infection.Kao [21] showed that novel pathogens may evolve towards a lower R0, even if this results in pathogen extinction. This is because the presence of exploitable heterogeneities, such as high variance in the number of potentially infectious contacts, increases R0; thus, pathogens that can exploit heterogeneities in the contact structure have an advantage over those that do not. The exploitation of heterogeneities results in a more rapid depletion of the potentially susceptible neighbourhood for an infected host. While the low R0 strategy is never evolutionarily stable, invading strains with higher R0 will also converge to the low R0 strategy if not sufficiently different from the resident strain. This is in contrast to the conventional belief that the emergence of novel pathogens is driven by maximisation of R0.In a randomly mixed epidemiological network,R0 can be approximated by(28)R0=〈linlout〉〈lout〉,
where lin and lout are, respectively, the number of inward and outward “truly infectious” links per node; the angled brackets represent the expected value of the relevant quantity [22].Meyers [23] showed that in a contact network framework, (29)R0=T(〈k2〉-〈k〉〈k〉),
where T is the mean probability of transmission between individuals and 〈k〉 and 〈k2〉 are the mean degree and mean square degree of the network. Here, R0 depends explicitly on the structure of the network (i.e., on 〈k〉 and 〈k2〉). A single pathogen may, therefore, have very different transmission dynamics depending on the population through which it spreads. If two networks have the same mean degree, k, then the one with the larger variance in degree, 〈k2〉-〈k〉2, will be more vulnerable to the spread of disease.In compartment models, infected hosts are assumed to have potentially disease-causing contacts with random individuals from the population according to a Poisson process that yields an average contact rate ofβ per unit time. The mass-action assumption of compartmental models is tantamount to assuming that the underlying contact patterns form a random graph with a Poisson degree distribution. Estimates of R0 that assume a mass-action model may, therefore, be invalid for populations with non-Poisson contact patterns and, in particular, will underestimate the actual growth rate of the disease in highly heterogeneous networks.


## 5.2. Individual-Level Models

Rahmandad and Sterman [24] noted that if R0>1 in an agent-based model then, due to the stochastic nature of interactions, it is possible that no epidemic occurs or that it ends early if, by chance, the few initially contagious individuals recover before generating new cases.Schimit and Monteiro [25] showed that in an individual-level model, R0 cannot be uniquely determined from some features of transient behaviour of the infective group. The value of R0 can be unambiguously determined from the asymptotical stable stationary concentrations, but this relies on waiting for the system to reach its permanent regime, which is not feasible in practice. The same value of R0 can be associated to networks with distinct values of clustering coefficients and average shortest path length. This result can affect the evaluation of the effectiveness concerning different strategies employed for controlling a disease. Because distinct values of topological properties can produce the same value of R0 in a model considering the spatial structure of the contact network, it is difficult to evaluate the effective contribution of each control measure. This is because the correspondences among R0 and the topological properties of the contact network are not one-to-one.


## 5.3. Metapopulation Models

Cross et al. [26] showed that when R0 is based on data collected from simulated epidemics mimicking epidemiological contact-tracing data, R0 can be substantially greater than one and yet not cause a pandemic. In populations with social or spatial structure, a chronic disease is more likely to invade than an acute disease with the same R0, because it persists longer within each group and allows for more host movement between groups.Under the settings where the rate of host population turnover was negligible relative to the rate of disease processes of infection and recovery, they showed thatR0>1 was insufficient for disease invasion when the product of the average group size and the expected number of between-group movements made by each individual while infectious was less than 1.Smith? et al. [15] examined a metapopulation model with travel between two regions, with reproductive ratios R0,1(0) and R0,2(0) for each region in the absence of travel, and R¯0,1 and R¯0,2 when only susceptibles travel, but infectives do not. They showed that if both R0,1(0)<1 and R0,2(0)<1, then there are conditions on the travel of susceptible such that R¯0,1>1 and R¯0,2<1. Thus, a disease which would otherwise be eradicated in both regions could be sustained in one of the regions if there were sufficient travel of the susceptibles (not infectives). Furthermore, if R0,1(0)<1 and R0,2(0)>1, then there are conditions on the travel of susceptibles such that R¯0,1>1 and R¯0,2>1. Thus, if one region sustains the disease on its own, while the other does not, then sufficient travel of susceptibles (not infectives) could sustain the disease in both regions.


## 5.4. Partial Differential Equation Models

Althaus et al. [27] examined an age-dependent partial differential equation model of in-host HIV infection. They showed that (30)R0=1∫0∞e-rag(a)da,
where the denominator is the Laplace transform of the generation time distribution g(a) and r is the growth rate. They found that estimates for R0 were generally smaller than those derived from the standard model when the generation time was taken into account.


## 6. Stochastic Effects

When stochastic effects, such as those inevitably found in nature, are included, the threshold atR0=1 may be disturbed. This includes assumptions about the distribution of transition times (assumed to be exponential in most models) as well as variations in individual parameters.Heffernan and Wahl [28] derived improved estimates of R0 for situations in which information about the dispersal of transition times is available to the clinical or epidemiological practitioner. Rather than rederiving R0 for a number of models (SIR, SEIR, etc.), they introduce a “correction factor”, ϕ: the ratio of R0 when the lifetimes are nonexponentially distributed to the value of R0 that would be calculated assuming exponential lifetimes. They were able to derive limiting values of ϕ and used this to gauge the sensitivity of R0 to dispersion in the underlying distributions.By combining the movement of hosts, transmission with groups, recovery from infection and the recruitment of new susceptibles, Cross et al. [29] expanded the earlier analysis of Cross et al. [26] to a much more broader set of disease-host relationships, exploring settings where the duration of immunity ranges from transient to lifelong or where the demographic processes occur on comparable (or faster) timescales to disease processes. The focus of this study was to investigate how recruitment of susceptibles affects disease invasion and how population structure can affect the frequency of superspreading events (SSEs). They found that the frequency of SSEs may decrease with the reduced movement and the group sizes due to the limited number of susceptibles available.The hierarchical nature of disease invasion in host metapopulations is illustrated by the classification tree analysis of the model results. Firstly, the pathogen must effectively transit within a group (R0>1), and then, the pathogen must persist within a group long enough to allow for movement between different groups. Hence, the infectious period, group size and recruitment of new susceptibles are as important as the local transmission rates in predicting the spread of pathogens across a metapopulation. It should be noted that in 35% of simulations when R0 was greater than one, the disease failed to invade.Tildesley and Keeling [30] examined whether R0 was a good predictor of the 2001 UK Foot and Mouth disease. They concluded that R0 explained just 29.3% of the standard deviation of the epidemic impact.  They also noted that R0=1 did not act as a threshold: at the value of R0=1, only 20% of initial seedings generated epidemics; this probability increased to around 50% for the largest reasonable R0 values. When heterogeneities exist in the population, infection is most likely to become focused within the high-risk individuals who are both more susceptible and more infectious. This highlights the stochastic nature of the disease in its early stages and the dependence of the ensuing epidemic on favourable local conditions in the neighbourhood of the initial infection.The probability of extinction, assuming exponentially distributed infectious periods, was(31)pextexp⁡=1R0,
when an epidemic began with a single infected individual. Thus, if R0≤1, then extinction occurs with probability 1, but if R0>1 then extinction occurs with some probability. See Chiang [31, Chapter 4] for more discussion.


## 7.R0 Failures

In this section, we note a variety of problems withR0 that are not covered in the previous sections. These include problems with the underlying structure of compartment models, the mismatch between an individual-based parameter and a population-level compartment model, and the failure of R0 to accurately measure an outbreak of a new disease.Breban et al. [32] argued that in order to associate an R0 to a model of ODEs, an individual-level model (ILM) which is compatible to the ODE model must be developed; only then can the R0 of the ILM be unambiguously calculated. These ILMs are growing (not static) network models, with individuals added to a network of who infected whom based on global or local network rules. Then, R0 is computed as the limit of the average number of outgoing links of individuals in a node that no longer accepts new links, as time goes to infinity. They showed that a broad range of R0 values were compatible with a given ODE model.For example, consider the basic model(32)dSdt=-βI,dIdt=βI-μI,
where β is the transmissibility and μ is the disease death rate; note that the transmission term is equivalent to the standard term βSI/(S+I) when the depletion of susceptibles is negligible so that S/(S+I)≈1.The expectedR0 value from this model using other methods (such as the next-generation method) is (33)R0=βμ.
The corresponding ILM consists of an infection rule, where an individual joining the infectious pool is infected by an infectious individual who is uniformly randomly selected, and a removal rule, where a uniformly randomly selected individual leaves the infectious pool. The flow of newly infected individuals is βI(t). Thus, the flow per already-infected individual is β. Since the removed individuals are randomly sampled from the infectious individuals, the average length of the infectious period equals the time expectation of the infectious period, which is 1/β  (rather than 1/μ, as calculated from the next-generation method).  It follows that R0=1, which is independent of the transmission and death rates from the corresponding ODE.  Thus, in this case, R0 does not signal epidemic growth as anticipated from other methods.Roberts [3] noted three fundamental properties commonly attributed to R0: (i) that an endemic infection can persist only if R0>1, (ii) R0 provides a direct measure of the control effort required to eliminate the infection, and (iii) pathogens evolve to maximise their R0 value. He demonstrated that all three statements can be false. The first, as we have noted, can fail due to the presence of backward bifurcations. The second can fail when control efforts are applied unevenly across different host types (such as a high-risk and a low-risk group), since R0 is determined by averaging over all host types and does not directly determine the control effort required to eliminate infection.The third can fail when two pathogens coexist at a steady state that exists and is stable whenever both single-pathogen steady states exist but are unstable. In this case, the order in which the pathogens are established in the host population matters. The established parasite has a role in determining a modified carrying capacity and the pathogen with the largest basic reproductive ratio does not necessarily exclude the other.Breban et al. [33] showed that two individual-level models having exactly the same expectations of the corresponding population-level variables may yield different R0 values. They showed that obtaining R0 from empirical contact-tracing data collected by epidemiologists and using this R0 as a threshold parameter for a population-level model could produce misleading estimates of the infectiousness of the pathogen, the severity of an outbreak, and the strength of the medical and/or behavioural interventions necessary for control. Thus, measuring R0 through contact tracing (as generally occurs during an outbreak investigation) may not be a useful measure for determining the strength of the necessary control interventions.Many different individual-level processes can generate the same incidence and prevalence patterns. Thus, assigning a meaningfulR0 value to an ODE model without knowledge of the underlying disease transmission network may be impossible. Only an epidemic threshold parameter can be used to design control strategies. Since R0 fails to possess this threshold quality, its usefulness may be vastly overstated.Meyers [23] compared the theoretical calculation of R0 with observed SARS data from China and showed that estimates of R0 seemed incompatible. The basic reproductive rate has two critical inputs:  (1)
intrinsic properties of the pathogen that determine the transmission efficiency per contact and the duration of the infectious period,(2)
the patterns of contacts between infected and susceptible hosts in the population.While the first factor may be fairly uniform across outbreaks, the second may depend significantly on context, varying both within and among populations. The problem with the SARS estimates stems from the mass-action assumption of compartmental models; that is, that all susceptible individuals are equally likely to become infected. When this assumption does not hold, the models may yield inaccurate estimates or estimates that do not apply to all populations.R0 estimates for SARS in the field were based largely on outbreak data from a hospital and a crowded apartment building, with anomalously high rates of close contacts among individuals.The author suggested that it might be inappropriate to extrapolate estimates forR0 from specific settings such as these to the population at large. Conversely, since the population at large is also unlikely to satisfy mass-action requirements, it may also be concluded that R0 is not a meaningful estimate of disease spread.


## 8. Alternatives toR0

A number of alternatives have been offered over the years, due to recognised problems withR0. Older examples include the critical community size  [5], the group-level reproductive number, R* [34], the type reproduction number [35] and the basic depression ratio [36]. See Heffernan et al. [2] for a summary of these methods. In this section, we review some of the more recent alternatives to R0.The actual reproduction number,Ra, is defined as a product of the mean duration of infectiousness and the ratio of incidence to prevalence [37–39]. R0 coincides with Ra when the transmission probability is constant but accounts for the more general situation when the transmission probability varies as a function of infection age (as happens in diseases such as HIV/AIDS).The effective reproduction numberR(t) measures the number of secondary cases generated by an infectious case once an epidemic is underway. In the absence of control measures, R(t)=R0S(t)/N, where S(t)/N is the proportion of the population susceptible [40, 41]. The estimation of reproductive numbers is typically an indirect process because some of the parameters on which these numbers depend are difficult or impossible to quantify directly. The effective reproductive number satisfies R(t)≤R0, with equality only when the entire population is susceptible [40].The effective reproduction number is of practical interest, since it is time dependent and can account for the degree of cross-immunity from earlier outbreaks. However, since it is based on the basic reproduction number, the effective reproduction number inherits many of the issues fromR0.Breban et al. [32] proposed Q0, the average number of secondary infections over the infectious population. The average number of secondary infections of actively infectious individuals, Q0(t), is computed as the average number of outgoing links of a node in the infected compartment at time t. Then, (34)R0≡lim⁡t→∞R0(t),Q0≡lim⁡t→∞Q0(t),
when the limits exist. Unfortunately, R0(t) is never defined in the paper, limiting the usefulness of this formulation of R0. However, by analogy with Q0(t), R0(t) is the average number of outgoing links of a removed compartment at time t (Romulus Breban, personal communication).For the SI model (32), under the assumption that every infection is uniquely assigned as a secondary infection for either a removed or an infected individual, (35)Ni(t)=I(t)Q0(t)+Nr(t)R0(t),
where Ni(t)=∫0tβI(u)du is the cumulative number of infected individuals that occur in the time interval (0,t] and Nr(t)=∫0tμI(u)du is the cumulative number of removed individuals.Their definition ofR0 evaluates the average number of secondary cases over removed individuals as the distribution of secondary cases becomes stationary. This definition does not imply a particular individual-level model; it depends exclusively on the structure of the disease-transmission network [42]. However, R0 is not measured at the start of an epidemic, which may limit its usefulness during an initial outbreak.Grassly and Fraser [43] demonstrated that standard epidemiological theory and concepts such as R0 do not apply when infectious diseases are affected by seasonal changes. They instead define (36)R¯0=D∫01β(t)dt,
where β(t) is the transmission parameter at time t and D is the average duration of infection. Thus, R¯0 can be interpreted as the average number of secondary cases arising from the introduction of a single infective into a wholly susceptible population at a random time of the year.They noted that the conditionR¯0<1 is not sufficient to prevent an outbreak, since chains of transmission can be established during high-infectious seasons if Dβ(t)>1. However, R¯0<1 is both necessary and sufficient for long-term disease extinction.Meyers [23] noted that in the case of networks, estimating the average transmissibility T may be more valuable than R0. This means reporting not only the number of new infections per case, but also the total estimated number of contacts during the infectious period of that case. Given the primary role of contact tracing in infectious disease control, these data are often collected. Unlike R0, T can be justifiably extrapolated from one location to another even if the contact patterns are quite disparate.Instead ofR0, the author offered a number of alternatives to determine whether an outbreak will occur, based on network modelling. These include the probability that an individual will spark an epidemic and the probability that a disease cluster will spark an epidemic.The probability that an individual with degreek will spark an epidemic is equal to the probability that transmission along at least one of the k edges emanating from that vertex will lead to an epidemic. For any of its k edges, the probability that the disease does not get transmitted along the edge is 1-T. The probability that disease is transmitted to the attached vertex but does not proceed into a full-blown epidemic is Tu, where u is the probability that a secondary infection does not spark an epidemic. Thus, the probability that an individual will spark an epidemic is 1-(1-T+Tu)k.The probability that an outbreak of sizeN sparks an epidemic is (37)1-(∑k=1∞kpk(1-T+Tu)k-1∑k=1∞kpk),
where pk is the relative frequency of vertices of degree k in the network.Kao et al. [22] defined an epidemiological network contact matrix M whose elements mij are either 1 or 0, depending on whether an infectious contact between nodes i and j is possible. The spectral radius of M is an alternative approximation for R0, which can be calculated via a weighted version of (28).This explicitly accounts for the full contact structure of the network, but the evaluation of extremely large, reasonably dense matrices (some highly active nodes may have hundreds of potentially infectious links) is difficult and time consuming, particularly when this evaluation process must be repeated multiple times. However, comparisons between the two approximations for subsets of a sheep network with several thousand nodes show little difference betweenR0 and the spectral radius of M (typically less than 5%).Kamgang and Sallet [44] used the special structure of Metzler matrices (real, square matrices with nonnegative off-diagonal entries) to define  𝒯0, an analytical threshold condition.  𝒯0 is a function of the parameters of the system such that if  𝒯0<1, the disease-free equilibrium is locally asymptotically stable, and if  𝒯0>1, the disease-free equilibrium is unstable.  𝒯0 has an association with R0, although it is a stronger condition; however, it has no direct biological interpretation. The algorithm for deriving  𝒯0, although highly mathematical in nature, allows computation of a threshold for high-dimensional epidemic models.Huang [45] defined four reproductive numbers associated with four types of transmission patterns, each depending on z, the ratio of the mean infectious period to the mean latent period. These four reproduction numbers are the following:  (1)
R0I, the minimal reproductive number associated with the slowest latency process and the fastest recovery process,(2)
R0II, the middle reproductive number associated with mean latency and recovery processes,(3)
R0III, the maximal reproductive number associated with the fastest latency process and the slow recovery process,(4)
R0IV, the largest reproductive number associated with the fastest latency process and the extremely slow recovery process.All four reproduction numbers are strictly increasing functions ofz and satisfy (38)R0I<R0II<R0III<R0IV.
These numbers allow a disease to be classified as mild (R0I<R0<R0II) or severe (R0III<R0<R0IV).Hosack et al. [46] noted that R0 does not necessarily address the dynamics of epidemics in a model that has an endemic equilibrium. They used the concept of reactivity to derive a threshold index for epidemicity, E0, which gives the maximum number of new infections produced by an infective individual at a disease-free equilibrium. They also showed that the relative influence of parameters on E0 and R0 may differ and lead to different strategies for control.IfR0 is derived from the next-generation operator so that (39)R0=ρ(FV-1),
then the threshold for endemicity is (40)E0=ρ(F+FT2V-1),
such that the system is reactive when E0>1 and nonreactive when E0<1. E0 describes the transitory behaviour of disease following a temporary perturbation in prevalence. If the threshold for epidemicity is surpassed, then disease prevalence can increase further even when the disease is not endemic. This suggests that epidemics can occur even in areas where long-term transmission cannot be maintained.Reluga et al. [47] defined the discounted reproductive number Rd. The discounted reproductive number is a measure of reproductive success that is an individual's expected lifetime offspring production discounted by the background population growth rate. It is calculated as (41)Rd=ρ(F(δI-V)-1),
where δ is the discount rate, I is the identity matrix, and F and V are from the next-generation matrix decomposition. Rd combines properties of both the basic reproductive number and the ultimate proliferation rate although it also inherits the nonuniqueness problems from the next-generation method.Nishiura  [4] developed a likelihood-based method for estimating R0 without assuming exponential growth of cases and offers a corrected value of the actual reproduction number. The author noted that R0 is extremely sensitive to dispersal of the progression of a disease or variations in the underlying epidemiological assumptions.


## 9. Discussion

Despite being a crucial concept in disease modelling, with a long history and frequent application,R0 is deeply flawed. It is not a measure of the number of secondary infections, it is never calculated consistently, and it fails to satisfy the threshold property. Rarely has an idea so erroneous enjoyed such popular appeal. Why, then, are we so attached to the concept?The answer is thatR0, despite all its flaws, is all that we have. No other concept has so effectively transcended mathematics, biology, epidemiology, and immunology. No other concept is so general that it can be understood in terms of compartment models, network models, partial differential equations, and metapopulation models. “The number of secondary infections” has an intuitive appeal that outlasts even the inaccuracy of that statement when applied to the concept.The threshold nature ofR0 is used to monitor and control severe real-time epidemics; control measures are often concluded if R0<1 [48], making the problems with R0 more than just theoretical. Due to the inconsistencies in calculation, different diseases cannot be compared unless the same method was used to calculate R0; if HIV has an R0 of 3 and swine flu has an R0 of 4, we cannot conclude that swine flu is worse than HIV if different methods were used to determine these values. All we can conclude is that both diseases have R0>1; however, as we have demonstrated, this does not necessarily guarantee disease persistence.Of the many different methods used to calculateR0, only the survival function reliably calculates the average number of secondary infections; however, this method is too cumbersome to use in most practical settings. The next-generation method is probably the most popular, yet it suffers from uniqueness problems and does not cope well with more than one disease state. Since R0 is rarely measured in the field, it instead relies upon after-the-fact calculations to determine the strength of disease spread [2]. This limits its usefulness even further.Policy decisions are being based upon the concept, with limited understanding of the complexity and errors that exist in the very structure of the concept. Funding decisions about where money should be best spent are based on estimates ofR0, resources are directed towards one disease over another, and monitoring programmes are abandoned, their objectives only half-realised, because of R0. Lives may be saved or lost, based on this imperfect and inconsistent measure.R0 is a quantity that relates to the initial phase of an epidemic. This makes practical sense in terms of disease prevention. However, it is also used to guide eradication efforts when a disease is endemic. Some methods derive an eradication threshold from an equilibrium value that may not be attained for a very long time. This suggests that different measures are needed in different stages of an epidemic in order to characterise transmissibility and to guide intervention strategies. When a new pathogen emerges, a quantity describing the initial spread is useful. When a disease is endemic, a quantity that applies to the long-term equilibrium is more appropriate.We note that although the definition ofR0 is broad, it is not a universal quantity that applies to all settings. In different settings, one should use a quantity that satisfies the following properties: that an endemic equilibrium only persist if R0>1 and that R0 provides a direct measure of the control effort required for eradication. The contact structure of the population, the variation of risk factors and the order of establishing parasites (if applicable) should accompany the identification of a meaningful R0.What is urgently needed is a simple, but accurate, measure of disease spread that has a consistent threshold property and which can be understood by nonmathematicians. IfR0 is to be used, it must be accompanied by a declaration of which method was used, which assumptions are underlying the model (e.g. mass-action transmission) and evidence that it is actually a threshold, with no backward bifurcation.  Without such caveats, the concept of R0 will continue to fail.



---

*Source: 527610-2011-08-16.xml*