EconBase
← Back to paper

Matching Theory and Evidence on Covid-19 using a Stochastic Network SIR Model

Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.

119,797 characters · 17 sections · 0 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Matching Theory and Evidence on Covid-19 using a Stochastic Network SIR Model

abstractThis paper develops an individual-based stochastic network SIR model for the empirical analysis of the Covid-19 pandemic. It derives moment conditions for the number of infected and active cases for single as well as multigroup epidemic models. These moment conditions are used to investigate the identification and estimation of the transmission rates. The paper then proposes a method that jointly estimates the transmission rate and the magnitude of under-reporting of infected cases. Empirical evidence on six European countries matches the simulated outcomes once the under-reporting of infected cases is addressed. It is estimated that the number of actual cases could be between 4 to 10 times higher than the reported numbers in October 2020 and declined to 2 to 3 times in April 2021. The calibrated models are used in the counterfactual analyses of the impact of social distancing and vaccination on the epidemic evolution, and the timing of early interventions in the UK and Germany. Keywords: Covid-19, multigroup SIR model, basic and effective reproduction numbers, transmission rates, vaccination, calibration and counterfactual analysis. JEL Classifications: C13, C15, C31, D85, I18, J18

\pagenumbering{gobble}

\baselineskip0.245in

\pagenumbering{arabic} \setcounter{page}{1} \doublespacing

Introduction

Since the outbreak of Covid-19, many researchers in epidemiology, behavioral sciences, and economics have applied various forms of compartmental models to study the disease transmission and potential outcomes under different intervention policies. The compartmental models are a major group of epidemiological models that categorize a population into several types or groups, such as susceptible (S), infected (I), and removed (recovered or deceased, R). Compartmental models owe their origin to the well-known SIR model pioneered by Kermack and McKendrick (1927) and have been developed in a number of important directions, allowing for multi-category (multi-location), under a variety of contact networks and transmission channels.\footnote{See Section (ref) of the online supplement for a review of the related literature. Comprehensive reviews can be found in Hethcote (2000) and Thieme (2013).}

In this paper, we develop a new stochastic network SIR model in which individual-specific infection and recovery processes are modelled, allowing for group heterogeneity and latent individual characteristics that distinguish individuals in terms of their degrees of resilience to becoming infected. The model is then used to derive individual-specific conditional probabilities of infection and recovery. In this respect, our modelling approach is to be distinguished from the individual-based models in epidemiology that specify the transition probabilities of individuals from one state to another.\footnote{See, for example, Rocha and Masuda (2016), Willem et al. (2017), and Nepomuceno, Resende, and Lacerda (2018).} In modelling the infection process, we consider an individual's contact pattern with others in the network, plus an individual-specific latent factor assumed to be exponentially distributed. The time from infection to recovery (or death) is assumed to be geometrically distributed. The individual processes are shown to aggregate up to the familiar multigroup SIR model. We allow for group heterogeneity and, in line with the literature, assume contact probabilities are homogeneous within groups but could differ across groups.

We then derive the probabilities of individuals within a given group being in a particular state at a given time, conditional on contact patterns, exposure intensities, and unobserved characteristics. These conditional probabilities are aggregated up to form a set of moment conditions that can be taken to the data on the number of infected and active cases both at the aggregate and group (or regional) levels. We make use of the moment conditions to investigate the identification of the underlying structural parameters. Most importantly, we show that whilst one cannot distinguish between average contact numbers and the degree of exposure to the virus upon contact, it is nevertheless possible to identify the basic and effective reproduction numbers from relatively short time series observations on infections and recoveries. Using Monte Carlo simulations, the small sample properties of the proposed estimator are shown to be satisfactory, with a high degree of precision even when using two and three weeks of rolling observations.

However, in practice, estimation of the transmission rate must take account of the well-known measurement problem where the number of infected cases is often grossly under-reported. This problem is further complicated since the degree of under-reporting varies over time and tends to fall as society becomes more familiar with the disease and testing becomes more widespread. To deal with this mismeasurement problem, we propose a new method that jointly estimates the transmission rates and the multiplication factor that measures the degree of under-reporting. Equipped with daily estimates of the transmission rates, we are then able to calibrate our epidemic model and investigate its properties under different network topology, group numbers, and different social distancing and vaccination strategies.

We apply the proposed joint estimation approach to examine how well the outcome of the proposed epidemic model matches the Covid-19 evidence in the case of six European countries (Austria, France, Germany, Italy, Spain, and the UK) from March 2020 to April 2021. We provide rolling estimates of the transmission rates and related effective reproduction numbers, as well as estimates of multiple factors. We then use the estimated transmission rates to calibrate the model parameters across the six countries. The stochastically simulated outcomes are shown to be reasonably close to the reported cases once the under-reporting issue has been addressed. We estimate that the degree of under-reporting declined from a multiple of $4$--$10$ to $2$--$3$ times during the study period across the countries considered.

Finally, we illustrate the use of our model for two different counterfactual exercises. First, we consider the effects of vaccination on the evolution of the epidemic using a multigroup setup, where we also evaluate the implications of age-based vaccine prioritization on the outcomes. Our model allows each individual to have their own degree of immunity, with vaccination increasing this individual-specific immunity by a factor of $20$ in the case of Moderna or Pfizer-BioNTech that are shown to be $94$--$95\%$ effective (Oliver et al., 2020, 2021). The multigroup model is particularly helpful in examining the implications of different vaccine prioritization strategies. Second, we investigate the potential outcomes if the first lockdown in Germany had been delayed for one or two weeks; and if the lockdown in the UK had started one or two weeks earlier. Such counterfactual analyses can be achieved by shifting the estimated transmission rates forward or backward. We show that early intervention is critical in managing the infection and controlling the total number of infected and active cases.

The problem of how to balance the public health risks from the spread of the epidemic with the economic costs associated with lockdowns and other mandatory social-distancing regulations will not be addressed in this paper. However, the proposed network SIR model with its individual-based architecture is eminently suited to this purpose. The proposed model can be combined with behavioral assumptions about how individuals trade off infection risk and economic well-being, thus generalizing the aggregate framework proposed in Chudik, Pesaran, and Rebucci (2021) to individual-based SIR models.

The rest of the paper is organized as follows. Section (ref) introduces the basic concepts and the classical multigroup SIR\ model. Section (ref) lays out our individual based stochastic model on a network. Section (ref) explains the calibration of our model to basic and effective reproduction numbers. Section (ref) documents the properties of the model. Sections (ref) discusses the estimation of the transmission rate. Section (ref) presents the calibration of the model to Covid-19 evidence. Section (ref) concentrates on the counterfactual analyses, and Section (ref) concludes.

To save space, a detailed review of the related literature is given in an online supplement, where we also provide supplementary theoretical derivations, additional estimation results, counterfactual outcomes, and data sources.

Basic concepts and the multigroup SIR model

We consider a population of $n$ individuals susceptible to the spread of a disease with some initially infected individuals. We suppose that the susceptible population can be categorized into $L$ groups of size $n_{\ell}$, $\ell=1,2,\ldots,L$, with $L$ fixed such that $n=\sum_{\ell=1}^{L}n_{\ell}$. It is further assumed that the group shares, $w_{\ell}=n_{\ell}/n>0,$ for all $n$ and as $n\rightarrow\infty$. The grouping could be based on demographic factors (age and/or gender), or other observed characteristics such as contact locations and/or schedules, mode of transmission, genetic susceptibility, group-specific vaccination coverage, as well as socioeconomic factors (see, e.g., Hethcote, 2000). Individual $i$ in group $\ell$ will be referred to as individual $\left( i,\ell\right) $, with $i=1,2,\ldots,n_{\ell}$ and $\ell=1,2,\ldots,L$. It is assumed that $n_{\ell}$ is relatively large but remains fixed over the course of the epidemic measured in days.

Suppose that individual $\left( i,\ell\right) $ becomes infected on day $t=t_{i\ell}^{\ast}$, and let $x_{i\ell,t}$ take the value of unity for all $t\geq t_{i\ell}^{\ast}$, and zero otherwise. In this way, we follow the convention that once an individual becomes infected, he/she is considered as infected thereafter, irrespective of whether that individual recovers or dies. Specifically, we set

equation[equation omitted — 152 chars of source]

The event of recovery or death of individual $\left( i,\ell\right) $ at time $t$ will be represented by $y_{i\ell,t}$, which will be equal to zero unless the individual is "removed" (recovered or dead). An individual $(i,\ell)$ is considered to be "active" if he/she is infected and not yet recovered. We denote the active indicator by $z_{i\ell,t},$ which is formally defined by

equation[equation omitted — 85 chars of source]

$z_{i\ell,t}$ takes the value of $0$ if individual $(i,\ell)$ has not been infected, or has been infected but recovered/dead. It takes the value of $1$ if he/she is infected and not yet recovered. Any individual $(i,\ell)$ who has not been infected is viewed "susceptible" and indicated by $s_{i\ell,t}=1$, where

equation[equation omitted — 74 chars of source]

It then readily follows that the total (cumulative) number of those "infected" in group $\ell$ at the end of day $t$ is given by

equation[equation omitted — 97 chars of source]

where the summation is over all individuals in group $\ell$. The total number of "recovered" in group $\ell$ in day $t$ is given by

equation[equation omitted — 97 chars of source]

The total number of "active" cases (individuals who are infected and not yet removed) in group $\ell$ in day $t$ is

equation[equation omitted — 171 chars of source]

The number of "susceptible" individuals in group $\ell$ in day $t$ is

equation[equation omitted — 109 chars of source]

Our model does not distinguish between recovery and death. Once an individual is removed (recovered or dead), following the SIR literature, we assume that he/she cannot be infected again. Under this assumption, recovery and death have the same effects on the evolution of the epidemic, and accordingly in what follows we shall not differentiate between recovery and death and refer to their total as "removed".

The classic multigroup SIR\ model in discretized form can be written as\footnote{See, e.g., Guo, Li, and Shuai (2006), Zhang et al. (2020), and the references therein.}

align[align omitted — 368 chars of source]

for $\ell=1,2,\ldots,L$ and $t=1,2,\ldots,T$, where $S_{\ell t}$, $I_{\ell t}$ and $R_{\ell t}$ are defined as above, $\gamma_{\ell}$ is the recovery rate which is assumed to be time-invariant and the same for all people in group $\ell$, and $\beta_{\ell\ell^{\prime}}$ is the transmission coefficient between $S_{\ell t}$ and $I_{\ell^{^{\prime}}t}$. Note that individuals in group $\ell^{^{\prime}}$ may transmit the disease to individuals in group $\ell$, with the new infections in group $\ell$ given by $S_{\ell t}\sum _{\ell^{^{\prime}}=1}^{L}\beta_{\ell\ell^{\prime}}I_{\ell^{^{\prime}}t}$.

An individual based stochastic epidemic model on a network

We now depart from the literature by explicitly modelling the individual indicators, $x_{i\ell,t}$ and $y_{i\ell,t}$, (and hence $z_{i\ell,t}$) and then simulate and aggregate up to match the theoretical predictions with realized aggregated outcomes. In this section, we first describe the infection and recovery processes at the individual level, we then show how they lead to the moment conditions at group levels, and finally derive the relation between aggregated outcomes from our model and the multigroup SIR\ deterministic model.

Modelling the infection and recovery processes

As an attempt to better integrate individual decisions to mitigate their health risk within the epidemic model, we propose to directly model $x_{i\ell,t}$ for each individual $(i,\ell)$, as compared to modelling the group aggregates $S_{\ell t}$, $R_{\ell t}$, and $I_{\ell t}$. We follow the micro-econometric literature and model the infection process using the latent variable, $x_{i\ell,t+1}^{\ast},$ which determines whether individual $(i,\ell)$ becomes infected. Specifically, we begin with the following Markov switching process for individual $(i,\ell)$:

equation[equation omitted — 127 chars of source]

where $I\left( \mathcal{A}\right) $ is the indicator function that takes the value of unity if $\mathcal{A}$ holds and zero otherwise. We suppose that $x_{i\ell,t+1}^{\ast}$ is composed of two different components. The first component relates to the contact pattern of individual $(i,\ell)$ with all other individuals in the active set, denoted by $z_{j\ell^{\prime},t},$ both within (when $\ell^{\prime}=\ell$) and outside of his/her group (when $\ell^{\prime}\neq\ell$). The second component is an unobserved individual-specific infection threshold variable, $\xi_{i\ell,t+1}>0$. Formally, we set

equation[equation omitted — 207 chars of source]

where the first component depends on the pattern of contacts, $d_{i\ell ,j\ell^{^{\prime}}}\left( t\right) $, whether the contacted individuals are infectious, $z_{j\ell^{^{\prime}},t}$, and the exposure intensity parameter, denoted by $\tau_{i\ell}$. $\mathbf{D}(t)=\left[ d_{i\ell,j\ell^{^{\prime}} }\left( t\right) \right] $ is the contact network matrix, such that $d_{i\ell,j\ell^{^{\prime}}}\left( t\right) =1$ if individual $\left( i,\ell\right) $ is in contact with individual $\left( j,\ell^{^{\prime} }\right) $ at time $t.$ $z_{j\ell^{\prime},t}=\left( 1-y_{j\ell^{\prime} ,t}\right) x_{j\ell^{\prime},t}$ is an infectious indicator, already defined by ((ref)), and takes the value of unity if individual $\left( j,\ell^{\prime}\right) $ is infected and not yet recovered, zero otherwise. The exposure intensity parameter, $\tau_{i\ell}$, is group-specific and depends on the average duration of contacts, whether the contacting individuals are wearing facemasks, and if they follow other recommended precautions.

The multigroup structure of the first component of ((ref)) covers a wide range of observable characteristics, and can be extended to allow for differences in age, location, and medical pre-conditions. There are also many unobservable characteristics that lead to different probabilities of infection, even for individuals with the same contact patterns and exposure intensities. To allow for such latent factors, we have introduced the individual-specific positive random variable ($\xi_{i\ell,t+1}>0$) which represents the individual's degree of resilience to becoming infected and varies across $(i,\ell)$ and $t.$ Ceteris paribus, an individual with a low value of $\xi_{i\ell,t+1}$ is more likely to become infected. $\xi_{i\ell,t+1}$ is assumed to be independently distributed over $i$, $\ell$ and $t$, and follows an exponential distribution with the cumulative distribution function given by

equation[equation omitted — 123 chars of source]

where $\mu_{i\ell}=E\left( \xi_{i\ell,t+1}\right) $.

To complete the specification of the infection process we also need to model $y_{i\ell,t}$, namely the recovery indicator. We assume that recovery depends on the number of days since infection. Specifically, the recovery process for individual $i$ is defined by

equation[equation omitted — 124 chars of source]

where $z_{i\ell,t}=\left( 1-y_{i\ell,t}\right) x_{i\ell,t}$, $\zeta _{i\ell,t+1}\left( t_{i\ell}^{\ast}\right) =1$ if individual $\left( i,\ell\right) $ recovers at time $t+1$, having been infected exactly at time $t_{i\ell}^{\ast}$ and not before, and $\zeta_{i\ell,t+1}\left( t_{i\ell }^{\ast}\right) =0,$ otherwise. The analysis of recovery simplifies considerably if we assume time to removal, denoted by $T_{i\ell,t}^{\ast }=t-t_{i\ell}^{\ast},$ follows the geometric distribution (for $t-t_{i\ell }^{\ast}=1,2,\ldots$)

equation[equation omitted — 233 chars of source]

Then the probability of recovery at time $t+1$\ having remained infected for $t-t_{i\ell}^{\ast}-1$ days (also known as the "hazard function") is given by

equation[equation omitted — 329 chars of source]

which is the same across all individuals within a given group and, most importantly, does not depend on the number of days since infection.\footnote{A more general specification that allows the recovery probability to depend on the number of days being infected is considered in Section (ref) of the online supplement.} Therefore, using ((ref)) in ((ref)), the recovery micro-moment condition simplifies to

equation[equation omitted — 150 chars of source]

which implies that

equation[equation omitted — 169 chars of source]

where $\mathbf{y}_{\ell,t}=(y_{1\ell,t},y_{2\ell,t},....,y_{n_{\ell}\ell ,t})^{\prime}$ and $\mathbf{z}_{\ell,t}=(z_{1\ell,t},z_{2\ell,t} ,....,z_{n_{\ell}\ell,t})^{\prime}.$

We assume that $d_{i\ell,j\ell^{^{\prime}}}\left( t\right) $, the elements of the $n\times n$ network matrix $\mathbf{D}\left( t\right) =\left[ d_{i\ell,j\ell^{^{\prime}}}\left( t\right) \right] ,$ are independent draws with $E\left[ d_{i\ell,j\ell^{^{\prime}}}\left( t\right) \right] =p_{\ell\ell^{^{\prime}}}$, namely, the probability of contacts is homogeneous within groups but differs across groups. Let $\boldsymbol{d}_{i}^{^{\prime} }\left( t\right) $ be the $i^{th}$ row of $\mathbf{D}\left( t\right) $. Also let $\mathbf{z}_{t}$ be a column vector consisting of $z_{j\ell,t}$, for $j=1,2,\ldots,n_{\ell}$ and $\ell=1,2,\ldots,L$. Then using ((ref)) we have

align*[align* omitted — 753 chars of source]

where $\theta_{i\ell}=\tau_{i\ell}/\mu_{i\ell}$, which represents the net exposure effect. Since in general individual contact patterns are not observed, we also need to derive $E\left( x_{i\ell,t+1}|x_{i\ell ,t},\mathbf{z}_{t},\theta_{i\ell}\right) $. To this end, we first note that

align*[align* omitted — 510 chars of source]

and since by assumption $d_{i\ell,j\ell^{^{\prime}}}\left( t\right) $\ are independently distributed, we then have

align*[align* omitted — 562 chars of source]

However, recall that $z_{j\ell^{\prime},t}=1$ if individual $\left( j,\ell^{\prime}\right) $ is currently infected (namely if at time $t$ he/she is a member of the active set, $\mathcal{I}_{t}$), otherwise $z_{j\ell ^{\prime},t}=0$. In the latter case $1-p_{\ell\ell^{^{\prime}}}+\exp\left( -\theta_{i\ell}z_{j\ell^{^{\prime}},t}\right) p_{\ell\ell^{^{\prime}}}=1$, and hence

align[align omitted — 510 chars of source]

where $I_{\ell^{\prime}t}=\sum_{j=1}^{n_{\ell^{\prime}}}z_{j\ell^{\prime} ,t}=C_{\ell^{\prime}t}-R_{\ell^{\prime}t}$. See also ((ref)).

Moment conditions at group and aggregate levels

We will first derive the moment conditions at the group level. Denote the per capita infected and recovered in group $\ell$ by $c_{\ell t}=C_{\ell t}/n_{\ell}$ and $r_{\ell t}=R_{\ell t}/n_{\ell}$, respectively, and note that $i_{\ell t}=c_{\ell t}-r_{\ell t}=I_{\ell t}/n_{\ell}$. Let $p_{\ell \ell^{^{\prime}}}=k_{\ell\ell^{\prime}}/n_{\ell^{\prime}},$ where $k_{\ell \ell^{\prime}}$ is the mean daily contacts from group $\ell^{\prime}$\ for an individual in group $\ell$.\footnote{The parameters of the model, including $\tau_{\ell}$, $k_{\ell\ell^{\prime}}$, $p_{\ell\ell^{^{\prime}}}$, $\beta_{\ell\ell^{\prime}},$ and $\mu_{i\ell}$ can be time-varying due to behavioral changes, vaccination, or other reasons. However, we suppress the time subscript $t$ to simplify the exposition.}\ To preserve the symmetry of contact probabilities, the mean contact numbers must satisfy the so-called reciprocity condition, $n_{\ell}k_{\ell\ell^{\prime}}=n_{\ell^{\prime}} k_{\ell^{\prime}\ell}$ (see, e.g., Willem et al., 2020). That is, the total number of contacts that people in group $\ell$ have with people in group $\ell^{\prime}$ must be the same as the number of contacts that people in group $\ell^{\prime}$ have with those in group $\ell$. In practice, $n_{\ell}$ is often quite large, with $k_{\ell\ell^{\prime}}$ relatively small (often less than $30$). Therefore, it is reasonable to assume that $k_{\ell \ell^{\prime}}$ is fixed in $n_{\ell^{\prime}}$ and hence $p_{\ell \ell^{^{\prime}}}=O\left( n_{\ell^{\prime}}^{-1}\right) $. Then we have

align[align omitted — 603 chars of source]

Suppose that $\theta_{i\ell}$ is small enough such that $1-e^{-\theta_{i\ell} }\approx\theta_{i\ell}$ is sufficiently accurate. Also recall that $w_{\ell }=n_{\ell}/n>0,$ for all $\ell$, then $n_{\ell}$ rises at the same rate as $n$. It follows that

equation[equation omitted — 265 chars of source]

and hence

equation[equation omitted — 319 chars of source]

Let $\chi_{\ell t}=\sum_{\ell^{\prime}=1}^{L}i_{\ell^{\prime}t}k_{\ell \ell^{\prime}}=\mathbf{i}_{t}^{\prime}\mathbf{k}_{\ell\circ}$, where $\mathbf{i}_{t}=(i_{1t},i_{2t},\ldots,i_{Lt})^{\prime}$ and $\mathbf{k} _{\ell\circ}=\left( k_{\ell1},k_{\ell2},\ldots,k_{\ell L}\right) ^{\prime}.$ Let $H\left( \chi_{\ell t},\boldsymbol{\varphi}_{\ell}\right) =E\left( e^{\chi_{\ell t}\theta_{i\ell}}|\chi_{\ell t}\right) $, where the expectations are taken with respect to the distribution of $\theta_{i\ell}$ for a given $\ell$, and $\boldsymbol{\varphi}_{\ell}$ refers to the parameters of the distribution of $\theta_{i\ell}$ over $i$ in group $\ell$. Note that $H\left( \chi_{\ell t},\boldsymbol{\varphi}_{\ell}\right) $ is the moment generating function of $\theta_{i\ell}$, assumed to be the same across all individuals. Then using the above result in the micro infection moment conditions, ((ref)), gives

equation[equation omitted — 205 chars of source]

Let $\mathbf{x}_{\ell t}=(x_{1\ell,t},x_{2\ell,t},...,x_{n_{\ell}\ell ,t})^{\prime}$ and note that since $x_{i\ell,t}$ is a subset of $\mathbf{x} _{\ell t}$, then \[ E\left[ E\left( x_{i\ell,t+1}|\mathbf{x}_{\ell t},\mathbf{i}_{t}\right) \left\vert x_{i\ell,t}\right. \right] =E\left( x_{i\ell,t+1}|x_{i\ell ,t},\mathbf{i}_{t}\right) , \] and the moment condition ((ref)) also implies that

equation[equation omitted — 213 chars of source]

for $i=1,2,...,n$. Averaging the above conditions over $i$ for a given group $\ell$, and recalling that $c_{\ell,t+1}=\sum_{i=1}^{n}x_{i\ell,t+1}/n_{\ell} $, we obtain\

equation[equation omitted — 217 chars of source]

We will return to the heterogeneous $\theta_{i\ell}$ in the counterfactual analysis of vaccination to be discussed in Section (ref), where $\theta_{i\ell}$ is associated with the vaccine effectiveness for individual $(i,\ell)$. In order to derive analytical results and achieve identification in estimation, in what follows, we assume $\theta_{i\ell }=\theta_{\ell}=\tau_{\ell}/\mu_{\ell}$ for all $i$ in group $\ell$. Also note that $\tau_{\ell}$ and $\mu_{\ell}$ are not separately identified. Without loss of generality, we normalize $\mu_{\ell}=1$. Under these conditions, $\theta_{i\ell}=\tau_{\ell}$, the group-level infection moment condition can be written as \[ E\left( \left. \sum_{i=1}^{n_{\ell}}x_{i\ell,t+1}\right\vert C_{\ell t},\mathbf{I}_{t}\right) =E\left( C_{\ell,t+1}|\text{ }\mathbf{z} _{t}\right) =n_{\ell}-\left( n_{\ell}-C_{\ell t}\right) \prod _{\ell^{^{\prime}}=1}^{L}\left( 1-p_{\ell\ell^{^{\prime}}}+p_{\ell \ell^{^{\prime}}}e^{-\tau_{\ell}}\right) ^{I_{\ell^{\prime}t}}, \] which can be written equivalently as (recall that $I_{\ell t}=C_{\ell t}-R_{\ell t}$)

equation[equation omitted — 348 chars of source]

Also aggregating the micro recovery moment conditions, ((ref)), we have

equation[equation omitted — 183 chars of source]

To sum up, in per capita terms, we obtain the following $2L\,$dimensional system of moment conditions (for $\ell=1,2,\ldots,L$)

align[align omitted — 436 chars of source]

Given time series data on $\mathbf{c}_{t}=(c_{1t},c_{2t},\ldots,c_{Lt} )^{\prime}$ and $\mathbf{r}_{t}=(r_{1t},r_{2t},\ldots,r_{Lt})^{\prime}$, the above moment conditions can be used to estimate the structural parameters, $\gamma_{\ell}$, $\tau_{\ell}$ and $p_{\ell\ell^{\prime}}=p_{\ell^{\prime} \ell}$.

In relating the theory to the data, one may need to further aggregate across groups to the population level if group-level data are unavailable or unreliable. It is interesting to note that the multigroup model does not lead to a model for the aggregates, $C_{t}=\sum_{\ell=1}^{L}C_{\ell t}$ and $I_{t}=\sum_{\ell=1}^{L}I_{\ell t}$, without additional restrictions. To see this, using ((ref)) in ((ref)) and under the assumption that $\theta_{i\ell}=\tau_{\ell}$, we obtain

equation[equation omitted — 241 chars of source]

where $\beta_{\ell\ell^{\prime}}=\left( 1-e^{-\tau_{\ell}}\right) k_{\ell\ell^{\prime}}\approx\tau_{\ell}k_{\ell\ell^{\prime}}$. The approximation holds since $\tau_{\ell}$ is small. Notice that $\sum_{\ell =1}^{L}w_{\ell}c_{\ell}=C_{t}/n=c_{t}$ and $\sum_{\ell=1}^{L}w_{\ell}i_{\ell }=I_{t}/n=i_{t}$. If we multiply both sides of ((ref)) by $w_{\ell }$ and sum across $\ell=1,2,\ldots,L$, we obtain

equation[equation omitted — 255 chars of source]

It is now clear that the group moment condition for infected cases, ((ref)), does not aggregate up to the moment conditions in terms of $c_{t}$ and $i_{t}$, unless $\beta_{\ell\ell^{\prime}}/n_{\ell^{\prime}}$ is the same across all $\ell$ and $\ell^{^{\prime}}$. It is also straightforward to see that the group moment condition for recovery, ((ref)), does not aggregate up either unless $\gamma_{\ell}=\gamma$ for all $\ell$.

In the case of a single group, we have $\tau_{\ell}=\tau$, $k_{\ell \ell^{\prime}}=k$, and $\beta_{\ell\ell^{\prime}}=\beta\approx\tau k,$ for all $\ell$ and $\ell^{^{\prime}}$. Then ((ref)) simplifies to

equation[equation omitted — 150 chars of source]

Also, if $\gamma_{\ell}=\gamma$ for all $\ell$, the recovery moment condition, ((ref)), becomes

equation[equation omitted — 114 chars of source]

Given aggregate data on $c_{t}$, $i_{t},$ and $r_{t}$, one can estimate $\beta$ and $\gamma$ using the moment conditions ((ref)) and ((ref)), respectively. Interestingly, it can be shown that the multigroup SIR model given by ((ref))--((ref)) is a linearized-deterministic version of the above moment conditions. The relationship between our model and the classical SIR model is set out in Section (ref) of the online supplement.

Basic and effective reproduction numbers

In this section, we consider the calibration of our model to a given basic reproduction number assuming no intervention, and derive the effective reproduction numbers in terms of mean contact patterns, exposure intensities, and the recovery rate. We also consider the problem of identifying contact patterns from the exposure rates in single and multigroup contexts.

Basic reproduction number

The basic reproduction number, denoted by $\mathcal{R}_{0}$, is defined as "the average number of secondary cases produced by one infected individual during the infected individual's entire infectious period assuming a fully susceptible population" (Del Valle, Hyman, and Chitnis, 2013). By construction, $\mathcal{R}_{0}$ measures the degree to which an infectious disease spreads when left unchecked. The infection spreads if $\mathcal{R} _{0}>1$ and abates if $\mathcal{R}_{0}<1$.

In order to derive $\mathcal{R}_{0}$ for our multigroup model, we suppose that on day $1$ a fraction $w_{\ell}=n_{\ell}/n$ of each group $\ell$ becomes infected, which represents the equivalent of one individual becoming infected as required by the definition of $\mathcal{R}_{0}$. That is, on day $1$, $R_{\ell1}=0$, $I_{\ell1}=C_{\ell1}=w_{\ell}$, for $\ell=1,2,\ldots,L$. Also, in view of our model of the recovery, the probability that an individual infected on day $1$ remains infected on day $s\geq1$ is given by $\Pr\left( T_{i\ell,t}^{\ast}\geq s\right) =\left( 1-\gamma_{\ell}\right) ^{s-1}$, for $s=1,2,\ldots$.\footnote{Using ((ref)), note that $\Pr\left( T_{i\ell,t}^{\ast}\geq s\right) =1-\Pr\left( T_{i\ell,t}^{\ast}<s\right) =1-\sum_{j=1}^{s-1}\Pr\left( T_{i\ell,t}^{\ast}=j\right) =1-\sum_{j=1} ^{s-1}\gamma_{\ell}\left( 1-\gamma_{\ell}\right) ^{j-1}=\left( 1-\gamma_{\ell}\right) ^{s-1}.$} Hence, we have

equation[equation omitted — 163 chars of source]

where $\mathbf{w=(}w_{1},w_{2},\ldots,w_{L})^{\prime}$. Now using ((ref)) for $s=2$ we have

equation[equation omitted — 235 chars of source]

Due to the large number of possibilities that follow after the second day of the epidemic, it is not possible to derive similar analytical expressions for $E\left( C_{_{\ell,s}}|\mathbf{w}\right) $, $s=3,4,\ldots$. But since the weights of these future expected values decay geometrically, and at the start of the epidemic the number of infected is likely to be very small relative to the susceptible population, we think it is reasonable to follow the literature (Farrington and Whitaker, 2003; Elliott and Gourieroux, 2020) and assume that $E\left( C_{_{\ell,s+1}}|\mathbf{w}\right) \approx E\left( C_{\ell 2}|\mathbf{w}\right) $, for $s\geq1$.\footnote{Elliott and Gourieroux (2020) make a similar assumption that the expected number of susceptibles is constant over time in deducing the reproduction numbers (p. 7). The constant recovery intensity is a standard assumption in the SIR\ literature. In contrast, our calibration to the reproduction numbers differs significantly from Elliott and Gourieroux (2020) in that we directly consider individual $i$ in group $\ell$ and his/her contact network, rather than modelling population groups classified by their S, I, or R\ status.} Under this assumption, the following approximate expression for $\mathcal{R}_{0}$ obtains: \[ \mathcal{R}_{0}\approx\sum_{\ell=1}^{L}\sum_{s=1}^{\infty}\left( 1-\gamma_{\ell}\right) ^{s-1}E\left( C_{\ell2}|\mathbf{w}\right) =\sum_{\ell=1}^{L}\gamma_{\ell}^{-1}E\left( C_{\ell2}|\mathbf{w}\right) . \] In the case where the recovery rates are the same across the groups ($\gamma_{\ell}=\gamma$), the above expression simplifies further and we have $\mathcal{R}_{0}\approx\gamma^{-1}\sum_{\ell=1}^{L}E\left( C_{\ell 2}|\mathbf{w}\right) $.\footnote{In the case of Covid-19, it is universally assumed that $\gamma_{\ell}=\gamma=1/14$, which we also adopt in our empirical analysis.} Now using ((ref)) gives

equation[equation omitted — 250 chars of source]

To see how the above result relates to the well-known expression $\mathcal{R}_{0}=\beta/\gamma$, consider the case of a single group with $p=k/(n-1)\thickapprox k/n$ as $n$ is large. Then the expression in ((ref)) reduces to

equation[equation omitted — 150 chars of source]

where the last result follows by $1-e^{-\tau}\approx\tau$. Hence the model can be calibrated to any choice of $\mathcal{R}_{0}$ and $\gamma$ by setting the average number of contacts, $k$, and/or the exposure intensity parameter, $\tau$. It is clear that $\tau$ and $k\,$are not separately identified---only their product is identified. In addition, we would obtain the standard result $\gamma\mathcal{R}_{0}=\beta$ for SIR\ models if we set $\beta=np(1-e^{-\tau })\thickapprox\tau k$.

Returning to the multigroup case, expression ((ref)) continues to apply if the population is homogeneous in the sense that $p_{\ell \ell^{^{\prime}}}=p$, $\tau_{\ell}=\tau$, for all $\ell$ and $\ell^{^{\prime} }$. But in the more realistic case of group heterogeneity, we can use ((ref)) to calibrate $\tau_{\ell}$ and/or $p_{\ell\ell^{^{\prime}}}$ for given choices of $\mathcal{R}_{0}$ and $\gamma$. Since we have assumed that $n_{\ell}$ is large and $L$ is fixed, ((ref)) can be further simplified with a linear approximation derived as follows. Let $A_{n,\ell \ell^{^{\prime}}}=\left( 1-p_{\ell\ell^{^{\prime}}}+p_{\ell\ell^{^{\prime}} }e^{-\tau_{\ell}}\right) ^{w_{\ell^{^{\prime}}}}$, and use a similar argument as in ((ref)) to obtain (recall that $p_{\ell\ell^{\prime}} =O(n^{-1})$)

align*[align* omitted — 378 chars of source]

Then $A_{n,\ell\ell^{^{\prime}}}=\exp\left( -\tau_{\ell}w_{\ell^{^{\prime}} }p_{\ell\ell^{^{\prime}}}\right) +O\left( n^{-2}\right) $. Using this in ((ref)) gives

align[align omitted — 538 chars of source]

As before setting $p_{\ell\ell^{^{\prime}}}=k_{\ell\ell^{\prime}} /n_{\ell^{\prime}}$, for $n$ sufficiently large, the above expression can be written equivalently as

equation[equation omitted — 202 chars of source]

where $\beta$ and $\beta_{\ell}$ are the aggregate and group-specific transmission rates, respectively.

Similar to the case of a single group, equation ((ref)) implies that $\tau_{\ell}$ and $p_{\ell\ell^{^{\prime}}}$ are not separately identified; only their products are identified (or equivalently, $\beta_{\ell\ell^{\prime }}\approx\tau_{\ell}k_{\ell\ell^{\prime}}$ are identified). To see this more formally, consider the simple case of two groups $(L=2)$. Then for sufficiently large $n$, using ((ref)) with $L=2$ we have

align[align omitted — 345 chars of source]

where the last line follows by the symmetry of contact probabilities: $p_{\ell\ell^{^{\prime}}}=p_{\ell^{^{\prime}}\ell}$. It is clear from ((ref)) that only $\tau_{1}p_{11}$, $\left( \tau_{1}+\tau_{2}\right) p_{12}$, and $\tau_{2}p_{22}$ can be identified given $w_{1}$, $w_{2}$, $n$ and $\gamma$. More generally, for finite $L\geq2$, $\tau_{\ell}p_{\ell \ell^{^{\prime}}}$ are identified for any $\ell$ and $\ell^{^{\prime} }=1,2,\ldots,L$.

Effective reproduction numbers and mitigation policies

In reality, the average number of secondary cases will vary over time as a result of the decline in the number of susceptible individuals (due to immunity or death) and/or changes in behavior (due to mitigation strategies such as social distancing, quarantine measures, travel restrictions and wearing of facemasks). The effective reproduction number, which we denote by $\mathcal{R}_{et}$,\footnote{We use this notation in order to clearly distinguish the effective reproduction number from the number of removed cases, $R_{t}$.} is the expected number of secondary cases produced by one infected individual in a population that includes both susceptible and non-susceptible individuals at time $t.$ In a multigroup setting, we represent "one infected individual" by the vector of population proportions, $\mathbf{w=(}w_{1},w_{2},\ldots,w_{L})^{\prime}$. The evolution of $\mathcal{R}_{et}$ is determined by the remaining number of susceptibles by groups, $S_{\ell t}=n_{\ell}-C_{\ell t},$ for $\ell=1,2,\ldots,L$. Formally, $\mathcal{R}_{et}$ is defined by

equation[equation omitted — 147 chars of source]

In the absence of any interventions, using ((ref)) we have

equation[equation omitted — 262 chars of source]

Recalling that $\left( 1-p_{\ell\ell^{^{\prime}}}+p_{\ell\ell^{^{\prime}} }e^{-\tau_{\ell}}\right) ^{w_{\ell^{^{\prime}}}}=\exp\left( -\tau_{\ell }w_{\ell^{^{\prime}}}p_{\ell\ell^{^{\prime}}}\right) +O\left( n^{-2}\right) $, then for $n$ sufficiently large we have the following approximate expression for $\mathcal{R}_{et}$: \[ \gamma\mathcal{R}_{et}=\sum_{\ell=1}^{L}S_{\ell t}\left( \sum_{\ell^{\prime }=1}^{L}\tau_{\ell}w_{\ell^{^{\prime}}}p_{\ell\ell^{^{\prime}}}\right) +O\left( n^{-1}\right) , \] Setting $p_{\ell\ell^{^{\prime}}}=k_{\ell\ell^{\prime}}/n_{\ell^{\prime}}$ we can alternatively write $\gamma\mathcal{R}_{et}$ as (recall that $w_{\ell^{\prime}}=n_{\ell^{\prime}}/n$ and $s_{\ell t}=S_{\ell t}/n_{\ell}$)

equation[equation omitted — 123 chars of source]

where $\beta_{\ell}$ is already defined by ((ref)).

In the case of a single group or when $\beta_{\ell}=\beta$ is homogeneous across groups, the above expression simplifies to $\gamma\mathcal{R} _{et}=\beta\left( \sum_{\ell=1}^{L}w_{\ell}s_{\ell t}\right) =\beta s_{t}$, which can be written equivalently as $\mathcal{R}_{et}=(1-c_{t})\mathcal{R} _{0}$. In the absence of any interventions $\mathcal{R}_{et}$ declines as $c_{t}$ rises, and $\mathcal{R}_{et}$ falls below $1$ when $c_{t} >(\mathcal{R}_{0}-1)/\mathcal{R}_{0}$. The value $(\mathcal{R}_{0} -1)/\mathcal{R}_{0}$ is often referred to as the herd immunity threshold. For the multigroup case, using ((ref)) and ((ref)), the condition for herd immunity is more complicated and is given by (for $n$ sufficiently large) \[ \frac{\sum_{\ell=1}^{L}w_{\ell}\beta_{\ell}s_{\ell t}}{\gamma}=\left[ \frac{\sum_{\ell=1}^{L}w_{\ell}\beta_{\ell}\left( 1-c_{\ell t}\right) } {\sum_{\ell=1}^{L}w_{\ell}\beta_{\ell}}\right] \mathcal{R}_{0}<1, \] and the herd immunity threshold becomes \[ \frac{\sum_{\ell=1}^{L}w_{\ell}\beta_{\ell}c_{\ell t}}{\sum_{\ell=1} ^{L}w_{\ell}\beta_{\ell}}>\frac{\mathcal{R}_{0}-1}{\mathcal{R}_{0}}. \] This formula clearly shows that for herd immunity to apply, the group-specific infection rate, $c_{\ell t}$, must be sufficiently large -- shielding one group requires higher infection rates in other groups with larger population weights. To see this, let us consider a simple example of two groups ($L=2$) with a homogeneous transmission rate across the two groups ($\beta_{1} =\beta_{2}=\beta$), and note that $0<w_{1},w_{2}<1$. Suppose that policymakers want to shield Group 1, which may comprise elderly people, from infection. In the extreme case where all individuals in Group 1 are protected, namely, $c_{1t}=0$, then herd immunity requires $c_{2t}>(\mathcal{R}_{0}-1)/\left( \mathcal{R}_{0}w_{2}\right) $, which is higher than the threshold value of $(\mathcal{R}_{0}-1)/\mathcal{R}_{0}$ where the population groups are treated symmetrically.

Social intervention might be necessary if the herd immunity threshold is too high and could lead to significant hospitalization and deaths. In such cases, intervention becomes necessary to reduce the transmission rates $\beta_{\ell} $, thus introducing independent policy-induced reductions in the transmission rates. In the presence of social policy interventions, the effective reproduction number for the multigroup can be written as

equation[equation omitted — 148 chars of source]

where (using ((ref))) $\beta_{\ell t}=\sum_{\ell^{^{\prime}} =1}^{L}\tau_{\ell t}k_{\ell\ell^{\prime},t}.$ Reductions in $\beta_{\ell t}$ can come about either by reducing the average number of contacts within and across groups, $k_{\ell\ell^{\prime},t}$, or by reducing the group-specific exposure intensity parameter, $\tau_{\ell t}$, or both. Since only the product of $\tau_{\ell t}$ and $k_{\ell\ell^{\prime},t}$ is identified, in our simulations we fix the contact patterns and calibrate the desired value of $\beta_{\ell t}$ by setting the value of $\tau_{\ell t}$ for each $\ell$ to achieve a desired $\mathcal{R}$ number. Of course, one would obtain equivalent results if the average number of contacts is assumed to be time-varying and the exposure intensity parameter is assumed constant. In the case of a single group or when $\beta_{\ell t}=\beta_{t}$ for all $\ell$, we have

equation[equation omitted — 88 chars of source]

where $(1-c_{t})$ is the herding component. It is also worth bearing in mind that at the outset of epidemic outbreaks the value of $c_{t}$ is close to zero which ensures that $\mathcal{R}_{e0}=\beta_{0}/\gamma=\mathcal{R}_{0}$.

Calibration and simulation of the model

Although it is difficult to obtain an analytical solution to the individual-based stochastic epidemic model, we can study its properties by simulations. This section focuses on the baseline scenario of no containment measures or mutation of the virus so that the transmission rate is constant. We will discuss simulation results with time-varying transmission rate under social distancing and vaccination in Section (ref). In light of the recent studies on the value of $\mathcal{R}_{0}$ for Covid-19, we set $\mathcal{R}_{0}=3$.\footnote{A summary of published estimates of $\mathcal{R}_{0}$ is provided in Table 1 of D'Arienzo and Coniglio (2020).} For the recovery rate, in view of the World Health Organization guidelines of two weeks self-isolation, we set $\gamma=1/14$.\footnote{ Similar guidelines issued by the US and the UK can be found at \url{https://www.cdc.gov/coronavirus/2019-ncov/if-you-are-sick/quarantine.html} and \url{https://www.nhs.uk/conditions/coronavirus-covid-19/self-isolation-and-treatment/how-long-to-self-isolate/}, respectively (last accessed October 2020).} It follows that $\beta =\gamma\mathcal{R}_{0}=3/14.$

We consider dividing the population into $L=5$ age groups: $[0,15)$, $[15,30)$, $[30,50)$, $[50,65)$, and $65$+ years old, and, of course, one can readily consider a different number of groups based on other characteristics if such data are available. We use the data on Germany as an illustration. The social contact surveys by Mossong et al. (2008) provide rich data on the contact patterns in Germany. We update the contact matrix by age with the most recent population data such that the reciprocity condition, $n_{\ell} k_{\ell\ell^{\prime}}=n_{\ell^{\prime}}k_{\ell^{\prime}\ell}$, is satisfied. The population shares for the five age groups are $\mathbf{w}=(0.13,$ $0.17,$ $0.28,$ $0.20,$ $0.21)^{\prime}$, and the resulting (pre-pandemic) contact matrix is

footnotesize\begin{equation} \mathbf{K}=\left( k_{\ell\ell^{^{\prime}}}\right) =\left( \begin{tabular} [c]{rrrrr} 3.43 & 1.10 & 2.34 & 0.67 & 0.47\\ 0.87 & 4.55 & 2.72 & 1.14 & 0.41\\ 1.11 & 1.64 & 3.74 & 1.42 & 0.78\\ 0.45 & 0.96 & 1.99 & 2.30 & 0.92\\ 0.31 & 0.34 & 1.08 & 0.91 & 1.70 \end{tabular} \right) , \end{equation}

where the element, $k_{\ell\ell^{^{\prime}}}$, represents the average number of daily contacts reported by participants in group $\ell$ with someone of group $\ell^{\prime}$. The larger diagonal values in ((ref)) indicate that people tend to mix more with others of the same age group---a phenomenon well documented by contact surveys across different countries. In order to calibrate $\tau_{\ell}$ across groups, we match the ratio of infection probabilities of groups with the ratio of reported cases. Specifically, denoting the reference group by $\ell_{0}$ and the ratio of reported infections of group $\ell$ to group $\ell_{0}$ by $\lambda_{\ell\ell_{0}}$, then $\lambda_{\ell\ell_{0}}$ should match the ratio of the related probabilities, namely, \[ \lambda_{\ell\ell_{0}}=\frac{E\left( x_{i\ell,t+1}|\mathbf{z}_{t}\right) }{E\left( x_{i\ell_{0},t+1}|\mathbf{z}_{t}\right) }=\frac{1-\prod _{\ell^{^{\prime}}=1}^{L}\left( 1-p_{\ell\ell^{^{\prime}}}+p_{\ell \ell^{^{\prime}}}e^{-\tau_{\ell}}\right) ^{n_{\ell^{\prime}}i_{\ell^{\prime }t}}}{1-\prod_{\ell^{^{\prime}}=1}^{L}\left( 1-p_{\ell_{0}\ell^{^{\prime}} }+p_{\ell_{0}\ell^{^{\prime}}}e^{-\tau_{\ell_{0}}}\right) ^{n_{\ell^{\prime} }i_{\ell^{\prime}t}}}. \] Using ((ref)), we now have \[ \lambda_{\ell\ell_{0}}\approx\frac{1-\exp\left( -\tau_{\ell}\sum _{\ell^{\prime}=1}^{L}i_{\ell^{\prime}t}k_{\ell\ell^{\prime}}\right) } {1-\exp\left( -\tau_{\ell_{0}}\sum_{\ell^{\prime}=1}^{L}i_{\ell^{\prime} t}k_{\ell_{0}\ell^{\prime}}\right) }\approx\frac{\tau_{\ell}}{\tau_{\ell_{0} }}\frac{\sum_{\ell^{\prime}=1}^{L}i_{\ell^{\prime}t}k_{\ell\ell^{\prime}} }{\sum_{\ell^{\prime}=1}^{L}i_{\ell^{\prime}t}k_{\ell_{0}\ell^{\prime}}}. \] For the purpose of calibration, we further assume that $i_{\ell t}=w_{\ell }f_{t}$, where $f_{t}$ can be viewed as the latent common driver of the epidemic at time $t$. It then follows that

equation[equation omitted — 242 chars of source]

We now use data on infected cases in Germany by the five age groups at the end of 2020 (before the rollout of Covid-19 vaccines) to calibrate the relative transmission rates by groups. Setting the first age group as the reference group ($\ell_{0}=1$), we obtain $\boldsymbol{\lambda}=\left( \lambda _{\ell\ell_{0}}\right) =\left( 1,\text{ }2.83,\text{ }3.81,\text{ }2.94,\text{ }2.39\right) ^{^{\prime}}$, which in conjunction with ((ref)) yields $\boldsymbol{\tau}=\left( \tau_{\ell}\right) =\tau_{1}\left( 1,\text{ }2.21,\text{ }3.04,\text{ }3.15\text{ },3.93\right) ^{\prime}$. To calibrate $\tau_{1}$, we use ((ref)) and obtain $\tau_{1}=0.011$ setting $\mathcal{R}_{0}=3$ and $\gamma=1/14$.

For each replication, the simulation begins with $1/1000$ of total population randomly infected on day $t=1$, that is, $c_{1}^{\left( b\right) } =i_{1}^{\left( b\right) }=0.001$ and $r_{1}^{\left( b\right) }=0$, where $b$ denotes the $b^{th}$ replication, for $b=1,2,\ldots,B$.\footnote{We find that there will not be outbreaks in many replications if the simulation begins with less than $1/1000$ of the population initially infected.} Then from $t=2$ onwards, the infection and recovery processes follow ((ref)) and ((ref)), respectively, for $\ell=1,2,\ldots,L.$ The proportion of infections for each age group is computed as $c_{\ell t}^{\left( b\right) }=C_{\ell t}^{\left( b\right) }/n_{\ell}=\sum_{i=1}^{n_{\ell}}x_{i\ell ,t}^{\left( b\right) }/n_{\ell}$, and the daily new cases are computed by $\Delta c_{\ell t}^{\left( b\right) }=c_{\ell t}^{\left( b\right) }-c_{\ell,t-1}^{\left( b\right) }$. The aggregate infections and new cases are computed as $c_{t}^{\left( b\right) }=\sum_{\ell=1}^{L}w_{\ell}c_{\ell t}^{\left( b\right) }$ and $\Delta c_{t}^{\left( b\right) }=c_{t}^{\left( b\right) }-c_{t-1}^{\left( b\right) }$, respectively. Details on the generation of random networks are given in Section (ref) of the online supplement. Note that the contact network randomly changes every day (and also across replications). This feature captures the random nature of many encounters an individual has on a daily basis. We consider $B=1,000$ replications and set the population size to $n=10,000$. We also tried larger population sizes, but, as will be seen below, the interquartile range of the simulated new cases is already very tight when $n=10,000$. Some simulation results for $n=50,000$ and $n=100,000$ are provided in Figure (ref) of the online supplement.

Figure (ref) displays the simulated proportion of group-specific and aggregate new cases in fan chart style with the $10^{th}$, $25^{th}$, $50^{th}$, $75^{th},$ and $90^{th}$ percentiles over $1,000$ replications. The mean values are very close to the median and not shown. We also report the maximum proportion of infected for each group averaged across replications, i.e., $c_{\ell}^{\ast}=B^{-1}\sum_{b=1}^{B}\max_{t}c_{\ell t}^{(b)}$, and the maximum proportion of aggregate infected, $c^{\ast} =B^{-1}\sum_{b=1}^{B}\max_{t}c_{t}^{(b)}$. The duration of the epidemic, denoted by $T^{\ast}$, is computed as the number of days to reach zero active cases averaged across replications.\footnote{Note that the model implies that the disease will not spread again once $i_{t}$ becomes zero.}

Figure (ref) shows that if the disease transmits at a fixed $\mathcal{R}_{0}=3$, the youngest age group will have the lowest maximum proportion of infections, ending up with $62$ percent infected in comparison to over $90$ percent infected in the other groups. The uncontrolled epidemic is expected to end about $215$ days after the outbreak. The daily new cases for the five groups peak around the same time (about $50$ days on average), with the highest daily infection ranging from $2.0$ percent in the youngest group to $3.7$ percent in the middle-aged group (Group 3). As a whole, the maximum aggregate infection rate will reach $90$ percent, with daily new cases peaking at $2.9$ percent of the population.

figure[figure omitted — 1,489 chars of source]

In order to examine whether the number of groups affects the aggregate outcomes, we carry out simulations using a single group model and compare the results with the aggregate outcomes using the multigroup model. When there is only one group, the contact network reduces to the Erd\H{o}s-R\'{e}nyi random network (simply referred to as the random network below), where each pair of the nodes (or individuals) are connected at random with a uniform probability $p=k/\left( n-1\right) \approx k/n,$ where $k$ is the mean degree of the network (or the mean number of contacts per individual).\footnote{The degree of a node in a network is the number of connections it has (or the number of edges attached to it).} We set the average number of contacts to $k=10$ based on the literature on social contacts in the pre-Covid period, and then set the exposure intensity parameter to $\tau=\beta/k=\gamma\mathcal{R}_{0}/k$, where, as before, $\gamma=1/14$ and $\mathcal{R}_{0}=3$. Figure (ref) of the online supplement shows that the simulated aggregate outcomes are very close under the single- and multi-group models. This finding is reasonable because the simulations were performed with the same fixed $\beta=\gamma\mathcal{R}_{0}$. The heterogeneity in $\tau_{\ell}$ and $p_{\ell\ell^{^{\prime}}}$ affects the epidemic curves for each group but does not seem to impact the aggregate outcomes. This result suggests that an aggregate analysis may be justified if the primary focus is on the spread of the infection across the population as a whole rather than on particular age/type groups.

Lastly, to investigate whether the results are robust to different network topologies, we considered another widely used contact network---the power law random network, in which a small number of nodes (individuals) may have a relatively high number of links (contacts). Figure (ref) of the online supplement shows that the simulation outcomes obtained by the Erd\H{o}s-R\'{e}nyi and the power law networks with the same average number of contacts are very similar.

Estimation of transmission rates

The previous section investigates the properties of the model, assuming the transmission rates are given. This section turns to detailing how to estimate the transmission rate using data on infected cases. We first derive the method of moments estimation of the transmission rate when there are no measurement errors and present the finite sample properties of the estimators using Monte Carlo techniques. We then allow for under-reporting of infected cases and propose a recursive joint estimation of the transmission rate and the degree of under-reporting by a simulated method of moments.

Estimation without measurement errors

Let us first consider the case of a single group, and recall that the moment condition for this case is given by ((ref)), which is replicated here for convenience

equation[equation omitted — 155 chars of source]

We can estimate the transmission rate, $\beta$, using ((ref)) by nonlinear least squares (NLS) given time series data on $\left\{ c_{t} ,i_{t}\right\} $. The recovery rate, $\gamma$, can be estimated using the recovery equation, ((ref)), $E\left( r_{t+1}|r_{t},c_{t}\right) =\left( 1-\gamma\right) r_{t}+\gamma c_{t}$. Nevertheless, in reality, $r_{t}$ is often not recorded in a timely manner and $\gamma$ is estimated from the hospitalization data. We therefore set $\gamma=1/14$ in our estimation and calibration exercises,\footnote{For the rationale behind setting $\gamma=1/14$, see Footnote (ref).} and discuss the properties of the moment estimator of $\gamma$ in the online supplement.

In the absence of any interventions (voluntary or mandatory), we have $\beta=\gamma\mathcal{R}_{0}$, where, as before, $\mathcal{R}_{0}$ is the basic reproduction number. It follows that $\mathcal{R}_{0}$ can be estimated by $\mathcal{\hat{R}}_{0}=\hat{\beta}/\gamma$, where $\hat{\beta}$ is the NLS estimate of $\beta$ using ((ref)). Under social interventions, the recovery equation holds (since $\gamma$ is unaffected), but the moment condition for $c_{t+1}$ now depends on the time-varying transmission rate, $\beta_{t}$. For $\gamma$ in the range of $1/14$ to $1/21$, it is reasonable to use two or three weeks rolling windows when estimating $\beta_{t}$. For a window of size $W$, we have

equation[equation omitted — 193 chars of source]

Note that even though\ the time series $\left\{ c_{t},i_{t}\right\} $ over the course of the epidemic are non-stationary, the rolling estimation is based on short-$T$ series ($T=14$ or $21$).

To examine the finite sample performance of $\hat{\beta}_{t}\left( W\right) $, we estimate $\beta$ using the simulated data generated from the single group stochastic SIR\ model on a random network with mean contact $k=10$ and assuming $1/1000$ of the population is randomly infected on day $1$. The true value of the transmission rate is set to $\beta=3/14$ such that $\mathcal{R} _{0}=\beta/\gamma=3$. We consider population sizes $n=10,000$, $50,000,$ and $100,000$, and set the the number of replications to $B=1,000$. Recall that $T$ is fixed and $n\rightarrow\infty$. To alleviate noise induced by zero and near-zero observations at the start and final stages of the epidemic, the rolling estimation of $\beta$ is carried out over the $4^{th}-15^{th}$ weeks after the outbreak.

Since the value of $\beta$ is quite small, we present the estimation results in terms of $\mathcal{\hat{R}}_{0}\left( W\right) =\hat{\beta}_{t}\left( W\right) /\gamma$. Table (ref) summarizes the bias and root mean square error (RMSE) of the 2-weekly rolling estimates of $\mathcal{R}_{0},$ averaged over the four non-overlapping 3-weekly sub-samples, for different population sizes. The bias is computed as $B^{-1}\sum_{b=1}^{B}\left[ \mathcal{\hat{R}}_{0}^{\left( b\right) }\left( W\right) -\mathcal{R} _{0}\right] $, and the RMSE is computed by $\sqrt{B^{-1}\sum_{b=1}^{B}\left[ \mathcal{\hat{R}}_{0}^{\left( b\right) }\left( W\right) -\mathcal{R} _{0}\right] ^{2}}$, where $\mathcal{\hat{R}}_{0}^{\left( b\right) }\left( W\right) =\hat{\beta}_{t}^{(b)}\left( W\right) /\gamma$ and $\hat{\beta }^{(b)}\left( W\right) $ is the estimate of $\beta$ in the $b^{th}$ simulated sample. As can be seen from Table (ref), although $\mathcal{\hat{R}}_{0}$ tends to slightly underestimate $\mathcal{R}_{0}$, its bias and RMSE\ are quite small in all experiments and sub-samples. The RMSE declines as the population size $n$ increases, but $\mathcal{R}_{0}$ can be estimated reasonably well even with $n=10,000$. Comparing the results over different epidemic stages, the RMSE is relatively larger at the early and late stages of the epidemic. This finding is not surprising since it is difficult to obtain precise estimates when $c_{t}$ and $i_{t}$ are near zero. We also considered the 3-weekly rolling estimates, reported in Table (ref) of the online supplement, and as can be seen are very close to the 2-weekly estimates, with slightly better performance in the early and late stages of the epidemic. We will hereafter mainly focus on the 2-weekly rolling estimation.

table[table omitted — 1,558 chars of source]

Similar moment conditions can also be used to estimate the parameters of the multigroup model. If time series data on $\left\{ c_{\ell t},i_{\ell t}\right\} $, for $\ell=1,2,\ldots,L$ ($L$ is finite) are available, we can estimate $\beta_{\ell\ell^{\prime}}$ using the moment conditions ((ref)), namely, \[ E\left( \left. \frac{1-c_{\ell,t+1}}{1-c_{\ell t}}\right\vert \text{ }\mathbf{i}_{t}\right) =\exp\left( -\sum_{\ell^{\prime}=1}^{L}\beta _{\ell\ell^{\prime}}i_{\ell^{\prime}t}\right) +O\left( n^{-1}\right) . \] Then, as we have discussed in Section (ref), $\beta _{\ell\ell^{\prime}}$ is identifiable from ((ref)) for given values of $\gamma$ and $\mathbf{w=(}w_{1},w_{2},\ldots,w_{L})^{\prime}.$

Estimation allowing for measurement errors

It is widely recognized that in practice $c_{t}$ and $r_{t}$ are under-reported. The magnitude of under-reporting is measured by the multiplication factor (MF) in the literature (see, e.g., Gibbons et al., 2014). It is expected that the MF will decline over time since data quality will improve as more testing is conducted, but, in any case, MF is certainly greater than one. Denoting the multiplication factor by $m_{t},$ and denoting the observed values of $c_{t}$ and $i_{t}$ by $\tilde{c}_{t}$ and $\tilde{\imath}_{t}$, respectively, we have $c_{t}=m_{t}\tilde{c}_{t}$ and $i_{t}=c_{t}-r_{t}=m_{t}\tilde{\imath}_{t}$ (assuming that $r_{t} =R_{t}/n=m_{t}\tilde{r}_{t}$, where $\tilde{R}_{t}$ is the observed value of $R_{t}$). Then the moment condition in terms of the observed values ($\tilde{c}_{t}$ and $\tilde{\imath}_{t}$) can be written as

equation[equation omitted — 225 chars of source]

It can be seen from ((ref)) that $m_{t}$ is not identified when $\tilde {c}_{t}$ and $\tilde{\imath}_{t}$ are very small in the early stage of the epidemic. When $\tilde{c}_{t}$ becomes large enough, we can estimate $m_{t}$ by the simulated method of moments based on ((ref)). In practice, $m_{t}$ varies slowly, and it is reasonable to assume $m_{t}=m_{t-1}$ within a short time interval (two or three weeks). Then we have

equation[equation omitted — 153 chars of source]

where $i_{t}^{\left( b\right) }$ denotes the simulated value of $i_{t}$ in the $b^{th}$ replication and $B$ is the total number of replications. Solving ((ref)) for $m_{t}$ yields

equation[equation omitted — 224 chars of source]

It is now clear that one can estimate $m_{t}$ by ((ref)) for given values of $\beta_{t}$, and estimate $\beta_{t}$ by ((ref)) if $m_{t}$ is known. Accordingly, we propose a method that estimates $\beta_{t}$ and $m_{t}$ jointly. The algorithm is described in detail in Section (ref) of the online supplement. We apply the procedure recursively using 2- and 3-weekly rolling windows in the next section to examine how the transmission rates and under-reporting of cases changed over time in a number of countries and evaluate how our model matches the Covid-19 evidence.

Matching the model with evidence from a number of European countries

We now assess how our model matches with the recorded cases in six European countries: Austria, France, Germany, Italy, Spain, and the United Kingdom, while taking account of under-reporting of infections.\footnote{We also examined how our model matches the Covid-19 evidence in the US. The estimates of $\mathcal{R}_{et}$ for the country as a whole and for each of the $48$ mainland states and the District of Columbia are displayed in Section (ref) of the online supplement. The estimates of the multiplication factor and the comparison between the realized and calibrated new cases are presented in Figure (ref) of the online supplement.} The Covid-19 outbreak in Europe began with Italy in early February 2020, with the recorded number of infections accelerating rapidly from February 21 onward. A rapid rise in infections took place about one week later in Spain, France, and Germany, followed by UK and Austria at the end of February.

To estimate $\beta_{t}$, we need observations on per capita infected and active cases, $c_{t}$ and $i_{t}$. Using the recorded number of infected cases, $\tilde{C}_{t}$, and population data, the per capita cases, $\tilde {c}_{t}$, are readily available, where as before we use the tilde symbol to indicate observed values. Since $\tilde{I}_{t}=\tilde{C}_{t}-\tilde{R}_{t}$, we can obtain $\tilde{I}_{t}$ if the number of removed (the sum of those who recovered and the deceased) cases, $\tilde{R}_{t}$, is available. Unfortunately, the recovery data is either not reported or is subject to severe measurement error/reporting issues in many countries. For all six countries, we therefore estimate the number of removed using the recursion $\tilde{R}_{t}=\left( 1-\gamma\right) \tilde{R}_{t-1}+\gamma\tilde{C}_{t-1} $, for $t=2,3,\ldots,T$, where the recovery rate $\gamma$ is set to $1/14,$ and values of $\tilde{R}_{t}$ are generated starting with $\tilde{R} _{1}=\tilde{C}_{1}=0$. We then compute $\tilde{I}_{t}$ by subtracting the estimated $\tilde{R}_{t}$ from the recorded $\tilde{C}_{t}$ .\footnote{Specifically, among the six countries, the recorded data on recovery are unavailable for Spain and UK; they are of poor quality for France and Italy; they are relatively close to our estimated recovery for Austria and Germany. We have also calibrated our model using the recorded recovery data for Austria and Germany as a robustness check and obtained similar results.} That is, in this empirical exercise, we only need data on Covid-19 cases per capita, $\tilde{c}_{t}$. To alleviate the wide fluctuations in the data due to irregular update schedules and reporting/recording delays, we smooth the series by taking the 7-day moving average before they are used in the estimation and calibration.

We adopt the joint estimation approach proposed in the last section to calibrate and evaluate our stochastic network model. The procedure is applied recursively using 2- and 3-weekly rolling windows. Recall that in the early stage of an epidemic, $\tilde{c}_{t}$ is small, and MF\ is not identified. We choose $\tilde{c}^{0}=0.01$ as the threshold value.\footnote{For the countries we considered, there is virtually no difference in the estimates of $\beta _{t}$ when $c_{t}\leq0.01$ if MF takes the value from $2$ to $7$.} When $\tilde{c}_{t}\leq\tilde{c}^{0}$, we use an initial guess of the multiplication factor, MF $=5$, in the estimation of $\beta_{t}$. Since the early Covid data are quite noisy, we start the rolling estimation when the daily new cases exceed one per $100,000$ people and use the estimates below $\mathcal{R}_{0}=3$ in the calibration. Specifically, the simulations begin with $1/1000$ of the population randomly infected on day $1$. To render the calibrations comparable across the countries in our sample, during the first week after the outbreak we set the value of $\beta$ such that $\mathcal{R} _{0}$ equals its first estimate, $\hat{\beta}_{t}/\gamma$, that is less than $3$. Then, from the second week onwards, we set $\beta_{t}$ to the rolling estimates computed from the realized data (with MF $=5$) until $\tilde{c}_{t}$ reaches $\tilde{c}^{0}$ on day $t^{0}$. As shown in Section (ref), it makes little difference to the aggregate outcomes whether we carry out the simulations using single- or multi-group models. Since we are interested in comparing the calibrated outcomes with realized cases, we conduct simulations using the single group model with the Erd\H{o}s-R\'{e}nyi random network in this exercise. The first estimate of MF is computed as the ratio of the average calibrated cases to realized cases on day $t^{0}$. When $\tilde{c}_{t}>\tilde{c}^{0}$, we perform the joint rolling estimation of $\beta_{t}$ and $m_{t}$ using ((ref)) and ((ref)). We present the 2-weekly estimation and calibration results in the main paper. The results using the 3-weekly rolling windows are very close and are given in the online supplement.\footnote{Figures (ref) and (ref) of the online supplement compares the 2- and 3-weekly estimates of $\mathcal{R}_{et}$ and MF, respectively. We find that the 2- and 3-weekly estimates of $\mathcal{R} _{et}$ are quite similar. The 3-weekly estimates of MF tend to be slightly higher than the 2-weekly estimates, but overall they are very close.} Since the moment condition, ((ref)), used in the joint estimation was derived assuming no vaccination, we end the joint estimation when the recorded share of the population fully vaccinated reaches $10$ percent. The population size in simulations is set to $n=50,000$. To ease the computational burden, the number of replications is set to $B=500$.

figure[figure omitted — 1,805 chars of source]

Figure (ref) shows the evolution of $\mathcal{\hat{R}}_{et}$ over the period March 2020--April 2021 for the six countries, where $\mathcal{\hat {R}}_{et}=\left( 1-\hat{m}_{t}\tilde{c}_{t}\right) \hat{\beta}_{t}/\gamma$ and $\hat{m}_{t}$ is displayed in Figure (ref).\footnote{In the early stages of the epidemic when $c_{t}$ is small, estimates of $\mathcal{R}_{et}$ and $\beta_{t}/\gamma$ are very close, even if we set MF to $10$. See Figure (ref) in the online supplement where $\mathcal{\hat{R}}_{et}$ are compared with $\hat{\beta}_{t}/\gamma$ for the six countries.} It also marks the start and end dates of the respective lockdowns.\footnote{The lockdown dates across countries can be found at \url{https://en.wikipedia.org/wiki/COVID-19_pandemic_lockdowns}. One could consider more accurate measures of the strictness of the lockdowns. For example, Chudik, Pesaran, and Rebucci (2021) studied the impact of the OxCGRT's stringency index on the estimated $\mathcal{R}_{et}$.} It should be noted that the epidemic tends to expand (contract) if $\mathcal{\hat{R}}_{et}$ is above (below) unity.\footnote{See also Figure (ref) of the online supplement for graphs of recorded daily new cases for these countries.} Among the six countries, Italy started the first national lockdown on March 9, 2020, followed by Spain, Austria, and France about a week later, and Germany and the UK two weeks later (on March 23, 2020). As can be seen from Figure (ref), $\mathcal{\hat{R}}_{et}$ fell below one in mid to late April 2020 in all these countries except for the UK, which took a bit longer before falling below unity in early May. On average, it took $36$ days to bring $\mathcal{\hat{R}}_{et}$ down below one from the start of the lockdown, with Germany being the fastest ($27$ days) and the UK being the slowest ($47$ days). By the end of May 2020, $\mathcal{\hat{R}}_{et}$ were brought down below $0.5$ in all six countries except for the UK, where the lowest value of $\mathcal{\hat{R}}_{et}$ occurred in early July. As lockdowns eased, not surprisingly, the transmission rates started to rise. This new surge in estimates of $\mathcal{R}_{et}$ led some of the countries to announce their second lockdowns in early November 2020. The estimates of $\mathcal{R} _{et}$ fell below one again in December 2020, but the second trough in $\mathcal{\hat{R}}_{et}$ is higher than the first in all countries except France. As the pandemic progressed, $\mathcal{\hat{R}}_{et}$ displayed different patterns (timing and magnitudes of peaks and troughs) across the six countries due to different containment measures. By late April 2021 (the end of our sample), $\mathcal{\hat{R}}_{et}$ is estimated to be close to one in all these countries, but they appear to be rising in Spain and the UK and falling in the other countries.

figure[figure omitted — 1,009 chars of source]

Figure (ref) plots the estimated MF for the six countries. The results offer evidence of substantial under-reporting in the pandemic's early stages, with the magnitude of under-reporting falling over time in all these countries. Closer inspection of the figure shows that the number of cases in Austria, Germany, and Italy in late October-mid November 2020 was underestimated by $5$--$6$ times, which declined to $2$--$3$ times in late April 2021. About a quarter of actual infections were recorded in Spain in September 2020, compared to about a half being detected in late April 2021. France and the UK have a greater level of under-reporting in the early stages---the number of cases was underestimated by as much as a factor of $9$--$10$, which fell to $2$--$4$ during the study period. Overall, the magnitude of these estimated MF seems reasonable and comparable to the estimates obtained by other approaches in the literature.\footnote{See, for example, Jagodnik et al. (2020), Li et al. (2020), Havers et al. (2020), Kalish et al. (2021), and Rahmandad, Lim, and Sterman (2021). See also Section (ref) of the online supplement.}

figure[figure omitted — 1,103 chars of source]

Figure (ref) presents the calibrated new cases and the 7-day moving average of the reported new cases multiplied by the estimated MF. The fan charts depict the $10^{th}$ through the $90^{th}$ percentiles of the calibrated data. It can be seen that once we have taken account of under-reporting, the calibrated cases match with the recorded cases fairly well. It is noteworthy that our model is able to catch the multiple waves\ of Covid-19 cases over the course of the epidemic. Lastly, it is interesting to see how the total cases per capita compare across the six countries, with and without adjustments for under-reporting. Figure (ref) of the online supplement displays the reported total cases and the case numbers after adjusting for under-reporting using the MF estimates. The results show that the number of total cases could have been underestimated three to five times in these countries as of early August 2021. These comparisons clearly show the importance of adjusting the number of infected cases due to under-reporting, which can be reasonably estimated using our joint estimation procedure.

Counterfactual exercises

Having shown that the outcomes of the calibrated model closely match the evidence, we now demonstrate how the model can be used for two counterfactual analyses of interest. First, we investigate the impact of social distancing and vaccination on the evolution of the epidemic. To simplify the exposition, we consider an epidemic with two waves and investigate if the second wave can be avoided by vaccination. Second, we consider counterfactual outcomes that could have resulted from different timing of the early interventions in Germany and the UK.

Social distancing and vaccination

In order to understand the impact of non-pharmaceutical interventions and vaccination on controlling the epidemic, we perform counterfactual analyses using the multigroup model with the five age groups introduced in Section (ref). The different age groups in the model also allow us to consider the implications of prioritizing the Covid-19 vaccine by age.

In reality, the transmission rate varies over time due to both voluntary and mandatory social distancing as well as other mitigation measures such as vaccination. Here we use social distancing to refer broadly to all types of non-pharmaceutical interventions (including lockdown measures). We assume that, in the absence of a vaccination program, the (scaled) transmission rate, $\beta_{t}/\gamma$, equals $3$ in the first two weeks after the epidemic outbreak, falls to $0.9$ linearly over the next three weeks, and remains at $0.9$ for eight weeks. When social distancing is relaxed, the transmission rate increases to $1.5$ linearly over the next three weeks and remains at $1.5$ thereafter. Note that the effective reproduction number, $\mathcal{R} _{t}$, could still fall below unity due to herding.\footnote{Figure (ref) of the online supplement displays the time profile of the transmission rate under this social distancing policy.} Also note that $\beta_{t}=\sum_{\ell=1}^{L}w_{\ell}\beta_{\ell t}$, where $\beta_{\ell t}=\sum_{\ell^{\prime}=1}^{L}\tau_{\ell t}k_{\ell\ell^{\prime},t}$. We assume that $\beta_{\ell t}$ has the same rate of change as $\beta_{t}$, for all $\ell$.

To model vaccination, we assume, for simplicity, that a single-dose shot vaccine with efficacy of $\epsilon_{v}$ becomes available when $i_{t}=i^{0}$. The vaccine takes full effect immediately, and its immune protection does not wane over time.\footnote{On August 5, 2021, Moderna reported that its vaccine efficacy remained almost the same through six months after the second shot. \url{https://time.com/6087722/moderna-vaccine-long-term-efficacy/}} The effectiveness of vaccination can be measured by the parameter $\mu_{i\ell}$, which is the mean of the infection threshold variable, $\xi_{it}$, defined in ((ref)). We assume that $\mu_{i\ell}$ takes the value $\mu^{0}$ if individual $(i,\ell)$ is not vaccinated and takes the value $\mu^{1}$ after vaccination. In the case of a single group, the probability of an individual getting infected when the proportion of active cases is $i^{0}$, for any given value of $\mu_{i}$, is

equation[equation omitted — 173 chars of source]

which declines with $\mu_{i}$. By definition, the vaccine efficacy should equal the percentage reduction in the probability of infection. Then the value of $\mu_{i}$ associated with efficacy $\epsilon_{v}$ is given by

align[align omitted — 315 chars of source]

Using the result in ((ref)), we have

equation[equation omitted — 186 chars of source]

Combining ((ref)), ((ref)), and ((ref)), we obtain $\tau ki^{0}/\mu^{1}\approx\left( 1-\epsilon_{v}\right) \tau ki^{0}/\mu^{0}$, which simplifies to\footnote{The exact solution is very close to its approximation given by ((ref)).}

equation[equation omitted — 72 chars of source]

This result also holds in the multigroup model.\footnote{A\ proof is given in Section (ref) of the online supplement.} Intuitively, ((ref)) states that an individual becomes $1/\left( 1-\epsilon_{v}\right) $ times more immune relative to his/her level of immunity after vaccination.

In simulations, an individual's degree of resilience, $\xi_{it}$, are i.i.d. draws from an exponential distribution with mean $\mu^{0}$ $(\mu^{1})$ before (after) vaccination. Without loss of generality, we normalize $\mu^{0}=1.$ The Pfizer-BioNTech and Moderna vaccines have been shown to have $95\%$ and $94.1\%$ efficacy rates in preventing symptomatic laboratory-confirmed Covid-19 infection, respectively (Oliver et al., 2020, 2021). Accordingly, using ((ref)) we have $\mu^{1}=20$ for $\epsilon_{v}=0.95$, namely Pfizer and Moderna vaccines increase the level of immunity by a factor of $20$.\footnote{We also considered $\epsilon_{v}=0.66$, which is in line with a $66.3\%$ efficacy rate of the Johnson & Johnson vaccine (Oliver et al., 2021). The results are presented in Section (ref) of the online supplement.}

We suppose that $75$ percent of the population is vaccinated over $12$ weeks, with an equal number of people vaccinated each day.\footnote{We chose to consider constant daily vaccination rate and a relatively short period as an example. Of course, one could consider increasing daily rate over a more extended period if desired.} We consider two vaccination schemes---random vaccination and vaccination in decreasing age order. Under the former, people are randomly selected without replacement for vaccination irrespective of their age. In the latter, older people are vaccinated first. Individuals within an age group are randomly selected for vaccination on each day when their group is eligible. After all people in the oldest group have been vaccinated, the second-oldest group becomes eligible. This process continues until $75$ percent of the population is vaccinated. In both schemes, we assume that vaccination eligibility does not depend on whether an individual is susceptible, infected, or recovered.

Let us first consider the random vaccination scheme. Figure (ref) compares the aggregate outcomes when there are (a) no containment measures, (b) social distancing only, (c) vaccination only, and (d) combined social distancing and vaccination. Specifically, the transmission rate, $\beta_{t}/\gamma$, is fixed at $3$ in the absence of social distancing (i.e., in cases (a) and (c)). In case (c), the vaccination starts from the $4^{th}$ week after the outbreak$.$ In case (d), the vaccination starts during the last month of social distancing (i.e., the $10^{th}$ week after the outbreak). In cases (c) and (d), $75$ percent of the population is randomly vaccinated over $12$ weeks, and the vaccine efficacy is set at $\epsilon _{v}=0.95$.

figure[figure omitted — 2,336 chars of source]

The results show that social distancing can quickly bring down the daily case rate and thus reduce the demands on the healthcare system. However, as social distancing restrictions are relaxed, a second wave is expected to emerge. The second wave may have a higher peak than the first wave if the transmission rates rise too fast due to increasing contact rates or the exposure intensity. Comparing graphs (a) and (b) reveal that social distancing alone can reduce the maximum proportion of cases from $90$ percent in an uncontrolled epidemic to $50$ percent. Nonetheless, the duration of the epidemic could more than double, and an enduring epidemic may entail high social and economic costs. If vaccination is the only containment tool, it must be implemented soon enough to curb the spread of the disease, especially for a highly contagious disease with $\mathcal{R}_{0}$ about $3$ (or even higher as evidenced by the new variants). Graph (c) shows that even if (in a very unlikely scenario) a highly effective vaccine becomes available during the $4^{th}$ week after the outbreak, $59$ percent of the population could end up getting infected. In reality, developing a new vaccine takes considerable time. Therefore, vaccination is not a substitute for non-pharmaceutical interventions, which are necessary to slow the spread of the disease, allowing more time for vaccine development. Vaccination can effectively prevent the second wave if it is rolled out during the last month of social distancing, as shown in Graph (d). Under the assumption that $75$ percent of the population end up getting vaccinated with efficacy of $\epsilon_{v}=0.95$, the maximum proportion of infected could reduce to $12$ percent, and the number of highest daily new infections could be more than $10$ and $7$ times lower than that in cases (a) and (c), respectively. We also examined the implications of different vaccination coverages, start times, speeds of delivery, and vaccine efficacies. The results are summarized in Section (ref) of the online supplement, where we provide counterfactual outcomes assuming (i) 50 percent versus 75 percent vaccination coverage, (ii) vaccination starts at the end of social distancing versus during the last month of social distancing, (iii) 75 percent of the population getting vaccinated over 8 versus 12 weeks, and (iv) $\epsilon_{v}=0.95$ versus $0.66$.

figure[figure omitted — 2,498 chars of source]

Figure (ref) compares the simulation outcomes under random vaccination and vaccination in decreasing age order for each age group and the entire population, assuming the same social distancing policy as described above. In this experiment, the vaccination starts during the last month of social distancing. $75$ percent of the population is vaccinated over $12$ weeks, and the vaccine efficacy is $\epsilon_{v}=0.95$.\footnote{We also examined the case if $50$ percent of the population is vaccinated over eight weeks. See Figure (ref) of the online supplement.} The results show that the maximum proportion of infected in the oldest group is reduced by $2$ percentage points if the old gets vaccinated first, compared to random vaccination. Not surprisingly, the cost of protecting the elderly is reflected in the higher infection rates in the younger age groups, increasing the maximum cases by $1$, $4$, and $2$ percentage points for age groups 1 to 3, respectively. The vaccine effectively curbs the spread of the disease and prevents the second surge of cases in the two senior groups. A comparison of the aggregate outcomes reveals that prioritizing the old results in a higher level of overall infections and a longer duration of the epidemic, owning to higher contact rates of the younger population. Of course, how to prioritize vaccines is a complex decision requiring further information on the rates of hospitalization and death in each age group. It also requires evaluating the social and economic costs of high infection rates and lockdown measures among young people.

Counterfactual outcomes of early interventions in UK and Germany

We now turn to different counterfactual outcomes that could have resulted from different timing of the first lockdowns in Germany and the UK, focusing on the first wave\ of Covid-19 that leveled off at the end of June 2020 in both countries.\footnote{For example, Neil Ferguson, once an advisor to the UK government, stated on June 10, 2021, that "Had we introduced lockdown measures a week earlier, we would have reduced the final death toll by at least a half". See \url{https://www.politico.com/news/2020/06/10/boris-johnson-britain-coronavirus-response-312668}.} In particular, we investigate the quantitative effect of bringing forward the lockdown in the UK on the number of infected cases, as compared to the effect of delaying the lockdown in Germany. To this end, we shift the estimated $\beta_{t}$ values backward or forward for one or two weeks. As shown in Figure (ref), if the German lockdown had been delayed by one week, the maximum proportion of infected cases would have increased from $2.2$ to $5.0$ percent, and the maximum proportion of active cases would have risen from $0.6$ to $1.5$ percent. In contrast, if the UK lockdown had been brought forward by one week, the model predicts that the maximum proportion of infected cases would have reduced from $5.3$ to $2.3$ percent, and the maximum number of active cases would have reduced from $1.2$ to $0.5$ percent. These results suggest that the UK could have achieved a similarly low level of infected cases per capita as Germany if it had implemented social distancing sooner. The maximum proportion of infected (active cases) is estimated to rise further to $10.8$ ($3.2$) percent if the German lockdown was delayed by two weeks, and the maximum proportion of infected (active cases) is estimated to decrease further to $1.1$ ($0.3$) percent if the UK lockdown was brought forward by two weeks.\footnote{See Figure (ref) of the online supplement.} In summary, this counterfactual exercise shows that it is critical to take measures to lower the effective reproduction number as early as possible if a policymaker aims to control the number of infected and active cases.

figure[figure omitted — 2,522 chars of source]

Concluding remarks

This paper has developed a stochastic network SIR model for empirical analyses of the Covid-19 pandemic across countries or regions. Moment conditions are derived for the number of infected and active cases for the single group as well as multigroup models. It is shown how these moment conditions can be used to identify the structural parameters and provide rolling estimates of the transmission rate in different phases of the epidemic. To allow for time-varying under-reporting of cases, it proposes a method that jointly estimates the transmission rate and the multiplication factor using a simulated method of moments approach. In empirical applications to six European countries, the estimates of the transmission rate are used to calibrate the proposed epidemic model. It is shown that the simulated outcomes are reasonably close to the reported cases once the under-reporting of cases is addressed. The multiplication factors are found to be declining over the course of the pandemic. It is estimated that the actual number of infections could be between $4$--$10$ times higher than the number of reported cases around October 2020, whereas only $2$--$3$ times higher in April 2021. The multigroup model is used for counterfactual analyses of the impact of social distancing and vaccination on the evolution of the epidemic. It is shown that lockdown measures are needed to slow down the spread of a highly contagious disease such as Covid-19, buying time for the development of vaccines and treatments. Vaccination can prevent additional waves of epidemics as social distancing is eased after lockdowns if it is introduced early enough. The calibrated model is also used for empirically-based counterfactual analyses of the first lockdowns in Germany and the UK. It is shown that the UK could have achieved an outcome similar to that experienced by Germany during the first wave if she had started the lookdown just one week earlier. Almost symmetrically, Germany would have experienced much higher infection rates (similar to the UK's experience) if she had started the lockdown one week later.

References

singlespace\begin{footnotesize} \begin{list}{{8pt}{0.1in} {-0.1in}{-0.1in} {0.6em}} • Chudik, A., M. H. Pesaran, and A. Rebucci (2021). COVID-19 time-varying reproduction numbers worldwide: An empirical analysis of mandatory and voluntary social distancing. NBER working paper No. 28629. • D'Arienzo, M. and A. Coniglio (2020). Assessment of the SARS-CoV-2 basic reproduction number, R0, based on the early phase of COVID-19 outbreak in Italy. Biosafety and Health 2(2), 57--59. • Del Valle, S. Y., J. M. Hyman, and N. Chitnis (2013). Mathematical models of contact patterns between age groups for predicting the spread of infectious diseases. Mathematical Biosciences and Engineering 10, 1475. • Elliott, S. and C. Gourieroux (2020). Uncertainty on the reproduction ratio in the SIR model. arXiv preprint: 2012.11542. • Farrington, P., and H., Whitaker (2003). Estimation of Effective Reproduction Numbers for Infectious Diseases Using Serological Survey Data. Biostatistics, 4, 621--632. • Gibbons, C. L., M.-J. J. Mangen, D. Plass, A. H. Havelaar, R. J. Brooke, P. Kramarz, ... M. E. Kretzschmar (2014). Measuring underreporting and underascertainment in infectious disease datasets: A comparison of methods. BMC Public Health 14(1), 147. • Guo, H., M. Y. Li, and Z. Shuai (2006). Global stability of the endemic equilibrium of multigroup SIR epidemic models. Canadian Applied Mathematics Quarterly 14(3), 259--284. • Havers, F. P., C. Reed, T. Lim, J. M. Montgomery, J. D. Klena, A. J. Hall, ... N. J. Thornburg (2020). Seroprevalence of antibodies to SARS-CoV-2 in 10 sites in the United States, March 23-May 12, 2020. JAMA Internal Medicine, \textit{180}(12), 1576--1586. • Hethcote, H. W. (2000). The mathematics of infectious diseases. \textit{SIAM Review 42}(4), 599--653. • Jagodnik, K., F. Ray, F. M. Giorgi, and A. Lachmann (2020). Correcting under-reported COVID-19 case numbers: estimating the true scale of the pandemic. medRxiv preprint doi: 10.1101/2020.03.14.20036178. • Kalish, H., C. Klumpp-Thomas, S. Hunsberger, H. A. Baus, M. P. Fay, N. Siripong, ... K. Sadtler (2021). Undiagnosed SARS-CoV-2 seropositivity during the first six months of the COVID-19 pandemic in the United States. \textit{Science Translational Medicine} \textit{13}(601), 1--11. • Kermack, W. and A. McKendrick (1927). A contribution to the mathematical theory of epidemics. \textit{Proceedings of the Royal Society of London. Series A. 115}(772), 700--721. • Li, R., S. Pei, B. Chen, Y. Song, T. Zhang, W. Yang, and J. Shaman (2020). Substantial undocumented infection facilitates the rapid dissemination of novel coronavirus (SARS-CoV-2). \textit{Science 368}(6490), 489--493. • Mossong, J., N. Hens, M. Jit, P. Beutels, K. Auranen, R. Mikolajczyk, ... W. J. Edmunds (2008). Social contacts and mixing patterns relevant to the spread of infectious diseases. \textit{PLoS Med 5}(3), e74. • Nepomuceno, E., D. F. Resende, and M. J. Lacerda (2018). A survey of the individual-based model applied in biomedical and epidemiology. \textit{Journal of Biomedical Research and Reviews,} \textit{1}(1): 11--24. • Oliver, S., J. Gargano, M. Marin, M. Wallace, K. G. Curran, M. Chamberland, ... K. Dooling (2020). The advisory committee on immunization practices' interim recommendation for use of Pfizer-BioNTech COVID-19 vaccine --- United States, December 2020. \textit{MMWR. Morbidity and Mortality Weekly Report 69}(50), 1922--1924. • Oliver, S., J. Gargano, M. Marin, M. Wallace, K. G. Curran, M. Chamberland, ... K. Dooling (2021). The advisory committee on immunization practices' interim recommendation for use of Moderna COVID-19 vaccine --- United States, December 2020. \textit{MMWR. Morbidity and Mortality Weekly Report 69}(5152), 1653--1656. • Oliver, S. E., J. W. Gargano, H. Scobie, M. Wallace, S. C. Hadler, J. Leung, ... K. Dooling (2021). The advisory committee on immunization practices' interim recommendation for use of Janssen COVID-19 vaccine --- United States, February 2021. \textit{MMWR. Morbidity and Mortality Weekly Report 70}(9), 329--332. • Rahmandad, H., T. Y. Lim, and J. Sterman (2021). Behavioral dynamics of COVID-19: estimating underreporting, multiple waves, and adherence fatigue across 92 nations. \textit{System Dynamics Review 37}(1), 5--31. • Rocha, L. E. and N. Masuda (2016). Individual-based approach to epidemic processes on arbitrary dynamic contact networks. \textit{Scientific Reports 6}, 31456. • Thieme, H. R. (2013). \textit{Mathematics in population biology}, Volume 12 of \textit{Princeton Series in Theoretical and Computational Biology}. Princeton University Press. • Willem, L., F. Verelst, J. Bilcke, N. Hens, and P. Beutels (2017). Lessons from a decade of individual-based models for infectious disease transmission: A systematic review (2006--2015). \textit{BMC infectious diseases} \textit{17}(1), 612. • Willem, L., T. Van Hoang, S. Funk, P. Coletti, P. Beutels, and N. Hens (2020). SOCRATES: An online tool leveraging a social contact data sharing initiative to assess mitigation strategies for COVID-19. \textit{BMC Research Notes 13}(1), 1--8. • Zhang, J., M. Litvinova, Y. Liang, Y. Wang, W. Wang, S. Zhao, ... H. Yu (2020). Supplementary materials for "Changes in contact patterns shape the dynamics of the COVID-19 outbreak in China". Available at: science.sciencemag.org/content/368/6498/1481/suppl/DC1. \end{list} \end{footnotesize}