2.1. Algorithm details and remarks

A few remarks are in order to clarify Algorithm 1.

• We abuse notation by writing $\xi_t^i \in u$ or $\hat{\xi}_t^j \in u$ to indicate that $\xi_t^i$ or $\hat{\xi}_t^j$ is in bin u, even though the bins form a partition of the particles, and not (necessarily) a partition of the state space.

• A child is simply a copy of its parent: if $\hat{\xi}_t^j$ is a child of $\xi_t^i$ , then $\hat{\xi}_t^j = \xi_t^i$ . The indices of the children are not important, so specifying the number of children of each parent is enough to define them.

• Since the weights do not change in the mutation step and the the selection step preserves the total weight, the total weight is constant in time: $\omega_t^1 + \cdots + \omega_t^N = 1$ for all $t \ge 0$ . We discuss the importance of this in Remark 1.

• We assume that the bins and particle allocation at time t are included in ${\mathcal F}_t$ , the $\sigma$ -algebra generated by the information from Algorithm 1 just before the tth selection step. We also write $\hat{\mathcal F}_t$ for the $\sigma$ -algebra generated by the information from Algorithm 1 just after the tth selection step. See Section 5.1 for details.

• We assume multinomial sampling in the bins because it leads to simple explicit variance expressions in terms of intrabin variances. In Remark 3 we comment on residual sampling, which performs much better than multinomial resampling and still admits nice variance formulas.

For our ergodic theorem, the bins and particle allocation can be arbitrary. To actually do better than direct Monte Carlo, they must be judiciously chosen. The most common strategy is to define bins based on a carefully constructed partition of state space – particles occupy the same bin when they belong to the same element of the partition – and then allocate children approximately uniformly among these bins. Some knowledge about the underlying problem is needed to choose the bins, but this strategy has had considerable success, as the references in the introduction attest (see [Reference Zuckerman and Chong54] for a mostly current list of application papers). We propose a different strategy in our companion paper [Reference Aristoff and Zuckerman3] that uses our variance analysis below. We summarize that strategy in Section 4 below.

3.1. Unbiased property

We begin with the unbiased property of weighted ensemble. This property was previously noted in [Reference Zhang, Jasnow and Zuckerman52], and proved in a slightly different setting in [Reference Aristoff2].

Theorem 2. (Unbiased property.) For each $t\ge 0$ and all bounded measurable g,

\begin{equation*}{\mathbb E}\biggl[\sum_{i=1}^N \omega_t^i g(\xi_t^i)\biggr] = \int K^t g\,{\mathrm{d}}\nu,\end{equation*}

where $\nu$ is the weighted ensemble initial distribution,

\[\int g\,{\mathrm{d}}\nu \;:\!=\; {\mathbb E}\biggl[\sum_{i=1}^N \omega_0^i g(\xi_0^i)\biggr].\]

We could interpret Theorem 2 as follows. If $(X_t)$ is a Markov chain with kernel K and initial distribution $\nu$ , then

\[ {\mathbb E}\biggl[\sum_{i=1}^N \omega_t^i g(\xi_t^i)\biggr] = {\mathbb E}[g(X_t)]. \]

In this sense, weighted ensemble gives unbiased estimates of the law of the underlying Markov chain.

3.2. Ergodic theorem

To ensure that weighted ensemble is ergodic, the underlying Markov kernel K must be ergodic in some sense. We assume that K is uniformly ergodic [Reference Douc, Moulines and Stoffer27].

Assumption 1. There is $c>0$ , $\lambda \in [0,1)$ and a probability measure $\mu$ such that

\[\|K^t(\xi,\cdot) - \mu(\!\cdot\!)\|_{TV} \le c \lambda^t \quad \textit{for all $\xi$ and all $t\ge 0$.}\]

Here and below, f is a fixed bounded measurable function.

Theorem 3. (Ergodic theorem.) If Assumption 1 holds, then with probability 1,

(6) \begin{equation}\lim_{T \to\infty} \dfrac{1}{T}\sum_{t=0}^{T-1}\sum_{i=1}^N \omega_t^i f(\xi_t^i) =\int f\,{\mathrm{d}}\mu.\end{equation}

Convergence of the mean of the time average in (6), at the same rate as direct Monte Carlo, follows from Assumption 1 and the unbiased property (Theorem 2). For the ergodic theorem to hold, and for weighted ensemble to beat direct Monte Carlo, the variance of the time average should be sufficiently small. Well-behaved variance is not automatic for unbiased methods; see Remark 1 below.

3.3. Variance formulas

Here, we give exact, finite N formulas for the variance of weighted ensemble, based on a martingale decomposition. To get nice concise formulas, we need some notation. Define the intrabin distributions

\begin{equation*}\eta_t^u = \sum_{i\colon \xi_t^i \in u} \dfrac{\omega_t^i}{\omega_t(u)} \delta_{\xi_t^i},\end{equation*}

where $\delta_{\xi}$ is the Dirac delta distribution centered at $\xi$ . Define also

\begin{equation*}h_{t,T} = \sum_{s=0}^{T-t-1} K^sf.\end{equation*}

For a probability measure $\eta$ and bounded measurable function g, define

\begin{equation*}\eta(g) = \int g\,{\mathrm{d}}\eta, {\rm{Var}}_\eta(g) = \eta(g^2)-\eta(g)^2,\end{equation*}

and in particular, let ${\rm{Var}}_K g(\xi) \;:\!=\; {\rm{Var}}_{K(\xi,\cdot)}(g) = Kg^2(\xi) - (Kg)^2(\xi)$ .

Theorem 4. (Variance formulas.) For each time $T >0$ ,

(7) \begin{align} {\rm{Var}}\biggl(\dfrac{1}{T}\sum_{t=0}^{T-1}\sum_{i=1}^N \omega_t^if(\xi_t^i) \biggr) \end{align}

(8) \begin{align} &\!\!\!\!\!\!\!\!\!\!\!= \dfrac{1}{T^2}{\rm{Var}}\biggl(\sum_{i=1}^N \omega_0^i h_{0,T}(\xi_0^i)\biggr) \quad {\rm{(initial \, condition \, variance)}} \end{align}

(9) \begin{align} & \quad\quad + \dfrac{1}{T^2}\sum_{t=0}^{T-2}{\mathbb E}\biggl[\sum_u \dfrac{\omega_t(u)^2}{N_t(u)}{\rm{Var}}_{\eta_t^u}(Kh_{t+1,T})\biggr]\quad \text{(selection variance)} \end{align}

(10) \begin{align} & \quad\quad + \dfrac{1}{T^2}\sum_{t=0}^{T-2} {\mathbb E}\biggl[\sum_u\dfrac{\omega_t(u)^2}{N_t(u)}\eta_t^u({\rm{Var}}_K h_{t+1,T})\biggr] \quad \text{(mutation variance)}. \end{align}

The expression in (8) can be interpreted as the variance coming from the initial condition, while the expressions in (9) and (10) can be understood as the variances arising from each selection and mutation step, respectively.

Using Theorem 4, the proof of the ergodic theorem is straightforward. Under Assumption 1, we can show that variances of $h_{t,T}$ and $Kh_{t,T}$ are uniformly bounded in t and T. This makes the weighted ensemble variance

\[{\rm{Var}}\biggl(\dfrac{1}{T}\sum_{t=0}^{T-1}\sum_{i=1}^N \omega_t^if(\xi_t^i) \biggr) = {\mathrm{O}}(1/T).\]

The ergodic theorem then follows from standard arguments. Note that the variance expressions (7)–(10) by themselves do not require Assumption 1.

Beyond the ergodic theorem, these variance formulas are interesting in their own right, since they could be used to design binning and particle allocation schemes that minimize the weighted ensemble variance. Indeed, that is what we have done in our companion paper [Reference Aristoff and Zuckerman3] (see also [Reference Aristoff2]). We discuss this more in the next section.

4.1. Example

Below is a simple example illustrating the variance reduction in weighted ensemble, compared to direct Monte Carlo, in the context of our variance analysis above. Consider a Markov chain on three states, with transition matrix

\[K(1,2) = K(2,3) = \delta,\quad K(1,1) = K(2,1) = 1-\delta, \quad K(3,1) = 1,\]

where $\delta>0$ is small, and we take $f(1) = f(2) = 0$ and $f(3) = 1$ .

For weighted ensemble, we will assign particles to the same bin if and only if they occupy the same point in space. Following the usual method in the literature, we allocate an approximately equal number of children to each bin.

In this case the variance of weighted ensemble can be estimated as

(15) \begin{equation} \dfrac{1}{T^2}\sum_{t=0}^{T-2}{\mathbb E}\biggl[\sum_{i=1}^3\dfrac{\omega_t(i)^2}{N_t(i)}{\rm{Var}}_K h_{t+1,T}(i)\biggr]\approx \dfrac{1}{T}\sum_{i=1}^3 \dfrac{\mu(i)^2}{N/3}{\rm{Var}}_K h_{t+1,T}(i)\approx \dfrac{6\delta^3}{NT},\end{equation}

where the first approximation in (15) holds for large enough N and T, and the second approximation uses direct calculations, dropping terms of higher order than $\delta^3$ . The variance of direct Monte Carlo (see Remark 2) can be estimated by

(16) \begin{equation} \dfrac{1}{NT^2}\sum_{t=0}^{T-2}\nu(K^t({\rm{Var}}_Kh_{t+1,T}))\approx \dfrac{1}{NT}\mu({\rm{Var}}_K h_{t+1,T})\approx \dfrac{\delta^2-\delta^3}{NT}\end{equation}

for large T. Figure 1 shows numerical confirmation of these estimates. Note the smaller variance (higher order in $\delta$ ) for weighted ensemble, compared to direct Monte Carlo.

Figure 1. Comparison of weighted ensemble with direct Monte Carlo when $\delta = 0.001$ and $T = 500$ . (a) Average values of $ ({{1}/{T}})\sum_{t=0}^{T-1}\sum_{i=1}^N \omega_t^i f(\xi_t^i) $ versus N, computed from $10^4$ independent trials, for weighted ensemble and direct Monte Carlo. Error bars widths are $\sigma_T/10^2$ , where $\sigma_T^2$ are the empirical variances. (b) Weighted ensemble empirical standard deviation compared with (15). (c) Direct Monte Carlo empirical standard deviation compared with (16).

5.1. Notation

Below, ${\mathcal F}_t$ is the $\sigma$ -algebra generated by the parents and their weights from times $0 \le s \le t$ , the children and their weights from times $0 \le s \le t-1$ , and the bins and particle allocations from times $0 \le s \le t$ . Meanwhile, $\hat{\mathcal F}_t$ is the $\sigma$ -algebra generated by ${\mathcal F}_t$ together with the children and their weights at time t. Throughout, g denotes a bounded measurable function, and c a positive constant whose value can change between different equations. We will use the notation

(17) \begin{equation}{\rm{par}}(\hat{\xi}_t^j) = \xi_t^i \Longleftrightarrow \xi_t^i \text{ is the parent of }\hat{\xi}_t^j.\end{equation}

5.2. One-step means

Lemma 1. For each $i=1,\ldots,N$ and $t \ge 0$ ,

(18) \begin{equation}{\mathbb E}[\omega_{t+1}^i g(\xi_{t+1}^i)\mid \hat{\mathcal F}_t] = \hat{\omega}_t^i Kg(\hat{\xi}_t^i).\end{equation}

At each time $t \ge 0$ , for each bin u,

(19) \begin{equation}{\mathbb E}\biggl[\sum_{i\colon \hat{\xi}_t^i \in u}\hat{\omega}_t^i g(\hat{\xi}_t^i)\biggm|{\mathcal F}_t\biggr] =\sum_{i\colon \xi_t^i \in u} \omega_t^i g(\xi_t^i).\end{equation}

Proof. By (4), ${\mathbb E}[g(\xi_{t+1}^i)\mid \hat{\mathcal F}_t] = Kg(\hat{\xi}_t^i)$ . Now (18) follows from this and (5). Meanwhile, by (2), ${\mathbb E}[N_t^i\mid {\mathcal F}_t] = N_t(u)\omega_t^i/\omega_t(u)$ if $\xi_t^i \in u$ . Thus, by (3) and (17),

\begin{align*}{\mathbb E}\biggl[\sum_{i\colon \hat{\xi}_t^i \in u}\hat{\omega}_t^i g(\hat{\xi}_t^i)\biggm|{\mathcal F}_t\biggr]&= \sum_{i\colon \xi_t^i \in u}{\mathbb E}\biggl[\sum_{j\colon {\rm{par}}(\hat{\xi}_t^j) = \xi_t^i}\hat{\omega}_t^j g(\hat{\xi}_t^j)\biggm|{\mathcal F}_t\biggr] \\&= \sum_{i\colon \xi_t^i \in u}\dfrac{\omega_t(u)}{N_t(u)}g(\xi_t^i){\mathbb E}[N_t^i \mid {\mathcal F}_t] \\&= \sum_{i\colon \xi_t^i \in u}\omega_t^i g(\xi_t^i).\end{align*}

Lemma 2. (One-step means.) For each time $t \ge 0$ ,

(20) \begin{align} &{\mathbb E}\biggl[\sum_{i=1}^N \omega_{t+1}^i g(\xi_{t+1}^i)\biggm|\hat{\mathcal F}_t\biggr] = \sum_{i=1}^N \hat{\omega}_t^i Kg(\hat{\xi}_t^i), \end{align}

(21) \begin{align} &{\mathbb E}\biggl[\sum_{i=1}^N \hat{\omega}_t^i g(\hat{\xi}_t^i)\biggm|{\mathcal F}_t\biggr] = \sum_{i=1}^N {\omega}_t^i g({\xi}_t^i). \end{align}

Proof. This follows immediately from Lemma 1, by summing over the particles in (18) to get (20), and summing over the bins in (19) to get (21).

5.3. Proof of the unbiased property

Proof of Theorem 2. By repeated application of Lemma 2 with the tower property,

(22) \begin{equation}{\mathbb E}\biggl[\sum_{i=1}^N \omega_{t}^i g(\xi_{t}^i)\biggm|{\mathcal F}_0\biggr] = \sum_{i=1}^N \omega_0^i K^t g(\xi_0^i).\end{equation}

Taking expectations in (22) gives the result.

5.4. Doob martingale and variance decomposition

Below, define the Doob martingale

\begin{equation*}D_t = {\mathbb E}\biggl[\sum_{s=0}^{T-1}\sum_{i=1}^N \omega_s^i f(\xi_s^i)\biggm|{\mathcal F}_t\biggr], \quad \hat{D}_t = {\mathbb E}\biggl[\sum_{s=0}^{T-1}\sum_{i=1}^N \omega_s^i f(\xi_s^i)\biggm|\hat{\mathcal F}_t\biggr].\end{equation*}

Proposition 1. For $0 \le t \le T-1$ , we have

(23) \begin{align}D_t = \sum_{s=0}^{t} \sum_{i=1}^N \omega_s^i f(\xi_s^i) + \sum_{i=1}^N \omega_{t}^i Kh_{t+1,T}(\xi_{t}^i), \quad \hat{D}_t = \sum_{s=0}^{t} \sum_{i=1}^N \omega_s^i f(\xi_s^i) + \sum_{i=1}^N \hat{\omega}_{t}^i Kh_{t+1,T}(\hat{\xi}_t^i).\end{align}

Proof. This comes from Lemma 2 by repeated application of the tower property.

Proposition 2. (Martingale variance decomposition.) For each $T>0$ ,

(24) \begin{align}&{\rm{Var}}\biggl(\dfrac{1}{T}\sum_{t=0}^{T-1}\sum_{i=1}^N \omega_t^i f(\xi_t^i)\biggr) \notag \\*&\quad = \underbrace{\dfrac{1}{T^2}{\rm{Var}}(D_0)}_{{initial \, condition \, variance}} +\underbrace{\dfrac{1}{T^2}\sum_{t=0}^{T-2}{\mathbb E}[(\hat{D}_t-D_t)^2]}_{{selection \, variance}} + \underbrace{\dfrac{1}{T^2}\sum_{t=0}^{T-2}{\mathbb E}[(D_{t+1}-\hat{D}_t)^2]}_{{mutation \, variance}}.\end{align}

Proof. It is straightforward to check that all the martingale differences $D_{t+1}-\hat{D}_t$ and $\hat{D}_t-D_t$ are uncorrelated with each other and with $D_0$ . The proof is finished by writing

\[D_{T-1} = \sum_{t=0}^{T-1}\sum_{i=1}^N \omega_t^i f(\xi_t^i)\]

as a telescoping sum of the martingale differences, computing ${\mathbb E}[D_{T-1}^2]$ in terms of the martingale differences, and subtracting ${\mathbb E}[D_{T-1}]^2 = {\mathbb E}[D_0]^2$ from the resulting expression.

In the proof of Theorem 4 below, we will see that the initial condition variance, selection variance, and mutation variance in (24) are the same as those in (8)–(10).

5.5. Proof of the variance formulas

Proof of Theorem 4. Using the formula (23) for the Doob martingale, together with the one-step mean formula (20), the weight update formula (5), and the conditional independence of particle evolution in (4),

(25) \begin{align} {\mathbb E}[(D_{t+1} - \hat{D}_{t})^2\mid \hat{\mathcal F}_{t}] &= {\rm{Var}}\biggl(\sum_{i=1}^N \omega_{t+1}^i h_{t+1,T}(\xi_{t+1}^i)\biggm|\hat{\mathcal F}_t\biggr) \notag \\ &= \sum_{i=1}^N (\hat{\omega}_t^i)^2{\rm{Var}}_{K}h_{t+1,T}(\hat{\xi}_t^i).\end{align}

By (25), the weight update formula (3), and the bin mean formula (19),

\begin{align*} {\mathbb E}[(D_{t+1}-\hat{D}_t)^2] &= {\mathbb E}\biggl[\sum_{u} \dfrac{\omega_t(u)}{N_t(u)} \sum_{i\colon \hat{\xi}_t^i \in u}\hat{\omega}_t^i {\rm{Var}}_{K}h_{t+1,T}(\hat{\xi}_t^i)\biggr] \\[-2pt] &= {\mathbb E}\biggl[\sum_{u} \dfrac{\omega_t(u)}{N_t(u)} \sum_{i\colon {\xi}_t^i \in u}\omega_t^i{\rm{Var}}_{K}h_{t+1,T}(\xi_t^i)\biggr].\end{align*}

In light of Proposition 2, this gives (10). Now we turn to (9). Using the formula (23) for the Doob martingale, the one-step mean formula (21), and the fact that selections in distinct bins are conditionally independent,

(26) \begin{align} {\mathbb E}[(\hat{D}_t - D_t)^2\mid {\mathcal F}_t] &= {\rm{Var}}\biggl(\sum_{i=1}^N \hat{\omega}_{t}^i Kh_{t+1,T}(\hat{\xi}_{t}^i)\biggm|{\mathcal F}_t\biggr)\notag \\[-2pt] &= \sum_u {\rm{Var}}\biggl(\sum_{i\colon \hat{\xi}_t^i\in u} \hat{\omega}_t^i Kh_{t+1,T}(\hat{\xi}_t^i)\biggm|{\mathcal F}_t\biggr).\end{align}

Using (26), the weight update formula (3), and the conditional independence property of multinomial resampling (see e.g. [Reference Douc, Cappe and Moulines26, equation (6)),

\begin{align*} {\mathbb E}[(\hat{D}_t - D_t)^2\mid {\mathcal F}_t] &= \sum_u \omega_t(u)^2 {\rm{Var}}\biggl(\dfrac{1}{N_t(u)}\sum_{i\colon \hat{\xi}_t^i\in u} Kh_{t+1,T}(\hat{\xi}_t^i)\biggm|{\mathcal F}_t\biggr) \\[-2pt] &= \sum_u \dfrac{\omega_t(u)^2}{N_t(u)}{\rm{Var}}_{\eta_t^u}(Kh_{t+1,T}).\end{align*}

Taking expectations in this expression and appealing to Proposition 2 gives (9).

5.6. Remarks on the variance formulas

Remark 1. We briefly comment on the variance for other methods of the selection and mutation type (1). Consider a method with the same mutation step, but a different selection step that is still unbiased in the sense of (21).

In this case the same variance decomposition (24) applies, with selection variance

(27) \begin{equation}\dfrac{1}{T^2}\sum_{t=0}^{T-2}{\mathbb E}[(\hat{D}_t-D_t)^2] =\dfrac{1}{T^2}\sum_{t=0}^{T-2}{\mathbb E}\biggl[ {\rm{Var}}\biggl(\sum_{i=1}^N \hat{\omega}_t^i Kh_{t+1,T}(\hat{\xi}_t^i)\biggm|{\mathcal F}_t\biggr)\biggr].\end{equation}

For a method like weighted ensemble in which the total weight is always 1, $Kh_{t+1,T}$ can be replaced with $Kh_{t+1,T} - (T-t-1)\int f\,{\mathrm{d}}\mu$ in (27) without otherwise changing the equation. Assumption 1 shows that $Kh_{t+1,T} - (T-t-1)\int f\,{\mathrm{d}}\mu$ is uniformly bounded in t and T. As a result, the selection variance (27) is ${\mathrm{O}}(1/T)$ , and the ergodic theorem remains valid.

We observed numerically that if the total weight varies at each time, then the weights tend to approach zero, and the variance (7) is of order T as $T \to \infty$ . Indeed, $Kh_{t+1,T}$ is typically of order T as $T \to \infty$ , which suggests that in this case the selection variance (27) is of the order of

\[\sum_{t=0}^{T-2}{\mathbb E}\biggl[{\rm{Var}}\biggl(\sum_{i=1}^N \hat{w}_t^i\biggm|{\mathcal F}_t\biggr)\biggr]\quad\text{as $T \to \infty$.}\]

Of course, the ergodic theorem fails if the variance (7) goes to infinity as $T \to \infty$ .

Remark 2. Note that direct Monte Carlo – which we define as independent, equally weighted particles – is a special case of weighted ensemble in which each particle always has weight $1/N$ , every parent always gets its own bin u, and $N_t(u)$ is always 1. In this case the selection variance is zero, while the mutation variance is

(28) \begin{equation}\dfrac{1}{NT^2} \sum_{t=0}^{T-2}{\mathbb E}\biggl[\dfrac{1}{N}\sum_{i=1}^N {\rm{Var}}_Kh_{t+1,T}(\xi_t^i)\biggr] = \dfrac{1}{NT^2} \sum_{t=0}^{T-2}\nu(K^t({\rm{Var}}_Kh_{t+1,T})).\end{equation}

Remark 3. If the selections in the bins use residual multinomial resampling [Reference Douc, Cappe and Moulines26] instead of multinomial resampling, then the selection variance (9) becomes

\begin{equation*}\dfrac{1}{T^2}\sum_{t=0}^{T-2} {\mathbb E}\biggl[\sum_u\dfrac{\omega_t(u)^2}{N_t(u)^2}r_t(u){\rm{Var}}_{\gamma_t^u} (Kh_{t+1,T})\biggr],\end{equation*}

where $r_t(u)$ and $\gamma_t^u$ are defined from the residuals

\[r_t^i = \dfrac{N_t(u)\omega_t^i}{\omega_t(u)} - \biggl\lfloor \dfrac{N_t(u)\omega_t^i}{\omega_t(u)}\biggr\rfloor\]

by

\[r_t(u) = \sum_{i\colon \xi_t^i \in u} r_t^i, \quad \gamma_t^u = \sum_{i\colon \xi_t^i\in u} \dfrac{r_t^i}{r_t(u)}\delta_{\xi_t^i}.\]

We omit proof, but include this formula in case it is useful for optimizations, since residual multinomial resampling performs much better than multinomial resampling yet still admits simple explicit variance expressions.

5.7. Proof of the ergodic theorem

Lemma 3. If Assumption 1 holds, then as $T \to \infty$ ,

\begin{equation*}{\rm{Var}}\biggl(\dfrac{1}{T}\sum_{t=0}^{T-1}\sum_{i=1}^N \omega_t^i f(\xi_t^i)\biggr) = {\mathrm{O}}(1/T).\end{equation*}

Proof. By Assumption 1, we have

(29) \begin{equation}|K^tf(x)-K^tf(y)| \le c\lambda^t \quad \text{for all {x},{y} and all $t\ge 0$,}\end{equation}

where now $c>0$ is a different constant. Thus

\begin{equation*}|h_{t+1,T}(x)-h_{t+1,T}(y)| \le \sum_{s=0}^{T-t-2}|K^sf(x)-K^sf(y)| \le \dfrac{c}{1-\lambda} \;:\!=\; C.\end{equation*}

This shows that ${\rm{Var}}_\eta h_{t+1,T} \le C^2$ , and similarly ${\rm{Var}}_\eta Kh_{t+1,T} \le C^2$ , for any probability distribution $\eta$ . As a result, the selection and mutation variances (9) and (10) are both ${\mathrm{O}}(1/T)$ as $T \to \infty$ . By similar arguments,

\[{\rm{Var}}\biggl(\sum_{i=1}^N \omega_0^i h_{0,T}(\xi_0^i)\biggr) \le C^2,\]

which makes the initialization variance (8) ${\mathrm{O}}(1/T^2)$ . Thus the variance in (7) is ${\mathrm{O}}(1/T)$ .

Proof of Theorem 3 . Define

\[ \theta_T = \dfrac{1}{T}\sum_{t=0}^{T-1}\sum_{i=1}^N \omega_t^i f(\xi_t^i). \]

Then, for $0 < S \le T$ ,

(30) \begin{equation}|\theta_T - \theta_S| = \biggl|\biggl(\dfrac{S}{T}-1\biggr)\theta_S + \dfrac{1}{T} \sum_{t=S}^{T-1}\sum_{i=1}^N\omega_t^i f(\xi_t^i)\biggr| \\\le 2\sup|f|\biggl(1-\dfrac{S}{T}\biggr).\end{equation}

By Lemma 3, ${\rm{Var}}(\theta_T) = {\mathrm{O}}(1/T)$ . So by Chebyshev’s inequality, there is $c>0$ so that

\begin{equation*}{\mathbb P}(|\theta_T - {\mathbb E}[\theta_T]| \ge cT^{1/3-1/2}) \le \dfrac{1}{T^{2/3}}\end{equation*}

for large enough T. With $T_n = n^2$ , by the Borel–Cantelli lemma, there is $n_0$ such that

(31) \begin{equation}|\theta_{T_n} - {\mathbb E}[\theta_{T_n}]| < cT_n^{1/3-1/2} \quad\text{a.s. for all $ n \ge n_0$.}\end{equation}

By the unbiased property and Assumption 1,

(32) \begin{equation}\lim_{T \to \infty} {\mathbb E}[\theta_T] = \int f\,{\mathrm{d}}\mu.\end{equation}

Now given $S>0$ , we can choose $T_n$ so that $T_n \le S \le T_{n+1}$ and write

(33) \begin{equation}|\theta_S - \int f\,{\mathrm{d}}\mu| \le |\theta_S - \theta_{T_n}| + |\theta_{T_n} - {\mathbb E}[\theta_{T_n}]|+ |{\mathbb E}[\theta_{T_n}] - \int f\,{\mathrm{d}}\mu|.\end{equation}

By (30)–(32), with probability 1, the right-hand side of (33) vanishes as $S \to \infty$ .