Commentary on the use of the reproduction number R during the COVID-19 pandemic

In the early stages of a new outbreak of an infectious disease, we can define an initial $R$ value, known as the basic reproduction number $R_0$, that is the average number of individuals infected by each infectious individual in a fully susceptible population^3–5. An outbreak resulting from one infected individual may die out within a few infection generations by chance^6,7. Otherwise, if $R_0>1$, the incidence of cases will grow exponentially, with on average $R_0^n$ cases in the $n\text{th}$ generation. Already, this simple description introduces a number of concepts and assumptions. An individual’s infection generation specifies their position in the chain of infections, the $(n-1)\text{th}$ generation infects the $n^{\text{th}}$ generation, and so on. It also assumes an underlying scenario (model) in which the average number of susceptibles infected by each infective stays the same over successive infection generations, and ignores the depletion of susceptibles. (We refer to those members of the population who are uninfected and susceptible to infection as susceptibles, and those that are infected and infectious as infectives.) The potential importance of these assumptions depends on the contact structure of the population, to which we return below.

Thus, $R_0$ (and other $R$ values to be defined later) is not just a property of the infectious agent (pathogen). It depends on demography, and whatever human behaviour is associated with the possibility of infectious contact (an effective contact is one that results in transmission if made with a susceptible, while a contact in the common sense of the word has a certain probability of transmission). For the simplest models, $R_0>1$ implies that an introduction of infection will result in an epidemic. Furthermore, if there were no interventions or changes in behaviour, then the proportion of the population infected during the entire course of an epidemic would be approximately the non-zero solution of the equation $P=1-{\mathrm{e}}^{-R_{0}P}$ (e.g. if $R_0=2$, then $P\approx 0.8$). This result is referred to as the final size equation and underscores the fact that during an epidemic it is not generally true that everybody will be infected at some point.

Individuals may vary considerably in their susceptibility to infection and in their propensity to pass it on through their biology or behaviour. Age is often an important determinant. If the population is grouped in some way, so that for instance some groups have higher $R$ values than others, then the overall outbreak is expected to grow as described by an $R_0$ that depends on all of these values, and also depends on how each group infects the others, i.e. on the $R$ values between groups as well as within them ( $R_0$ is then the dominant eigenvalue of the matrix of $R$ values^3,5,8). The first few stages of the outbreak may be atypical, depending on which group is first infected.

For the simplest mathematical model of the beginning of an outbreak, it is assumed that because only a small fraction of the population has been infected, all potential contacts are with susceptibles. This may be an unrealistic assumption because human interaction networks tend to be clustered (e.g. through households, workplaces or schools). Growth through successive generations of infection, which is the basis for defining $R_0$, does not translate simply into time because the generation interval of an infection (the time interval back from the instant when a susceptible is infected, to that when their infector was infected) is variable, and infection generations may overlap temporally. Typically, growth in the early stages is faster than the simple assumption of a fixed average generation time would suggest and this is a major problem in estimating $R_0$ from early outbreak data. In addition, the implicit assumption is that all infectives are identifiable as such. If there is a significant proportion of asymptomatic cases, an estimate of $R_0$ may be affected by the time from when an asymptomatic infective has become infected to when he/she is expected to infect susceptibles. If this timing is the same for asymptomatic and symptomatic cases, then the estimate for $R_0$ will be unaffected.

When the pandemic is well-established in a country (or region), with large numbers of cases most of which are internal to the country, an ‘effective reproduction number’ at time $t$, $R_t$ (sometimes denoted $R_e$ or $R_{eff}$), is a useful descriptor of the progress of the outbreak (Figure 1). Again, the concept is of an average of how many new cases each infectious case causes. The value of $R_t$ may be affected by interventions: typically the aim is to reduce $R_t$ below one and to as small a value as possible. For models including detailed, and therefore complex, contact networks there may be more than one way of defining $R_t$; however, definitions should always agree that the value of $R_t$ is 1 when the expected number of new infections is constant.

The relevance of the assumptions here (large numbers of cases, mostly internal to the region) is that in such circumstances we expect $R_t$ to have a fairly stable value that changes substantially over time only when interventions are introduced or cease. The definition of $R_t$ here is in terms of actual new infectious cases, i.e. excluding potentially infectious contacts with individuals who have been infected and are immune to reinfection. As the number of immune individuals grows large compared to the entire population, the spread of infections will gradually slow, because many contacts will be with immune individuals, and hence the value of $R_t$ will be reduced. The level of immunity at which $R_t=1$ is the herd immunity threshold (see the ‘Methods’ section on vaccination and herd immunity below).

If the outbreak is at a low level either because it has run its course or because of successful interventions, the definition and the use of an $R$ value are problematic (Figure 1). At low levels of prevalence, there will (as in the early stages of the outbreak) be greater statistical variability. Additionally, there are likely to be heterogeneities associated with the infection being unevenly spread among different subgroups of the population (possibly dependent on age, behaviour or geographical location⁹), with some parts of the population having had more exposure than others. There may also be local variability in interventions, and it may not be easy to allow for the effect of some cases being introductions from outside the population under consideration. If the outbreak is fragmented, particularly when close to elimination, it will make more sense to think of it as composed of separate local outbreaks, which can be modelled separately, rather than trying to specify an average $R$ value overall.

If the population is heterogeneous or structured, defining a reproduction number needs care, as the number of new cases, an infective is expected to cause will depend on both their infectiousness and how well connected they are. It has been shown that in the early stages of an epidemic, when the relevant contact structures of a population are not known and interventions are not targeted, assuming a homogeneous contact structure results in conservative estimates of $R_0$ and the required control effort. However, designing targeted intervention strategies requires reliable information on infectious contact structures¹⁰. There are several basic ways to use structured population models to capture departures from the simplest epidemic models. The four most common are (i) household models, (ii) multi-type models, (iii) network models, and (iv) spatial models.

In a household model, every person in the population is assumed to be part of a single household, which is typically small, and may even be of size one. Those in the same household have a higher probability of infecting each other than is the case for two people chosen randomly from the population. In this model, reproduction numbers can still be defined^11,12. The most commonly used is the household reproduction number $R_{*}$, which is the expected number of members of other households that are infected by people from a primary infected household. It is still possible to consider the average number of susceptibles infected by a single infectious person. However, for this to be useful, the average has to be computed in a sophisticated way, because the number of people a person can infect will depend on how many members of the same household are still susceptible when s/he becomes infectious¹³.

A second way of modelling heterogeneity in the population is to assume that the population can be subdivided into groups. The groups may be defined through age bands, social activity levels, health status, type of job, place of residence and so on. Characteristics such as susceptibility, infectivity and frequency of contact may depend on an individual’s group, but all those in a single group have the same characteristics. It is often assumed that all these groups are large. If there are regular inter-group contacts then the largest eigenvalue of the so-called next-generation matrix^5,8 has many similar properties to those of $R_0$ for an epidemic spreading in a homogeneously mixing population, although the final size equation is generally not satisfied.

A third way of introducing heterogeneity is to represent the population by a network, where transmission is only possible between people sharing a link in the network. For many network models, it is still possible to define a reproduction number¹⁴. It is important to note that the person initially infected in a population is often atypical and should be ignored in computing or estimating the reproduction number. A useful extension is a mixture of a network model and a homogeneous mixing model, in which both regular and casual contacts are captured. In this extension, a reproduction number with the desired threshold properties can be defined¹⁵.

Sometimes most transmission is restricted to people living close to each other, and spatial models are useful when physical location should be incorporated. For these, it is often difficult to define a reproduction number because there is no phase in which the number of infected is growing exponentially^16,17. If standard estimation methods are used where there is a considerable spatial component then the estimates will be close to one, even when the spread is highly supercritical and transmission needs to be much reduced to control the epidemic.

In this section, we discuss two frequently used simple statistical methods to estimate $R$ and common issues associated with them. The Cori and Wallinga–Teunis methods estimate subtly different versions of $R_t$; the Cori method generates estimates of the instantaneous reproduction number and the Wallinga–Teunis method generates estimates of the case reproduction number^33,37. The key difference is that the instantaneous reproduction number gives an average $R_t$ for a homogeneous population at a single point in time, whereas the case reproduction number can accommodate individual heterogeneity, but blurs over several dates of transmission. Furthermore, the case reproduction number is a leading estimator of the instantaneous reproduction number, i.e. it depends on data from after the time for which the reproduction number is to be estimated, and must be adjusted accurately to infer the impact of time-specific interventions³⁸.

The instantaneous reproduction number represents the expected number of infections generated at time $t$ by currently infectious individuals³³. For real-time analysis, one of the benefits of estimating the instantaneous reproduction number is that it does not require information about future changes in transmissibility, and it reflects the effectiveness of control measures in place at time $t$. But as an aggregate measure of transmission by all individuals infected in the past (who may now be shedding virus), it does not easily consider heterogeneity in transmission. In contrast, the case reproduction number represents the expected number of infections generated by an individual who is first infected at time $t$ and has yet to progress through the full course of viral shedding. This leads to ‘right censoring’ when the case reproduction number is estimated in real-time; if all infections generated by individuals who were infected at time $t$ have not yet been observed, then the data must be adjusted^39–41 or the case reproduction number will be underestimated.

The Cori method and the Wallinga–Teunis method involve inferring the values of $R_t$ that are most consistent with observed incidence data (for a review, see Gostic et al.³⁸). In the Cori method, typically this inference is carried out by assuming that $R_t$ is constant over fixed time windows. Smoothing windows are used to avoid spurious fluctuations in estimates of $R_t$. These can occur if imperfect observation and reporting effects, rather than actual bursts in transmission, are the main source of noise in the data. Cross-validation and proper scoring rules can be used to avoid under- or oversmoothing $R_t$ estimates⁴².

An important concept, basic to both methods, is the intrinsic generation time also referred to as the infectiousness profile. The intrinsic generation interval is a theoretical quantity derived from the renewal equation of Lotka and Euler^30,43. It describes the time distribution of potentially infectious contacts made by an index case and is independent of population susceptibility⁴⁴. In practice, the intrinsic generation interval is not observable, and it must be estimated carefully from observed serial intervals within contact tracing or household data^44–47. The serial interval is generally defined as the duration of time between the onset of symptoms in an index case and in a secondary case⁴⁸. In the early stages of an outbreak, accurate estimation should adjust for right truncation of observations, for changes over time in population susceptibility, and for interventions such as case isolation, which may shorten the generation interval by limiting transmission events late in the course of infectiousness^44,45,49.

Both the Cori and Wallinga–Teunis methods are conceptually based on separating the infectiousness of an infective into two components, total amount and timing. The timing is expressed by the generation time distribution while the total amount is expressed by $R_t$. The variation of (average) infectivity over time is ascribed, at least in practical implementations of the methods, to changes in $R_t$, while the intrinsic generation time is assumed to remain fixed. This is a simplification that may lead to inaccurate estimation of $R_t$, since, in reality, the observed generation time distribution varies over time, both because of the epidemic dynamics^48,50,51, because of the epidemic affecting different subgroups of the population, with possibly different generation time distributions over time^52,53, and, more importantly, because of interventions that affect the length or efficacy of the infectious period⁴⁹. An additional complication is that the ‘intrinsic’ generation interval of the Cori and Wallinga–Teunis estimators includes potentially infectious contacts with both susceptible and immune individuals, whereas only contacts with susceptible individuals cause new infections, and are observed in contact tracing^44,45. Even when using an accurately estimated fixed generation time distribution, both $R_t$ estimators are numerically sensitive to the specified mean and variance of the intrinsic generation interval⁵⁴.

Fundamentally, $R_t$ is a measure of transmission. Ideally, it would be estimated from data on the total number of incident infections (i.e. transmission events) occurring each day. But in practice, only a small fraction of infections are observed, and notifications do not occur until days or weeks after the moment of infection. Temporally accurate $R_t$ estimation requires adjusting for lags to observation, which can be estimated as the sum of the incubation period and delays from symptom onset to case observation^54,55. Delays not only shift observations into the future, but they also blur infections incident on a particular date across many dates of observation. This blurring can be particularly problematic when working with long and variable delays (e.g. from infection to death), and when $R_t$ is changing. Deconvolution^56–59, or $R_t$ estimation models that include forward delays⁶⁰ can be used to adjust lagged observations. Simpler approaches may be justifiable under some circumstances. If observation delays are relatively short and not highly variable, and if $R_t$ is not rapidly changing, simply shifting unadjusted $R_t$ estimates back in time by the mean delay can provide a reasonable approximation to the true value (see Challen et al.,⁵⁴ in this volume, for an in-depth discussion). The advantages and disadvantages of each approach are reviewed in Gostic et al.³⁸. Changes over time in case ascertainment can also bias $R_t$ estimates, so ideally data should be drawn from structured surveillance (see, e.g., the REACT study⁶¹) or adjusted for known changes in testing or reporting effort^61,62.

In practice, $R_t$ can be estimated from a time series of new symptom onset reports, cases, hospitalisations or deaths. Choosing an appropriate data stream involves weighing representativeness, timeliness of reporting, consistency of ascertainment, and length of lag. For example, reported deaths may be reasonably unaffected by changes over time in ascertainment, but adjusting for long lags to observation can be challenging, and deaths may not be representative of overall transmission (e.g. if the epidemic shifts towards younger age groups)^63,64. Extensions of existing statistical models for $R_t$ estimation could potentially integrate multiple kinds of data, by assuming that, for example, cases, hospitalisations and deaths, arise from a shared, latent infection process, with different delays³⁸. A mechanistic model can also pull multiple data streams together by modelling the different processes underlying each data stream. Problems can arise if different data streams disagree on the progress of the pandemic. However, if the disagreement is caused by a shift in delays from events to reporting in different data streams, a mechanistic model can highlight these changes. Sometimes different data streams can be used for model validation.

All methods used to estimate $R_t$ must decide on the length of the time window over which it is to be estimated. All data used to estimate $R_t$ are noisy. The shorter the time window used for estimation, the higher will be the noise-to-signal ratio and, therefore, the uncertainty in the estimate of $R_t$. In contrast, longer time windows will produce estimates with lower uncertainty, but sudden changes in transmission may not be detected if the time window is too long.