0

1. Introduction

The basic reproductive ratio—also known as the basic reproductive number, the basic reproduction number, the control reproduction number, or R

0—is one of the foremost concepts in epidemiology [1–3]. R

0 is the most widely used epidemiological measurement of the transmission potential in a given population. It is a measure of initial disease spread, such that if R

0 > 1, then the disease can invade an otherwise susceptible population and hence persist, whereas if R

0 < 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.

Despite its place at the forefront of mathematical epidemiology, the concept of R

0 is deeply flawed. Defining R

0 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 R

0 for a specific disease relies on secondary methods. There are many methods to calculate R

0 from mathematical models, few of which agree with each other and few of which produce the average number of secondary infections. Methods to calculate R

0 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. R

0 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..

Furthermore, there are many diseases that can persist with R

0 < 1, while diseases with R

0 > 1 can die out, reducing the utility of the concept as a threshold. R

0 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 R

0. 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 with R

0 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 R

0 in the literature. The number of articles that use R

0 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 R

0.

2.1. The Survival Function

&!start&0&!end&

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 time a) ×(the probability that anewly infected individualremains infectiousfor at least time a)da

The 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 Figure 1, where one human infects two mosquitoes, who each infect three humans, the survival function produces R

0 = 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.

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 time t, 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. 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 correct R

0, 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

&!start&0&!end&

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 R

0,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

&!start&0&!end&

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 a

0 = 0 is a threshold if a

j ≥ 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 a

j ≥ 0 under the constraint a

0 = 0 (so that the largest eigenvalue at a

0 = 0 is 0).

This method is significantly easier to use than finding the largest eigenvalue, although verifying that a

0 = 0 necessarily corresponds to the largest eigenvalue can become complicated for some models. However, similar to the above, rearranging a

0 = 0 to produce R

0,C = 1 is not a unique process and does not always produce the average number of secondary infections.

2.4. The Next-Generation Method

&!start&0&!end&

The next-generation method, developed by Diekmann et al. and Diekmann and Heesterbeek, and popularised by van den Driessche and Watmough, 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 matrices F 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 R

0,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 R

0-like threshold at one. However, as we have seen, this does not produce a unique result.

van den Driessche and Watmough note that other decompositions of F and V can be chosen, which lead to different values for R

0,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 R

0 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 Figure 1. For convenience, label these R

H = 2 and R

M = 3. Then the number of humans infected from a primary human (via mosquitoes) is R

0 = 2 × 3 = 6. (This is also the value calculated by the survival function.)

However, the next-generation method would calculate R0,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 (so R

H = 2), each tick infects four animals (so R

T = 4) and each animal infects three humans (so R

A = 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, R

0,N becomes a progressively higher surd.

The next-generation method is likely the most frequently used method to calculate R

0. It has been used extensively to calculate R

0-like values from host-vector models (see Wonham et al., Gaff et al. or Gubbins et al.). 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

&!start&0&!end&

In de Camino-Beck et al., a graph-theoretic method for calculating R

0 is given. Starting from the definition of R

0 = ρ(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 λ = R

0. The rules are as follows.

Rule 1To reduce the loop −a

ii < 0 to −1 at node i, every arc entering i has weight divided by a

ii.

Rule 2For a trivial node i 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 of V

−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. The R

0 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 Figure 2. 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.

2.6. Existence of the Endemic Equilibrium

&!start&0&!end&

R

0 can also be calculated from the endemic equilibrium, in the case where there is a bifurcation at R

0 = 1 such that the endemic equilibrium does not exist for R

0 < 1. The existence of the endemic equilibrium is, thus, a threshold that can be rearranged to produce an R

0-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 for R

0 < 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 R

0 and does not always produce the average number of secondary infections.

van den Bosch et al. 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>1 or αμK+q>1,

and note that either would suffice as an R

0 value, but were unable to resolve the question of which was appropriate.

2.7. Summary of R

0 Methods

&!start&0&!end&

In summary, there are many methods available for calculating R

0, 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 calculating R

0 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 R

0 = 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 (H

S and H

I, resp.), while mosquitoes may be susceptible or infected (M

S and M

I, 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.

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[−μ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 R

0,9 > 1 and is eradicated if R

0,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. called T

0.

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, R

0,N, R

0,J, R

0,C, R

0,9, and R

0,e are all distinct. Although R

0,E = R

0,C in this case, this does not hold in general. Note that R

0,C and R

0,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 R

0-like expressions vary with σ, when all other parameters are held constant.

Anecdotally, we have found that most biologists believe that there is one and only one R

0 value for a given disease and, in order to verify an R

0 is correct, all that is required is to check that it has a threshold at 1. We have demonstrated here that R

0 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 for R

0 < 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 at R

0 = 1 may result in persistence of the disease when R

0 < 1. In this case, the disease will always persist for R

0 > 1. However, there is a point R

a < 1 such that the endemic equilibrium exists for R

a < R

0 < 1 and a third, unstable, equilibrium also exists between the two stable equilibria. Hence, an endemic disease is only eradicated if 0 < R

0 < R

a. For R

a < R

0 < 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 R

0 < 1 is reached. See Figure 5(a).

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 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 R

0 > R

0(σ

0), where R

0(σ

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. 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.

For r > 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. 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 R

0 < 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. 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 R

0 = 1. See Figure 5(b).

5.1. Networks

&!start&0&!end& in Spatial Contexts

Green et al. 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 length u is given implicitly by

(27)R(m,u)=∫0u(m−R)CIS(t)ψ(t)dt,

where C

IS(t) denotes the contact per susceptible neighbour and ψ(t) denotes the infectiousness at time t after infection.

Kao showed that novel pathogens may evolve towards a lower R

0, 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 R

0; 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 R

0 strategy is never evolutionarily stable, invading strains with higher R

0 will also converge to the low R

0 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 R

0.

In a randomly mixed epidemiological network, R

0 can be approximated by

(28)R0=〈linlout〉〈lout〉,

where l

in and l

out are, respectively, the number of inward and outward “truly infectious” links per node; the angled brackets represent the expected value of the relevant quantity.

Meyers 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 〈k

2〉 are the mean degree and mean square degree of the network. Here, R

0 depends explicitly on the structure of the network (i.e., on 〈k〉 and 〈k

2〉). 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, 〈k

2〉−〈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 R

0 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

&!start&0&!end& in Spatial Contexts

Rahmandad and Sterman noted that if R

0 > 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 showed that in an individual-level model, R

0 cannot be uniquely determined from some features of transient behaviour of the infective group. The value of R

0 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 R

0 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 R

0 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 R

0 and the topological properties of the contact network are not one-to-one.

5.3. Metapopulation Models

&!start&0&!end& in Spatial Contexts

Cross et al. showed that when R

0 is based on data collected from simulated epidemics mimicking epidemiological contact-tracing data, R

0 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 R

0, 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 that R

0 > 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. examined a metapopulation model with travel between two regions, with reproductive ratios R

0,1

(0) and R

0,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 R

0,1

(0) < 1 and R

0,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 R

0,1

(0) < 1 and R

0,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

&!start&0&!end& in Spatial Contexts

Althaus et al. 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 R

0 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 at R

0 = 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 derived improved estimates of R

0 for situations in which information about the dispersal of transition times is available to the clinical or epidemiological practitioner. Rather than rederiving R

0 for a number of models (SIR, SEIR, etc.), they introduce a “correction factor”, ϕ: the ratio of R

0 when the lifetimes are nonexponentially distributed to the value of R

0 that would be calculated assuming exponential lifetimes. They were able to derive limiting values of ϕ and used this to gauge the sensitivity of R

0 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. expanded the earlier analysis of Cross et al. 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 (R

0 > 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 R

0 was greater than one, the disease failed to invade.

Tildesley and Keeling examined whether R

0 was a good predictor of the 2001 UK Foot and Mouth disease. They concluded that R

0 explained just 29.3% of the standard deviation of the epidemic impact. They also noted that R

0 = 1 did not act as a threshold: at the value of R

0 = 1, only 20% of initial seedings generated epidemics; this probability increased to around 50% for the largest reasonable R

0 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 R

0 ≤ 1, then extinction occurs with probability 1, but if R

0 > 1 then extinction occurs with some probability. See Chiang [31, Chapter 4] for more discussion.

7. R

&!start&0&!end& Failures

In this section, we note a variety of problems with R

0 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 R

0 to accurately measure an outbreak of a new disease.

Breban et al. argued that in order to associate an R

0 to a model of ODEs, an individual-level model (ILM) which is compatible to the ODE model must be developed; only then can the R

0 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, R

0 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 R

0 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 βS

I/(S + I) when the depletion of susceptibles is negligible so that S/(S + I) ≈ 1.

The expected R

0 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 R

0 = 1, which is independent of the transmission and death rates from the corresponding ODE. Thus, in this case, R

0 does not signal epidemic growth as anticipated from other methods.

Roberts noted three fundamental properties commonly attributed to R

0: (i) that an endemic infection can persist only if R

0 > 1, (ii) R

0 provides a direct measure of the control effort required to eliminate the infection, and (iii) pathogens evolve to maximise their R

0 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 R

0 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. showed that two individual-level models having exactly the same expectations of the corresponding population-level variables may yield different R

0 values. They showed that obtaining R

0 from empirical contact-tracing data collected by epidemiologists and using this R

0 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 R

0 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 meaningful R

0 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 R

0 fails to possess this threshold quality, its usefulness may be vastly overstated.

Meyers compared the theoretical calculation of R

0 with observed SARS data from China and showed that estimates of R

0 seemed incompatible. The basic reproductive rate has two critical inputs:

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. R

0 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 for R

0 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 R

0 is not a meaningful estimate of disease spread.

8. Alternatives to R

&!start&0&!end&

A number of alternatives have been offered over the years, due to recognised problems with R

0. Older examples include the critical community size, the group-level reproductive number, R*, the type reproduction number and the basic depression ratio. See Heffernan et al. for a summary of these methods. In this section, we review some of the more recent alternatives to R

0.

The actual reproduction number, R

a, is defined as a product of the mean duration of infectiousness and the ratio of incidence to prevalence [37–39]. R

0 coincides with R

a 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 number R(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) = R

0

S(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) ≤ R

0, with equality only when the entire population is susceptible.

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 from R

0.

Breban et al. proposed Q

0, the average number of secondary infections over the infectious population. The average number of secondary infections of actively infectious individuals, Q

0(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, R

0(t) is never defined in the paper, limiting the usefulness of this formulation of R

0. However, by analogy with Q

0(t), R

0(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 N

i(t) = ∫0

t

βI(u)du is the cumulative number of infected individuals that occur in the time interval (0, t] and N

r(t) = ∫0

t

μI(u)du is the cumulative number of removed individuals.

Their definition of R

0 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. However, R

0 is not measured at the start of an epidemic, which may limit its usefulness during an initial outbreak.

Grassly and Fraser demonstrated that standard epidemiological theory and concepts such as R

0 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 condition R¯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 noted that in the case of networks, estimating the average transmissibility T may be more valuable than R

0. 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 R

0, T can be justifiably extrapolated from one location to another even if the contact patterns are quite disparate.

Instead of R

0, 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 degree k 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 size N sparks an epidemic is

(37)1−(∑k=1∞kpk(1−T+Tu)k−1∑k=1∞kpk),

where p

k is the relative frequency of vertices of degree k in the network.

Kao et al. defined an epidemiological network contact matrix M whose elements m

ij 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 R

0, 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 between R

0 and the spectral radius of M (typically less than 5%).

Kamgang and Sallet 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 R

0, 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 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:

All four reproduction numbers are strictly increasing functions of z and satisfy

(38)R0I

These numbers allow a disease to be classified as mild (R

0

I < R

0 < R

0

II) or severe (R

0

III < R

0 < R

0

IV).

Hosack et al. noted that R

0 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, E

0, 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 E

0 and R

0 may differ and lead to different strategies for control.

If R

0 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 E

0 > 1 and nonreactive when E

0 < 1. E

0 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. defined the discounted reproductive number R

d. 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. R

d 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 developed a likelihood-based method for estimating R

0 without assuming exponential growth of cases and offers a corrected value of the actual reproduction number. The author noted that R

0 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, R

0 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 that R

0, 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 of R

0 is used to monitor and control severe real-time epidemics; control measures are often concluded if R

0 < 1, making the problems with R

0 more than just theoretical. Due to the inconsistencies in calculation, different diseases cannot be compared unless the same method was used to calculate R

0; if HIV has an R

0 of 3 and swine flu has an R

0 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 R

0 > 1; however, as we have demonstrated, this does not necessarily guarantee disease persistence.

Of the many different methods used to calculate R

0, 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 R

0 is rarely measured in the field, it instead relies upon after-the-fact calculations to determine the strength of disease spread. 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 of R

0, resources are directed towards one disease over another, and monitoring programmes are abandoned, their objectives only half-realised, because of R

0. Lives may be saved or lost, based on this imperfect and inconsistent measure.

R

0 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 of R

0 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 R

0 > 1 and that R

0 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 R

0.

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. If R

0 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 R

0 will continue to fail.