Central limit theorems for SIR epidemics and percolation on configuration model random graphs

TL;DR

This paper establishes central limit theorems for the final size of SIR epidemics on configuration model graphs, providing explicit variance formulas, and connects these results to percolation theory on such graphs.

Contribution

It introduces new CLTs for epidemic final sizes on configuration models, including explicit variance expressions, covering various initial infective scenarios and linking to percolation properties.

Findings

CLTs for epidemic final size with explicit variance formulas

Results apply to both Molloy-Reed and Newman-Strogatz-Watts graphs

Connections established between epidemic sizes and percolation cluster sizes

Abstract

We consider a stochastic SIR (susceptible $\to$ infective $\to$ recovered) epidemic defined on a configuration model random graph, in which infective individuals can infect only their neighbours in the graph during an infectious period which has an arbitrary but specified distribution. Central limit theorems for the final size (number of initial susceptibles that become infected) of such an epidemic as the population size $n$ tends to infinity, with explicit, easy to compute expressions for the asymptotic variance, are proved assuming that the degrees are bounded. The results are obtained for both the Molloy-Reed random graph, in which the degrees of individuals are deterministic, and the Newman-Strogatz-Watts random graph, in which the degrees are independent and identically distributed. The central limit theorems cover the cases when the number of initial infectives either (a) tends to infinity or (b) is held fixed as $n \to \infty$. In (a) it is assumed that the fraction of the population that is initially infected converges to a limit (which may be $0$) as $n \to \infty$, while in (b) the central limit theorems are conditional upon the occurrence of a large outbreak (more precisely one of size at least $\log n$). Central limit theorems for the size of the largest cluster in bond percolation on Molloy-Reed and Newman-Strogatz-Watts random graphs follow immediately from our results, as do central limit theorems for the size of the giant component of those graphs. Corresponding central limit theorems for site percolation on those graphs are also proved.

Tables

Table 1. Table 1: Simulation results against theoretical (asymptotic) calculations for final size of epidemics with ϵ = 0.05 italic-ϵ 0.05 \epsilon=0.05 and p I = 0.3 subscript 𝑝 𝐼 0.3 p_{I}=0.3 on NSW networks.

Equations

$$\mathrm{P}(X=k)=\sum _{d=k}^{\infty }\mathrm{P}(D=d)\mathrm{E}\left[\left(\genfrac{}{}{0ex}{}{d}{k}\right)(1-{\mathrm{e}}^{-\lambda I}{)}^k{\mathrm{e}}^{-(d-k)\lambda I}\right](k=0,1,\dots ).$$

$$\underset{n\to \infty }{\lim }\sqrt{n}\left(n^{-1}{v}_i^{(n)}-p_i\right)=0(i=0,1,\dots ,d_{\max }).$$

$$z-q_I={\mu }_D^{-1}{p}_I{f}_{D_{ϵ}}^{(1)}(z)$$

$$\rho =1-ϵ-f_{D_{ϵ}}(z).$$

$$\sqrt{n}\left(n^{-1}{T}_{\mathrm{MR}}^{(n)}-\rho \right)\stackrel{D}{⟶}\mathrm{N}(0,{\sigma }_{\mathrm{MR}}^2)$$

$$\sqrt{n}\left(n^{-1}{T}_{\mathrm{NSW}}^{(n)}-\rho \right)\stackrel{D}{⟶}\mathrm{N}(0,{\sigma }_{\mathrm{NSW}}^2),$$

$${\sigma }_{\mathrm{MR}}^2$$

$$+h(z)\left[2p_{I}z{f}_{D_{ϵ}}^{(1)}(z^2)-(z-q_I){\mu }_D\right]+1-ϵ-\rho -f_{D_{ϵ}}(z^2)$$

$$+\left(q_I^{(2)}-q_I^2\right)h(z{)}^2\left[f_D^{(2)}(1)-f_{D_{ϵ}}^{(2)}(z)\right],$$

$$h(z)=\frac{p_I^{-1}(q_I-z)}{1-p_I{\mu }_D^{-1}{f}_{D_{ϵ}}^{(2)}(z)},$$

$${\sigma }_{\mathrm{NSW}}^2$$

$$+h(z{)}^2\left\{\left[p_I{q}_I-2z(z-q_I)\right]{\mu }_D+(z-q_I{)}^2\left({\sigma }_D^2+\frac{1-2ϵ}{1-ϵ}{\mu }_D^2\right)\right\}$$

$$+\left(q_I^{(2)}-q_I^2\right)h(z{)}^2\left[f_D^{(2)}(1)-(1-ϵ)f_D^{(2)}(z)\right],$$

$$h(z)=\frac{p_I^{-1}(q_I-z)}{1-p_I{\mu }_D^{-1}(1-ϵ)f_D^{(2)}(z)}.$$

$$R_0={\mu }_{\tilde{D}-1}{p}_I=\left({\mu }_D+{\mu }_D^{-1}{\sigma }_D^2-1\right)p_I,$$

$$z-q_I={\mu }_D^{-1}{p}_I{f}_D^{(1)}(z)$$

$$\sqrt{n}\left(n^{-1}{T}_{\mathrm{MR}}^{(n)}-\rho \right)|G^{(n)}\stackrel{D}{⟶}\mathrm{N}(0,{\tilde{\sigma }}_{\mathrm{MR}}^2)$$

$$\sqrt{n}\left(n^{-1}{T}_{\mathrm{NSW}}^{(n)}-\rho \right)|G^{(n)}\stackrel{D}{⟶}\mathrm{N}(0,{\tilde{\sigma }}_{\mathrm{NSW}}^2),$$

$$\sqrt{n}\left(n^{-1}{T}_{\mathrm{MR}}^{(n)}-\rho \right)\stackrel{D}{⟶}\mathrm{N}(0,{\tilde{\sigma }}_{\mathrm{MR}}^2)$$

$$\sqrt{n}\left(n^{-1}{T}_{\mathrm{NSW}}^{(n)}-\rho \right)\stackrel{D}{⟶}\mathrm{N}(0,{\tilde{\sigma }}_{\mathrm{NSW}}^2),$$

$$(z-1+\pi )={\mu }_D^{-1}\pi f_D^{(1)}(z)$$

$$\sqrt{n}\left(n^{-1}{C}_{\mathrm{MR},\mathrm{B}}^{(n)}-\rho \right)\stackrel{D}{⟶}\mathrm{N}(0,{\sigma }_{\mathrm{MR},\mathrm{B}}^2),$$

$$\sqrt{n}\left(n^{-1}{C}_{\mathrm{NSW},\mathrm{B}}^{(n)}-\rho \right)\stackrel{D}{⟶}\mathrm{N}(0,{\sigma }_{\mathrm{NSW},\mathrm{B}}^2),$$

$$\sqrt{n}\left(n^{-1}{C}_{\mathrm{MR},\mathrm{S}}^{(n)}-\pi \rho \right)\stackrel{D}{⟶}\mathrm{N}(0,{\sigma }_{\mathrm{MR},\mathrm{S}}^2),$$

$$\sqrt{n}\left(n^{-1}{C}_{\mathrm{NSW},\mathrm{S}}^{(n)}-\pi \rho \right)\stackrel{D}{⟶}\mathrm{N}(0,{\sigma }_{\mathrm{NSW},\mathrm{S}}^2),$$

$${\sigma }_{\mathrm{MR},\mathrm{B}}^2$$

$$+h(z)\left[2\pi z{f}_D^{(1)}(z^2)-(z-1+\pi ){\mu }_D\right]+1-\rho -f_D(z^2)$$

$${\sigma }_{\mathrm{NSW},\mathrm{B}}^2$$

$$+h(z{)}^2\left\{\left[\pi (1-\pi )-2z(z-1+\pi )\right]{\mu }_D+(z-1+\pi {)}^2\left({\sigma }_D^2+{\mu }_D^2\right)\right\},$$

$$h(z)=\frac{{\pi }^{-1}(1-\pi -z)}{1-\pi {\mu }_D^{-1}{f}_D^{(2)}(z)};$$

$${\sigma }_{\mathrm{MR},\mathrm{S}}^2$$

$$+\pi (1-\pi )[\rho -2h(z)(1-z){\mu }_D]$$

$$q^{(n)}(\mathbf{i},\mathbf{i}+\mathbf{l})=n{\beta }_{\mathbf{l}}^{(n)}(n^{-1}\mathbf{i})(\mathbf{i}\in S^{(n)},\mathbf{l}\in \Delta ),$$

$$q^{(n)}(\mathbf{i},\mathbf{i}+\mathbf{l})=n\left[{\beta }_{\mathbf{l}}(n^{-1}\mathbf{i})+O(n^{-1})\right](\mathbf{i}\in S^{(n)},\mathbf{l}\in \Delta ).$$

Full text

1 Introduction

There has been considerable work in the past two decades on models for the spread of epidemics on random networks; see, for example, the recent book [Kiss et al.(2017)]. The usual paradigm is that individuals in a population are represented by nodes in a random graph and infected individuals are able to transmit infection only to their neighbours in the graph. The graph is often constructed using the configuration model (see, for example, [van der Hofstad(2016)], Chapter 7), which allows for an arbitrary but specified degree distribution. The most-studied type of epidemic model is the SIR (susceptible $\to$ infective $\to$ recovered) model. In this model individuals are classified into three types: susceptibles, infectives and recovered. If a susceptible individual is contacted by an infective then it too becomes an infective and remains so for a time, called its infectious period, that is distributed according to a non-negative random variable $I$ having an arbitrary but specified distribution. An infective individual recovers at the end of its infectious period and is then immune to further infection. During its infectious period, an infective contacts its susceptible neighbours in the graph independently at the points of Poisson processes having rate $\lambda$. The graph is assumed to be static and the population closed (i.e. there are no births or deaths), so eventually the epidemic process terminates. The final size of the epidemic is the number of initial susceptibles that are infected during its course. The final size is a key epidemic statistic, not only as a measure of the impact of an epidemic but also in an inferential setting, since often it can be observed more reliably than the precise temporal spread. The main aim of this paper is to develop central limit theorems for the final size of an SIR epidemic on configuration model graphs as the population size $n\to\infty$.

The configuration model, which was introduced by [Bollobás(1980)], is a model for random graphs with a given degree sequence. There are two distinct approaches for constructing configuration model graphs with a given degree distribution as $n\to\infty$. In both approaches, individuals are assigned a number of half-edges, corresponding to their degree, and then these half-edges are paired uniformly at random. In [Molloy and Reed(1995)], the degrees of indiviudals are prescribed deterministically whilst in [Newman et al.(2001)] they are i.i.d. (independent and identically distributed) copies of a random variable $D$, that describes the limiting degree distribution. We refer to the former as the MR random graph and to the latter as the NSW random graph. Subject to suitable conditions on the degree sequences in the MR model and $D$ in the NSW model, law of large number limits for SIR epidemics on the two graphs are the same. That is not the case for central limit theorems, as for finite $n$, there is greater variability in the degrees of individuals in the NSW model than in the MR model; indeed, in the NSW model, such variability is of the same order of magnitude as the variability in the epidemic. Thus, though the asymptotic means are the same, the asymptotic variances in the central limit theorems for final size are greater for the epidemic on the NSW random graph.

There have been numerous studies, some fully rigorous and some heuristic, of SIR epidemics on configuration model networks making various assumptions concerning the infectious period random variable $I$. For example, assuming $I$ is constant, [Andersson(1998)] derives a law of large numbers for the final size of an epidemic on an NSW random graph when a strictly positive fraction of the population is initially infected in the limit as $n\to\infty$ and [Britton et al.(2007)] obtain a similar result for epidemics on an MR random graph initiated by a single infective. In the latter case, a large outbreak is possible only if the basic reproduction number $R_{0}>1$; see (2.8) in Section 2.4. In a highly influential paper, [Newman(2002)] uses heuristic percolation arguments to obtain a number of results, including the fraction of the population infected by a large outbreak, for SIR epidemics on NSW random graphs with $I$ having an arbitrary but specified distribution. Several authors have studied the case when $I$ has an exponential distribution, so the model becomes Markovian. [Decreusefond et al.(2012)] obtain a law of large numbers type result for the epidemic process on an NSW random graph with a strictly positive fraction initially infected, which yields a rigorous proof of the deterministic approximation of [Volz(2008)] (see also [Miller(2011)] and [Miller et al.(2012)]). [Bohman and Picollelli(2012)] obtain law of large numbers results for both the process and final size of an epidemic with one intial infective on an MR random graph with bounded degrees. [Janson et al.(2012)] obtain similar results under weaker conditions on the degree sequences considering the cases when the fraction initially infected, in the limit as $n\to\infty$, is either strictly positive or zero (assuming of course there is at least one initial infective). In the latter case, the limiting “deterministic” process involves a random time translation reflecting the time taken for the number of infectives to reach order $n$; a similar result is obtained by [Barbour and Reinert(2013)] assuming a bounded degree sequence and an arbitrary but specified distribution for $I$.

There has been very little work to date on central limit theorems for SIR epidemics on configuration model networks. A functional central limit theorem for the SI epidemic (in which ${\rm P}(I=\infty)=1$, so infectives remain infectious forever) on an MR random graph with unbounded degrees is obtained by [KhudaBukhsh et al.(2017)], who note that their method is not straightforward to extend to an SIR model. Assuming that $I$ follows an exponential distribution and bounded degrees, [Ball et al.(2018)] use an effective degree approach ([Ball and Neal(2008)]) and density dependent population processes ([Ethier and Kurtz(1986)], Chapter 11) to obtain functional central limit theorems for SIR epidemics on MR and NSW random graphs, in which susceptible individuals can also drop their edges to infective neighbours. They also conjecture central limit theorems for the final size of such epidemics (and hence as a special case for the final size of standard SIR epidemics, without dropping of edges), assuming either a strictly positive fraction or a constant number of initial infectives (in which case the central limit theorem is conditional on the occurrence of a large outbreak). However, the arguments are not fully rigorous and the result for a constant number of initial infectives is based purely on the existence of equivalent results for other (non-network) SIR epidemic models. Another limitation of [Ball et al.(2018)] is the assumption that $I$ is exponentially distributed, which is unrealistic for most real-life diseases.

In the present paper, we address these shortcomings and derive fully rigorous central limit theorems for the final size of SIR epidemics on MR and NSW random graphs having bounded degrees, when the infectious period $I$ follows an arbitrary but specified distribution. We consider the cases when the limiting fraction of the population is (i) strictly positive and (ii) zero. For the latter we treat the situtations where the number of initial infectives either (i) is held fixed independent of $n$ or (ii) tends to $\infty$ as $n\to\infty$. The mean parameter $\rho$ in the central limit theorems, which coincides with the corresponding law of large numbers limit, depends on the solution $z$ of a non-linear equation (see (2.2) and (2.9) in Section 2.4). Given $z$, the variance parameter in the central limit theorems is fully explicit and hence easy to compute.

If $I$ is constant, say ${\rm P}(I=1)=1$ then the above SIR model is essentially bond percolation with probability $\pi=1-{\rm e}^{-\lambda}$ and if ${\rm P}(I=\infty)=\pi=1-{\rm P}(I=\infty)=0$ then it is closely related to site percolation with probability $\pi$; see, for example, [Durrett(2007)], page 15, and [Janson(2009a)]. Central limit theorems for the size of the giant component (largest cluster) of bond percolation on the MR and NSW random graphs follow immediately from our results. (Corresponding theorems for site percolation are also obtained using our methodology.) Further, setting $\pi=1$ yields central limit theorems for the giant component of those graphs (Remark 2.6); cf. [Barbour and Röllin(2017)] who obtain a central limit theorem for the giant component of the MR random graph and [Ball and Neal(2017)] who derive the asymptotic variance of the giant component of MR and NSW random graphs, all allowing for unbounded degrees.

The proofs involve constructing the random graph and epidemic on it simultanaeously, modifying the infection mechanism so that when a susceptible is infected it decides which of its half-edges it will try to infect along with (its remaining half-edges becoming recovered half-edges), with the times of those infection attempts (relative to the time of infection of the susceptible) being realisations of i.i.d. exponential random variables. The distribution of the final outcome of the epidemic, and hence also its final size, is invariant to this modification. The process describing the evolution of the numbers of susceptibles of different degrees, infective half-edges and recovered half-edges is an asymptotically density dependent population process ([Ethier and Kurtz(1986)], Chapter 11, and [Pollett(1990)]). The asymptotic distribution of the final outcome of the epidemic is studied by considering a boundary crossing problem for a random time-scale transformation of that process. The proofs extend, at least in principle, to SIR epidemics and percolation on extensions of the configuration model that include fully-connected cliques ([Trapman(2007)], [Gleeson(2009)], [Ball et al.(2010)] and [Coupechoux and Lelarge(2014)]), though explicit calculation of the asymptotic variances may be difficult.

The remainder of the paper is organised as follows. The MR and NSW random graphs are defined in Section 2.1 and the SIR epidemic model is described in Section 2.2. The main central limit theorems (Theorems 2.1-2.3 for SIR epidemics and Theorem 2.7 for percolation) are stated in Section 2.4, together with some remarks giving comparisons of their variance parameters and applications to giant components indicated above. Some numerical illustrations, which show that the central limit theorems can yield good approximations even for relatively small graphs, are given in Section 3. The proofs are given in Section 5. They make extensive use of asymptotically density dependent population processes and in particular require a version of the functional central limit theorem for such processes to include asymptotically random initial conditions. For ease of reference, the required results for such processes are collected together in Section 4. Some brief concluding comments are given in Section 6. Calculation of the asymptotic variances for the central limit theorems is lengthy, though straightforward, so this and a few other details are deferred to an appendix.

1.1 Notation

All vectors are row vectors and ⊤ denotes transpose. With the dimension being obvious from the context, $I$ denotes an identity matrix and $\boldsymbol{0}$ and $\boldsymbol{1}$ denotes vectors all of whose elements are [math] and $1$, respectively. For $x\in\mathbb{R}$, the usual floor and ceiling functions are denoted by $\left\lfloor{x}\right\rfloor$ and $\left\lceil{x}\right\rceil$, respectively. Thus $\left\lfloor{x}\right\rfloor$ is the greatest integer $\leq x$ and $\left\lceil{x}\right\rceil$ is the smallest integer $\geq x$. For a positive integer $k$, the $k$th derivative of a real-valued function $f$ is denoted by $f^{(k)}$. The cardinality of a set $A$ is denoted by $|A|$. Sums are zero if vaccuous. We use $\stackrel{{\scriptstyle p}}{{\longrightarrow}}$, $\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}$ and $\stackrel{{\scriptstyle D}}{{\longrightarrow}}$ to denote convergence in probability, convergence almost sure and convergence in distribution, respectively.

Further, ${\rm U}(0,1)$ denotes a uniform random variable on $(0,1)$; ${\rm Exp}(1)$ denotes an exponential random variable with mean $1$; ${\rm N}(0,\sigma^{2})$ denotes a univariate normal random variable with mean [math] and variance $\sigma^{2}$; and ${\rm N}(\boldsymbol{0},\Sigma)$ denotes a multivariate normal random variable with mean vector $\boldsymbol{0}$ and variance matrix $\Sigma$, whose dimension again is obvious from the context. For a positive integer $n$ and $p\in[0,1]$, ${\rm Bin}(n,p)$ denotes a binomial random variable with $n$ trials and success probability $p$. Also, if $I$ is a non-negative random variable and $\lambda\in(0,\infty)$ then ${\rm Bin}(n,1-{\rm e}^{-\lambda I})$ denotes a mixed-Binomial random variable obtained by first sampling $I_{1}$ from the distribution of $I$ and then, given $I_{1}$, sampling independently from ${\rm Bin}(n,1-{\rm e}^{-\lambda I_{1}})$. Similarly, if $I$ is a non-negative random variable and $D$ is a non-negative integer-valued random variable, then ${\rm Bin}(D,1-{\rm e}^{-\lambda I})$ denotes a mixed-Binomial random variable, where the realisations of $D$ and $I$ are independent. Thus, if $X\sim{\rm Bin}(D,1-{\rm e}^{-\lambda I})$, then

[TABLE]

Note that we allow the possibility $D=0$.

2 Model and main results

2.1 Random graph

Consider a population of $n$ indivdiuals labelled $1,2,\dots,n$. For $i=1,2,\dots,n$, let $D^{(n)}_{i}$ denote the degree of individual $i$. We assume that $0\leq D^{(n)}_{i}\leq d_{\rm max}$ for all $i$, i.e. that there is a maximum degree $d_{\rm max}$. In the MR random graph the degrees are prescribed, while in the NSW random graph $D^{(n)}_{1},D^{(n)}_{2},\dots,D^{(n)}_{n}$ are i.i.d. copies of a random variable $D$ having probability mass function given by ${\rm P}(D=i)=p_{i}$ $(i=0,1,\dots,d_{\rm max})$. Under both models, the network (random graph) is formed by attaching $D^{(n)}_{i}$ half-edges to individual $i$, for $i=1,2,\dots,n$, and then pairing up the $D^{(n)}_{1}+D^{(n)}_{2}+\dots+D^{(n)}_{n}$ half-edges uniformly at random to give the edges in the random graph, which we denote by $\mathcal{G}^{(n)}$. In the NSW model, if $D^{(n)}_{1}+D^{(n)}_{2}+\dots+D^{(n)}_{n}$ is odd there is a left-over stub, which is ignored. (Of course in the MR model the prescribed degrees can be chosen so that $D^{(n)}_{1}+D^{(n)}_{2}+\dots+D^{(n)}_{n}$ is even.)

We are interested in asymptotic results as the number of individuals $n\to\infty$. In the MR random graph, for $i=0,1,\dots,d_{\rm max}$, let $v^{(n)}_{i}=\sum_{k=1}^{n}1_{\{D^{(n)}_{k}=i\}}$ be the number of individuals having degree $i$. We assume that

[TABLE]

Note that (2.1) implies $\lim_{n\to\infty}n^{-1}v^{(n)}_{i}=p_{i}$ $(i=0,1,\dots,d_{\rm max})$.

In both models the random graph may have some imperfections, specifically self-loops and multiple edges, but they are sparse in the network as $n\to\infty$; more precisely, the number of such imperfections converges in distrubution to a Poisson random variable as $n\to\infty$ ([Durrett(2007)], Theorem 3.1.2). Moreover, [Janson(2009b)] implies that the probability that the random graph is simple is bounded away from [math], as $n\to\infty$, and law of large numbers results continue to hold if the graph is conditioned on being simple. However, that is not necessarily the case for convergence in distribution; see [Janson(2010)], Remark 1.4, and [Barbour and Röllin(2017)], Remark 2.5. Thus whether or not our central limit theorems continue to hold when the random graph is conditioned on being simple is an open question.

2.2 SIR epidemic

An SIR epidemic, denoted by $\mathcal{E}^{(n)}$, is constructed on the above network as follows. Initially, at time $t=0$, a number of individuals are infective and the remaining individuals are susceptible. (Precise statements concerning the initial infectives are made later.) Distinct infectives behave independently of each other and of the construction of the network. Each infective remains infectious for a period of time that is distributed according to a random variable $I$, having an arbitrary but specified distribution, after which it becomes recovered. During its infectious period, an infective contacts its neighbours in the network independently at the points of Poisson processes each having rate $\lambda$, so the probability that a given neighbour is contacted is $p_{I}=1-\phi(\lambda)$, where $\phi(\theta)=\rm E[\exp(-\theta I)]$ $(\theta\geq 0)$ is the Laplace transform of $I$. If a contacted individual is susceptible then it becomes an infective, otherwise nothing happens. The epidemic ends when there is no infective individual in the population.

For $t\geq 0$, let $X^{(n)}_{i}(t)$ be the number of degree-$i$ susceptible individuals at time $t$ ($i=0,1,\dots,d_{\rm max}$) and let $Y^{(n)}(t)$ be the total number of infectives at time $t$. Let $\tau^{(n)}=\inf\{t\geq 0:Y^{(n)}(t)=0\}$ be the time of the end of the epidemic. Then $T^{(n)}_{i}=X^{(n)}_{i}(0)-X^{(n)}_{i}(\tau^{(n)})$ is the total number of degree-$i$ susceptibles that are infected by the epidemic. Let $T^{(n)}=\sum_{i=0}^{d_{\rm max}}T^{(n)}_{i}$ be the total number of susceptibles infected by the epidemic, i.e. the final size of the epdiemic. We are primarily interested in the asymptotic distribution of $T^{(n)}$ as $n\to\infty$.

The proofs allow for the possibility that $p_{I}=1$, i.e. ${\rm P}(I=\infty)=1$. In that case the set of individuals that are infected during the epidemic $\mathcal{E}^{(n)}$ comprises all individuals in the components of $\mathcal{G}^{(n)}$ that contain at least one initial infective. Thus central limit theorems for the size of the giant component in MR and NSW random graphs follow immediately from our results; see Remark 2.6 below.

2.3 Bond and site percolation

In bond percolation on $\mathcal{G}^{(n)}$, each edge in $\mathcal{G}^{(n)}$ is deleted indpendently with probability $1-\pi$, while in site percolation $\mathcal{G}^{(n)}$, each vertex (together with all incident edges) is deleted independently with probability $1-\pi$ ([Janson(2009a)]). Interest is often focused on the size, $C^{(n)}$ say, of the largest connected component in the resulting graph.

2.4 Main results

For $i=0,1,\dots,d_{\rm max}$, let $a^{(n)}_{i}$ be the number of degree-$i$ initial infectives in the epidemic $\mathcal{E}^{(n)}$ and let $a^{(n)}=\sum_{i=0}^{d_{\rm max}}a^{(n)}_{i}$ denote the total number of initial infectives. In the epidemic on the MR random graph, we assume that $a^{(n)}_{0},a^{(n)}_{1},\dots,a^{(n)}_{d_{\rm max}}$ are prescribed. In the epidemic on the NSW random graph, we assume that $a^{(n)}$ is prescribed and that the $a^{(n)}$ initial infectives are chosen by sampling uniformly at random without replacement from the $n$ individuals in the population.

Let $\epsilon^{(n)}=n^{-1}a^{(n)}$ and $\epsilon^{(n)}_{i}=n^{-1}a^{(n)}_{i}$ $(i=0,1,\dots,d_{\rm max})$. Suppose that $\epsilon^{(n)}\to\epsilon$ as $n\to\infty$ and that $\lim_{n\to\infty}\sqrt{n}(\epsilon^{(n)}-\epsilon)=0$. For the epidemic on the MR random graph, suppose further that, for $i=0,1,\dots,d_{\rm max}$, there exists $\epsilon_{i}$ such that $\lim_{n\to\infty}\sqrt{n}(\epsilon^{(n)}_{i}-\epsilon_{i})=0$. For the epidemic on the NSW random graph, let $\epsilon_{i}=\epsilon p_{i}$ ($i=0,1,\dots,d_{\rm max}$). Let $\mu_{D}={\rm E}[D]$, $\sigma_{D}^{2}=\rm var(D)$ and, for $s\in[0,1]$, $f_{D}(s)=\sum_{i=0}^{d_{\rm max}}p_{i}s^{i}$ and $f_{D_{\epsilon}}(s)=\sum_{i=0}^{d_{\rm max}}(p_{i}-\epsilon_{i})s^{i}$. Let $q_{I}=1-p_{I}=\phi(\lambda)$ be the probability that an infective fails to contact a given neighbour and $q_{I}^{(2)}=\phi(2\lambda)$ be the probability that a given infective fails to contact two given neighbours.

The first theorem concerns the case $\epsilon>0$, so in the limit as $n\to\infty$ a strictly positive fraction of the population is initially infective. Let $T^{(n)}_{\rm MR}$ and $T^{(n)}_{\rm NSW}$ denote the final size of the epidemic $\mathcal{E}^{(n)}$ on the MR and NSW random graphs, respectively. Let $z$ be the unique solution in $[0,1)$ of

[TABLE]

and

[TABLE]

Theorem 2.1

Suppose that $p_{i}\epsilon_{i}>0$ for at least one $i>0$ and, if $p_{I}=1$ then $p_{1}-\epsilon_{1}>0$. Then, as $n\to\infty$,

[TABLE]

and

[TABLE]

where

[TABLE]

with

[TABLE]

and

[TABLE]

with

[TABLE]

The second theorem concerns the case when the number of initial infectives is held fixed as $n\to\infty$, so $\epsilon=0$. More specifically, in the epidemic on the MR random graph, we assume that $a^{(n)}_{i}=a_{i}$ $(i=0,1,\dots,d_{\rm max})$ for all $n\geq a=\sum_{i=1}^{d_{\rm max}}a_{i}$ and in the epidemic on the NSW random graph, we assume that $a^{(n)}=a$ for all $n\geq a$. It is well known that, for large $n$, the process of infectives in the early stages of such an epidemic can be approximated by a Galton-Watson branching process, $\mathcal{B}$ say, in which, except for the initial generation, the offspring distribution is ${\rm Bin}(\tilde{D}-1,1-{\rm e}^{-\lambda I})$, where $\tilde{D}$ and $I$ are independent and $\tilde{D}$ has the size-biased degree distribution ${\rm P}(\tilde{D}=k)=\mu_{D}^{-1}kp_{k}$ $(k=1,2,\dots,d_{\rm max})$; see, for example, [Ball and Sirl(2013)]. This offspring distribution has mean

[TABLE]

where $\mu_{\tilde{D}-1}={\rm E}[\tilde{D}-1]$. The quantity $R_{0}$ is called the basic reproduction number of the epidemic.

For $\mathcal{E}^{(n)}$, we say that a major outbreak occurs if and only if the event $G^{(n)}=\{T^{(n)}\geq\log n\}$ occurs. Now $\lim_{n\to\infty}{\rm P}(G^{(n)})={\rm P}(B=\infty)$, where $B$ is the total progeny (not including the initial generation) of the branching process $\mathcal{B}$ (cf. Theorem 5.1). Thus, in the limit as $n\to\infty$, a major outbreak occurs with non-zero probability if and only if $R_{0}>1$.

Suppose that that $R_{0}>1$. Now let $z$ be the unique solution in $[0,1)$ of

[TABLE]

and $\rho=1-f_{D}(z)$. Note that if $p_{I}\in(0,1)$, or $p_{I}=1$ and $p_{1}>0$, then $z>0$.

Theorem 2.2

Suppose that $R_{0}>1$ and, if $p_{I}=1$ then $p_{1}>0$. Then, as $n\to\infty$,

[TABLE]

and

[TABLE]

where $\tilde{\sigma}^{2}_{\rm MR}$ is given by (2.4) and (2.5), with $f_{D_{\epsilon}}$ replaced by $f_{D}$, and $\tilde{\sigma}^{2}_{\rm NSW}$ is obtained by setting $\epsilon=0$ in (2.6) and (2.7).

The next theorem concerns the case when $\epsilon=0$ but the number of initial infectives $a^{(n)}\to\infty$ as $n\to\infty$.

Theorem 2.3

Suppose that $R_{0}>1$, $\epsilon=0$, $\sum_{i=1}^{d_{\rm max}}a^{(n)}_{i}\to\infty$ as $n\to\infty$ and, if $p_{I}=1$ then $p_{1}>0$. Then, as $n\to\infty$,

[TABLE]

and

[TABLE]

where $\tilde{\sigma}^{2}_{\rm MR}$ and $\tilde{\sigma}^{2}_{\rm NSW}$ are as in Theorem 2.2.

Remark 2.4

Note that $q_{I}^{(2)}=q_{I}^{2}$ if $I$ is almost surely constant, otherwise $q_{I}^{(2)}>q_{I}^{2}$ by Jensen’s inequality. Also $f_{D}^{(2)}(1)-f_{D_{\epsilon}}^{(2)}(z)>0$, so as one would expect on intuitive grounds, if $p_{I}$ is held fixed, the asymptotic variance $\sigma^{2}_{\rm MR}$ is smallest when the infectious period is constant. A similar comment holds for $\sigma^{2}_{\rm NSW},\tilde{\sigma}^{2}_{\rm MR}$ and $\tilde{\sigma}^{2}_{\rm NSW}$.

Remark 2.5

Although it is not transparent from (2.4) and (2.6), it is seen easily from the proof that, again as one would expect on intuitive grounds, $\sigma^{2}_{\rm NSW}\geq\sigma^{2}_{\rm MR}$ and $\tilde{\sigma}^{2}_{\rm NSW}\geq\tilde{\sigma}^{2}_{\rm MR}$, with strict inequalities unless the support of the degree random variable $D$ is concentrated on a single point (when the two models are identical); see Appendix A.3.4.

Remark 2.6

In the setting of Theorem 2.2, if $a^{(n)}=1$ and $p_{I}=1$ then with probability tending to $1$ as $n\to\infty$, the event $G^{(n)}$ occurs if and only if the initial infective belongs to the giant component of $\mathcal{G}^{(n)}$. Thus setting $p_{I}=1$ in Theorem 2.2 yields central limit theorems for giant components of MR and NSW random graphs.

The final theorem is concerned with percolation. Let $R_{0}$ be given by (2.8) with $p_{I}=\pi$. Let $C^{(n)}_{\rm MR,B}$ and $C^{(n)}_{\rm MR,S}$ denote the size of the largest connected component after bond and site percolation, respectively, on an MR graph; define $C^{(n)}_{\rm NSW,B}$ and $C^{(n)}_{\rm NSW,S}$ analogously for the NSW graph. Then, if $R_{0}>1$, for each of these four choices for $C^{(n)}$, there exists $\epsilon>0$ such that $\lim_{n\to\infty}{\rm P}(C^{(n)}\geq\epsilon n)=1$; see [Janson(2009a)], Theorems 3.5 and 3.9, which also give law of large number limits for the $C^{(n)}$ under weaker conditions than here. The following theorem gives associated central limit theorems.

Theorem 2.7

Suppose that $\pi\in(0,1)$ and $R_{0}>1$. Let $z$ be the unique solutiuon in $(0,1)$ of

[TABLE]

and $\rho=1-f_{D}(z)$. Then, as $n\to\infty$,

[TABLE]

[TABLE]

[TABLE]

and

[TABLE]

where

[TABLE]

and

[TABLE]

with

[TABLE]

[TABLE]

*and $\sigma^{2}_{\rm NSW,S}$ is given by (2.17), with $\sigma^{2}_{\rm MR,B}$ replaced by $\sigma^{2}_{\rm NSW,B}$. *

3 Illustrations

In this section, we illustrate the central limit theorems in Section 2.4 by using simulations to explore their applicability to graphs with finite $n$. We also investigate briefly the dependence of the limiting variances in the central limit theorems on the degree distribution, graph type and infectious period distribution. We consider four degree distributions:

1. (i)

$D\equiv d$, i.e. $D$ is constant with $p_{d}=1$; 2. (ii)

$D\sim{\rm Po}(\mu_{D})$, i.e. $D$ is Poisson with mean $\mu_{D}$; 3. (iii)

$D\sim{\rm Geom}(p)$, i.e. $p_{k}=(1-p)^{k}p$ $(k=0,1,\dots)$; 4. (iv)

$D\sim{\rm Power}(\alpha,\kappa)$, i.e. $p_{k}=ck^{-\alpha}{\rm e}^{-\frac{k}{\kappa}}$ $(k=1,2,\dots)$, where $\alpha,\kappa\in(0,\infty)$ and the normalising constant $c={\rm Li}_{\alpha}({\rm e}^{-\frac{1}{\kappa}})$, with ${\rm Li}_{\alpha}$ being the polylogarithm function.

The fourth distribution is a power law with exponential cut-off (see, for example, [Newman(2002)]) that has been used extensively in the physics literature. Note that, with $\theta={\rm e}^{-\frac{1}{\kappa}}$ and $\beta={\rm Li}_{\alpha}(\theta)$, $f_{D}(s)=\beta^{-1}{\rm Li}_{\alpha}(\theta s)$, $f_{D}^{(1)}(s)=\frac{1}{\beta s}{\rm Li}_{\alpha-1}(\theta s)$ and $f_{D}^{(2)}(s)=\frac{1}{\beta s^{2}}[{\rm Li}_{\alpha-2}(\theta s)-{\rm Li}_{\alpha-1}(\theta s)]$, which enables $R_{0}$ and the asymptotic means and variances in the central limit theorems to be calculated. Distributions (ii)-(iv) have unbounded support, so do not satisfy the conditions of our theorems. The asymptotic ditributions in this section are calculated under the assumption that the theorems still apply.

Each simulation consists of first simulating one graph, by simulating $D^{(n)}_{1},D^{(n)}_{2},\dots,D^{(n)}_{n}$ and pairing up the half-edges, and then simulating one epidemic, or percolation process, on it. For an MR graph with (limiting) degree random variable $D$ on $n$ nodes, the degrees are given by $D^{(n)}_{i}=\inf\{d\in\mathbb{Z}_{+}:F_{D}(d)>i/(n+1)\}$ $(i=1,2,\dots,n)$, where $F_{D}$ is the distribution function of $D$. For heavy-tailed $D$, in particular, this choice of MR degrees converges faster with $n$ to the intended $D$ than one based on rounding $np_{i}$ $(i=0,1,\dots)$ to nearest integers. Two choices of infectious period distribution are used in the simulations: (i) $I\equiv 1$ (i.e. $\rm P(I=1)=1$) and (ii) $I$ is [math] or $\infty$, matched to have a common $p_{I}$. Note that these yield the minimum and maximum asymptotic variances for fixed $p_{I}$. The parameters of the degree distributions are chosen so that $\mu_{D}=5$.

We first consider epidemics on an NSW network in which a fraction $\epsilon=0.05$ of the population is initially infective. Table 1 shows estimates of $\rho_{n}=n^{-1}\rm E[T^{(n)}_{\rm NSW}]$ and $\sigma_{n}^{2}=n^{-1}\rm var[T^{(n)}_{\rm NSW}]$ for epidemics with $p_{I}=0.3$ and various population size $n$, based on 100,000 simulations for each set of parameters, together with the corresponding asymptotic ($n=\infty$) values given by Theorem 2.1. The table indicates that the asymptotic approximations are useful for even moderate $n$. The approximations are better for the model with $I\equiv 1$ than that with $I=0$ or $\infty$, as one might expect since there is less variability in the process with $I\equiv 1$, and improve with increasing $\rm var(D)$. The approximations are noticeably worse when $n=200$ than for the other values of $n$. Histograms of the final size of 100,000 simulated epidemics on an NSW network with $n=200$ and $D\sim{\rm Geom}(1/6)$, together with the corresponding density of ${\rm N}(n\rho,n\sigma^{2}_{\rm NSW})$ with $\rho$ and $\sigma^{2}_{\rm NSW}$ given by Theorem 2.1, are shown in Figure 1. For the epidemic with $I\equiv 1$, the asymptotic normal distrbution gives an excellent approximation, even though $n$ is rather small. The approximation is markedly worse for model with $I=0$ or $\infty$, owing to the increased likelihood of small outbreaks. If $\epsilon$ is held fixed, the probability of a small outbreak decreases approximately exponentially with $n$, as the number of initial infectives is proportional to $n$, and the approximation improves significantly, particularly in the $I=0$ or $\infty$ model.

Turning to Theorem 2.2, Figure 2 shows simulations of the final size of epidemics on an NSW network with one initial infective and $I$ constant. Note that with a single initial infective there is always a non-negligible probability of a minor outbreak, even if $n$ is large. As $n\to\infty$, the probability of a major epidemic converges to the survival probability, $p_{\rm maj}$ say, of a branching process; see, for example, [Ball and Sirl(2013)], Section 2.1.1. Theorem 5.1. Moreover, as $I$ is constant, $p_{\rm maj}=\rho$; see, for example, [Kenah and Robins(2007)]. Superimposed on each histogram in Figure 2 is the density of ${\rm N}(n\rho,n\tilde{\sigma}^{2}_{\rm NSW})$ multiplied by $p_{\rm maj}$, which approximates the component of the distribution of $T^{(n)}_{\rm NSW}$ corresponding to a major outbreak. Even with the very small $n$, the approximation is very good when $D\sim{\rm Geom}(1/6)$ or $D\sim{\rm Power}(1,13.796)$. For both of these degree distributions there is a clear distiction between major and minor outbreaks. That is not the case for the other two degree distributions, though the approximation is still useful.

The upper panels of Figure 3 shows the dependence of $\rho$ (left panel) and $\tilde{\sigma}^{2}_{\rm MR}$ and $\tilde{\sigma}^{2}_{\rm NSW}$ (right panel) on $p_{I}$. The latter two are for the model with $I$ constant. (Recall that, given $p_{I}$, the scaled asymptotic mean $\rho$ is independent of the distribution of $I$.) The asymptotic scaled variances all decrease with $p_{I}$ and converge to the asymptotic scaled variance of the giant component on the relevant graph as $p_{I}\uparrow 1$ (cf. Remark 2.6). The asymptotic scaled variances tend to $\infty$ as $p_{I}\downarrow p_{C}$, where $p_{C}=\mu_{\tilde{D}-1}^{-1}$ is the critical value of $p_{I}$ so that $R_{0}=1$. Note that $\tilde{\sigma}^{2}_{\rm NSW}\geq\tilde{\sigma}^{2}_{\rm MR}$ for all choices of $D$, cf. Remark 2.5. The lower panels of Figure 3 show plots of $(\tilde{\sigma}^{2}_{\rm NSW}-\tilde{\sigma}^{2}_{\rm MR})/\tilde{\sigma}^{2}_{\rm MR}$ against $p_{I}$. Note the plots for the Poisson and geometric degree distributions are both increasing with $p_{I}$, while that for the ${\rm Power}(1,13.796)$ distribution is unimodal.

The final illustrations are concerned with percolation. Figure 4 shows plots of estimates $\hat{\sigma}_{n}^{2}$ of the scaled variance $n^{-1}\rm var(C^{(n)})$ of the size of the largest component, based on $n_{\rm sim}=10,000$ simulations for each choice of parameters, together with $95$% equal-tailed confidence intervals given by $\left[(n_{\rm sim}-1)/q_{2},(n_{\rm sim}-1)/q_{1}\right]$, where $q_{1}$ and $q_{2}$ are respectively the 2.5% and 97.5% quantiles of the $\chi^{2}_{n_{{\rm sim}-1}}$ distribution. In all cases $\pi=0.3$. The filled and unfilled markers correspond to percolation on NSW and MR networks, respectively. The $n\to\infty$ scaled variances, given by Theorem 2.7, are shown by horizontal dashed and solid lines for NSW and MR networks, respectively. The figure suggests that the asymptotic approximations are again generally good, even for moderately-sized networks. For fixed $n$, the approximation is better for ${\rm Power}(1,13.796)$ distribution than for the ${\rm Po}(5)$ distribution. The plot for site percolation when $D\sim{\rm Po}(5)$ is a bit odd, as $\hat{\sigma}_{n}^{2}$ is not monotonic in $n$. This is explored further in Figure 5, which is for percolation on NSW random graphs with $D\sim{\rm Po}(5)$. Note that the distribution of $C^{(n)}$ is clearly bimodal for site percolation with $n=200$ and $\pi=0.3$. The lower plots in Figure 4 demonstrate that increasing $n$ or $\pi$ alleviates the issue of small largest components.

Finally, Figure 6 shows histograms of the size of the largest component in 100,000 simulated bond and site percolations on an MR random graph with $n=200,\pi=0.3$ and $D\sim{\rm Power}(1,13.796)$. Two asymptotic normal approximations are superimposed. The solid lines are the densities of ${\rm N}(n\rho,n\sigma^{2}_{\rm MR,B})$ (bond percolation) and ${\rm N}(n\pi\rho,n\sigma^{2}_{\rm MR,S})$ (site percolation), with $\rho,\sigma^{2}_{\rm MR,B}$ and $\sigma^{2}_{\rm MR,S}$ obtained by setting $D\sim{\rm Power}(1,13.796)$ in Theorem 2.7. An improved approximation (dashed lines) is obtained by instead setting $D$ to be the empirical ditribution of $D^{(n)}_{1},D^{(n)}_{2},\dots,D^{(n)}_{n}$. The difference between the approximations is more noticeable for bond percolation. (The support of the histogram has been truncated slightly to make the difference clearer.) The difference decreases with $n$ and is appreciably greater with heavy-tailed degree distributions. Observe from Figures 5 and 6 that the asymptotic normal approximation underestimates the left tail and overestimated the right tail of the distribution of $C^{(n)}$. This phenomenon is present also in the asymptotic normal approximation of the epidemic final size $T^{(n)}$.

4 Density dependent population processes

This section collects together some results for density dependent population processes that are required in the paper. It is based on [Ethier and Kurtz(1986)], Chapter 11, and [Pollett(1990)], though the statement of the functional central limit theorem is slightly more general than that in those references. The notation is local to this section.

For $n=1,2,\dots$, let $\{\boldsymbol{X}^{(n)}(t):t\geq 0\}=\{(X^{(n)}_{1}(t),X^{(n)}_{2}(t),\dots,X^{(n)}_{p}(t)):t\geq 0\}$ be a continuous-time Markov chain with state space $S^{(n)}\subseteq\mathbb{Z}^{p}$ and transition intensities of the form

[TABLE]

where $\Delta$ is the set of possible jumps from a typical state $\boldsymbol{i}=(i_{1},i_{2},\dots,i_{p})$ and the $\beta^{(n)}_{\boldsymbol{l}}:E\to\mathbb{R}$ are nonnegative functions defined on an open set $E\subseteq\mathbb{R}^{p}$. We assume that $\Delta$ is finite. The theory in [Ethier and Kurtz(1986)], Chapter 11, and [Pollett(1990)] allows $\Delta$ to be infinite but only the finite case is required in our application and the results are easier to state in that setting. Suppose that $\beta_{\boldsymbol{l}}(\boldsymbol{x})=\lim_{n\to\infty}\beta^{(n)}_{\boldsymbol{l}}(\boldsymbol{x})$ exists for all $\boldsymbol{l}\in\Delta$ and all $\boldsymbol{x}\in E$; the corresponding family of processes is then called asymptotically density dependent by [Pollett(1990)]. In [Ethier and Kurtz(1986)], Chapter 11, a family of processes which satisfies (4.1) with $\beta^{(n)}_{\boldsymbol{l}}$ replaced by $\beta_{\boldsymbol{l}}$ is called a density dependent family, and it is noted that the results usually carry over with little additional effort to the more general form where

[TABLE]

Let

[TABLE]

The following weak law of large numbers follows from [Ethier and Kurtz(1986)], Theorem 11.2.1, allowing for random initial conditions and asymptotic density dependence; see also [Kurtz(1970)], Theorem 3.1 and [Pollett(1990)], Theorem 3.1.

Theorem 4.1

Suppose that for each compact $K\in H$,

[TABLE]

and there exists $M_{K}>0$ such that

[TABLE]

Suppose also that $n^{-1}\boldsymbol{X}^{(n)}(0)\stackrel{{\scriptstyle p}}{{\longrightarrow}}\boldsymbol{x}_{0}$ as $n\to\infty$. Let $\boldsymbol{x}(t)$ $(t\geq 0)$ be given by

[TABLE]

Then, for all $t\geq 0$,

[TABLE]

Write $F(\boldsymbol{x})=\left(F_{1}(\boldsymbol{x}),F_{2}(\boldsymbol{x}),\dots,F_{p}(\boldsymbol{x})\right)$ and $\boldsymbol{x}=(x_{1},x_{2},\dots,x_{p})$. Let $\partial F(\boldsymbol{x})$ and $G(\boldsymbol{x})$ be the $p\times p$ matrix functions defined by

[TABLE]

Further, let $\Phi(t,u)=[\phi_{ij}(t,u)]$ $(0\leq u\leq t<\infty)$ be the solution of the matrix differential equation

[TABLE]

Note that

[TABLE]

The following functional central limit theorem follows from [Ethier and Kurtz(1986)], Theorem 11.2.3, (see also [Kurtz(1970)], Theorem  3.5, and [Pollett(1990)], Theorem 3.2) allowing for asymptotically random initial conditions. Let $\Rightarrow$ denote weak convergence in the space of right-continuous functions $f:[0,\infty)\to\mathbb{R}^{p}$ having limits from the left (i.e. càdlàg functions), endowed with the Skorohod metric.

Theorem 4.2

Suppose that for each compact $K\in H$,

[TABLE]

and that $\beta_{\boldsymbol{l}}$ $(\boldsymbol{l}\in\Delta)$ and $\partial F$ are continuous. Suppose also that

[TABLE]

where $\boldsymbol{V}(0)\sim N(\boldsymbol{0},\Sigma_{0})$. Then, as $n\to\infty$,

[TABLE]

where $\{\boldsymbol{V}(t):t\geq 0\}$ is a zero-mean Gaussian process with covariance function given, for $t_{1},t_{2}\geq 0$, by

[TABLE]

We also require the following theorem concerned with the first exit of $\{\boldsymbol{X}^{(n)}(t):t\geq 0\}$ from a region of its state space. The theorem is a subset of [Ethier and Kurtz(1986)], Theorem 11.4.1, expressed now for asymptotically density dependent processes having random initial conditions.

Theorem 4.3

Suppose that the conditions of Theorem 4.2 are satisfied. Let $\varphi:\mathbb{R}^{p}\to\mathbb{R}$ be continuously differentiable. Suppose that $\varphi(\boldsymbol{x}(0))>0$. Let $\tau^{(n)}=\inf\{t\geq 0:\varphi(n^{-1}\boldsymbol{X}^{(n)}(t))\leq 0\}$ and $\tau=\inf\{t\geq 0:\varphi(\boldsymbol{x}(t))\leq 0\}$. Suppose that $\tau<\infty$ and $\nabla\varphi(\boldsymbol{x}(\tau))\cdot F(\boldsymbol{x}(\tau))<0$, where $\cdot$ denotes inner vector product. Then, as $n\to\infty$,

[TABLE]

Corollary 4.4

Suppose that the conditions of Theorem 4.3 are satisfied. Let

[TABLE]

where $\bigotimes$ denotes outer vector product. Then,

[TABLE]

where

[TABLE]

• Proof.

Corollary 4.4 follows immediately from Theorems 4.2 and 4.3 on noting that

[TABLE]

$\Box$

5 Proofs

5.1 Alternative construction of final size $T^{(n)}$

We describe first the well-known construction of the final outcome of the epidemic $\mathcal{E}^{(n)}$ using a random directed graph. We then use that construction to show that $T^{(n)}$ can be realised using the location of the first exit of an asymptotically density dependent population process from a given region.

Given a realisation of $\mathcal{G}^{(n)}$, construct a directed random graph $\tilde{\mathcal{G}}^{(n)}$, having vertex set $\mathcal{N}^{(n)}=\{1,2,\dots,n\}$, as follows. For each $i=1,2,\dots,n$, by sampling from its infectious period distribution $I$ and then the relevant Poisson processes, draw up a list of individuals $i$ would make contact with if $i$ were to become infected. Then, for each ordered pair of individuals, $(i,j)$ say, a directed edge from $i$ to $j$ is present in $\tilde{\mathcal{G}}^{(n)}$ if and only if $j$ is in $i$’s list. Let $\mathcal{I}^{(n)}$ denote the set of initial infectives in $\mathcal{E}^{(n)}$. For distinct $i,j\in\mathcal{N}^{(n)}$, write $i\leadsto j$ if and only if there is a chain of directed edges from $i$ to $j$ in $\tilde{\mathcal{G}}^{(n)}$. Let $\mathcal{T}^{(n)}$ be the set of initial susceptibles that are infected in $\mathcal{E}^{(n)}$. Then $\mathcal{T}^{(n)}=\{j\in\mathcal{N}^{(n)}\setminus\mathcal{I}^{(n)}:i\leadsto j\mbox{ for some i }\in\mathcal{I}^{(n)}\}$ and the final size $T^{(n)}$ of $\mathcal{E}^{(n)}$ is given by the cardinality of $\mathcal{T}^{(n)}$.

Note that $\mathcal{T}^{(n)}$ is determined purely by the presence/absence of directed edges in $\tilde{\mathcal{G}}^{(n)}$ and does not depend on the times of the corresponding infections in $\mathcal{E}^{(n)}$. (This implies that the distribution of $\mathcal{T}^{(n)}$, and hence also $T^{(n)}$ is invariant to the introduction of a latent/exposed period into the model, i.e. the time elapsing after infection of an individual before it is able to infect other individuals.) It follows that the final outcome $\mathcal{T}^{(n)}$ has the same distribution as that of a related epidemic $\tilde{\mathcal{E}}^{(n)}$, with set of initial infectives $\mathcal{I}^{(n)}$, in which for any infective, $i$ say, it is determined upon infection which, if any, of its neighbours $i$ will contact and those contacts take place at the first points of independent Poisson processes, each having rate $1$. More precisely, suppose $i$ is infected at time $t_{0}$ and $i$ has $d$ neighbours, $i_{1},i_{2},\dots,i_{d}$ say. Let $I_{i}$ be a realisation of $I$ and, given $I_{i}$, let $\chi_{i1},\chi_{i2},\dots,\chi_{id}$ be i.i.d. Bernoulli random variables with success probability $1-\exp(-\lambda I_{i})$. Let $W_{i1},W_{i2},\dots,W_{id}$ be an independent set of i.i.d. unit-mean exponential random random variables. Then $i$ contacts $i_{j}$ if and only if $\chi_{ij}=1$, and in that case the contact occurs at time $t_{0}+W_{ij}$. Of course, the $I$ and $W$ random variables for any set of distinct infectives are mutually independent.

Given the degrees $D^{(n)}_{1},D^{(n)}_{2},\dots,D^{(n)}_{n}$ and the set of initial infectives $\mathcal{I}^{(n)}$, the random graph $\mathcal{G}^{(n)}$ and the epidemic on it $\tilde{\mathcal{E}}^{(n)}$ can be constructed simultaneously as follows (cf. [Ball and Neal(2008)]). The process starts with no half-edge paired. The individuals in $\mathcal{I}^{(n)}$ become infected at time $t=0$ and all other individuals are susceptible. When an individual is infected, it determines immediately which, if any, of its half-edges it will transmit infection along and when it will make those contacts, according to the probabilistic law described above; the half-edges that the indivdiual will infect along and not infect along are then called infective and recovered half-edges, respectively. When infection is transmitted along a half-edge that half-edge is paired with a half-edge chosen uniformly at random from all unpaired half-edges, forming an edge in the network. If the chosen half-edge is attached to a susceptible individual then that individual is infected and determines immediately along which, if any, of its remaining half-edges it will transmit infection. If the chosen half-edge is infective or recovered then nothing happens, apart from the two half-edges being replaced by an edge. The process continues until there is no infective half-edge remaining. (In the epidemic on the NSW random graph, if $D^{(n)}_{1}+D^{(n)}_{2}+\dots+D^{(n)}_{n}$ is odd then it is possible for the process to end with one unpaired half-edge, which is infective, but we can ignore that possibility because under the conditions of the theorems its probability tends to [math] as $n\to\infty$.)

Note that as we are interested only in the final outcome of the epidemic, it is not necessary to keep track of the degrees of individuals to which infective and recovered half-edges are attached; we need to know just the total numbers of such half-edges. For $t\geq 0$, let $X^{(n)}_{i}(t)$ be the number of degree-$i$ susceptible individuals at time $t$ $(i=0,1,\dots,d_{\rm max})$, and let $Y^{(n)}_{E}(t)$ and $Z^{(n)}_{E}(t)$ be the total numbers of infective and recovered half-edges, respectively, at time $t$. Let $\{\boldsymbol{W}^{(n)}(t)\}=\{\boldsymbol{W}^{(n)}(t):t\geq 0\}$, where $\boldsymbol{W}^{(n)}(t)=(\boldsymbol{X}^{(n)}(t),Y^{(n)}_{E}(t),Z^{(n)}_{E}(t))$ and $\boldsymbol{X}^{(n)}(t)=(X^{(n)}_{0}(t),X^{(n)}_{1}(t),\dots,X^{(n)}_{d_{\rm max}}(t))$.

The process $\{\boldsymbol{W}^{(n)}(t)\}$ is a continuous-time Markov chain, whose initial state $\boldsymbol{W}^{(n)}(0)$ is random, even for the epidemic on the MR random graph as the numbers of infective and recovered half-edges created by the initial infectives are random. In the epidemic on the NSW random graph, $X^{(n)}_{i}(t)$ $(i=0,1,\dots,d_{\rm max})$ are also random. Before giving the possible transition intensities of $\{\boldsymbol{W}^{(n)}(t)\}$ some more notation is required.

Let $H^{(n)}$ denote the state space of $\{\boldsymbol{W}^{(n)}(t)\}$ and $\boldsymbol{n}=(n^{X}_{0},n^{X}_{1},\dots,n^{X}_{d_{\rm max}},n^{Y}_{E},n^{Z}_{E})$ denote a typical element of $H^{(n)}$. Let $n^{X}_{E}=\sum_{i=1}^{d_{\rm max}}in^{X}_{i}$. Define the unit vectors $\boldsymbol{e}^{\rm S}_{i}$ $(i=0,1,\dots,d_{\rm max})$, $\boldsymbol{e}^{\rm I}$ and $\boldsymbol{e}^{\rm R}$ on $H^{(n)}$, where, for example, $\boldsymbol{e}^{\rm S}_{i}$ has a one in the $i$th suscpetible component and zeros elsewhere. For $i=1,2,\dots,d_{\rm max}$ and $k=0,1,\dots,i-1$, let $p_{i,k}$ be the probability that if a degree-$i$ susceptible is infected it subsequently transmits infection along $k$ of its remaining $i-1$ half-edges. (Note that a degree-[math] susceptible cannot be infected.) Thus $p_{i,k}={\rm P}(X=k)$, where $X\sim{\rm Bin}(i-1,1-\exp(-\lambda I))$. The transition intensities of $\{\boldsymbol{W}^{(n)}(t)\}$ are as follows.

1. (i)

For $i=1,2,\dots,d_{\rm max}$ and $k=0,1,\dots,i-1$, an infective half-edge is paired with a degree-$i$ susceptible yielding $k$ infective half-edges and $i-1-k$ recovered half-edges

[TABLE] 2. (ii)

an infective half-edge is paired with an infective half-edge

[TABLE] 3. (iii)

an infective half-edge is paired with a recovered half-edge

[TABLE]

Note that these transition intensities are independent of the population size $n$. We index then by $n$ to connect with theory of density dependent population processes.

Let $\tau^{(n)}=\inf\{t\geq 0:Y^{(n)}_{E}(t)=0\}$. Then $X^{(n)}_{i}(\tau^{(n)})$ $(i=0,1,\dots,d_{\rm max})$ give the numbers of susceptibles of the different degrees at the end of the epidemic and $T^{(n)}=\sum_{i=1}^{d_{\rm max}}(X^{(n)}_{i}(0)-X^{(n)}_{i}(\tau^{(n)})$.

5.2 Random time-scale transformation

We wish to apply Corollary 4.4 to obtain a central limit theorem for $T^{(n)}$ but that corollary cannot be applied directly to $\{\boldsymbol{W}^{(n)}(t)\}$ as $\tau^{(n)}\stackrel{{\scriptstyle p}}{{\longrightarrow}}\infty$ as $n\to\infty$. Thus we consider the following random time-scale transformation of $\{\boldsymbol{W}^{(n)}(t)\}$; cf. [Ethier and Kurtz(1986)], page 467, [Watson(1980)] and [Janson et al.(2012)].

For $t\in[0,\tau^{(n)}]$, let

[TABLE]

where $X^{(n)}_{E}(u)=\sum_{i=1}^{d_{\rm max}}iX^{(n)}_{i}(u)$. Let $\tilde{\tau}^{(n)}=A^{(n)}(\tau^{(n)})$. For $t\in[0,\tilde{\tau}^{(n)}]$, let $U^{(n)}(t)=\inf\{u\geq 0:A^{(n)}(u)=t\}$ and $\tilde{\boldsymbol{W}}^{(n)}(t)=\boldsymbol{W}^{(n)}\left(U^{(n)}(t)\right)$. Then $\{\tilde{\boldsymbol{W}}^{(n)}(t)\}=\{\tilde{\boldsymbol{W}}^{(n)}(t):0\leq t\leq\tilde{\tau}^{(n)}\}$ is a continuous-time Markov chain with transition intensities:

1. (i)

for $i=1,2,\dots,d_{\rm max}$ and $k=0,1,\dots,i-1$,

[TABLE] 2. (ii)

[TABLE] 3. (iii)

[TABLE]

For $t\geq 0$, let $\tilde{\boldsymbol{X}}^{(n)}(t)=(\tilde{X}^{(n)}_{0}(t),\tilde{X}^{(n)}_{1}(t),\dots,\tilde{X}^{(n)}_{d_{\rm max}}(t))$, so $\tilde{\boldsymbol{W}}^{(n)}(t)=(\tilde{\boldsymbol{X}}^{(n)}(t),\tilde{Y}^{(n)}_{E}(t),\tilde{Z}^{(n)}_{E}(t))$.

For reasons that will become clear later, we need to extend $\{\tilde{\boldsymbol{W}}^{(n)}(t)\}$ so that it is defined beyond time $\tilde{\tau}^{(n)}$ and allow $\tilde{Y}^{(n)}_{E}(t)$ and $\tilde{Z}^{(n)}_{E}(t)$ to take negative values. Thus we add two new transition intensities:

1. (ii*′*)

for $n^{Y}_{E}<0$,

[TABLE] 2. (iii*′*)

for $n^{Z}_{E}<0$,

[TABLE]

Note that the transition intensities (ii) and (iii) are defined only for $n^{Y}_{E}>0$ and $n^{Z}_{E}>0$, respectively.

The state space of $\{\tilde{\boldsymbol{W}}^{(n)}(t)\}$ is a subset of

[TABLE]

Let $\boldsymbol{l}_{ik}^{(1)}=-\boldsymbol{e}^{\rm S}_{i}+(k-1)\boldsymbol{e}^{\rm I}+(i-k-1)\boldsymbol{e}^{\rm R}$ $(i=1,2,\dots,d_{\rm max};k=0,1,\dots,i-1)$, $\boldsymbol{l}^{(2)}_{+}=-2\boldsymbol{e}^{\rm I}$, $\boldsymbol{l}^{(2)}_{-}=2\boldsymbol{e}^{\rm I}$, $\boldsymbol{l}^{(3)}_{+}=-\boldsymbol{e}^{\rm I}-\boldsymbol{e}^{\rm R}$ and $\boldsymbol{l}^{(3)}_{+}=\boldsymbol{e}^{\rm I}+\boldsymbol{e}^{\rm R}$. The set of possible jumps of $\{\tilde{\boldsymbol{W}}^{(n)}(t)\}$ from a typical state $\boldsymbol{n}$ is $\Delta=\Delta_{1}\cup\Delta_{2}\cup\Delta_{3}$, where $\Delta_{1}=\{\boldsymbol{l}^{(1)}_{ik}:i=1,2,\dots,d_{\rm max};k=0,1,\dots,i-1\}$, $\Delta_{2}=\{\boldsymbol{l}^{(2)}_{+},\boldsymbol{l}^{(2)}_{-}\}$ and $\Delta_{3}=\{\boldsymbol{l}^{(3)}_{+},\boldsymbol{l}^{(3)}_{-}\}$. Let $\boldsymbol{w}=(\boldsymbol{x},y_{E},z_{E})$, where $\boldsymbol{x}=(x_{0},x_{1},\dots,x_{d_{\rm max}})$, and $H=[0,1]^{d_{\rm max}+1}\times(-\infty,d_{\rm max}]^{2}$. The intensities of the jumps of $\{\tilde{\boldsymbol{W}}^{(n)}(t)\}$ admit the form

[TABLE]

where the functions $\tilde{\beta}_{\boldsymbol{l}}^{(n)}:H\to\mathbb{R}_{+}$ $(\boldsymbol{l}\in\Delta)$ are given by

[TABLE]

with

[TABLE]

The family of processes $\{\tilde{\boldsymbol{W}}^{(n)}(t)\}$ is asymptotically density dependent (see Section 4).

5.3 Proof of Theorem 2.1

Note that $\boldsymbol{W}^{(n)}(\tau^{(n)})=\tilde{\boldsymbol{W}}^{(n)}(\tilde{\tau}^{(n)})$, so the final size of the epidemic is given by

[TABLE]

Note also that $\tilde{\tau}^{(n)}=\inf\{t\geq 0:\tilde{Y}^{(n)}_{E}(t)=0\}$. We use Corollary 4.4 to obtain a central limit theorem for $\tilde{\boldsymbol{W}}^{(n)}(\tilde{\tau}^{(n)})$, and hence for $T^{(n)}$. The asymptotic variance matrix in the central limit theorem for $\tilde{\boldsymbol{W}}^{(n)}(\tilde{\tau}^{(n)})$ is not in closed form. However, we derive a closed-form expression for the asymptotic variance of $T^{(n)}$. The main concepts of the proof are presented here, with some detailed but elementary calculations deferred to Appendix A.

5.3.1 Deterministic model

As at (4.2), define the drift function $\tilde{F}(\boldsymbol{w})=\sum_{\boldsymbol{l}\in\Delta}\boldsymbol{l}\tilde{\beta}_{\boldsymbol{l}}(\boldsymbol{w})$. Using (5.3),

[TABLE]

For $t\geq 0$, let $\tilde{\boldsymbol{w}}(t)=(\tilde{x}_{0}(t),\tilde{x}_{1}(t),\dots,\tilde{x}_{d_{\rm max}}(t),\tilde{y}_{E}(t),\tilde{z}_{E}(t))$ be defined by

[TABLE]

Thus $\tilde{\boldsymbol{w}}(t)$ satisfies the differential equations

[TABLE]

where $\tilde{x}_{E}=\sum_{i=1}^{d_{\rm max}}i\tilde{x}_{i}$, having solution (see Appendix A.1)

[TABLE]

Let $\tilde{\eta}_{E}(t)=\tilde{x}_{E}(t)+\tilde{y}_{E}(t)+\tilde{z}_{E}(t)$. Note that, for $t\geq 0$,

[TABLE]

5.3.2 Initial conditions

Consider first the epidemic on the MR random graph. Recall that in the epidemic $\mathcal{E}^{(n)}$, for $i=0,1,\dots,d_{\rm max}$, there are $v^{(n)}_{i}$ individuals of degree $i$, of whom $a^{(n)}_{i}$ are initially infected. Also, $\left(\tilde{Y}^{(n)}_{E}(0),\tilde{Z}^{(n)}_{E}(0)\right)$ is given by the total numbers of infective and recovered half-edges created by the initial infectives. Thus $\tilde{X}^{(n)}_{i}(0)=v^{(n)}_{i}-a^{(n)}_{i}$ ($i=0,1,\dots,d_{\rm max}$) and

[TABLE]

where $Y_{ij}$ $(i=1,2,\dots,d_{\rm max};j=1,2,\dots)$ are independent and $Y_{ij}\sim{\rm Bin}(i,1-\exp(-\lambda I))$. Let $\epsilon_{E}=\sum_{i=1}^{d_{\rm max}}i\epsilon_{i}$. Then, recalling that $\epsilon^{(n)}_{i}=n^{-1}a^{(n)}_{i}$,

[TABLE]

For $i=1,2,\dots,d_{\rm max}$, let $\sigma_{Y,i}^{2}={\rm var}(Y_{i1})$. Now $\lim_{n\to\infty}\sqrt{n}(\epsilon^{(n)}_{i}-\epsilon_{i})=0$ $(i=0,1,\dots,d_{\rm max})$, by assumption, so the central limit theorem and Slutsky’s theorem imply that

[TABLE]

as $n\to\infty$, where

[TABLE]

(A closed-form expression for $\sigma_{Y}^{2}$ is given by (A.27) in Appendix A.3.2.) Since $\tilde{\boldsymbol{X}}^{(n)}(0)$ is non-random, it follows using (2.1) that

[TABLE]

where

[TABLE]

and

[TABLE]

Turning to the epidemic on the NSW random graph, recall that in $\mathcal{E}^{(n)}$ the number of initial infectives $a^{(n)}$ is prescribed and the $a^{(n)}$ initial infectives are chosen by sampling uniformly at random without replacement from the population. Thus the network can be constructed using two independent sets of i.i.d. copies of $D$, viz. $D_{1}^{\prime},D_{2}^{\prime},\dots,D_{n-a^{(n)}}^{\prime}$ for the initial susceptibles and $D_{1},D_{2},\dots,D_{a^{(n)}}$ for the initial infectives. Let $(Y_{E},Z_{E})$ be the bivariate random variable obtained by first sampling $D$ and then letting $(Y_{E},Z_{E})=(Y_{E},D-Y_{E})$, where $Y_{E}|D~{}\sim{\rm Bin}(D,1-\exp(-\lambda I))$. Let $\sigma_{Y_{E}}^{2}={\rm var}(Y_{E})$, $\sigma_{Y_{E},Z_{E}}={\rm cov}(Y_{E},Z_{E})$ and $\sigma_{Y_{E}}^{2}={\rm var}(Y_{E})$. (Closed-form expressions for $\sigma_{Y_{E}}^{2}$, $\sigma_{Z_{E}}^{2}$ and $\sigma_{Y_{E},Z_{E}}$ are given in (A.3.3)-(A.40) in Appendix A.3.3.) Let $\boldsymbol{p}=(p_{0},p_{1},\dots,p_{d_{\rm max}})$ and $\Sigma_{XX}$ be the $(d_{\rm max}+1)\times(d_{\rm max}+1)$ matrix with elements

[TABLE]

Recalling that $\lim_{n\to\infty}\sqrt{n}(\epsilon^{(n)}-\epsilon)=0$, where $\epsilon^{(n)}=n^{-1}a^{(n)}$, a similar argument to the above shows that

[TABLE]

where

[TABLE]

and

[TABLE]

5.3.3 Central limit theorem

Noting from (5.2) that

[TABLE]

it is easily checked that $\{\tilde{\boldsymbol{W}}^{(n)}(t)\}$ satisfies the conditions of Theorem 4.3. Now $\tilde{\tau}^{(n)}=\inf\{t\geq 0:\varphi(n^{-1}\tilde{\boldsymbol{W}}^{(n)}(t))\leq 0\}$, where $\varphi(\boldsymbol{w})=\varphi(\boldsymbol{x},y_{E},z_{E})=y_{E}$.

Suppose that $\tilde{x}_{i}(0)=p_{i}-\epsilon_{i}$ $(i=0,1,\dots,d_{\rm max})$, $\tilde{y}_{E}(0)=p_{I}\sum_{i=1}^{d_{\rm max}}i\epsilon_{i}$ and $\tilde{z}_{E}(0)=q_{I}\sum_{i=1}^{d_{\rm max}}i\epsilon_{i}$, so $\tilde{\eta}_{E}(0)=\mu_{D}$. Let $\tilde{\tau}=\inf\{t\geq 0:\varphi(\tilde{\boldsymbol{w}}(t))\leq 0\}=\inf\{t\geq 0:\tilde{y}_{E}(t)=0\}$. Then it follows from (5.11) that $\tilde{\tau}$ satisfies the equation

[TABLE]

We show in Appendix A.2 that, under the conditions of Theorem 2.1, the equation (5.22) has a unique solution in $(0,\infty)$. Note that $z={\rm e}^{-\tilde{\tau}}$, where $z$ is defined at (2.2). Also, using (5.10) the deterministic final fraction of the population that is susceptible is given by

[TABLE]

The corresponding deterministic final size is $\rho=1-\epsilon-f_{D_{\epsilon}}({\rm e}^{-\tilde{\tau}})$, agreeing with (2.3).

Let $a(\tilde{\tau})=\nabla\varphi(\boldsymbol{w}(\tilde{\tau}))\cdot\tilde{F}(\boldsymbol{w}(\tilde{\tau}))$. Then, using (5.6), and noting that $\tilde{y}_{E}(\tilde{\tau})=0$,

[TABLE]

Thus, using (5.10) and (5.13),

[TABLE]

since $\tilde{\eta}_{E}(0)=\mu_{D}$. We show in Appendix A.2 that $a(\tilde{\tau})<0$, so we may apply Corollary 4.4.

Writing $\tilde{\boldsymbol{w}}=(\tilde{w}_{0},\tilde{w}_{1},\dots,\tilde{w}_{p})$, where $p=d_{\rm max}+2$, let $\tilde{\Phi}(t,u)=[\tilde{\phi}_{ij}(t,u)]$ $(0\leq u\leq t<\infty)$, where

[TABLE]

Also, let

[TABLE]

where

[TABLE]

and $\Sigma_{0}=\Sigma^{\rm MR}_{0}$ or $\Sigma^{\rm NSW}_{0}$ dependening on whether the epidemic is on an MR or an NSW random graph. Then application of Corollary 4.4 yields

[TABLE]

where

[TABLE]

Thus, recalling (5.4),

[TABLE]

where

[TABLE]

Now $\nabla\varphi(\tilde{\boldsymbol{w}}(\tilde{\tau}))=(\boldsymbol{0},1,0)$ and, using (5.6), $(\boldsymbol{1},0,0)[\tilde{F}(\tilde{\boldsymbol{w}}(\tilde{\tau}))]^{\top}=-\tilde{x}_{E}(\tilde{\tau})$, so $(\boldsymbol{1},0,0)B=(\boldsymbol{1},b(\tilde{\tau}),0)$, where

[TABLE]

Further, it follows from (5.10)-(5.12) and (5.24) that

[TABLE]

and

[TABLE]

Further, let

[TABLE]

where

[TABLE]

and, for $i=0,1,\dots,d_{\rm max}$,

[TABLE]

Noting that $\boldsymbol{c}(\tilde{\tau},u)\boldsymbol{l}^{\top}$ is a scalar, it follows from (5.25) and (5.28) that

[TABLE]

The asymptotic variances, $\sigma^{2}_{\rm MR}$ and $\sigma^{2}_{\rm NSW}$ in Theorem 2.1 can be obtained by substituting $\Sigma_{0}=\Sigma^{\rm MR}_{0}$ and $\Sigma_{0}=\Sigma^{\rm NSW}_{0}$, respectively, in (5.36) and using (5.31)-(5.35), together with (5.3) and (5.10)-(5.12), to evaluate the second term on the right-hand side of (5.36). The details are lengthy and may be found in Appendix A.3.

5.4 Proof of Theorem 2.2

We prove Theorem 2.2 for the epidemic on an NSW random graph. The proof for the epidemic on an MR random graph is similar but simpler, as there is no randomness in the degrees of individuals, and is thus omitted. The proof proceeds in two stages. First, in Section 5.4.1, we couple the early stages of the epidemic $\tilde{\mathcal{E}}^{(n)}$, defined in Section 5.1, to a two-type branching process $\tilde{\mathcal{B}}^{(n)}$ which assumes that all infective half-edges in $\tilde{\mathcal{E}}^{(n)}$ are paired with susceptible half-edges. The branching processes $\tilde{\mathcal{B}}^{(n)}$ $(n=1,2,\dots)$ are coupled to a limiting branching process $\tilde{\mathcal{B}}$. The couplings and standard properties of the limiting branching process $\tilde{\mathcal{B}}$ show that, with probability tending to $1$ as $n\to\infty$, a major outbreak occurs if and only if the branching process $\tilde{\mathcal{B}}$ does not go extinct, and yield weak convergence results concerning the composition of the population in $\tilde{\mathcal{E}}^{(n)}$ when the number of infective half-edges first reaches $\log n$ in the event of a major outbreak (see Theorem 5.1). Then, in Section 5.4.2, we use the random time-scale transformation introduced in Section 5.2 to determine the asymptotic distribution of the final size of a major outbreak. The argument proceeds as in the proof of Theorem 2.1 but the equation defining $\tilde{\tau}$ now has a solution at [math] and one at $\tilde{\tau}>0$ (see the discussion following (5.49)) and a lower bounding branching process for the epidemic $\tilde{\mathcal{E}}^{(n)}$ is used to show that $\tilde{\tau}>0$ is the relevant asymptotic hitting time.

For ease of presentation we assume that, for $n=1,2,\dots$, there is one initial infective in $\tilde{\mathcal{E}}^{(n)}$ (i.e. that $a=1$), who is chosen by sampling a half-edge uniformly at random from all $D^{(n)}_{1}+D^{(n)}_{2}+\dots+D^{(n)}_{n}$ half-edges and infecting the individual who owns that half-edge. The proofs are easily extended to $a>1$ and other ways of choosing the initial infective(s) but the details are more complicated.

5.4.1 Coupling of epidemic and branching processes

Let $(\Omega,\mathcal{F},{\rm P})$ be a probability space on which are defined the following independent sets of random varaibles:

1. (i)

$D_{1},D_{2},\dots$ i.i.d. $\sim D$; 2. (ii)

$U_{0},U_{1},\dots$ i.i.d. $\sim{\rm U}(0,1)$; 3. (iii)

$L_{1},L_{2},\dots$ i.i.d. $\sim{\rm Exp}(1)$; 4. (iv)

$Y_{01},Y_{02},\dots Y_{0d_{\rm max}}$, where $Y_{0i}\sim{\rm Bin}(i,1-{\rm e}^{-\lambda I})$; 5. (v)

for $i=1,2,\dots,d_{\rm max}$, $Y_{i1},Y_{i2},\dots$ i.i.d. $\sim{\rm Bin}(i-1,1-{\rm e}^{-\lambda I})$.

For $n=1,2,\dots$, let $p^{(n)}_{i}=n^{-1}\sum_{k=1}^{n}1_{\{D_{k}=n\}}$ $(i=0,1,\dots,d_{\rm max})$ and $\tilde{p}^{(n)}_{i}=ip^{(n)}_{i}/\mu_{D}^{(n)}$ $(i=1,2,\dots,d_{\rm max})$, where $\mu_{D}^{(n)}=n^{-1}\sum_{k=1}^{n}D_{k}$. Note that by the strong law of large numbers $\mu_{D}^{(n)}>0$ (and $\tilde{p}^{(n)}_{i}$ is well defined) for all sufficiently large $n$ almost surely. Let $\tilde{c}^{(n)}_{i}=\sum_{j=1}^{i}\tilde{p}^{(n)}_{j}$ $(i=1,2,\dots,d_{\rm max})$. For $x\in(0,1)$, let $\tilde{d}^{(n)}(x)=\min\{i:x\leq\tilde{c}^{(n)}_{i}\}$. Similarly, let $\tilde{p}_{i}=ip_{i}/\mu_{D}$ and $\tilde{c}_{i}=\sum_{j=1}^{i}\tilde{p}_{j}$ $(i=1,2,\dots,d_{\rm max})$, and let $\tilde{d}(x)=\min\{i:x\leq\tilde{c}_{i}\}$ $(0<x<1)$.

For $n=1,2,\dots$, construct on $(\Omega,\mathcal{F},{\rm P})$ a realisation of a two-type continuous-time Markov branching process $\tilde{\mathcal{B}}^{(n)}$, which approximates the process of infected and recovered half-edges in the epidemic $\tilde{\mathcal{E}}^{(n)}$, as follows. The types are denoted $I$ and $R$ depending on whether the individual corresponds to an infective or recovered half-edge. Only type-$I$ individuals have offspring and they do so at their moment of death. Type-$R$ individuals live forever. For $t\geq 0$, let $\hat{Y}^{(n)}_{E}(t)$ and $\hat{Z}^{(n)}_{E}(t)$ denote respectively the number of type-$I$ and type-$R$ individuals alive in $\tilde{\mathcal{B}}^{(n)}$ at time $t$. The initial state $(\hat{Y}^{(n)}_{E}(0),\hat{Z}^{(n)}_{E}(0))$ is determined as follows. Let $d^{(n)}_{0}=\tilde{d}^{(n)}(U_{0})$, which corresponds to the initial infective in $\tilde{\mathcal{E}}^{(n)}$ having degree $d^{(n)}_{0}$. If $d^{(n)}_{0}=0$ then $(\hat{Y}^{(n)}_{E}(0),\hat{Z}^{(n)}_{E}(0))=(0,0)$ and $\tilde{\mathcal{B}}^{(n)}$ goes extinct immediately. Alternatively, if $d^{(n)}_{0}>0$ then $(\hat{Y}^{(n)}_{E}(0),\hat{Z}^{(n)}_{E}(0))=(Y_{0d^{(n)}_{0}},d_{0}-Y_{0d^{(n)}_{0}})$. In that case, for $k=1,2,\dots$, the $k$th type-$I$ individual born in $\tilde{\mathcal{B}}^{(n)}$ (including the initial individuals) has degree $d^{(n)}_{k}=\tilde{d}^{(n)}(U_{k})$ and lives until age $L_{k}$, when it dies. Denote this individual by $i^{*}$. Suppose that $d^{(n)}_{k}=i$ and $i^{*}$ is the $l$th degree-$i$ individual (excluding the initial individuals) born in $\tilde{\mathcal{B}}^{(n)}$. Then, when $i^{*}$ dies, it leaves $Y_{il}$ type-$I$ and $i-1-Y_{il}$ type-$R$ offspring. Of course, reproduction stops in $\tilde{\mathcal{B}}^{(n)}$ if $\hat{Y}^{(n)}_{E}(t)=0$. Construct also on $(\Omega,\mathcal{F},{\rm P})$ a realisation of a two-type continuous-time branching process $\tilde{\mathcal{B}}$, defined analogously to $\tilde{\mathcal{B}}^{(n)}$ but using the function $\tilde{d}$ instead of $\tilde{d}^{(n)}$. For $t\geq 0$, let $Y_{E}(t)$ and $Z_{E}(t)$ denote respectively the numbers of type-$I$ and type-$R$ individuals alive in $\tilde{\mathcal{B}}$ at time $t$.

For $n=1,2,\dots$, construct on $(\Omega,\mathcal{F},{\rm P})$ a realisation of the epidemic $\tilde{\mathcal{E}}^{(n)}$, defined in Section 5.1, as follows. Give the $n$ individuals in $\tilde{\mathcal{E}}^{(n)}$ the labels $1,2,\dots,n$ in increasing order of degree. Now label the $n\mu_{D}^{(n)}$ half-edges $1,2,\dots,n\mu_{D}^{(n)}$, starting with the half-edges (if any) attached to individual $1$, then the half-edges (if any) attached to individual $2$ and so on. Thus half-edges attached to the same individual have consecutive labels. As in Section 5.1, for $t\geq 0$, let $Y^{(n)}_{E}(t)$ and $Z^{(n)}_{E}(t)$ denote respectively the number of infective and recovered half-edges in $\tilde{\mathcal{E}}^{(n)}$ at time $t$. The initial infective in $\tilde{\mathcal{E}}^{(n)}$ is the individual who owns the half-edge labelled $\left\lfloor{n\mu_{D}^{(n)}U_{0}}\right\rfloor+1$. Note that this individual has degree $d^{(n)}_{0}=\tilde{d}^{(n)}(U_{0})$. If $d^{(n)}_{0}=0$ then the epidemic stops immediately. Alternatively, if $d^{(n)}_{0}>0$ then the initial infective infects along $k^{(n)}_{0}=Y_{0d^{(n)}_{0}}$ of its half-edges with its remaining $d^{(n)}_{0}-k^{(n)}_{0}$ half-edges becoming recovered half-edges. The epidemic stops if $k^{(n)}_{0}=0$, otherwise the initial infective transmits infection along its infective half-edges at times $L_{1},L_{2},\dots,L_{k_{0}}$.

For $k=1,2,\dots,$ let $l_{k}^{(n)}=\left\lfloor{n\mu_{D}^{(n)}U_{k}}\right\rfloor+1$ and note that the half-edge having label $l_{k}^{(n)}$ is attached to an individual having degree $d^{(n)}_{k}=\tilde{d}^{(n)}(U_{k})$. When infection is transmitted along a half-edge that half-edge, $l_{*}$ say, is attempted to be paired with the half-edge having label $l_{k}^{(n)}$, where $k$ is the number of the $U_{0},U_{1},\dots$ that have been used already in the construction of $\tilde{\mathcal{E}}^{(n)}$. (Thus, for example, the first half-edge emanating from the initial infective is attempted to be paired with the half-edge having label $l_{1}^{(n)}$.) If the half-edge $l_{k}^{(n)}$ has already been paired or $l_{k}^{(n)}=l_{*}$ then the attempt fails and $l_{*}$ is attempted to be paired with the half-edge $l_{k+1}^{(n)}$, and so on until a valid pairing is obtained and an edge is formed. Suppose that a valid pairing is made with the half-edge having label $l_{V}$. Let $i^{*}$ be the individual that owns the half-edge $l_{V}$ and $i$ be the degree of $i^{*}$. If $i^{*}$ is susceptible then it becomes an infective, otherwise nothing happens apart from the formation of the edge. Suppose that $i^{*}$ is susceptible and is the $l$th degree-$i$ susceptible to be infected in $\tilde{\mathcal{E}}^{(n)}$, excluding the initial infective. Then $i^{*}$ infects along $Y_{il}$ of its half-edges and its remaining $i-1-Y_{il}$ half-edges become recovered half-edges. (When $i^{*}$ was infected one of its $i$ half-edges was paired to its infector.) Let $k_{1}=Y_{il}$. The times of these $k_{1}$, infections, relative to the infection time of $i^{*}$, are given by $L_{k^{*}+1},L_{k^{*}+2},\dots,L_{k^{*}+k_{1}}$, where $k^{*}$ is the number of infective half-edges created in $\tilde{\mathcal{E}}^{(n)}$ prior to the infection of $i^{*}$. The epidemic terminates when $Y^{(n)}_{E}(t)=0$, or when $Y^{(n)}_{E}(t)=1$ and all other half-edges have been paired.

As in Section 5.1, for $t\geq 0$, let $\boldsymbol{X}^{(n)}(t)=(X^{(n)}_{0}(t),X^{(n)}_{1}(t),\dots,X^{(n)}_{d_{\rm max}}(t))$, where $X^{(n)}_{i}(t)$ is the number of degree-$i$ susceptible individuals at time $t$ in $\tilde{\mathcal{E}}^{(n)}$. For $n=1,2,\dots$, let $\hat{\tau}^{(n)}=\inf\{t\geq 0:Y^{(n)}_{E}(t)\geq\log n\}$, where $\hat{\tau}^{(n)}=\infty$ if $Y^{(n)}_{E}(t)<\log n$ for all $t\geq 0$. Let $A_{\rm ext}=\{\omega\in\Omega:\lim_{t\to\infty}Y_{E}(t)=0\}$ denote the set on which type-$I$ individuals become extinct in the branching process $\tilde{\mathcal{B}}$ and let $\alpha$ denote the Malthusian parameter of $\tilde{\mathcal{B}}$. Then $\alpha$ is given by the unique real solution of the equation

[TABLE]

so $\alpha=R_{0}-1$, where $R_{0}$ is defined at (2.8).

Theorem 5.1

1. (a)

$\lim_{n\to\infty}{\rm P}\left(\hat{\tau}^{(n)}=\infty|A_{\rm ext}\right)=1$. 2. (b)

Suppose that $R_{0}>1$, so ${\rm P}(A_{\rm ext}^{C})>0$. Then, as $n\to\infty$,

1. (i)

${\rm P}\left(\hat{\tau}^{(n)}<\infty|A_{\rm ext}^{C}\right)\to 1$; 2. (ii)

$\frac{1}{\log n}Y^{(n)}_{E}(\hat{\tau}^{(n)})\left|\right.A_{\rm ext}^{C}\stackrel{{\scriptstyle p}}{{\longrightarrow}}1$; 3. (iii)

$\frac{1}{\log n}Z^{(n)}_{E}(\hat{\tau}^{(n)})\left|\right.A_{\rm ext}^{C}\stackrel{{\scriptstyle p}}{{\longrightarrow}}\alpha^{-1}q_{I}\mu_{\tilde{D}-1}$; 4. (iv)

$\sqrt{n}\left(n^{-1}\boldsymbol{X}^{(n)}(\hat{\tau}^{(n)})-\boldsymbol{p}\right)|A_{\rm ext}^{C}\stackrel{{\scriptstyle D}}{{\longrightarrow}}{\rm N}(\boldsymbol{0},\Sigma_{XX})$, where $\Sigma_{XX}$ is defined at (5.18).

• Proof.

The key observations underlying the proof are that (i) the processes $\{(\hat{Y}^{(n)}_{E}(t),\hat{Z}^{(n)}_{E}(t)):t\geq 0\}$ and $\{(Y^{(n)}_{E}(t),Z^{(n)}_{E}(t)):t\geq 0\}$ coincide up until at least the first time that an attempt is made to pair a half-edge with a half-edge belonging to an individual previously used in the construction of $\tilde{\mathcal{E}}^{(n)}$ and (ii) the branching processes $\tilde{\mathcal{B}}^{(n)}$ and $\tilde{\mathcal{B}}$ coincide up until the first time that $\tilde{d}^{(n)}(U_{k})\neq\tilde{d}(U_{k})$. Thus we show that the probability that the processes $\{(\hat{Y}^{(n)}_{E}(t),\hat{Z}^{(n)}_{E}(t)):0\leq t\leq\hat{\tau}^{(n)}\}$ and $\{(Y_{E}(t),Z_{E}(t)):0\leq t\leq\hat{\tau}^{(n)}\}$ coincide converges to $1$ as $n\to\infty$. The theorem then follows using standard results concerning the asymptotic behaviour of the branching process $\tilde{\mathcal{B}}$.

For $n=1,2,\dots$ and $k=1,2,\dots,n\mu_{D}^{(n)}$, let $\mathcal{C}^{(n)}(k)$ be the set of half-edges attached to the individual owning half-edge $k$ in $\tilde{\mathcal{E}}^{(n)}$. Thus $k\in\mathcal{C}^{(n)}(k)$ and $|\mathcal{C}^{(n)}(k)|\leq d_{\rm max}$. For $n=1,2,\dots$, let

[TABLE]

be the number of pairings required in $\tilde{\mathcal{E}}^{(n)}$ until an attempt is made to pair a half-edge with a half-edge belonging to a previoulsy used individual. Let $\boldsymbol{D}=(D_{1},D_{2},\dots)$. Then, for $m=1,2,\dots$,

[TABLE]

Now $\mu_{D}^{(n)}\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}\mu_{D}$ as $n\to\infty$ by the strong law of large numbers. Hence, for any $\gamma\in(0,1/2)$, ${\rm P}(M^{(n)}\leq n^{\gamma}|\boldsymbol{D})\stackrel{{\scriptstyle p}}{{\longrightarrow}}0$ as $n\to\infty$. Thus,

[TABLE]

as ${\rm P}(M^{(n)}\leq n^{\gamma}|\boldsymbol{D})$ $(n=1,2,\dots)$ is uniformly bounded.

Turning to the coupling of $\tilde{\mathcal{B}}^{(n)}$ and $\tilde{\mathcal{B}}$, for $n=\,2,\dots$,

[TABLE]

The first inequality in (5.4.1) follows from the total variation distance between the distributions $\tilde{p}^{(n)}_{i}$ $(i=1,2,\dots,d_{\rm max})$ and $\tilde{p}_{i}$ $(i=1,2,\dots,d_{\rm max})$, and the second inequality follows using the triangle inequality and elementary manipulation ([Ball and Neal(2017)], equation (3.2)). For $n=1,2,\dots$, let $\hat{M}^{(n)}=\min\{k\geq 0:\tilde{d}^{(n)}(U_{k})\neq\tilde{d}(U_{k})\}$. Then, for $k=0,1,\dots$,

[TABLE]

Now, for $i=1,2,d_{\rm max}$, as $np^{(n)}_{i}\sim{\rm Bin}(n,p_{i})$,

[TABLE]

Thus, for any $\gamma\in(0,1/2)$,

[TABLE]

For $t\geq 0$, let $T_{E}(t)$ be the total number of individuals (of either type) alive in the branching process $\tilde{\mathcal{B}}$ during $[0,t]$ and let $T_{E}(\infty)=\lim_{t\to\infty}T_{E}(t)$. Thus $T_{E}(\infty,\omega)<\infty$ if and only if $\omega\in A_{\rm ext}$. For $n=1,2,\dots$ and $t\geq 0$, let

[TABLE]

Let $\tau=\min\{t\geq 0:Y_{E}(t)=0\}$, where $\tau(\omega)=\infty$ if $\omega\in A_{\rm ext}^{C}$. Then (5.37) and (5.39) imply that, for $k=0,1,\dots$,

[TABLE]

Part (a) of the theorem follows since $A_{\rm ext}=\bigcup_{k=0}^{\infty}\{T_{E}(\infty)=k\}$.

Note that although $\tilde{\mathcal{B}}$ is a two-type branching process, since only type-$I$ individuals produce offspring, it can be treated as a single-type branching process consisting only of type-$I$ individuals, in which attached to each (type-$I$) individual is a random characteristic (see [Nerman(1981)]) given by its number of type-$R$ offspring. Then it follows from [Nerman(1981)], Theorem 5.4, that there exists a random variable $W\geq 0$, where $W(\omega)=0$ if and only if $\omega\in A_{\rm ext}$ such that, as $t\to\infty$,

[TABLE]

[TABLE]

and

[TABLE]

(It is easily checked that (5.40)-(5.42) hold by calculating the the appropriate $m_{\infty}^{\phi}$ in  [Nerman(1981)], Theorem 5.4, and using Corollary 3.2 of that paper. For a heuristic argument note, for example, that if $Y_{E}(u)\approx W{\rm e}^{\alpha u}$ $(u\geq 0)$ then $T_{E}(t)\approx\int_{0}^{t}\mu_{\tilde{D}-1}{\rm e}^{\alpha u}W\,{\rm d}u\approx\alpha^{-1}\mu_{\tilde{D}-1}W{\rm e}^{\alpha t}$, since type-$I$ individuals die at rate $1$ and have on average $\mu_{\tilde{D}-1}$ offspring.)

For $n=1,2,\dots$, let $\bar{\tau}^{(n)}=\inf\{t\geq 0:Y_{E}(t)\geq\log n\}$, where $\bar{\tau}^{(n)}(\omega)=\infty$ if $\omega\in A_{\rm ext}$. Then it follows from (5.40) that $\bar{\tau}^{(n)}(\omega)<\infty$ and $\lim_{n\to\infty}\bar{\tau}^{(n)}(\omega)=\infty$, for ${\rm P}$-almost all $\omega\in A_{\rm ext}^{C}$. Thus, for ${\rm P}$-almost all $\omega\in A_{\rm ext}^{C}$,

[TABLE]

since an individual has at most $d_{\rm max}$ offspring, and using (5.40)-(5.42),

[TABLE]

and

[TABLE]

Now (5.37), (5.39) and (5.45) imply that

[TABLE]

and (i)-(iii) of part (b) of the theorem follow using (5.41) and (5.43).

To prove part (b)(iv), note by the multivariate central limit theorem that

[TABLE]

Also, (5.45) and (5.46) imply that, for $i=0,1,\dots,d_{\rm max}$,

[TABLE]

so $\frac{1}{\sqrt{n}}\left(X^{(n)}_{i}(\hat{\tau}^{(n)})-X^{(n)}_{i}(0)\right)|A_{\rm ext}^{C}\stackrel{{\scriptstyle p}}{{\longrightarrow}}0$ as $n\to\infty$. Part (b)(iv) now follows using Slutsky’s theorem. $\Box$

5.4.2 Epidemic starting with $\left\lceil{\log n}\right\rceil$ infectives

Turning to the proof of Theorem 2.2, note that by Theorem 5.1, $\lim_{n\to\infty}{\rm P}\left(G^{(n)}\triangle A_{\rm ext}\right)=0$, where $\triangle$ denotes symmetric difference. Hence, using the strong Markov property, we can determine the asymptotic distribution of $T^{(n)}_{\rm NSW}|G^{(n)}$ by considering the random time-scale transformed process $\{\tilde{\boldsymbol{W}}^{(n)}(t)\}$, defined in Section 5.2, with initial state $\tilde{\boldsymbol{W}}^{(n)}(0)=(\boldsymbol{X}^{(n)}(\hat{\tau}^{(n)}),Y^{(n)}_{E}(\hat{\tau}^{(n)}),Z^{(n)}_{E}(\hat{\tau}^{(n)}))$. Thus, by Theorem 5.1(b),

[TABLE]

where $\Sigma^{\rm NSW}_{0}$ is obtained by setting $\epsilon=0$ in (5.21).

By Theorem 4.1, for any $t>0$,

[TABLE]

where, using (5.11)

[TABLE]

Note that $\tilde{y}_{E}(t)=0$ if and only if

[TABLE]

and that $t=0$ is a solution of (5.49). If $R_{0}\leq 1$ then $t=0$ is the only solution in $[0,\infty)$, but if $R_{0}>1$ then we show in Appendix A.2 that, under the conditions of Theorem 2.2, there is a unique solution in $(0,\infty)$ which we denote by $\tilde{\tau}$.

As in Section 5.3.3, let $\tilde{\tau}^{(n)}=\inf\{t\geq 0:\varphi(n^{-1}\tilde{\boldsymbol{W}}^{(n)}(t))\leq 0\}$. Recall that $R_{0}>1$ is assumed in Theorem 2.2. It follows from (5.48) that, $\min\left(\tilde{\tau}^{(n)},|\tilde{\tau}^{(n)}-\tilde{\tau}|\right)\stackrel{{\scriptstyle p}}{{\longrightarrow}}0$ as $n\to\infty$. We show below that there exists $\epsilon_{0}>0$ such that

[TABLE]

Theorem 4.3 and Corollary 4.4, as stated, cannot be applied in the present setting as $\varphi(\tilde{\boldsymbol{w}}(0))=0$. However, in the terminology of Theorem 4.3, the proof of [Ethier and Kurtz(1986)], Theorem 11.4.1, extends easily to the situation when $\tau^{(n)}=\inf\{t\geq\epsilon_{0}:\varphi(n^{-1}\boldsymbol{X}^{(n)}(t))\leq 0\}$ and $\tau=\inf\{t\geq\epsilon_{0}:\varphi(\boldsymbol{x}(t))\leq 0\}$, for fixed $\epsilon_{0}>0$. Hence, so do Theorem 4.3 and Corollary 4.4.

Let

[TABLE]

Then using (5.50) and the above extension of Corollary 4.4 (we show in Appendix A.2 that $a(\tilde{\tau})=\nabla\varphi(\boldsymbol{w}(\tilde{\tau}))\cdot\tilde{F}(\boldsymbol{w}(\tilde{\tau}))<0$), setting $\epsilon=0$ in the proof of Theorem 2.1 yields

[TABLE]

Theorem 2.2 follows using Slutsky’s theorem since, by (5.45),

[TABLE]

To complete the proof of Theorem 2.2 we use a device, first introduced by [Whittle(1955)] in the setting of a homegeneously mixing epidemic, to show that there exists $\epsilon_{0}>0$ such that (5.50) holds. For $n=1,2,\dots$, construct on $(\Omega,\mathcal{F},{\rm P})$ an epidemic $\check{\mathcal{E}}^{(n)}$ analogously to $\tilde{\mathcal{E}}^{(n)}$ except having initial state given in an obvious notation by

[TABLE]

Let $\check{\mathcal{B}}^{(n)}$ be the corresponding branching process which assumes all attempted pairing in $\check{\mathcal{E}}^{(n)}$ are valid. Let

[TABLE]

be the final size of $\check{\mathcal{E}}^{(n)}$. Fix $\delta\in(0,1)$. Then, provided $\check{T}^{(n)}\leq n\delta$, for a given pairing, the probability that an invalid pairing is attempted is at most $p^{(n)}(\delta)=\min(\delta nd_{\rm max}/(n\mu_{D}^{(n)}),1)=\min(\delta d_{\rm max}/\mu_{D}^{(n)},1)$. It follows that

[TABLE]

where $\check{T}^{(n)}_{B}(\delta)$ is the total number of type-$I$ individuals (excluding the initial individuals) ever alive in the branching process $\check{\mathcal{B}}^{(n)}(\delta)$ that is obtained from $\check{\mathcal{B}}^{(n)}$ by aborting each type-$I$ individual independently with probability $p^{(n)}(\delta)$.

Let $R^{(n)}_{0}(\delta)$ be the mean number of type-$I$ individuals spawned by a typical type-$I$ individual in $\check{\mathcal{B}}^{(n)}(\delta)$ and let $\pi^{(n)}(\delta)$ be the extinction probability for type-$I$ individuals in $\check{\mathcal{B}}^{(n)}(\delta)$ if the process were to start with one type-$I$ individual. As $n\to\infty$, $p^{(n)}_{i}\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}p_{i}$ ($i=0,1,\dots,d_{\rm max}$), by the strong law of large numbers, so $\tilde{p}^{(n)}_{i}\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}\tilde{p}_{i}$ ($i=1,2,\dots,d_{\rm max}$) and $\mu^{(n)}_{\tilde{D}-1}\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}\mu_{\tilde{D}-1}$, where $\mu^{(n)}_{\tilde{D}-1}=\sum_{i=1}^{d_{\rm max}}(i-1)\tilde{p}^{(n)}_{i}$. Thus,

[TABLE]

where $p(\delta)=\min(\delta d_{\rm max}/\mu_{D},1)$. Further, as $n\to\infty$, the offspring distribution of $\check{\mathcal{B}}^{(n)}(\delta)$ converges almost surely to that of $\tilde{\mathcal{B}}(\delta)$, where $\tilde{\mathcal{B}}(\delta)$ is the branching process obtained from $\tilde{\mathcal{B}}$ by aborting each type-$I$ individual independently with probability $p(\delta)$. Thus, by a simple extension of [Britton et al.(2007)], Lemma 4.1, $\pi^{(n)}(\delta)\stackrel{{\scriptstyle a.s.}}{{\longrightarrow}}\pi(\delta)$ as $n\to\infty$, where $\pi(\delta)$ is the extinction probability for type-$I$ individuals in $\tilde{\mathcal{B}}(\delta)$ starting from one type-$I$ individual.

Recall that $R_{0}>1$. Thus $\delta\in(0,1)$ can be chosen sufficently small such that $(1-p(\delta))R_{0}>1$, whence $\pi(\delta)<1$.

Now, $\check{Y}^{(n)}_{E}(0)\geq\log n$ by construction, so using (5.52), for such $\delta$,

[TABLE]

for ${\rm P}$-almost all $\omega\in A_{\rm ext}^{C}$. Hence

[TABLE]

Finally, for $t\geq 0$, let $\tilde{x}(t)=\sum_{i=0}^{d_{\rm max}}\tilde{x}_{i}(t)=f_{D}({\rm e}^{-t})$, using (5.10) with $\tilde{x}_{i}(0)=p_{i}$ $(i=0,1,\dots,d_{\rm max})$. By construction, $\boldsymbol{W}^{(n)}(\tau^{(n)})=\tilde{\boldsymbol{W}}^{(n)}(\tilde{\tau}^{(n)})$, so (5.53) and Theorem 4.1 applied to the process $\{\tilde{\boldsymbol{W}}^{(n)}(t)\}$, with initial state satisfying (5.47), imply that (5.50) holds with $\epsilon_{0}$ being given by the unique solution in $(0,\infty)$ of $1-f_{D}({\rm e}^{-\epsilon_{0}})=\delta$.

5.5 Proof of Theorem 2.3

The proof of Theorem 2.3 for the epidemic on the NSW random graph parallels in an obvious fashion the argument in Section 5.4.2 above without the need to condition on $A_{\rm ext}^{C}$, so the details are omitted. Again, the proof for the epidemic on the MR random graph is similar but simpler.

5.6 Proof of Theorem 2.7

Consider first bond percolation and note that when $I\equiv 1$ and $\lambda=-\log(1-\pi)$ then in the directed random graph $\tilde{\mathcal{G}}^{(n)}$ defined at the start of Section 5.1, all the possible directed edges in $\tilde{\mathcal{G}}^{(n)}$ (i.e. between pairs of neighbours in $\mathcal{G}^{(n)}$) are present independently, each with probability $\pi$. Further, when constructing the final outcome of the corresponding epidemic $\mathcal{E}^{(n)}$ using $\tilde{\mathcal{G}}^{(n)}$, for any pair $(i,j)$ of distinct individuals use is made of at most one of the directed edges $i\to j$ and $j\to i$ (if $i$ infects $j$, whether or not $j$ tries to infect $i$ is immaterial). Thus, in this situation $\tilde{\mathcal{G}}^{(n)}$ can be replaced by $\mathcal{G}^{(n)}_{\rm bond}$, obtained using bond percolation on $\mathcal{G}^{(n)}$, and in $\mathcal{E}^{(n)}$, an initial susceptible is ultimately infected if and only if in $\mathcal{G}^{(n)}_{\rm bond}$ there is a chain of edges connecting it to an initial infective.

Suppose that there is one initial infective. Then the final size of the epidemic (including the initial infective) is given by the size of the connected component of $\mathcal{G}^{(n)}_{\rm bond}$ that contains the initial infective. With probability tending to $1$ as $n\to\infty$, a major outbreak occurs in $\mathcal{E}^{(n)}$ if and only if the initial infective belongs to the largest connected component of $\mathcal{G}^{(n)}_{\rm bond}$, and the final size of a major epidemic is given by the size $C^{(n)}$ of that connected component. Further $C^{(n)}$ has the same asymptotic distribution as $T^{(n)}|G^{(n)}$. Thus (2.11) and (2.12) follow immediately on setting $I\equiv 1$ and $\lambda=-\log(1-\pi)$ in Theorem 2.2. (Note that $p_{I}=\pi$ and, since $I$ is constant, $q_{I}^{(2)}=q_{I}^{2}$.)

Turning to site percolation, consider the epidemic $\mathcal{E}^{(n)}$ with ${\rm P}(I=\infty)=\pi=1-{\rm P}(I=0)$, so each infective infects all of its neighbours with probability $\pi$ and none of them otherwise. Suppose that there is one initial infective, $i^{*}$ say, and let $\mathcal{T}^{(n)}=\{j\in\mathcal{N}^{(n)}\setminus\{i^{*}\}:i^{*}\leadsto j\}$ using the directed random graph $\tilde{\mathcal{G}}^{(n)}$. Then $\mathcal{T}^{(n)}\cup\{i^{*}\}$ differs from the connected component containing $i^{*}$ in $\mathcal{G}^{(n)}_{\rm site}$ (site percolation on $\mathcal{G}^{(n)}$) in that $\mathcal{T}^{(n)}$ also includes infected indiviudals having $I=0$, which are deleted in $\mathcal{G}^{(n)}_{\rm site}$. Thus, to obtain a central limit theorem for $C^{(n)}$, we need one for $V^{(n)}=|\{j\in\mathcal{N}^{(n)}\setminus\{i^{*}\}:I_{j}=\infty\mbox{ and }i^{*}\leadsto j\}|$, the final size of $\mathcal{E}^{(n)}$ counting only individuals with $I=\infty$. This can be obtained by augmenting the process $\{\boldsymbol{W}^{(n)}(t)\}$ as we now outline.

Let $\{\hat{\boldsymbol{W}}^{(n)}(t)\}=\{\hat{\boldsymbol{W}}^{(n)}(t):t\geq 0\}$, where

[TABLE]

and $V^{(n)}(t)$ is the total number of initial susceptibles that are infected during $(0,t]$ and have $I=\infty$. A typical element of $\hat{H}^{(n)}$, the state space of $\{\hat{\boldsymbol{W}}^{(n)}(t)\}$, is now $\hat{\boldsymbol{n}}=(n^{X}_{0},n^{X}_{1},\dots,n^{X}_{d_{\rm max}},n^{Y}_{E},n^{Z}_{E},n^{V})$. The transition intensities of $\{\hat{\boldsymbol{W}}^{(n)}(t)\}$ are essentially the same of those of $\{\boldsymbol{W}^{(n)}(t)\}$, given at the end of Section 5.1, except now transitions of type (i) are partitioned according to whether or not the infected susceptible has $I=\infty$. For $i=1,2,\dots,d_{\rm max}$, we have

[TABLE]

corresponding to the degree-$i$ infected susceptible having $I=0$, and

[TABLE]

corresponding to the degree-$i$ infected susceptible having $I=\infty$.

We apply the same random time-scale transformation to $\{\hat{\boldsymbol{W}}^{(n)}(t)\}$ as done to $\{\boldsymbol{W}^{(n)}(t)\}$ in Section 5.2. Denote the time-transformed process by $\{\bar{\boldsymbol{W}}^{(n)}(t)\}$. Let $\hat{\boldsymbol{l}}^{(1)}_{i,0}=\boldsymbol{l}_{i,0}^{(1)}$ and $\hat{\boldsymbol{l}}^{(1)}_{i,i-1}=\boldsymbol{l}_{i,i-1}^{(1)}+\boldsymbol{e}^{\rm V}$ $(i=1,2,\dots,d_{\rm max})$, and $\hat{\boldsymbol{l}}^{(j)}_{+}=\boldsymbol{l}^{(j)}_{+}$ and $\hat{\boldsymbol{l}}^{(j)}_{-}=\boldsymbol{l}^{(j)}_{-}$ $(j=2,3)$. The set of possible jumps of $\{\bar{\boldsymbol{W}}^{(n)}(t)\}$ from a typical state $\hat{\boldsymbol{n}}$ is now $\hat{\Delta}=\hat{\Delta}_{1}\cup\hat{\Delta}_{2}\cup\hat{\Delta}_{3}$, where $\hat{\Delta}_{1}=\bigcup_{i=1}^{d_{\rm max}}\{\hat{\boldsymbol{l}}^{(1)}_{i,0}\}\cup\{\hat{\boldsymbol{l}}^{(1)}_{i,i-1}\}$, $\hat{\Delta}_{2}=\{\hat{\boldsymbol{l}}^{(2)}_{+},\hat{\boldsymbol{l}}^{(2)}_{-}\}$ and $\hat{\Delta}_{3}=\{\hat{\boldsymbol{l}}^{(3)}_{+},\hat{\boldsymbol{l}}^{(3)}_{-}\}$.

Let $\hat{\boldsymbol{w}}=(\boldsymbol{x},y_{E},z_{E},v)$. The family of processes $\{\bar{\boldsymbol{W}}^{(n)}(t)\}$ is again asymptotically density dependent with corresponding functions $\bar{\beta}_{\boldsymbol{l}}(\hat{\boldsymbol{w}})$ $(\boldsymbol{l}\in\hat{\Delta})$ given by (cf. (5.3))

[TABLE]

The associated drift function is (cf. (5.6))

[TABLE]

For $t\geq 0$, let $\bar{\boldsymbol{w}}(t)=(\bar{x}_{0}(t),\bar{x}_{1}(t),\dots,\bar{x}_{d_{\rm max}}(t),\bar{y}_{E}(t),\bar{z}_{E}(t),\bar{v}(t))$ be defined analogously to $\tilde{\boldsymbol{w}}(t)$ at (5.6). Noting that $p_{I}=\pi$ and $q_{I}=1-\pi$, the corresponding deterministic model for $\bar{\boldsymbol{w}}(t)$ is given by (5.7)-(5.9), with $\tilde{x}_{i}$ replaced by $\bar{x}_{i}$ etc., augmented with

[TABLE]

Thus, (5.10)-(5.12) still hold, with the above change of notation, and, using (5.10),

[TABLE]

The stopping time $\tilde{\tau}^{(n)}$ is unchanged, except now $\varphi(\hat{\boldsymbol{w}})=\varphi(\boldsymbol{x},y_{E},z_{E},v)=y_{E}$. Recall that site percolation corresponds to one initial infective, so under both the MR and NSW models, $\epsilon_{i}=0$ and $\tilde{x}_{i}(0)=p_{i}$ $(i=0,1,\dots,d_{\rm max})$. Thus $\tilde{\tau}=\inf\{t\geq 0:\bar{y}_{E}(t)=0\}$ ($=\inf\{t\geq 0:\tilde{y}_{E}(t)=0\}$) is given by the unique solution of (5.49) in $(0,\infty)$. Hence $z={\rm e}^{-\tilde{\tau}}$ is given by the unique solution in $(0,1)$ of (2.10). Now $\tilde{v}(0)=0$, since $V^{(n)}(0)=0$ for all $n=1,2,\dots$, so using (5.56) and recalling that $\rho=1-f_{D}(z)$,

[TABLE]

Let $\bar{V}^{(n)}=\bar{V}^{(n)}(\tilde{\tau}^{(n)})$, so in the site percolation setting, $C^{(n)}\stackrel{{\scriptstyle D}}{{=}}1+\bar{V}^{(n)}|G^{(n)}$. Now application of Corollary 4.4, as at (5.27), yields

[TABLE]

where

[TABLE]

and $\bar{\Sigma}(\tilde{\tau})$ is obtained by making obvious changes to (5.24)-(5.26) to account for the extra dimension. Further, $\bar{V}^{(n)}=\bar{\boldsymbol{W}}^{(n)}(\tilde{\tau}^{(n)})(\boldsymbol{0},0,0,1)^{\top}$, so

[TABLE]

where

[TABLE]

A simple calculation yields $(\boldsymbol{0},0,0,1)\bar{B}=(\boldsymbol{0},-\pi b(\tilde{\tau}),0,1)$, where $b(\tilde{\tau})$ is given by (5.29) with $a(\tilde{\tau})$ obtained by replacing $f_{D_{\epsilon}}$ by $f_{D}$ in (5.23). Define $\bar{\Phi}(t,u)$ analogously to $\tilde{\Phi}(t,u)$ at (5.24). Then it follows using (5.10)-(5.12) and (5.56) that

[TABLE]

and $\left[(\boldsymbol{0},1,0,0)\bar{\Phi}(\tilde{\tau},u)\right]_{i}$ is [math], if $i=d_{\rm max}+3$, and given by the right-hand side of (5.30) (with $p_{I}=\pi$ and $q_{I}=1-\pi$, if $i=0,1,\dots,d_{\rm max}+2$.

It follows that $\bar{\boldsymbol{c}}(\tilde{\tau},u)=(\boldsymbol{0},0,0,1)\bar{B}\bar{\Phi}(\tilde{\tau},u)$ is given by

[TABLE]

where $\boldsymbol{c}(\tilde{\tau},u)$ is defined at (5.31), and (cf. (5.36))

[TABLE]

where $\bar{\Sigma}_{0}=\bar{\Sigma}^{\rm MR}_{0}$ or $\bar{\Sigma}^{\rm NSW}_{0}$, depending on whether the random graph is MR or NSW. Here, $\bar{\Sigma}^{\rm MR}_{0}$ and $\bar{\Sigma}^{\rm NSW}_{0}$ are the asymptotic variance matrices of $n^{-\frac{1}{2}}\bar{\boldsymbol{W}}^{(n)}(0)$ for the MR and NSW random graphs, respectively; cf. (5.15) and (5.19).

We show in Appendix A.4 that

[TABLE]

and, for both the MR and NSW random graphs, that

[TABLE]

Equations (2.13) and (2.14) then follow using Theorem 2.2, noting that $q_{I}^{(2)}=\pi(1-\pi)$.

6 Concluding comments

A shortcoming of our results, from a mathematical though not a practical viewpoint, is the requirement of a maximal degree $d_{\rm max}$. It seems likely that Theorems 2.1-2.3 continue to hold when that requirement is relaxed, subject to appropriate conditions on the degree sequence (MR model) or degree distribution $D$ (NSW model). This conjecture is supported by the recent work of [Barbour and Röllin(2017)], who use Stein’s method to obtain a central limit theorem for local graph statistics in the MR model, allowing for unbounded degree, and as an application obtain a central limit theorem for the size of the giant component with asymptotic variance given by $\tilde{\sigma}^{2}_{\rm MR}$ with $p_{I}=1$. To extend the present proof to models with unbounded degree would require a functional central limit theorem for density dependent population processes with countable state spaces. [Barbour and Luczak(2012)] give such a theorem but it is not applicable in our setting as it would require a finite upper bound on the number of neighbours an individual can infect.

The central limit theorems can be extended, at least in principle, to allow for the infection rate $\lambda$ to depend on the degree of an infective, and also to more general infection processes in which the set of its neighbours that are contacted by a given infective is a symmetric sampling procedure ([Martin-Löf(1986)]). In both cases it is straightforward to determine the limiting deterministic model in Section 5.3.1 and the equation corresponding to (5.22), which governs $\tilde{\tau}$, but calculation of the asymptotic variances in the central limit theorems is likely to be prohibitive.

The configuration model does not display clustering in the limit as $n\to\infty$ and several authors have considered modifications of the configuration model that introduce clustering. In [Trapman(2007)] and [Coupechoux and Lelarge(2014)], in the configuration model construction, for $d=1,2,\dots$, some individuals having $d$ half-edges are replaced by fully connected cliques, each of size $d$, with each member of a clique having exactly one half-edge. The half-edges are then paired up in the usual fashion. In [Gleeson(2009)] and [Ball et al.(2010)], the network is formed as in the configuration model and the population is also partitioned into fully connected cliques. In both models, the set of edges in the network is the union of those in cliques and the paired half-edges. The methodology in this paper can be extended to this general class of models as follows.

As in Section 5.1, the network and epidemic are constructed simultaneously. The objects counted are now fully susceptible cliques, typed by their size and degree composition, and infective and recovered half-edges. When infection is transmitted down a half-edge that half-edge is paired with a uniformly chosen half-edge as before. If it is paired with a susceptible half-edge, then an epidemic is triggered within the corresponding clique and associated half-edges, leading to the creation of further infective and recovered half-edges and, unless the clique epidemic infects the entire clique, a new susceptible clique having reduced size. Central limit theorems for the final size of epidemics on MR and NSW versions of such random graphs should follow using similar arguments to before but again calculation of the asymptotic variances may be difficult. If the infectious period is constant, so the epidemic model is equivalent to bond percolation on the network, the analysis may perhaps be simplified by first splitting the cliques into components determined by bond percolation and then using similar methods to the present paper treating the components as super-individuals.

Appendix A Detailed derivations

A.1 Deterministic solution (5.10)-(5.12)

First note that (5.10) follows immediately from (5.7). Substituting (5.10) into (5.9) yields

[TABLE]

and (5.12) follows using the integrating factor ${\rm e}^{t}$.

Multiplying (5.7) by $i$ and adding over $i=1,2,\dots,d_{\rm max}$ yields

[TABLE]

Recall that $\tilde{\eta}_{E}(t)=\tilde{x}_{E}(t)+\tilde{y}_{E}(t)+\tilde{z}_{E}(t)$. Summing (A.1), (5.8) and (5.9) gives

[TABLE]

and (5.13) follows. Equation (5.11) then follows since $\tilde{y}_{E}(t)=\tilde{\eta}_{E}(t)-\tilde{x}_{E}(t)-\tilde{z}_{E}(t)$ and, using (5.10), $\tilde{x}_{E}(t)=\sum_{i=1}^{d_{\rm max}}i\tilde{x}_{i}(0){\rm e}^{-it}$.

A.2 Properties of $\tilde{\tau}$ and $a(\tilde{\tau})$

In this appendix we show that $\tilde{\tau}\in(0,\infty)$ and $a(\tilde{\tau})<0$, first in the setting of Theorem 2.1 and then in the setting of Theorems 2.2 and 2.3. Let $G(s)=p_{I}f_{D_{\epsilon}}^{(1)}(s)-\mu_{D}(s-q_{I})$ $(0\leq s\leq 1)$. Then, from (5.22), $z={\rm e}^{-\tilde{\tau}}$ satisfies $G(z)=0$. Now $G(0)=p_{I}f_{D_{\epsilon}}^{(1)}(0)+\mu_{D}q_{I}>0$, unless $q_{I}=0$ and $p_{1}-\epsilon_{1}=0$. Also $f_{D_{\epsilon}}^{(1)}(1)<\mu_{D}$, since $p_{i}\epsilon_{i}>0$ for at least one $i>0$, so $G(1)<0$. Thus, under the conditions of Theorem 2.1, $G(s)$ has at least one zero in $(0,1)$. Moreover it has precisely one zero, $z$ say, as $G(s)$ is convex on $[0,1]$, since $G^{(2)}(s)=p_{I}f_{D_{\epsilon}}^{(3)}(s)\geq 0$ for all $s\in[0,1]$. Hence $\tilde{\tau}=-\log z\in(0,\infty)$, as required, and it follows from (5.23) that $a(\tilde{\tau})<0$ if and only if $G^{(1)}(z)<0$. Suppose, for contradiction, that $G^{(1)}(z)\geq 0$. Then, since $G^{(2)}(s)\geq 0$ for all $s\in[z,1]$,

[TABLE]

However, $G(1)<0$, so $G^{(1)}(z)<0$, as required.

Turning to the setting of Theorems 2.2 and 2.3, now let $G(s)=p_{I}f_{D}^{(1)}(s)-\mu_{D}(s-q_{I})$ $(0\leq s\leq 1)$. Then, from (5.49), $z={\rm e}^{-\tilde{\tau}}$ satisfies $G(z)=0$. Now $G(1)=0$ and, under the conditions of Theorems 2.2 and 2.3, $G(0)>0$. Also, using (2.8), $G^{(1)}(1)=p_{I}f_{D}^{(2)}(1)-\mu_{D}>0$, since $R_{0}>1$. Thus, since $G(s)$ is convex on $[0,1]$, there exists a unique $z\in(0,1)$ such that $G(z)=0$, so $\tilde{\tau}\in(0,\infty)$. Moreover, $G^{(2)}(s)>0$ for $s\in[z,1]$ and a similar contradiction argument to before shows that $G^{(1)}(z)<0$, whence $a(\tilde{\tau})<0$.

A.3 Asymptotic variances $\sigma^{2}_{\rm MR}$ and $\sigma^{2}_{\rm NSW}$

In this appendix we derive the expressions for $\sigma^{2}_{\rm MR}$ and $\sigma^{2}_{\rm NSW}$ given in Theorem 2.1. Recalling the partition $\Delta=\Delta_{1}\cup\Delta_{2}\cup\Delta_{3}$ defined just before (5.1), it follows from (5.36) that

[TABLE]

where

[TABLE]

and, for $i=1,2,3$,

[TABLE]

We calculate $\sum_{i=1}^{3}\sigma_{i}^{2}$ in Appendix A.3.1, $\sigma_{0,MR}^{2}$ and then $\sigma^{2}_{\rm MR}$ in Appendix A.3.2, and $\sigma_{0,NSW}^{2}$ and then $\sigma^{2}_{\rm NSW}$ in Appendix A.3.3.

A.3.1 Calculation of $\sum_{i=1}^{3}\sigma_{i}^{2}$

Noting that $\tilde{y}_{E}(t)\geq 0$ and $\tilde{z}_{E}(t)\geq 0$ for $0\leq t\leq\tilde{\tau}$, it follows from (5.3), (5.33) and (5.34) that

[TABLE]

and

[TABLE]

where

[TABLE]

Now

[TABLE]

where, using (5.3),

[TABLE]

Recalling (5.35), it follows using (5.3) that, for $i=1,2,\dots,d_{\rm max}$ and $k=0,1,\dots,i-1$,

[TABLE]

Recall that $p_{i,k}={\rm P}(X=k)$, where $X\sim{\rm Bin}(i-1,1-\exp(-\lambda I))$; see Section 5.1. Elementary calculation yields, for $i=1,2,\dots,d_{\rm max}$, that

[TABLE]

and

[TABLE]

Using (A.10) and (A.11), it follows from (A.9) and some algebra that, for $i=1,2,\dots,d_{\rm max}$

[TABLE]

Recall from (5.16) and (5.20) that, under both the MR and NSW models, $\tilde{x}_{i}(0)=p_{i}-\epsilon_{i}$ $(i=0,1,\dots,d_{\rm max})$, so (5.10) yields

[TABLE]

For $i,k=0,1,\dots$, let $i_{[k]}=i(i-1)\dots(i-k+1)$ denote a falling factorial, with the convention that $i_{[0]}=1$. Then it follows from (A.13) that, for $\theta\in\mathbb{R}$ and $k=1,2,\dots$,

[TABLE]

Thus, for $k=1,2,\dots$,

[TABLE]

Substituting (A.12) into (A.8) and using (A.14)-(A.16) yields

[TABLE]

It then follows using (A.4), (A.3.1) and (A.7) that

[TABLE]

where

[TABLE]

with

[TABLE]

Noting that $\tilde{\eta}_{E}(0)=\mu_{D}$, it follows using (5.33) and (5.13) that

[TABLE]

Recall from (5.16) and (5.20) that under both the MR and NSW models, $\tilde{x}_{E}(0)=\sum_{i=1}^{d_{\rm max}}i(p_{i}-\epsilon_{i})$ and $\tilde{z}_{E}(0)=q_{I}\sum_{i=1}^{d_{\rm max}}i\epsilon_{i}$. It then follows from (5.12) and (A.13) that

[TABLE]

so

[TABLE]

Also, setting $k=1$ in (A.14) yields $\tilde{x}_{E}(u)={\rm e}^{-u}f_{D_{\epsilon}}^{(1)}({\rm e}^{-u})$. Thus,

[TABLE]

The first integral in (A.19) involves only exponential functions and is easily evaluated. Using (A.6), the integrand in the second integral in (A.19) can be expressed as

[TABLE]

so that integral is also easily evaluated. Hence, omitting the details,

[TABLE]

Turning to $I_{5}$, note that

[TABLE]

so

[TABLE]

The integrals $I_{6}$ and $I_{7}$ are easily evaluated yielding

[TABLE]

To evaluate $I_{8}$ to $I_{10}$, let $J_{k}=\int_{0}^{\tilde{\tau}}{\rm e}^{-k(2\tilde{\tau}-u)}f_{D_{\epsilon}}^{(k)}({\rm e}^{-(2\tilde{\tau}-u)})\,{\rm d}u$ $(k=1,2,\dots)$. Then a simple reduction formula yields

[TABLE]

for $k=2,3,\dots$, and

[TABLE]

Applying these formulae to $I_{8}$ to $I_{10}$ yields after some algebra that

[TABLE]

Letting $z={\rm e}^{-\tilde{\tau}}$, it follows from (A.17) and the above equations that

[TABLE]

A.3.2 Calculation of $\sigma_{0,MR}^{2}$ and $\sigma^{2}_{\rm MR}$

To determine $\sigma_{0,MR}^{2}$, note from (A.3), (5.17) and (5.31) that

[TABLE]

using (5.33), (5.34) and $z={\rm e}^{-\tilde{\tau}}$. Now $\sigma_{Y}^{2}$ is given by (5.14), where, for $i=1,2,\dots,d_{\rm max}$, $\sigma_{Y,i}^{2}={\rm var}(Y_{i1})$ with $Y_{i1}\sim{\rm Bin}(i,1-\exp(-\lambda I))$. A simple calculation, conditioning on $I$, yields

[TABLE]

Now $\sum_{i=1}^{d_{\rm max}}i\epsilon_{i}=\mu_{D}-f_{D_{\epsilon}}^{(1)}(1)$ and $\sum_{i=1}^{d_{\rm max}}i(i-1)\epsilon_{i}=f_{D}^{(2)}(1)-f_{D_{\epsilon}}^{(2)}(1)$, so

[TABLE]

Using (A.2), adding (A.24) and (A.3.2), after substituting from (A.27), yields

[TABLE]

Using (5.22) and recalling that $z={\rm e}^{-\tilde{\tau}}$ yields

[TABLE]

Also, using (5.23), (5.29) and (A.14), with $k=1$, gives

[TABLE]

so

[TABLE]

Substituting (A.29) and (A.31) into (A.28), recalling (2.3) and noting from (A.29) and (A.30) that $h(z)=b(\tilde{\tau})z$ yields (2.4) after some algebra.

A.3.3 Calculation of $\sigma_{0,NSW}^{2}$ and $\sigma^{2}_{\rm NSW}$

To determine $\sigma_{0,NSW}^{2}$, first note that (A.3), (5.21) and (5.31) imply that

[TABLE]

where

[TABLE]

and

[TABLE]

Recalling (5.18) and (5.32),

[TABLE]

where from (5.35) and recalling that $z={\rm e}^{-\tilde{\tau}}$,

[TABLE]

Thus,

[TABLE]

and

[TABLE]

Note that in the NSW model

[TABLE]

Hence, using (A.29),

[TABLE]

and, using (A.29) and (A.31),

[TABLE]

It then follows using (A.33) and (A.35) that

[TABLE]

Recall the definition of $(Y_{E},Z_{E})$ just before (5.18). Now ${\rm E}[Y_{E}|D]=p_{I}D$ and, using (A.26), ${\rm var}(Y_{E}|D)=D(D-1)q_{I}^{(2)}+Dq_{I}-D^{2}q_{I}^{2}$, so

[TABLE]

Similar arguments show that

[TABLE]

and

[TABLE]

Thus, using (A.34), (5.33), (5.34) and $z={\rm e}^{-\tilde{\tau}}$,

[TABLE]

Noting from (A.36) that $f_{D_{\epsilon}}^{(1)}(1)=(1-\epsilon)\mu_{D}$ and $f_{D_{\epsilon}}^{(2)}(1)=(1-\epsilon)f_{D}^{(2)}(1)$, it follows using (A.3.2) that

[TABLE]

Note from (A.2) that $\sigma^{2}_{\rm NSW}-\sigma^{2}_{\rm MR}=\sigma_{0,NSW}^{2}-\sigma_{0,MR}^{2}$. Hence, using (A.32), (A.37) and (A.41),

[TABLE]

Recall that $h(z)=b(\tilde{\tau})z$. Note from (2.3) and (A.36) that

$f_{D}(z)=(1-\epsilon-\rho)/(1-\epsilon)$ and also that $f_{D_{\epsilon}}^{(2)}(z)=(1-\epsilon)f_{D}^{(2)}(z)$. Adding (2.4) and (A.42) then yields (2.6).

A.3.4 Proof of Remark 2.5

To prove the assertions in Remark 2.5, suppose that the support of $D$ is not concentrated on a single point. Then it is shown easily using (5.18) that the matrix $\Sigma_{XX}$ is positive definite and it follows immediately from (A.32), (A.33) and (A.41) that $\sigma_{0,NSW}^{2}>\sigma_{0,MR}^{2}$, whence, using (A.2), $\sigma^{2}_{\rm NSW}>\sigma^{2}_{\rm MR}$. Further, setting $\epsilon=0$ yields $\tilde{\sigma}^{2}_{\rm NSW}>\tilde{\sigma}^{2}_{\rm MR}$.

A.4 Asymptotic variances $\sigma^{2}_{\rm MR,S}$ and $\sigma^{2}_{\rm NSW,S}$

For $i=1,2,3$, let

[TABLE]

Note from (5.57) that $\bar{\boldsymbol{c}}(\tilde{\tau},u)(\hat{\boldsymbol{l}}_{+}^{(2)})^{\top}=-\pi\boldsymbol{c}(\tilde{\tau},u)(\boldsymbol{l}_{+}^{(2)})^{\top}$, so $\bar{\sigma}_{2}^{2}=\pi^{2}\sigma_{2}^{2}$ since $\bar{\beta}_{\boldsymbol{l}}(\bar{\boldsymbol{w}}(u))$ is independent of $\bar{v}$ for $\boldsymbol{l}\in\hat{\Delta}_{2}$. Similarly, $\bar{\sigma}_{3}^{2}=\pi^{2}\sigma_{3}^{2}$

Turning to $\bar{\sigma}_{1}^{2}$, note that for $i=1,2,\dots,d_{\rm max}$,

[TABLE]

whence,

[TABLE]

Now $\boldsymbol{l}^{(1)}_{i,0}-\boldsymbol{l}^{(1)}_{i,i-1}=(i-1)(\boldsymbol{e}^{\rm R}-\boldsymbol{e}^{\rm I})$, so, using (5.31),(5.33) and (5.34),

[TABLE]

where $c_{J}(\tilde{\tau},u)=b(\tilde{\tau}){\rm e}^{-(\tilde{\tau}-u)}$ (see (A.6)). Thus, using (5.54), $\bar{\sigma}_{1}^{2}=\int_{0}^{\tilde{\tau}}\bar{g}(u)\,{\rm d}u$, where

[TABLE]

and $\tilde{g}(u)$ is given by (A.8) with $p_{i,i-1}=\pi=1-p_{i,0}$.

Now $\tilde{x}_{i}(0)=p_{i}$ ($i=0,1,\dots,d_{\rm max}$), so (5.10) yields $\sum_{i=1}^{d_{\rm max}}i_{[k]}\tilde{x}_{i}(u)={\rm e}^{-ku}f_{D}^{(k)}({\rm e}^{-u})$ $(k=1,2)$. Hence,

[TABLE]

Thus, recalling that $z={\rm e}^{-\tilde{\tau}}$ and $\rho=1-f_{D}(z)$,

[TABLE]

and

[TABLE]

Further, (2.10) implies $\pi[\mu_{D}-f_{D}^{(1)}(z)]=\mu_{D}(1-z)$, so

[TABLE]

and (5.58) follows since $h(z)=b(\tilde{\tau})z$.

For the MR random graph, setting $\epsilon_{i}=0$ ($i=1,2,\dots,d_{\rm max}$) in (5.14), using (5.17) and recalling that $V^{(n)}(0)=0$ for all $n$, shows that $\hat{\Sigma}_{0}=0$ and (5.59) follows. For the NSW random graph, a similar argument setting $\epsilon=0$ in (5.21) and using (5.57) yields

[TABLE]

A simple calculation using (5.18) (or noting that $\Sigma_{XX}$ is the variance matrix of a single multinomial trial) gives $\boldsymbol{1}\Sigma_{XX}\boldsymbol{1}^{\top}=\boldsymbol{c}(\tilde{\tau},0)\Sigma_{XX}\boldsymbol{1}^{\top}=0$, so (5.59) follows from (A.43).