Past exam of the mathematics course of the University of Cambridge / 2016 / iii / Paper 208 (Time Series and Monte Carlo Inference)

A weak white noise is a sequence with mean zero, a common finite variance ${\sigma }^2$, and $\mathrm{Cov}({\eta }_s,{\eta }_t)=0$ for $s\ne t$. A strong white noise is a sequence of independent and identically distributed random variables with mean zero and finite variance. Strong white noise is therefore weak white noise, but the converse need not hold. Neither definition, by itself, requires a Gaussian distribution.

For the spectral questions take the unique weakly stationary process solving the equation, as is customary for a stable autoregressive moving-average model. The equation alone also admits nonstationary solutions differing by $C{4}^{-t}$, for which a spectral density of a stationary process need not exist. Let $B$ be the backshift operator. Since the autoregressive root is $4$, the stationary solution is causal andUsing the convention $\gamma (h)={\int }_{-\pi }^{\pi }{e}^{ih\lambda }f(\lambda )d\lambda$, the spectral density of a stationary process isThe stationary mean is zero, since $\mu -\mu /4=0$.

Here the innovation process consists of linear innovations, the linear one-step prediction errors: ${\epsilon }_t=X_t-{\mathrm{proj}}_{{\mathcal{H}}_{t-1}}{X}_t$, where ${\mathcal{H}}_{t-1}$ is the closed linear span of the past $X$'s in $L^2$. For Gaussian processes this is also the conditional expectation prediction error. For general non-Gaussian strong white noise these two notions can differ.

The given ${\eta }_t$ is not the linear innovation process. The moving-average factor has its zero at $1/3$, inside the unit disk, and is noninvertible as a causal moving-average filter. The identityis the moving-average root reflection that places this zero outside the unit disk. Thus the causal invertible representation has innovation variance $9$, rather than the given variance $1$.

For an explicit verification, defineThe filter has constant squared modulus $9$, so $\epsilon$ is weak white noise with variance $9$. The new moving-average factor has root $3$ and is invertible, while its autoregressive factor is causal. Hence the past spans of $X$ and $\epsilon$ agree, and $X_t-{\epsilon }_t$ belongs to that past span. Orthogonality of ${\epsilon }_t$ to past $\epsilon$ therefore identifies it as the linear innovation process. The variance difference proves that it cannot be ${\eta }_t$. If the original noise is Gaussian, the new linear innovations are independent Gaussian variables; without Gaussianity they need only be uncorrelated random variables.

The root-reflected, invertible autoregressive moving-average model isThis follows by multiplying the filter identity for $\epsilon$ by $1-B/3$. Its moving-average root $3$ and autoregressive root $4$ are both outside the unit disk. Thus it is the representation in terms of the innovation process, rather than merely another noise representation with the same spectrum.

Expanding the stable autoregressive inverse gives the infinite moving-average representation in linear innovations:ThereforeFor comparison, the causal representation in the originally supplied noise isBoth converge in $L^2$, since their coefficients are square summable. The second is causal, but its driving noise is not the linear innovation process.

Use the original-noise coefficients ${\psi }_0=1$ and ${\psi }_j=-11\cdot 4^{-j}$ for $j\ge 1$. The autocovariance function of a causal linear process is ${\gamma }_X(h)={\sum }_{j\ge 0}{\psi }_j{\psi }_{j+|h|}$, since the original noise has variance $1$. ThusFor $h\ge 1$,By covariance symmetry,In particular, all nonzero-lag covariances are negative, despite the positive autoregressive coefficient.

The printed strong white noise assumption does not imply a normal distribution or even a finite fourth moment. Thus it does not, by itself, determine the covariance of a quadratic transform. We give the intended Gaussian calculation, and then the general finite-fourth-moment answer.

If $\eta$ is Gaussian, $X$ is a centered Gaussian process. Put $v={\gamma }_X(0)=136/15$ and $U_t=X_t/\sqrt{v}$. Since $H_1(u)=u$ and $H_2(u)=u^2-1$, the centered transform isThe supplied Hermite polynomial identity makes the cross terms vanish and givesEquivalently,The constant $a$ has no effect on covariance.

For a general iid noise with $E{\eta }^4<\infty$, put $m_3=E{\eta }^3$ and ${\kappa }_4=E{\eta }^4-3$, its fourth cumulant. Independence and expansion of third and fourth moments give, for $k=|h|$,whereSumming the geometric series explicitly givesThe covariance of quadratic transforms of a linear process follows from a fourth-moment expansion consisting of the three Isserlis theorem pairings, plus the fourth-cumulant contribution when all four noise indices coincide. The third-moment contribution similarly requires three coincident indices. This proves the general formula without assuming a normal distribution. For Gaussian white noise $m_3={\kappa }_4=0$, recovering the simpler answer. If $c\ne 0$ and the noise has infinite fourth moment, $Z_t$ need not have finite variance, so an autocovariance function may not exist.

For a causal autoregressive model of order $p$, the population partial autocorrelation function cuts off after lag $p$, whereas the autocorrelation function generally decays. For an invertible moving-average model of order $q$, the autocorrelation function cuts off after lag $q$, whereas the partial correlation coefficients generally decay. For a mixed autoregressive moving-average model, both generally decay. Exponential decay may alternate in sign or show damped oscillations. Sample correlation coefficients only approximate these patterns, so isolated crossings of the significance bands are not exact order tests.

In PDF Figure 1, the autocorrelation function alternates sign, with a large negative lag-one correlation coefficient and geometrically decreasing magnitude. The partial autocorrelation function has essentially one substantial spike, negative at lag one and about $-0.8$; later values mostly lie within the displayed bands. Thus a stationary AR(1) with a negative coefficient is the natural first model, with $\hat{\varphi }$ initially near $-0.8$. The trace fluctuates around an approximately constant mean and shows no evident deterministic trend.

Fit that model, compare nearby low-order alternatives using likelihood function and an Akaike information criterion or Bayesian information criterion, and inspect the residual autocorrelation function and residual variance. A few small later PACF spikes are expected in a plot with many lags. The roughly $\pm 1.96/\sqrt{T}$ white noise bands are a guide, rather than simultaneous guarantees for every lag.

The sign alternation also suggests concentration of spectral power near the high-frequency end: the AR(1) denominator $1+{\varphi }^2-2\varphi \cos \lambda$ is smallest near $\lambda =\pi$ when $\varphi <0$. From estimated $\gamma (0)$ and lag-one correlation coefficient one can estimate innovation variance as ${\hat{\sigma }}_{\epsilon }^2\approx \hat{\gamma }(0)(1-{\hat{\varphi }}^2)$. The persistence parameter controls the decay rate and long-run variance of a stationary process. These plots do not establish a normal distribution, independence, or the absence of nonlinear dependence; correlation coefficients of squared residuals can provide a separate variance diagnostic.

With angular frequency $\lambda$, the periodogram isAt Fourier frequencies $2\pi j/T$ it is the squared modulus of the normalized discrete Fourier transform. If a nonzero mean is unknown, subtract the sample mean first.

For a zero-mean stationary process, expanding the square givesUnder absolute summability of the covariances this converges to the spectral density of a stationary process. Thus the periodogram is generally biased at finite $T$, but asymptotically unbiased under this short-memory time series condition.

It is nevertheless not a pointwise estimator with statistical consistency unless it is smoothed. For Gaussian white noise, at a nonzero Fourier frequency other than the Nyquist frequency, the real and imaginary Fourier components are independent normal variables. Exactly,where $f={\sigma }^2/(2\pi )$. The variance does not decrease with $T$. Under usual short-memory time series assumptions the same exponential limit is asymptotic for general processes. Averaging nearby frequencies or using a lag-window estimator reduces variance; a frequency bandwidth tending to zero while $T$ times that bandwidth tends to infinity can give statistical consistency. The raw plot remains useful for detecting strong periodic peaks, but increasing the record length alone does not remove its pointwise noise.

First, the absolute values in the printed limit are an error: the correct long-run variance of a stationary process is the signed sum of autocovariances. Directly,soEach coefficient tends to $1$ and has absolute value at most $1$. Absolute summability and the dominated convergence theorem therefore giveThe sum of absolute values is an upper bound, not the general limit. For a concrete counterexample to the printed assertion, take $X_t={\eta }_t-{\eta }_{t-1}$ with unit-variance iid noise. Then $\gamma (0)=2$, $\gamma (\pm 1)=-1$, and other covariances vanish. The absolute sum is $4$, but ${\overline{X}}_T=({\eta }_T-{\eta }_0)/T$, so $T\mathrm{Var}({\overline{X}}_T)=2/T\to 0$.

For correlated data the central limit theorem can still hold under suitable strong mixing of a stationary process and moment conditions, but the limiting variance is the long-run variance of a stationary process, rather than the one-observation variance. When it is positive,Positive serial dependence usually increases the standard error, while negative dependence can reduce it. The approximate effective sample size of a stationary sample is $T{\gamma }_X(0)/{\sum }_h{\gamma }_X(h)$ when that denominator is positive. A long-memory time series can require a different normalization or a different limit law. Absolute covariance summability alone does not prove a CLT: if $X_t=V{\eta }_t$ with a common independent random scale $V$ taking values $1$ and $2$ with equal probabilities and iid ${\eta }_t$ with the standard normal distribution, the off-diagonal covariances are zero, yet $\sqrt{T}{\overline{X}}_T$ has the nonnormal scale-mixture law $VN(0,1)$. The example is not an ergodic stationary process. If the long-run variance of a stationary process is zero, the usual nondegenerate square-root-$T$ CLT is unavailable.

For the causal AR(1), $|\varphi |<1$ andThe conditional distribution follows because the current innovation is independent of the past.

For a fixed known $\varphi$, condition on the observed $X_1$ and use the $n=T-1$ transitions $t=2,\dots ,T$. Their conditional maximum likelihood criterion is, up to constants,Differentiating in $\mu$ yieldsSince $X_t-\varphi X_{t-1}=(1-\varphi )\mu +{\epsilon }_t$,Thus it is an unbiased estimator, even conditionally on $X_1$, andIt has statistical consistency with mean-square convergence, and the iid-noise strong law of large numbers also gives almost-sure statistical consistency. Its conditional distribution is exactly normal with the displayed mean and variance. If $X_0$ is also observed and all $T$ transitions are used, replace $T-1$ by $T$.

The fixed-$\varphi$ qualification is necessary for the exact finite-sample claims. If $\varphi$ is jointly estimated, conditional likelihood function is linear regression with an intercept: writing ${\bar{X}}_{-}=n^{-1}{\sum }_{t=2}^T{X}_{t-1}$ and ${\bar{X}}_{+}=n^{-1}{\sum }_{t=2}^T{X}_t$, the unconstrained estimators areThis ratio is not generally an unbiased estimator and does not have the preceding finite-sample variance. Under the usual stationary regression conditions it has statistical consistency; its asymptotic variance is ${\sigma }_{\epsilon }^2/(1-\varphi {)}^2$. Indeed, it differs from ${\bar{X}}_{+}$ by $\hat{\varphi }(X_T-X_1)/(n(1-\hat{\varphi }))$, an asymptotically negligible endpoint term. Profiling an unknown innovation variance does not change the fixed-$\varphi$ estimate of $\mu$.

The labels in PDF Figure 2 match as follows:Series $1$ has a wandering level and an overall upward drift. Its empirical autocorrelation function stays close to $1$ and declines very slowly, matching panel $c$; this is characteristic of a trend or a record from a unit-root autoregressive process. Such an empirical ACF is not a population stationary ACF if the underlying model is a random walk.

Series $2$ shows relatively smooth runs of adjacent observations on the same side of its level. Positive short-lag dependence emphasizes low frequencies, matching the decreasing spectrum $d$.

Series $3$ is rougher, with much more rapid fluctuation. Negative short-lag dependence emphasizes high frequencies, matching the increasing spectrum $a$. A negative-coefficient AR(1), for example, has this type of spectrum; the figure does not uniquely identify its exact order or parameters.

Series $4$ has conspicuous seasonality. Its autocorrelation function should alternate between positive peaks at full seasonal periods and negative troughs between them, matching $b$. Here “characteristic function” is used descriptively for an ACF or spectrum, rather than for the probability-theoretic characteristic function of a random variable.

A nonstationary process has a statistical law that changes with the time origin. In the weak sense, this includes a time-dependent mean or variance, or a covariance that depends on the two times separately rather than only on their lag. A strictly stationary process requires invariance of every finite-dimensional distribution under a common time shift.

In the preceding plots, series $1$ suggests a stochastic trend or changing level, and series $4$ suggests a deterministic seasonal mean. These are the intended nonstationary examples. A plot alone cannot prove either conclusion: in particular, a stationary random-phase sinusoid can have a seasonal-looking trace and an oscillatory ACF. The interpretation of series $4$ as nonstationary concerns a seasonal mean tied to calendar time.

Three standard responses are to remove a fitted deterministic trend or seasonal mean; to use regular or seasonal differencing for an appropriate trend or seasonal component; and to stabilize a changing variance by a transformation or explicit seasonal scale model. For example, logarithms or a Box–Cox transformation can address level-dependent variance. The operation should match the source of nonstationarity; differencing a varying variance does not generally make it stationary.

Put $c_t=\cos (2\pi t/S)$ and assume ${\sigma }^2>0$. ThenThus $Y$ generally has a periodic mean, and $Z$ generally has a periodic variance. The edge cases allowed by $S\in {\mathbb{N}}^{\ast }$ matter: $Y$ is stationary when $S=1$, because the mean is then constant, and $Z$ is stationary for $S=1$ or $S=2$. In the latter case multiplying iid centered Gaussian white noise by $(-1{)}^t$ leaves its iid distribution unchanged. For every $S\ge 3$, $c_0^2=1$ but $c_1^2<1$, so $Z$ is nonstationary. For $S>1$, $Y$ is nonstationary because $c_t$ is not constant.

Use the seasonal difference operator ${\Delta }_S=1-B^S$. Since $c_t=c_{t-S}$,This is a stationary moving-average model of iid Gaussian white noise. Its covariance is $2{\sigma }^2$ at lag $0$, $-{\sigma }^2$ at lags $\pm S$, and zero otherwise.

For the variance-modulated process the same operation givesIts variance is $2{\sigma }^2{c}_t^2$ and its covariance at lag $S$ is $-{\sigma }^2{c}_t^2$. Therefore it remains nonstationary for $S\ge 3$, and the operation doubles the marginal variance at nonzero seasons. Seasonal differencing does not remove periodic variance: it removes a periodic mean, but the periodic variance generally remains. Practical alternatives are a periodic model or seasonal variance standardization; at seasons with $c_t=0$ the observations are deterministically zero, so division by $c_t$ is not possible there. For the already stationary special cases $S=1,2$, filtering preserves stationarity.

A periodically correlated process with period $S$ has its mean and two-time covariance unchanged when both times are shifted by $S$. For $Z$, the mean is zero andHence $Z$ is a periodically correlated process. Its covariance period can be smaller than $S$: when $S$ is even, $c_{t+S/2}^2=c_t^2$, so $S/2$ is already a covariance period.

Use the usual normalization that ${\epsilon }_t$ is mean-zero unit-variance white noise orthogonal to the past, and extend $\varphi (\nu ),\sigma (\nu )$ periodically in $\nu$. A nondegenerate causal periodic autoregressive model of order one requires $|{\prod }_{\nu =1}^T\varphi (\nu )|<1$. Write $V_{\nu }={\gamma }_{\nu }(0)$, with season indices understood modulo $T$.

For $h\ge 1$, multiply the recursion by $X_{nT+\nu -h}$ and use orthogonality of the current innovation to the past. The periodic Yule-Walker equations areAt lag zero, expansion of the squared recursion givesIn particular, ${\gamma }_{\nu }(1)=\varphi (\nu )V_{\nu -1}$, so an equivalent variance equation is $V_{\nu }=\varphi (\nu ){\gamma }_{\nu }(1)+\sigma (\nu {)}^2$. Negative lags are obtained from the given symmetry ${\gamma }_{\nu }(-h)={\gamma }_{\nu +h}(h)$.

Replacing covariances by their sample versions givesTake the nonnegative scale and require a positive estimated predecessor variance. Using matched seasonal sample covariance matrices makes the quantity under the square root nonnegative by Cauchy-Schwarz inequality; with unmatched covariance estimates one can instead fit seasonal ordinary least squares and use the nonnegative residual variance. Statistical consistency of the covariance estimates gives statistical consistency of these plug-in estimators. If the driving noise variance were not normalized to $1$, the equations would identify $\sigma (\nu {)}^2\mathrm{Var}({\epsilon }_t)$, not the two scales separately.

Draw independent $U_1,U_2$ uniformly on $(0,1)$, and setThis is the Box-Muller transform. The endpoint $U_1=0$ has probability zero; an implementation should avoid evaluating its logarithm.

For $r\ge 0$, $P(R\le r)=1-e^{-r^2/2}$, so $R$ has density $r{e}^{-r^2/2}$. It is independent of the angle, which has a uniform distribution on $[0,2\pi )$, and their joint density is $(2\pi {)}^{-1}r{e}^{-r^2/2}$. Transforming to Cartesian coordinates divides by the polar-coordinate Jacobian $r$, yieldingThus the outputs are independent random variables with the standard normal distribution, proving the algorithm.

Let $p=\Phi (1)-\Phi (0)$, and write $\phi$ for the standard normal distribution density. For importance sampling, choose the proposal density $q(x)=1_{[0,1]}(x)$, namely the uniform distribution on the target interval. Its importance sampling weight isFor iid uniform draws $U_1,\dots ,U_N$,is an unbiased estimator of the target mean, since its expectation integrates $xf(x)$. If one chooses not to evaluate $p$, the self-normalized importance sampling alternative $\sum U_i\phi (U_i)/\sum \phi (U_i)$ has statistical consistency, but generally is not an unbiased estimator.

Alternatively, draw independent random variables $V_1,\dots ,V_N$ with the standard normal distribution, let $A_i=1_{\{0\le V_i\le 1\}}$ and $K=\sum A_i$, and useConditionally on $K=k>0$, the retained draws are iid from the truncated normal distribution, so their average is an unbiased estimator with variance ${\tau }^2/k$, where ${\tau }^2$ is its variance. If $K=0$ the ratio is undefined; one can continue drawing until a prescribed number of retained observations is obtained.

Here importance sampling uses every draw inside the interval, whereas only $p\approx 0.3413$ of the normal proposals are retained. The advantage can also be quantified. Integration using ${\phi }^{\prime }(x)=-x\phi (x)$ givesFor uniform importance sampling,For the retained-draw ratio with a fixed total proposal budget, its asymptotic variance coefficient is ${\tau }^2/p\approx 0.233347$. Thus this importance choice has substantially smaller error for the same proposal count; it also avoids wasting about two thirds of the standard normal distribution draws. This comparison is specific to this well-matched importance sampling proposal; arbitrary proposal choices do not share this guarantee.

Use the spectral theorem for real symmetric matrices to write $\Sigma =Q\Lambda Q^T$, with $Q$ orthogonal and $\Lambda =\mathrm{diag}({\lambda }_1,\dots ,{\lambda }_d)$, ${\lambda }_i\ge 0$. Set $L=Q{\Lambda }^{1/2}$ and generateThen $L{L}^T=\Sigma$. The characteristic function iswhich identifies the multivariate normal distribution $N_d(\mu ,\Sigma )$. Equivalently, every linear combination is normal, with the required mean and variance. A Cholesky decomposition can be used for a positive-definite matrix. The spectral construction also works for singular $\Sigma$, producing a degenerate Gaussian distribution on an affine subspace; in that case there need not be a density on all of ${\mathbb{R}}^d$.

Use $v={\sigma }^2$ as the variance parameter, with prior density relative to $d\mu dv$. Let $\bar{x}=n^{-1}\sum x_i$ and $S=\sum (x_i-\bar{x}{)}^2$. Multiplying the prior by the likelihood function gives the unnormalized Bayesian posteriorThe printed request for an “improper posterior” should not be read as a claim that this posterior is always improper. Although the prior is improper, integrating over $\mu$ gives a marginal kernel $v^{-(n+1)/2}{e}^{-S/(2v)}$, which is integrable exactly when $n>1$ and $S>0$. Under these usual nondegenerate-data conditions it is a proper posterior. If those conditions fail, the displayed kernel is improper and there is no posterior probability distribution for a sampler to target.

Completing the square and using the inverse-gamma distribution give the full conditionalsThe inverse-gamma distribution convention is density proportional to $v^{-\alpha -1}{e}^{-\beta /v}$. Equivalently the precision satisfiesusing shape and rate, as in the question.

A Gibbs sampler starts from $v^{(0)}>0$ and alternatesThis is the normal mean-variance posterior with a log-uniform variance prior. For a marginal check, $v\mid x\sim \mathrm{IG}((n-1)/2,S/2)$: integrating out $\mu$ reduces the shape by $1/2$. This marginal shape is different from the full-conditional shape $n/2$.

Let $Y_i$ be iid with density $f$, set $H_i=h(Y_i)$, and put $I=E_{f}h(Y)$ and ${\tau }^2={\mathrm{Var}}_f(h(Y))$. The finite-second-moment assumption gives, for ${\tau }^2>0$,The sample variance $s_N^2=(N-1{)}^{-1}{\sum }_i(H_i-{\hat{I}}_N{)}^2$ consistently estimates ${\tau }^2$. The central limit theorem and Slutsky theorem therefore give the asymptotic confidence intervalIts confidence probability is approximately $1-\alpha$, conventionally written $100(1-\alpha )\%$; the printed $(1-\alpha )\%$ omits the factor $100$. For ${\tau }^2=0$, $h(Y)=I$ almost surely and the estimator has zero error. This interval concerns iid Monte Carlo estimators; correlated simulation output requires a long-run variance of a stationary process estimate instead of the iid sample variance.

The support condition for importance sampling is missing from the printed assertion. Assume $M>1$, $E_f|h|<\infty$, and $Mg(z)>f(z)$ for $f$-almost every $z$. For a proposal $Z\sim g$, rejection has conditional probability $1-f(Z)/(Mg(Z))$. Its unconditional probability is $(M-1)/M$, so the conditional density of a rejected draw isThus the rejected-draw weight is precisely $f(z)/r(z)$, andConditionally on any positive rejection count $n$ among $N$ independent proposals, the rejected locations have this iid conditional density. Henceis an unbiased estimator conditionally on $n$, and therefore also conditionally on the event $0<n<N$. This is importance sampling from the rejected-proposal distribution.

The nonstrict assumption $f\le Mg$ alone does not suffice: the rejected distribution may completely miss part of the target support. For example, let $g$ be uniform on $[0,1]$, $f(z)=21_{[0,1/2]}(z)$, $M=2$, and $h=1$. All rejected proposals lie in $(1/2,1]$, where $f=0$, so the displayed estimator is always $0$, whereas $E_{f}h=1$. This remains a counterexample after conditioning on $0<n<N$. A strict envelope, such as choosing $M$ larger than the essential supremum of $f/g$, repairs the missing support condition. Near-equality can still give high-variance weights, even when it remains an unbiased estimator.

On $x>0$, the exponential proposal is the proper density $g(x)=e^{-x}{1}_{\{x>0\}}$; its integral is $1$. The word “improper” in the question is therefore not applicable on its stated positive support.

The target-to-proposal ratio isIts maximum occurs at $x=1$, so the smallest valid envelope constant isOne implementation draws $Z=-\log U_1$ and accepts it when $U_2\le e^{-(Z-1{)}^2/2}$, for independent random variables with the uniform distribution on $(0,1)$. Accepted observations have the positive half-normal density by the usual rejection sampling calculation.

The printed identity needs an expectation on its right-hand side. For fixed $r\ge 1$, condition on $X_2,\dots ,X_r$ and use independence of $X_1$:The expression inside this expectation is generally random and cannot equal the unconditional probability by itself.

For each independent Monte Carlo method replicate, draw $R_j$ from its given distribution. If $R_j=r\ge 1$, draw $r-1$ independent $F$-variables, form $W_j={\sum }_{i=2}^r{X}_{j,i}$ and set $H_j=F(x-W_j)$. For $r=1$ the empty sum is $0$. If $R$ may equal zero, use $H_j=1_{\{x\ge 0\}}$ when $R_j=0$, since $S_0=0$. The law of total expectation givesSince $0\le H_j\le 1$, the estimator has statistical consistency by the strong law of large numbers. It is the conditional Monte Carlo estimator obtained by Rao-Blackwellization of the direct indicator $1_{\{S_R\le x\}}$, and the law of total variance gives $\mathrm{Var}(H_j)\le l(1-l)$.

This version assumes that the cumulative distribution function $F$ can be evaluated. An easy sampler alone does not automatically supply an easy cumulative distribution function evaluation. If only sampling is available, use the direct indicator estimator, or replace each conditional cumulative distribution function by the average of several independent indicators $1_{\{X_{j,1,k}\le x-W_j\}}$; the nested version remains an unbiased estimator. With $L$ such inner draws its variance per outer replicate is $\mathrm{Var}(H_j)+L^{-1}E[H_j(1-H_j)]$.

Let $\tilde{\pi }$ be a nonnegative target kernel and $q(y\mid x)$ a proposal distribution. At state $x$, propose $Y\sim q(\cdot \mid x)$, draw an independent uniform $U$, and move to $Y$ ifOtherwise keep $x$. This is the Metropolis–Hastings algorithm. Start on the positive target support, with the usual zero-ratio conventions. The unknown normalizing constant cancels.

An unnormalized kernel is sufficient, provided $0<\int \tilde{\pi }<\infty$. A genuinely improper target with infinite integral is not a probability density and cannot supply a stationary target probability distribution. The acceptance formula can still be written formally, but cannot be said to generate samples from that nonexistent probability law.

Because rejected proposals leave the state unchanged, the Markov kernel includes an atom:Detailed balance with a probability measure $\pi$ means the measure identityWhere both sides have ordinary densities, this reads $\pi (x)K(x,y)=\pi (y)K(y,x)$. The measure formulation also includes the rejection atom. Integrating it shows that $\pi$ is invariant.