# The Importance Markov Chain

Charly Andral  CEREMADE, CNRS, UMR 7534, Université Paris-Dauphine, PSL
University, 75016, Paris, France  Randal Douc  SAMOVAR, Telecom Sudparis,
Institut Polytechnique de Paris, 9 rue Charles Fourier, 91820, Evry, France
Hugo Marival  SAMOVAR, Telecom Sudparis, Institut Polytechnique de Paris, 9
rue Charles Fourier, 91820, Evry, France  Christian P. Robert  CEREMADE, CNRS,
UMR 7534, Université Paris-Dauphine, PSL University, 75016, Paris, France
Department of Statistics, University of Warwick, Coventry, CV4 7AL, UK

###### Abstract

The Importance Markov chain is a novel algorithm bridging the gap between
rejection sampling and importance sampling, moving from one to the other
through a tuning parameter. Based on a modified sample of an instrumental
Markov chain targeting an instrumental distribution (typically via a MCMC
kernel), the Importance Markov chain produces an extended Markov chain where
the marginal distribution of the first component converges to the target
distribution. For example, when targeting a multimodal distribution, the
instrumental distribution can be chosen as a tempered version of the target
which allows the algorithm to explore its modes more efficiently. We obtain a
Law of Large Numbers and a Central Limit Theorem as well as geometric
ergodicity for this extended kernel under mild assumptions on the instrumental
kernel. Computationally, the algorithm is easy to implement and preexisting
librairies can be used to sample from the instrumental distribution.

Keywords— Markov chain Monte Carlo,importance sampling, Monte Carlo methods,
ergodicity, regeneration

##  1 Introduction

In Monte Carlo methods [1] and in particular in computational Bayesian
statistics, sampling is used to construct estimates for quantities depending
on problem-specific distributions. As a first approach, one can simulate
independently according to another distribution, called the instrumental
distribution, and use this sample to build an estimate of the quantity of
interest (see, e.g., [2] for a general introduction to Monte Carlo methods).
The most well-known example is the importance sampling (IS) technique [3],
which produces a weighted sample to approximate $\pi(f)=\int
f(x)\pi(x)\mathrm{d}x$ where $\pi$ is a given distribution (by an abuse of
notation, we also denote $\pi$ its density with respect to a dominating
measure $\mathrm{d}x$). Importance sampling is based on rewriting the quantity
of interest as $\pi(f)=\int
f(x)\frac{\pi(x)}{\tilde{\pi}(x)}\tilde{\pi}(x)\mathrm{d}x$ for any density
$\tilde{\pi}$ that dominates $\pi$. Then, $\pi(f)$ can be estimated by
sampling independently $X_{1},X_{2},\cdots$ from the instrumental distribution
$\tilde{\pi}$ and by returning the estimate
$\tilde{I}=n^{-1}\sum_{i=1}^{n}\frac{\pi(X_{i})}{\tilde{\pi}(X_{i})}f(X_{i})$.
It is fundamental to recall here that importance sampling does not deliver a
sample distributed from $\pi$. In contrast, rejection sampling allows to
construct a perfect sample according to $\pi$ but at the cost that a portion
of the sampled points are rejected. To be more specific, if we assume that
$\pi\leqslant M\tilde{\pi}$ for some constant $M$, then we sample
independently $X_{1},X_{2},\cdots\sim_{iid}\tilde{\pi}$ and
$U_{1},U_{2},\cdots\sim_{iid}\mathcal{U}(0,1)$ until the condition
$U_{i}<\frac{\pi(X_{i})}{M\tilde{\pi}(X_{i})}$ is met. For the exit index $i$,
setting $Y=X_{i}$, it turns out that the law of the accepted candidate $Y$ is
then exactly $\pi$ [4].

Another approach is proposed by Markov Chain Monte Carlo (MCMC) methods:
instead of constructing an independent and identically distributed (iid)
sample, an MCMC algorithm provides a Markov chain (thus a dependent sample),
that converges to the distribution of interest. The most common MCMC algorithm
is the Metropolis-Hastings algorithm [5, 6]. Note that MCMC and IS are not
incompatible, and the idea of using a Markov chain for the instrumental
distribution appeared as soon as 1963 [7], and is mentioned by Hasting in [6].
More recently, many algorithms combine IS and MCMC [8, 9, 10].

The Importance Markov Chain (IMC) algorithm uses those ideas in a novel way.
Indeed, while most MCMC algorithms try to adjust the proposals to explore the
support of the target distribution efficiently, the IMC algorithm allows to
target a more friendly instrumental distribution which is then transformed
into the initial target with IS. More specifically, the instrumental Markov
chain is transformed into an augmented Markov chain targetting the
distribution of interest on its first marginal. This is different from classic
subsampling or thinning of the chain that preserve the distribution [11, 12,
13]. The instrumental distribution is considered as a given. Indeed, our aim
in this paper is to establish properties that are preserved by our
transformation for a given instrumental distribution—namely a law of large
numbers (LLN), a Central Limit Theorem (CLT) and geometric ergodicity.

Of course, adding a resampling step to a classical importance sampling based
on a $\tilde{\pi}$-sample $(\tilde{X}_{1},\ldots,\tilde{X}_{n})$ may lead to a
random variable with distribution $\hat{\pi}_{n}$ close to the target
distribution $\pi$. But the total variation norm between the two distributions
$\hat{\pi}_{n}$ and $\pi$ is typically of order $O(1/n)$ whereas our
Importance Markov chain, under mild assumptions, is geometrically ergodic,
showing that the decrease in the total variation norm may be geometrically
fast with respect to $n$.

The Importance Markov chain, in the specific setting with independent
proposals, is related to previous works on Self Regenerative Markov chains
[14, 15] and on Dynamic Weighting Monte Carlo [16, 17]. It was developed in
the dependent case in [18] but the framework there was restrained to a semi-
Markov formulation.

The article is organized as follows:

  1. 1.

We first define the Rejection Markov chain, a generalization of the rejection
sampling in a context of MCMC sampling. This part allows us to define the
rejection kernel used further on, and provides some intuition for the
algorithm in the specific case where there exists a known constant $M$ such
that the density ratio $\frac{\pi}{M\tilde{\pi}}$ is upper-bounded by $1$.

  2. 2.

We then generalize the Rejection Markov chain using repetitions to allow the
density ratio $\frac{\pi}{M\tilde{\pi}}$ to be greater than $1$, thereby
relieving the constraint of the first part. The idea is similar to IS as the
number of repetitions is proportional to the density ratio, to the exceptions
that: (1) the number of repetitions is a random integer and the constraint is
simply that its expectation is proportional to the density ratio; (2) the
output is a true random sample and not a weighted one as in classical IS. We
use an extended space to construct an augmented Markov chain composed of
repetitions of the instrumental chain as its first component and an integer as
its second, keeping track of the remaining number of remaining repetitions. We
then proceed to establish some theoretical properties, under mild assumptions,
notably a law of large numbers, a Central Limit Theorem, a geometric
ergodicity property and some uniqueness results.

  3. 3.

Finally, we illustrate the IMC on two synthetic examples. The first is a
multidimensional mixture of Gaussian distributions, using as instrumental
distribution a tempered version of the target, and a NUTS kernel. The second
focuses on an i.i.d. sample from the instrumental distribution, defined by a
normalizing flow trained to approximate a multimodal target up to dimension
25.

##  2 Notations

Let us denote $\mathbb{R}^{+}=[0,\infty)$, $\mathbb{N}=\\{1,2,...\\}$ and
${\mathbb{N}_{0}}=\mathbb{N}\cup\\{0\\}$. We use the standard convention that
$\prod_{k=m}^{n}=1$ if $m>n$. For integers $k\leqslant\ell$, the notation
$[k:\ell]$ stands for the set $\\{k,\ldots,\ell\\}$ and in case where
$k>\ell$, $[k:\ell]$ is the empty set. Moreover,
$u_{k:\ell}=(u_{k},u_{k+1},\ldots,u_{\ell})$ for all $k\leqslant\ell$ and
$u_{k:\infty}=(u_{\ell})_{\ell\geqslant k}$. If a space ${\mathsf{X}}$ is
equipped with a $\sigma$-field $\mathcal{X}$, we denote by
$\mathsf{F}_{+}({\mathsf{X}})$ the set of all nonnegative measurable functions
with respect to $\mathcal{X}$, that is, we make implicit the dependence on the
$\sigma$-field $\mathcal{X}$ in the notation $\mathsf{F}_{+}({\mathsf{X}})$.
Similarly $\mathsf{F}_{b}({\mathsf{X}})$ is the set of all bounded measurable
functions on ${\mathsf{X}}$ and
$\mathsf{F}_{b+}({\mathsf{X}})=\mathsf{F}_{b}({\mathsf{X}})\cap\mathsf{F}_{+}({\mathsf{X}})$.
Moreover, we denote by $\mathsf{M}_{1}({\mathsf{X}})$ the set of probability
measures on $({\mathsf{X}},\mathcal{X})$. For a non-negative real number $x$,
we denote the floor function by $\left\lfloor x\right\rfloor$ and the
fractional part by $\left\langle x\right\rangle$ and hence $x=\left\lfloor
x\right\rfloor+\left\langle x\right\rangle$. The positive part of $x$ is
written $(x)^{+}$.

If $P,Q$ are Markov kernels on ${\mathsf{X}}\times\mathcal{X}$, $h$ a
measurable function on ${\mathsf{X}}$, $\nu$ a measure on
$({\mathsf{X}},\mathcal{X})$ we define:

  * •

$Ph(x):=\int_{\mathsf{X}}P(x,\mathrm{d}y)h(y)$, for $x\in{\mathsf{X}}$,

  * •

$\nu P(\mathsf{A}):=\int_{\mathsf{X}}\nu(\mathrm{d}x)P(x,\mathsf{A})$ for
$A\in\mathcal{X}$,

  * •

$PQ(x,\mathsf{A}):=\int_{\mathsf{X}}P(x,\mathrm{d}y)Q(y,\mathsf{A})$ for
$x\in{\mathsf{X}}$ and $\mathsf{A}\in\mathcal{X}$,

  * •

$\nu(h):=\int_{\mathsf{X}}h(x)\nu(\mathrm{d}x)$, also denoted $\nu h$ if the
context is clear.

Furthermore, we simply denote $P^{0}:=I$ and for $k\in\mathbb{N}$,
$P^{k}:=PP^{k-1}=P^{k-1}P$.

##  3 The Rejection Markov chain

###  3.1 Formal definition

Let $({\mathsf{X}},\mathcal{X})$ be a measurable space. For a given Markov
kernel $Q$ on ${\mathsf{X}}\times\mathcal{X}$, we denote by
$\mathbb{P}_{\xi}^{Q}$ the probability measure induced on
$({\mathsf{X}}^{\mathbb{N}_{0}},\mathcal{X}^{\otimes{\mathbb{N}_{0}}})$ by the
Markov kernel $Q$ and the initial distribution $\xi$, and by
$\mathbb{E}^{Q}_{\xi}$ the associated expectation operator. If
$\xi=\delta_{x}$ for some $x\in{\mathsf{X}}$, we simply use
$\mathbb{E}_{x}^{Q}:=\mathbb{E}^{Q}_{\delta_{x}}$. On
$({\mathsf{X}}^{\mathbb{N}_{0}},\mathcal{X}^{\otimes{\mathbb{N}_{0}}})$, we
define $X_{\ell}$ as the projection on $\ell^{th}$-component, i.e., for any
$w=(w_{\ell})_{\ell\in{\mathbb{N}_{0}}}\in{\mathsf{X}}^{\mathbb{N}_{0}}$,
$X_{\ell}(w)=w_{\ell}$, and $\theta$ the shift operator on
${\mathsf{X}}^{\mathbb{N}_{0}}$ such that
$\theta:(w_{0},w_{1},...)\mapsto(w_{1},w_{2},...)$. For any measurable
function $\rho:{\mathsf{X}}\to[0,1]$, the Rejection sub-Markovian kernel $S$
is defined as follows:

$Sh(x)=\sum_{k=1}^{\infty}\mathbb{E}_{x}^{Q}\left[h(X_{k})\rho(X_{k})\prod_{i=1}^{k-1}(1-\rho(X_{i}))\right],$
(1)

where $x\in{\mathsf{X}}$ and $h$ is any nonnegative measurable function on
$({\mathsf{X}},\mathcal{X})$. A transition according to $S$ is obtained by
generating a Markov chain $\\{X_{k}:\,k\in{\mathbb{N}_{0}}\\}$ according to
the kernel $Q$ and by selecting the first accepted candidate among
$\\{X_{k}:\,k\in{\mathbb{N}_{0}}\\}$ with the success probability sequence
$\\{\rho(X_{k}):\,k\in{\mathbb{N}_{0}}\\}$. More precisely, define
${\mathsf{Y}}={\mathsf{X}}\times[0,1]$ and
$\mathcal{Y}=\mathcal{X}\otimes\mathcal{B}([0,1])$ and let $G$ be the Markov
kernel on ${\mathsf{Y}}\times\mathcal{Y}$ such that for all
$y=(x,u)\in{\mathsf{Y}}$ and all $\mathsf{A}\in\mathcal{Y}$,

$G(y,\mathsf{A})=\int_{{\mathsf{Y}}}\mathbf{1}_{\mathsf{A}}(x^{\prime},u^{\prime})Q(x,\mathrm{d}x^{\prime})\mathrm{d}u^{\prime}\,.$
(2)

Therefore if $Y^{\prime}=(X^{\prime},U^{\prime})\sim G(y,\cdot)$, then
$(X^{\prime},U^{\prime})$ are independent and marginally, $X^{\prime}\sim
Q(x,\cdot)$ and $U^{\prime}\sim\mathcal{U}(0,1)$. For the Markov chain
$\\{Y_{k}=(X_{k},U_{k}):\,k\in{\mathbb{N}_{0}}\\}$ with Markov kernel $G$,
define the first return time to the set
$\mathsf{D}=\left\\{y=(x,u)\in{\mathsf{Y}}:\,u\leqslant\rho(x)\right\\}$ by

$\sigma_{\mathsf{D}}=\inf\left\\{k\geqslant
1:\,Y_{k}\in\mathsf{D}\right\\}\,.$

Then,
$Sh(x)=\mathbb{E}_{\delta_{x}\otimes\gamma}^{G}[\mathbf{1}_{\\{\sigma_{\mathsf{D}}<\infty\\}}h(X_{\sigma_{\mathsf{D}}})]$
where $\gamma$ is any probability measure on $[0,1]$, showing on the side that
$S\mathbf{1}_{{\mathsf{X}}}(x)=\mathbb{P}_{\delta_{x}\otimes\gamma}^{G}(\sigma_{\mathsf{D}}<\infty)\leqslant
1$, which indeed implies that the kernel $S$ is sub-Markovian.

######  Proposition 1.

Let $Q$ be a Markov kernel on ${\mathsf{X}}\times\mathcal{X}$ with invariant
probability measure $\mu\in\mathsf{M}_{1}({\mathsf{X}})$ and let
$\rho:{\mathsf{X}}\to[0,1]$ be a measurable function. Provided that
$\mu(\rho)>0$, the probability measure $\nu$ defined by

$\nu(\mathsf{A})=\frac{\int_{\mathsf{A}}\rho(x)\mu(\mathrm{d}x)}{\int_{\mathsf{X}}\rho(x)\mu(\mathrm{d}x)}\,,\quad\mathsf{A}\in\mathcal{X}\,,$
(3)

is invariant with respect to $S$, i.e. $\nu S=\nu$.

###### Proof.

Define $\bar{\rho}:=1-\rho$. For any bounded function
$h\in\mathsf{F}_{+}({\mathsf{X}})$ and any $x\in{\mathsf{X}}$,

$\displaystyle Sh(x)$
$\displaystyle=\mathbb{E}_{x}^{Q}\left[\sum_{k=1}^{\infty}\rho(X_{k})h(X_{k})\prod_{i=1}^{k-1}\bar{\rho}(X_{i})\right]$
$\displaystyle=\mathbb{E}_{x}^{Q}\left[\rho(X_{1})h(X_{1})\right]+\mathbb{E}_{x}^{Q}\left[\bar{\rho}(X_{1})\sum_{\ell=1}^{\infty}\rho(X_{\ell+1})h(X_{\ell+1})\prod_{j=1}^{\ell-1}\bar{\rho}(X_{j+1})\right]$
$\displaystyle=Q(\rho h)(x)+Q(\bar{\rho}Sh)(x)\,.$ (4)

Integrating with respect to $\mu$ yields

$\mu Sh=\mu Q(\rho h)+\mu Q(\bar{\rho}Sh)=\mu(\rho h)+\mu(\bar{\rho}Sh)\,,$

where we used $\mu Q=\mu$ in the last equality. Since $h$ is bounded, $\mu
Sh<\infty$. Retrieving $\mu Sh$ on both sides, we finally obtain $\mu(\rho
Sh)=\mu(\rho h)$. Hence $\nu(Sh)=\nu(h)$. ∎

###  3.2 Application to sampling

Let $\pi\in\mathsf{M}_{1}({\mathsf{X}})$ be the target distribution and denote
by $\tilde{\pi}\in\mathsf{M}_{1}({\mathsf{X}})$ an instrumental distribution.
As for rejection or importance sampling, the goal is to produce a sample
targeting $\pi$ using a sample of $\tilde{\pi}$, here obtained by using a
Markov kernel $Q$. We denote $\\{\tilde{X}_{i}:\,i\in{\mathbb{N}_{0}}\\}$ a
Markov chain with transition kernel $Q$ and make the following hypothesis on
$Q$ :

  1. (H1)

The Markov kernel $Q$ admits $\tilde{\pi}$ as invariant probability measure.

We also need the following domination assumption, which is compulsory for
rejection sampling.

  1. (H${}_{\mbox{\tiny rej}}$)

There exists $M>0$ such that $\pi\leqslant M\tilde{\pi}.$

Then, we can use (3) to define $\rho$ such that $\pi$ is the invariant
probability measure for the Markov kernel $S$. Indeed, if $\mu=\tilde{\pi}$,
we get $\nu=\pi$ in (3) by defining

$\rho\propto\frac{\mathrm{d}\pi}{\mathrm{d}\tilde{\pi}}\,.$

If in addition, we want $\rho$ to take values in $[0,1]$, we may pick

$\rho(x)=\frac{1}{M}\frac{\mathrm{d}\pi}{\mathrm{d}\tilde{\pi}}(x)\,,$ (5)

for $\tilde{\pi}$-almost all $x\in{\mathsf{X}}$. From Proposition 1, we deduce
immediately:

######  Theorem 2.

Assume (H1) and (H${}_{\mbox{\tiny rej}}$) and take $\rho$ as defined in (5).
Then $S$ is $\pi$-invariant.

As expected, the rejection kernel $S$ is related to the classical rejection
sampling method. More precisely, in the rejection sampling, we create samples
$\\{X_{k}:\,k\in{\mathbb{N}_{0}}\\}$ distributed according to $\pi$ by
subsampling among a batch of iid random variables
$\\{\tilde{X}_{k}:\,k\in{\mathbb{N}_{0}}\\}$ distributed according to
$\tilde{\pi}$.

Subsampling is done by accepting each candidate sample $\tilde{X}_{k}$ with
probability $\rho(\tilde{X}_{k})$ where $\rho$ is defined as in (5).
Therefore, rejection sampling consists in applying a transition according to
$S$ to the particular case where $Q(x,\cdot)=\tilde{\pi}(\cdot)$,
$\mu=\tilde{\pi}$, $\rho$ as in (5) and hence $\nu=\pi$.

Our Markov rejection algorithm closely resembles the rejection algorithm
except that we no longer need to subsample from an i.i.d. batch of random
variables exactly distributed according to $\tilde{\pi}$, which can be
restrictive. Instead, we rely on a Markov chain targeting $\tilde{\pi}$, i.e.
generated by a Markov kernel with invariant probability $\tilde{\pi}$, which
can be achieved for example via a Metropolis Hastings algorithm. Note that it
is sufficient to know $\tilde{\pi}$ up to a normalizing constant.

##  4 The Importance Markov chain

As for classical rejection sampling, the Markov chain rejection sampling
suffers from the drawback that (H${}_{\mbox{\tiny rej}}$) may be satisfied
only with a prohibitively large $M$ (or worse, (H${}_{\mbox{\tiny rej}}$) may
not even be satisfied). The sampling is in that case inefficient since the
average acceptance ratio is equal to $1/M$. Actually, $\rho$ can be
interpreted in the following way: given $\tilde{X}_{k-1}$, we draw a new point
$\tilde{X}_{k}\in{\mathsf{X}}$ according to $Q$ and insert
$\tilde{N}_{k}\sim\text{Ber}(\rho(\tilde{X}_{k}))$ replica in current sample.
Therefore $\rho(\tilde{X}_{k})=\mathbb{E}[\tilde{N}_{k}|\tilde{X}_{k}]$. This
can be viewed as the equivalent of the weight of $\tilde{X}_{k}$ in a
importance sampling context.

The idea with $\varrho:{\mathsf{X}}\to\mathbb{R}^{+}$ is similar but we now
offer to replicate $\tilde{X}_{k}$ a random number of times $\tilde{N}_{k}$
where the conditional expectation of the random integer
$\tilde{N}_{k}\in{\mathbb{N}_{0}}$ w.r.t. $\tilde{X}_{k}$ is equal to
$\varrho(\tilde{X}_{k})$. This relates to [15] where the author conditions the
weights of his estimate to be nonnegative integers.

###  4.1 The extended Markov chain

Let $\pi$ and $\tilde{\pi}$ be two probability measures on
$({\mathsf{X}},\mathcal{X})$ and let $\kappa$ be a positive real number.
Assume that $\pi$ is dominated by $\tilde{\pi}$ and let
$\varrho_{\kappa}:{\mathsf{X}}\to\mathbb{R}^{+}$ be a measurable function such
that

$\varrho_{\kappa}(x)=\kappa\frac{\mathrm{d}\pi}{\mathrm{d}\tilde{\pi}}(x)\,,$
(6)

for $\tilde{\pi}$-almost all $x\in{\mathsf{X}}$. For $\kappa=1$ let us simply
denote

$\varrho:=\varrho_{1}.$ (7)

Let $\tilde{R}$ be a Markov kernel on
${\mathsf{X}}\times{\mathcal{P}}({\mathbb{N}_{0}})$ and by an abuse of
notation, let us write $\tilde{R}(x,n)=\tilde{R}(x,\\{n\\})$ for any
$(x,n)\in{\mathsf{X}}\times{\mathbb{N}_{0}}$. The distribution
$\tilde{R}(\tilde{X}_{k},\cdot)$ will be used to draw the number
$\tilde{N}_{k}$ of replications of $\tilde{X}_{k}$ under the unbiasedness
assumption:

  1. (H2)

For all $x\in{\mathsf{X}}$,

$\sum_{n=0}^{\infty}\tilde{R}(x,n)n=\varrho_{\kappa}(x)\,.$

Algorithm 1 Importance Markov chain, semi-Markov version

1:$X=[\ ]$

2:Set an arbitrary $\tilde{X}_{0}$.

3:for $k\leftarrow 1$ to $n$ do

4: Draw $\tilde{X}_{k}\sim Q(\tilde{X}_{k-1},\cdot)$ and
$\tilde{N}_{k}\sim\tilde{R}(\tilde{X}_{k},\cdot)$

5: Append $\tilde{N}_{k}$ replicas of $\tilde{X}_{k}$ to $X$

6:end for

7:output: $X$

As seen below, if the Markov kernel $Q$ is $\tilde{\pi}$-invariant, then the
output sequence $X=(X_{0},X_{1},\ldots)$ of Algorithm 1 targets $\pi$.
However, $X$ is not a Markov chain per se, and in order to study its ergodic
properties, we need add a second component $N$ to the sequence $X$ so that the
augmented sequence $(X,N)$ becomes a Markov chain. This is done by rewriting
Algorithm 1 as Algorithm 2.

Algorithm 2 Importance Markov chain (IMC)

1:$\ell\leftarrow 0$

2:Set an arbitrary $\tilde{X}_{0}$.

3:for $k\leftarrow 1$ to $n$ do

4: Draw $\tilde{X}_{k}\sim Q(\tilde{X}_{k-1},\cdot)$ and
$\tilde{N}_{k}\sim\tilde{R}(\tilde{X}_{k},\cdot)$

5: Set $N_{\ell}=\tilde{N}_{k}$

6: while $N_{\ell}\geqslant 1$ do

7: Set $(X_{\ell},N_{\ell})\leftarrow(\tilde{X}_{k},N_{\ell}-1)$

8: Set $\ell\leftarrow\ell+1$

9: end while

10:end for

Let us describe the transition of the extended Markov chain
$\\{(X_{\ell},N_{\ell}):\,\ell\in{\mathbb{N}_{0}}\\}$. From Algorithm 2, we
can see that $(X_{\ell},N_{\ell})$ is updated according to two different
moves. Either we are already inside the while loop described in lines 6-9 of
Algorithm 2, in which case, $N_{\ell}\geqslant 1$ and the update is simply
$(X_{\ell+1},N_{\ell+1})=(X_{\ell},N_{\ell}-1)$, or we are outside the while
loop, in which case, $N_{\ell}=0$ and $X_{\ell}=\tilde{X}_{k}$ for some
$k\in{\mathbb{N}_{0}}$. Then, the update of $(X_{\ell},N_{\ell})$ happens when
we enter again the while loop, in which case
$(X_{\ell+1},N_{\ell+1})=(\tilde{X}_{T},\tilde{N}_{T}-1)$ where
$T=\inf\left\\{n>k:\,\tilde{N}_{n}\neq 0\right\\}$.

The associated Markov kernel $P$ on
$({\mathsf{X}}\times{\mathbb{N}_{0}})\times(\mathcal{X}\otimes\mathcal{P}({\mathbb{N}_{0}}))$
is then defined by: for all
$h\in\mathsf{F}_{+}({\mathsf{X}}\times{\mathbb{N}_{0}})$,

$Ph(x,n)=\mathbf{1}_{\\{n\geqslant
1\\}}h(x,n-1)+\mathbf{1}_{\\{n=0\\}}\cdot\sum_{k,\ell=1}^{\infty}\mathbb{E}_{x}^{Q}\left[h(X_{k},\ell-1)\tilde{R}(X_{k},\ell)\prod_{i=1}^{k-1}\tilde{R}(X_{i},0)\right].$
(8)

To simplify this expression, let us introduce additional notation. Write
$\rho_{\tilde{R}}(x)=\tilde{R}(x,[1:\infty))$ and let $S$ be the Markov kernel
on ${\mathsf{X}}\times\mathcal{X}$ defined by

$Sf(x)=\sum_{k=1}^{\infty}\mathbb{E}_{x}^{Q}\left[f(X_{k})\rho_{\tilde{R}}(X_{k})\prod_{i=1}^{k-1}(1-\rho_{\tilde{R}}(X_{i}))\right]\,,$
(9)

for $f\in\mathsf{F}_{+}({\mathsf{X}}).$

The kernel $S$ is of the same form as in (1) except that $\rho$ is now
replaced by $\rho_{\tilde{R}}$. Note that by construction, $\rho_{\tilde{R}}$
is $[0,1]$-valued and $\rho_{\tilde{R}}(x)$ can be interpreted as the
probability for an $\tilde{X}_{i}=x$ drawn from $Q$ to be accepted for the
chain $(X_{\ell})$, in which case, we keep at least one replica of
$\tilde{X}_{i}$.

Then, the extended Markov kernel $P$ writes, for
$h\in\mathsf{F}_{+}({\mathsf{X}}\times{\mathbb{N}_{0}})$:

$Ph(x,n)=\mathbf{1}_{\\{n\geqslant
1\\}}h(x,n-1)+\mathbf{1}_{\\{n=0\\}}\sum_{n^{\prime}=0}^{\infty}\int_{{\mathsf{X}}}S(x,\mathrm{d}x^{\prime})R(x^{\prime},n^{\prime})h(x^{\prime},n^{\prime})\,,$
(10)

where $R$ is the Markov kernel on
${\mathsf{X}}\times{\mathcal{P}}({\mathbb{N}_{0}})$ defined by

$R(x,n):=\tilde{R}(x,n+1)/\rho_{\tilde{R}}(x)\,,\quad(x,n)\in{\mathsf{X}}\times{\mathbb{N}_{0}}\,,$
(11)

and where we, again, make the abuse of notation $R(x,n):=R(x,\\{n\\})$.

Note that (8) and (10) give two different but equivalent decompositions of
$P$, for the sampling step (when $n=0$). In (8), we sample $\tilde{X}_{i}$
according to $Q$ and then use $\tilde{R}(\tilde{X}_{i},\cdot)$ to draw a
number of replicas $\tilde{N}_{i}$, until it is larger than $1$, in which case
we retain $\tilde{X}_{i}$ and $\tilde{N}_{i}-1$, the number of remaining
replicas. In (10), we bypass the rejection step by drawing directly a new
accepted point $X_{i}$ using $S$ and then the number of remaining replicas
from $R(X_{i},\cdot)$, which corresponds to the law of $\tilde{N}-1$
conditionnally on $\\{\tilde{N}\geqslant 1\\}$ when
$\tilde{N}\sim\tilde{R}(X_{i},\cdot)$.

######  Remark 1.

The unbiasedness assumption (H2) is closely related to the notion of correctly
weighted density developed in [16, 17]. Write
$\hat{\pi}(\mathrm{d}x\mathrm{d}n)=\tilde{\pi}(\mathrm{d}x)\tilde{R}(x,\mathrm{d}n)$
a joint distribution on ${\mathsf{X}}\times{\mathbb{N}_{0}}$. Then, under
(H2),
$\sum_{n\in\mathbb{N}}\tilde{\pi}(\mathrm{d}x)\tilde{R}(x,n)n=\kappa\pi(\mathrm{d}x)$
so $\hat{\pi}$ is correctly weighted. And by construction, the kernel
$Q(x,\mathrm{d}y)\tilde{R}(y,dn)$ that generates the samples
$(\tilde{X}_{i},\tilde{N}_{i})_{i\in\mathbb{N}}$ of Algorithm 1 admits
$\hat{\pi}$ as an invariant probability distribution (see Lemma 6).

######  Remark 2.

While importance sampling requires exact simulations from $\tilde{\pi}$, the
IMC method only relies on a Markov kernel $Q$ targetting $\tilde{\pi}$. This
allows us to extend the set of usable instrumental distributions.

######  Remark 3.

Perhaps surprisingly, Metropolis-Hastings (MH) algorithms can be cast into the
framework of importance Markov chains. Indeed, take a MH algorithm with
proposition kernel $A(x,dy)$ and acceptance rate $\alpha(x,y)$, targeting
$\pi$. Following the framework of [19], the accepted points
$\\{\tilde{X}_{i}:\,i\in{\mathbb{N}_{0}}\\}$ form a Markov chain with Markov
kernel $Q(x,dy)$ proportional to $\alpha(x,y)A(x,dy)$. Before moving to a new
accepted point, $\tilde{X}_{i}=x$ is repeated a random number of times
$\tilde{N}_{i}$ that follows (conditionally on $\tilde{X}_{i}=x$) a geometric
distribution with success probability
$p(x):=\int_{\mathsf{X}}\alpha(x,y)A(x,dy)$. Then it can be shown that $Q$ is
$\tilde{\pi}$-invariant where $\tilde{\pi}(dx)\propto p(x)\pi(dx)$, hence (H1)
holds. Define $\tilde{R}(x,\cdot)$ as the geometric distribution with
parameter $p(x)$. Choosing
$\kappa=1/\int_{\mathsf{X}}\rho(x)\pi(\mathrm{d}x)$, $\varrho_{\kappa}$
defined in (6) writes $\varrho_{\kappa}=1/p(x)$, hence (H2) holds. Then, with
these choices of $Q$ and $\tilde{R}$, $(\tilde{X}_{i},\tilde{N}_{i})$
corresponds to the output of the IMC algorithm defined in Algorithm 1.

###  4.2 Invariant probability measure

####  4.2.1 Existence

Let $\bar{\pi}$ be the measure on ${\mathsf{X}}\times{\mathbb{N}_{0}}$ defined
by: for any $h\in\mathsf{F}_{+}({\mathsf{X}}\times{\mathbb{N}_{0}})$,

$\displaystyle\bar{\pi}(h)$
$\displaystyle=\kappa^{-1}\sum_{n=1}^{\infty}\int_{{\mathsf{X}}}\tilde{\pi}(\mathrm{d}x)\tilde{R}(x,n)\sum_{k=0}^{n-1}h(x,k)$
$\displaystyle=\kappa^{-1}\sum_{\ell=0}^{\infty}\int_{{\mathsf{X}}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)R(x,\ell)\sum_{k=0}^{\ell}h(x,k)\,,$
(12)

where the last equality follows from (11) and the change of variable
$\ell=n-1$.

######  Proposition 3.

Assume (H1) and (H2). Let $P$ be the Markov kernel defined in (10) and let
$\bar{\pi}$ be the probability measure on ${\mathsf{X}}\times{\mathbb{N}_{0}}$
defined in (4.2.1). Then,

  1. (i)

the Markov kernel $P$ is $\bar{\pi}$-invariant,

  2. (ii)

the marginal of $\bar{\pi}$ on the first component is $\pi$.

###### Proof.

We start with (i). Let
$h\in\mathsf{F}_{+}({\mathsf{X}}\times{\mathbb{N}_{0}})$. Interchanging the
sum in $\ell$ and the sum in $k$ in (4.2.1) yields

$\bar{\pi}(h)=\kappa^{-1}\sum_{k=0}^{\infty}\int_{{\mathsf{X}}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)R(x,[k:\infty))h(x,k)\,.$
(13)

We now replace $h$ by $Ph$ and combine with the expression of $Ph$ given in
(10), we then obtain

$\displaystyle\bar{\pi}(Ph)$
$\displaystyle=\kappa^{-1}\sum_{k=1}^{\infty}\int_{{\mathsf{X}}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)R(x,[k:\infty))h(x,k-1)+\kappa^{-1}\sum_{n^{\prime}=0}^{\infty}\int_{{\mathsf{X}}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)\int_{{\mathsf{X}}}S(x,\mathrm{d}x^{\prime})R(x^{\prime},n^{\prime})h(x^{\prime},n^{\prime})$
$\displaystyle=\kappa^{-1}\sum_{n^{\prime}=0}^{\infty}\int_{{\mathsf{X}}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)R(x,[n^{\prime}+1:\infty))h(x,n^{\prime})+\kappa^{-1}\sum_{n^{\prime}=0}^{\infty}\int_{{\mathsf{X}}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)R(x,n^{\prime})h(x,n^{\prime})\,,$

where the last equality follows (a) from the change of variable
$n^{\prime}=k-1$ for the first term of the rhs and (b) from Proposition 1
applied, under (H1), to $\rho=\rho_{\tilde{R}}$ and $\mu=\tilde{\pi}$ for the
second term. Noting that
$R(x,[n^{\prime}+1:\infty))+R(x,n^{\prime})=R(x,[n^{\prime}:\infty))$, we
finally get

$\bar{\pi}(Ph)=\kappa^{-1}\sum_{n^{\prime}=0}^{\infty}\int_{{\mathsf{X}}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)R(x,[n^{\prime}:\infty))h(x,n^{\prime})=\bar{\pi}(h)\,,$

where (13) is used to obtain the last equality.

We now turn to (ii). For any $\mathsf{A}\in\mathcal{X}$, applying (4.2.1) with
$h(x,k)=\mathbf{1}_{\mathsf{A}}(x)$ yields under (H2)

$\bar{\pi}(\mathsf{A}\times{\mathbb{N}_{0}})=\kappa^{-1}\int_{\mathsf{X}}\tilde{\pi}(\mathrm{d}x)\left(\sum_{n=1}^{\infty}\tilde{R}(x,n)n\right)\mathbf{1}_{\mathsf{A}}(x)=\kappa^{-1}\tilde{\pi}(\rho_{\tilde{R}}\mathbf{1}_{\mathsf{A}})=\pi(\mathsf{A})\,.$

∎

###  4.3 Uniqueness

######  Proposition 4.

Assume (H1) and (H2). If any invariant measure for $Q$ is proportional to
$\tilde{\pi}$ (defined in (H1)), then $\bar{\pi}$ defined in (4.2.1) is the
unique invariant probability measure for $P$.

###### Proof.

See A. ∎

The uniqueness of the invariant probability measure $\bar{\pi}$ for $P$, as
stated in Proposition 4 allows to obtain the Birkhoff ergodic theorem ([20,
Theorem 5.2.9]): for any measurable function
$g:{\mathsf{X}}\times{\mathbb{N}_{0}}\rightarrow\mathbb{R}$ such that
$\bar{\pi}(\left\lvert g\right\rvert)<\infty$,

$\lim_{n\rightarrow\infty}n^{-1}\sum_{k=0}^{n-1}g(X_{k},N_{k})=\bar{\pi}(g),\quad\mathbb{P}_{\bar{\pi}}^{P}-a.s.$

Although reassuring, the law of large numbers holds
$\mathbb{P}_{\bar{\pi}}^{P}-a.s.$, i.e. the initial distribution is set to be
the invariant probability measure $\bar{\pi}$, which is not realistic from a
practical point of view. Consequently, we will now turn to conditions under
which the law of large numbers holds, irrespective to the initial
distribution.

###  4.4 Law of large numbers

To establish a strong law of large numbers for the kernel $P$, we rely on the
single hypothesis that the instrumental kernel $Q$ satisfies a law of large
numbers. More precisely if the instrumental kernel $Q$ satisfies a law of
large numbers for any initial distribution
$\xi\in\mathsf{M}_{1}({\mathsf{X}})$, Theorem 5 will show that it is also the
case for the importance Markov kernel $P$.

  1. (H${}_{\mbox{\tiny lln}}$)

For every $\xi\in\mathsf{M}_{1}({\mathsf{X}})$ and measurable function
$g:{\mathsf{X}}\rightarrow\mathbb{R}$ such that $\tilde{\pi}(\left\lvert
g\right\rvert)<\infty$,

$\lim_{n\rightarrow\infty}n^{-1}\sum_{k=0}^{n-1}g(\tilde{X}_{k})=\tilde{\pi}(g),\quad\mathbb{P}_{\xi}^{Q}-a.s.$

######  Theorem 5.

Assume (H1) and (H${}_{\mbox{\tiny lln}}$). Then, for every
$\xi\in\mathsf{M}_{1}({\mathsf{X}}\times{\mathbb{N}_{0}})$ and measurable
function $g:{\mathsf{X}}\times{\mathbb{N}_{0}}\rightarrow\mathbb{R}$ such that
$\bar{\pi}(\left\lvert g\right\rvert)<\infty$,

$\lim_{n\rightarrow\infty}n^{-1}\sum_{k=0}^{n-1}g(X_{k},N_{k})=\bar{\pi}(g),\quad\mathbb{P}_{\xi}^{P}-a.s.$

###### Proof.

The proof relies on [21, Proposition 3.5], which relates (H${}_{\mbox{\tiny
lln}}$) to a property on the harmonic functions for $Q$ (i.e. measurable
functions $h$ such that $Qh=h$). More precisely, it states that for any Markov
kernel $Q$ satisfying $\tilde{\pi}Q=\tilde{\pi}$ for some
$\tilde{\pi}\in\mathsf{M}_{1}({\mathsf{X}})$, (H${}_{\mbox{\tiny lln}}$) is
equivalent to (H${}_{\mbox{\tiny hrm}}$) defined as follows:

  1. (H${}_{\mbox{\tiny hrm}}$)

Any bounded harmonic function $h:{\mathsf{X}}\rightarrow\mathbb{R}$ for $Q$ is
constant.

Hence, proving Theorem 5 is equivalent to showing that any bounded harmonic
function for $P$ is constant.

Let $\bar{h}:{\mathsf{X}}\times{\mathbb{N}_{0}}\rightarrow\mathbb{R}$ be a
bounded harmonic function for $P$. Then for $n>0$,

$\bar{h}(x,n)=P\bar{h}(x,n)=\bar{h}(x,n-1),$

where the last equality comes from (10). Thus $\bar{h}$ does not depend on its
second argument and we can define a measurable function
$h:{\mathsf{X}}\rightarrow\mathbb{R}$ such that

$h(x)=\bar{h}(x,n)$ (14)

for all $(x,n)\in{\mathsf{X}}\times{\mathbb{N}_{0}}$. Now,

$h(x)=\bar{h}(x,0)=P\bar{h}(x,0)=\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}S(x,dx^{\prime})R(x^{\prime},dn)\bar{h}(x,n)=Sh(x)$

using the expression of $P$ in (10) as well as (14). From the recursive
expression for $S$ in (4) we have:

$\displaystyle h(x)$ $\displaystyle=Sh(x)$ $\displaystyle=Q(\rho
h)(x)+Q((1-\rho)Sh)(x)$ $\displaystyle=Q(\rho h)(x)+Q((1-\rho)h)(x)$
$\displaystyle=Qh(x).$

Therefore $h$ is harmonic for $Q$, and since it is also bounded,
(H${}_{\mbox{\tiny lln}}$) implies that it is constant. Then (14) shows that
$\bar{h}$ is constant which concludes the proof. ∎

###  4.5 Central Limit Theorem

Let us now establish a Central Limit Theorem (CLT) associated to the
Importance Markov chain for a particular function $h$, based on a similar
hypothesis on the instrumental kernel $Q$ with the function $\varrho h$. More
precisely, the Central Limit Theorem for $Q$ stems from the existence of a
solution to the Poisson equation for this kernel. This is a quite common
sufficient condition for a CLT, and although other conditions exist
(references are given for example in [20, Chap 21]), we choose this one for
its simplicity in the proofs. To be more specific, we are interested in
measurable functions $h:{\mathsf{X}}\to\mathbb{R}$ such that if we define

$h_{0}:=h-\pi h,$ (15)

the following condition holds:

  1. (H${}_{\mbox{\tiny Poiss}}$)

The Poisson equation associated to $\varrho h_{0}$ for the kernel $Q$ on
${\mathsf{X}}$ admits a $\tilde{\pi}$-square integrable solution $H$, i.e. for
all $x\in{\mathsf{X}}$,

$H(x)-QH(x)=\varrho h_{0}(x)\text{ and }\tilde{\pi}H^{2}<\infty.$

In addition,
$\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}n^{2}h(x)^{2}\tilde{\pi}(\mathrm{d}x)\tilde{R}(x,\mathrm{d}n)<\infty$.

######  Remark 4.

Note that $\tilde{\pi}(\varrho h_{0})=0$ since
$\varrho=\frac{\mathrm{d}\pi}{\mathrm{d}\tilde{\pi}}$, hence this term does
not appear in the Poisson equation.

Under (H${}_{\mbox{\tiny Poiss}}$), [20, Theorem 21.2.5] ensures that the
Markov chain $(\tilde{X}_{i})$ generated by the kernel $Q$ satisfies a Central
Limit Theorem for the function $\varrho h_{0}$ :

$\frac{1}{\sqrt{n}}\sum_{i=1}^{n}\varrho(\tilde{X}_{i})h_{0}(\tilde{X}_{i})\overset{\mathbb{P}^{Q}_{\tilde{\pi}}-law}{\leadsto}\mathcal{N}(0,\tilde{\sigma}^{2}(\varrho
h_{0})),$ (16)

where $\tilde{\sigma}^{2}(\varrho h_{0})=2\tilde{\pi}\left(\varrho
h_{0}H\right)-\tilde{\pi}((\varrho h_{0})^{2})$. Lemma 11 of A.3 combined with
(H${}_{\mbox{\tiny lln}}$) then extends the weak convergence under
${\mathbb{P}^{Q}_{\tilde{\pi}}}$ in the equation above to a weak convergence
under ${\mathbb{P}^{Q}_{\xi}}$ for any $\xi\in\mathsf{M}_{1}({\mathsf{X}})$.
We can now state the CLT for the Importance Markov chain in a formal manner.

######  Theorem 6.

Assume (H1), (H2), (H${}_{\mbox{\tiny lln}}$) and let
$h:{\mathsf{X}}\to\mathbb{R}$ be a measurable function that satisfies
(H${}_{\mbox{\tiny Poiss}}$). Then there exists a constant $\sigma^{2}(h)>0$
such that

$\frac{1}{\sqrt{n}}\sum_{i=1}^{n}\left(h(X_{i})-\pi
h\right)\overset{\mathbb{P}^{P}_{\chi}-law}{\leadsto}\mathcal{N}(0,\sigma^{2}(h)),$

where the distribution $\chi$ is defined by
$\chi(f)=\int\xi(\mathrm{d}x)S(x,\mathrm{d}x^{\prime})R(x^{\prime},\mathrm{d}n^{\prime})f(x^{\prime},n^{\prime})$.
Moreover, we have the following expression of $\sigma^{2}(h)$:

$\sigma^{2}(h)=\kappa\tilde{\sigma}^{2}(\varrho
h_{0})+\kappa^{-1}\hat{\sigma}^{2}(h_{0},\kappa),$ (17)

where

-

$\tilde{\sigma}^{2}(\varrho h_{0})$ is the variance obtained in (16),

-

$\hat{\sigma}^{2}(h_{0},\kappa):=\int_{\mathsf{X}}h_{0}^{2}(x)\mathbb{V}\mathrm{ar}^{\tilde{R}}_{x}[N]\tilde{\pi}(dx)$,

-

$\mathbb{V}\mathrm{ar}^{\tilde{R}}_{x}[N]:=\int_{\mathbb{N}_{0}}\tilde{R}(x,\mathrm{d}n)n^{2}-\left(\int_{\mathbb{N}_{0}}\tilde{R}(x,\mathrm{d}n)n\right)^{2}$.

###### Proof.

See A.3. ∎

######  Remark 5.

Note that the variance $\sigma^{2}(h)$ can be decomposed into two terms: $(1)$
$\tilde{\sigma}^{2}(\varrho h)$ is the variance coming from the instrumental
chain, while $(2)$ $\hat{\sigma}^{2}(h,\kappa)$ is the variance brought by the
random number of repetitions of the instrumental chain.

###  4.6 Minimizing the asymptotic variance

####  4.6.1 Optimal choice of the kernel $\tilde{R}$

Following Remark 5, one can notice that the expression
$\hat{\sigma}^{2}(h,\kappa)=\int_{\mathsf{X}}h^{2}(x)\mathbb{V}\mathrm{ar}^{\tilde{R}}_{x}[N]\tilde{\pi}(dx)$
directly depends on the variance of $N$ under $\tilde{R}(x,\cdot)$. Therefore,
minimizing the variance associated to $\tilde{R}(x,\cdot)$, for
$x\in{\mathsf{X}}$, leads to minimization of the asymptotic variance of the
chain as defined in Theorem 6. To help tuning $\tilde{R}$, we state the
following lemma:

######  Lemma 1.

Let $N$ be an integer-value random variable on some probability space
$(\Omega,\mathcal{F},\mathbb{P})$ such that $\mathbb{E}[N]=\rho$ for a fixed
$\rho\in\mathbb{R}^{+}$. Then,

$\mathbb{V}\mathrm{ar}(N)\geqslant\left\langle\rho\right\rangle(1-\left\langle\rho\right\rangle).$

This bound is reached for $N=\left\lfloor\rho\right\rfloor+S$, where $S\sim
Ber(\left\langle\rho\right\rangle)$

###### Proof.

Using $0=\mathbb{E}[N]-\rho=\mathbb{E}[(N-\rho)^{+}]-\mathbb{E}[(N-\rho)^{-}]$
and
$(N-\rho)^{2}\geqslant(1-\left\langle\rho\right\rangle)(N-\rho)^{+}+\left\langle\rho\right\rangle(N-\rho)^{-}$,
we get

$\mathbb{E}[(N-\rho)^{2}]\geqslant\mathbb{E}[(N-\rho)^{+}]=\mathbb{E}[(N-\rho)^{-}].$

  * •

If $\mathbb{P}(N>\rho)\geqslant\left\langle\rho\right\rangle$, then

$\mathbb{E}[(N-\rho)^{2}]\geqslant\mathbb{E}[(N-\rho)^{+}]\geqslant(1-\left\langle\rho\right\rangle)P(N>\rho)\geqslant\left\langle\rho\right\rangle(1-\left\langle\rho\right\rangle).$

  * •

If
$\mathbb{P}(N>\rho)<\left\langle\rho\right\rangle\Leftrightarrow\mathbb{P}(N\leqslant\rho)>1-\left\langle\rho\right\rangle$,
then

$\mathbb{E}[(N-\rho)^{2}]\geqslant\mathbb{E}[(N-\rho)^{-}]\geqslant\left\langle\rho\right\rangle
P(N\leqslant\rho)\geqslant\left\langle\rho\right\rangle(1-\left\langle\rho\right\rangle).$

∎

In the case where $\varrho_{\kappa}$ can be computed, we can use Lemma 1 to
define $\tilde{R}$ as:

$\tilde{R}_{\scriptscriptstyle\mathrm{opt}}=(1-\left\langle\varrho_{\kappa}(x)\right\rangle)\delta_{\left\lfloor\varrho_{\kappa}(x)\right\rfloor}+\left\langle\varrho_{\kappa}(x)\right\rangle\delta_{\left\lfloor\varrho_{\kappa}(x)\right\rfloor+1}.$
(18)

This $\tilde{R}_{\scriptscriptstyle\mathrm{opt}}$ implies the following
expression for $R_{\scriptscriptstyle\mathrm{opt}}$:

$R_{\scriptscriptstyle\mathrm{opt}}=(1-\left\langle\varrho_{\kappa}(x)\right\rangle)\delta_{\left\lfloor\varrho_{\kappa}(x)\right\rfloor}+\left\langle\varrho_{\kappa}(x)\right\rangle\delta_{(\left\lfloor\varrho_{\kappa}(x)\right\rfloor-1)^{+}}.$

####  4.6.2 Optimal upper bound

With the optimal choice of $\tilde{R}$ given in (18), we have
$\mathbb{V}\mathrm{ar}^{\tilde{R}_{\scriptscriptstyle\mathrm{opt}}}_{x}[N]=\left\langle\varrho_{\kappa}(x)\right\rangle(1-\left\langle\varrho_{\kappa}(x)\right\rangle)\leqslant
1/4$. Hence,

$\displaystyle\sigma^{2}(h)$ $\displaystyle=\kappa\tilde{\sigma}^{2}(\varrho
h_{0})+\kappa^{-1}\hat{\sigma}^{2}(h_{0},\kappa)$
$\displaystyle\leqslant\kappa\tilde{\sigma}^{2}(\varrho
h_{0})+\kappa^{-1}\tilde{\pi}(h_{0}^{2})/4.$

Therefore, for a given function $h$, optimizing the rhs of the above
inequality yields:

$\kappa=\frac{1}{2}\sqrt{\frac{\tilde{\pi}(h_{0}^{2})}{\tilde{\sigma}^{2}(\varrho
h_{0})}}$ (19)

from which we deduce the upper bound

$\sigma^{2}(h)\leqslant\sqrt{\tilde{\pi}(h_{0}^{2})\tilde{\sigma}^{2}(\varrho
h_{0})}\,.$

Note that the choice of $\kappa$ in (19) depends on the function $h$ which is
usually not given beforehand in practice. We propose another way of choosing
$\kappa$ in Section 6.1.2, more adapted to practical concerns.

###  4.7 Geometric ergodicity

For any set $\mathsf{A}\in\mathcal{X}$, we use the notation
$\bar{\mathsf{A}}={\mathsf{X}}\setminus\mathsf{A}$ to denote its complement.
Define the set
$\mathsf{C}_{\eta}:=\left\\{x\in{\mathsf{X}}:\,\rho_{\tilde{R}}(x)\geqslant\eta\right\\}$
for $\eta\in[0,1]$ and note that

$x\in\mathsf{C}_{\eta}\Longleftrightarrow 1-\rho_{\tilde{R}}(x)\leqslant
1-\eta\,.$ (20)

Finally we denote by $\sigma_{\mathsf{C}}=\inf\left\\{k\geqslant
1:\,X_{k}\in\mathsf{C}\right\\}$ the first return time of the set
$\mathsf{C}$.

######  Lemma 2.

Assume that for some $\eta\in(0,1)$, $\mathsf{C}_{\eta}$ is a
$(1,\varepsilon\nu)$-small set for the kernel $Q$. Then, there exists a
probability measure $\tilde{\nu}$ on ${\mathsf{X}}\times{\mathbb{N}_{0}}$
satisfying

  1. (i)

$\mathsf{C}_{\eta}\times\\{0\\}$ is a $(1,\varepsilon\tilde{\nu})$-small set
for the kernel $P$.

  2. (ii)

if $\nu(\mathsf{C}_{\eta}\cap\\{\tilde{R}(.,1)>0\\})>0$, then

$\tilde{\nu}(\mathsf{C}_{\eta}\times\\{0\\})>0.$

###### Proof.

See A.4.1. ∎

We now introduce the following assumption:

  1. (H3)

There exists $\beta_{0}>1$ such that

$\sup_{x\in{\mathsf{X}}}\int_{\mathbb{N}_{0}}\beta_{0}^{n}\tilde{R}(x,\mathrm{d}n)<\infty.$

In other words, the support of $\tilde{R}(x,\cdot)$ is uniformly bounded on
${\mathsf{X}}$. In particular, this assumption is satisfied whenever
$\frac{\mathrm{d}\pi}{\mathrm{d}\tilde{\pi}}$ is upper-bounded, and
$\tilde{R}(x,\cdot)$ is the distribution of
$\lfloor\varrho_{\kappa}(x)\rfloor+U$ where
$U\sim\mathrm{Ber}(\langle\varrho_{\kappa}(x)\rangle)$

######  Remark 6.

Actually, although condition (ii) of Lemma 2 is not verified for any kernel
$\tilde{R}$, it is always possible to transform it slightly into a new kernel
$\tilde{R}^{\prime}$ satisfying this condition while keeping its other
properties untouched, namely assumptions (H2) and (H3). Indeed, define
$\mathsf{C}_{\eta}^{-}=\\{x\in\mathsf{C}_{\eta}:\tilde{R}(x,0)>0\\}$ and
observe that
$\nu(\mathsf{C}_{\eta}\cap\\{\tilde{R}(.,1)>0\\})>\nu(\mathsf{C}_{\eta}^{-}\cap\\{\tilde{R}(.,1)>0\\})$
since $\mathsf{C}_{\eta}^{-}\subset\mathsf{C}_{\eta}$. We will now construct a
kernel $\tilde{R}^{\prime}$ such that for all $x\in\mathsf{C}_{\eta}^{-}$,
$\tilde{R}^{\prime}(x,1)>0$, which satisfies the desired condition as
$\nu(\mathsf{C}_{\eta}^{-}\cap\\{\tilde{R}^{\prime}(.,1)>0\\})=\nu(\mathsf{C}_{\eta}^{-})>0$.
Let $x\in\mathsf{C}_{\eta}^{-}$, implying that $\tilde{R}(x,0)>0$. Due to the
unbiasedness assumption, there exists $k\in{\mathbb{N}_{0}},k>1$ such that
$\tilde{R}(x,k)>0$. Now define $\tilde{R}^{\prime}$ such that:

  * •

$\tilde{R}^{\prime}(x,n)=\tilde{R}(x,n)$ if $n\notin\\{0,1,k\\}$,

  * •

$\tilde{R}^{\prime}(x,1)=\epsilon\tilde{R}(x,0)$,

  * •

$\tilde{R}^{\prime}(x,0)=\tilde{R}(x,0)-\frac{(k-1)\epsilon}{k}\tilde{R}(x,0)$,

  * •

$\tilde{R}^{\prime}(x,k)=\tilde{R}(x,k)-\frac{\epsilon}{k}\tilde{R}(x,0).$

Note that $\epsilon$ can be chosen small enough to guarantee
$\tilde{R}^{\prime}(x,k)>0$. One can easily check that

  * •

$\sum_{n\geqslant 0}\tilde{R}^{\prime}(x,n)=\sum_{n\geqslant
0}\tilde{R}(x,n)=1$,

  * •

$\sum_{n\geqslant 0}\tilde{R}^{\prime}(x,n)n=\sum_{n\geqslant
0}\tilde{R}(x,n)n=\varrho_{\kappa}(x)$.

Hence (H2) and (H3) hold for $\tilde{R}^{\prime}$.

Recall that $\mathsf{A}\in\mathcal{X}$ is accessible for the kernel $Q$ if and
only if for all $x\in{\mathsf{X}}$, there exists $n\in{\mathbb{N}_{0}}$ such
that $Q^{n}(x,\mathsf{A})>0$ (see for example Lemma 3.5.2 of [20]). We then
have the following lemma:

######  Lemma 3.

Assume (H3). Let $\mathsf{A}\in\mathcal{X}$ be an accessible set for $Q$ such
that $\inf_{x\in\mathsf{A}}\rho_{\tilde{R}}(x)>0$. Then $\mathsf{A}$ is
accessible for $S$ and $\mathsf{A}\times\\{0\\}$ is accessible for $P$.

###### Proof.

See A.4.1. ∎

Let us now assume a drift condition on $Q$ in order to deduce an upper bound
for $\sup_{x\in
C_{\eta}}\mathbb{E}_{x}^{S}\left[\beta^{\sigma_{\mathsf{C}_{\eta}}}\right]$.

  1. (H${}_{\mbox{\tiny dft}}$)

There exist constants $(\eta_{0},\lambda)\in(0,1)^{2}$ and a measurable
function $V:{\mathsf{X}}\to[1,\infty)$ such that

    1. (a)

for any $x\notin C_{\eta_{0}}$, we have $QV(x)\leqslant\lambda V(x)$,

    2. (b)

$b_{\infty}=\sup_{x\in C_{\eta_{0}}}\frac{QV}{V}(x)<\infty$,

    3. (c)

$\sup_{x\in C_{\eta_{0}}}V(x)<\infty$.

######  Remark 7.

Under (H${}_{\mbox{\tiny dft}}$), we have
$QV(x)\leqslant(b_{\infty}\vee\lambda)V(x)$ for any $x\in{\mathsf{X}}$. A
straightforward induction yields for any $n\in{\mathbb{N}_{0}}$, $1\leqslant
Q^{n}V(x)\leqslant(b_{\infty}\vee\lambda)^{n}V(x)$. This implies that
$b_{\infty}\vee\lambda\geqslant V(x)^{-1/n}$. Since $n$ is arbitrary,
$b_{\infty}\vee\lambda\geqslant 1$. Finally, $b_{\infty}\geqslant 1>\lambda$
and for all $x\in{\mathsf{X}}$,

$QV(x)\leqslant b_{\infty}V(x).$ (21)

######  Lemma 4.

Assume (H3) and (H${}_{\mbox{\tiny dft}}$). Then there exist constants
$(\eta,\beta_{r},\beta_{\star})\in(0,\eta_{0})\times(1,\infty)\times(0,\infty)$
such that for all $(x,n)\in{\mathsf{X}}\times{\mathbb{N}_{0}}$ :

$\mathbb{E}_{(x,n)}^{P}\left[\beta_{r}^{\sigma_{\mathsf{C}_{\eta}\times\\{0\\}}}\right]\leqslant\beta_{\star}\beta_{r}^{n}V(x)\,.$

###### Proof.

See A.4.2. ∎

Before stating the main theorem, let us introduce a smallness assumption:

  1. (H${}_{\mbox{\tiny sml}}$)

For any $\eta\in(0,1)$ there exist a probability measure
$\nu\in\mathsf{M}_{1}({\mathsf{X}})$ and a constant $\epsilon>0$ such that,
$\mathsf{C}_{\eta}$ is a $(1,\epsilon\nu)$-small set and
$\nu(\mathsf{C}_{\eta}\cap\\{\tilde{R}(\,\cdot\,,1)>0\\})>0$.

######  Theorem 7.

Assume (H1), (H2), (H3), (H${}_{\mbox{\tiny sml}}$) and (H${}_{\mbox{\tiny
dft}}$). Then $P$ has a unique invariant probability measure $\pi$ and there
exist constants $\delta,\beta_{r}>1$, $\zeta<\infty$, such that for all
$\xi\in\mathsf{M}_{1}({\mathsf{X}}\times{\mathbb{N}_{0}})$,

$\sum_{k=1}^{\infty}\delta^{k}d_{TV}(\xi
P^{k},\bar{\pi})\leqslant\zeta\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}\beta_{r}^{n}V(x)\,\xi(\mathrm{d}x\mathrm{d}n).$
(22)

###### Proof.

According to Theorem 11.4.2 of [20], there exists $\zeta_{0}<\infty$ such
that:

$\sum_{k=1}^{\infty}\delta^{k}d_{TV}(\xi
P^{k},\bar{\pi})\leqslant\zeta_{0}\mathbb{E}_{\xi}^{P}[\beta_{r}^{\sigma_{\mathsf{C}_{\eta}\times\\{0\\}}}],$
(23)

provided that:

  1. (i)

$\mathsf{C}_{\eta}\times\\{0\\}$ is an accessible
$(1,\varepsilon\tilde{\nu})$-small set for $P$ satisfying

$\tilde{\nu}(\mathsf{C}_{\eta}\times\\{0\\})>0,$

  2. (ii)

$\sup_{x\in\mathsf{C}_{\eta}}\mathbb{E}_{(x,0)}^{P}[\beta_{r}^{\sigma_{\mathsf{C}_{\eta}\times\\{0\\}}}]<\infty$
for some $\beta_{r}>1$.

Let us start by proving (i). By Lemma 2, there exists
$\tilde{\nu}\in\mathsf{M}_{1}({\mathsf{X}})$ such that
$\mathsf{C}_{\eta}\times\\{0\\}$ is a $(1,\varepsilon\tilde{\nu})$-small set
for $P$ satisfying $\tilde{\nu}(\mathsf{C}_{\eta}\times\\{0\\})>0$. Moreover
it is accessible for $P$, by Lemma 3 since $\mathsf{C}_{\eta}$ is accessible
for $Q$ and $\inf_{x\in\mathsf{C}_{\eta}}\rho_{\tilde{R}}(x)\geqslant\eta>0$.
Hence (i).

It remains to show (ii). We can apply Lemma 4 to get, for $\beta>1$ :

$\sup_{x\in\mathsf{C}_{\eta}}\mathbb{E}_{(x,0)}^{P}[\beta_{r}^{\sigma_{\mathsf{C}_{\eta}\times\\{0\\}}}]\leqslant\beta_{\star}\sup_{x\in\mathsf{C}_{\eta}}V(x)<\infty,$

which shows (23). Then (22) is obtained from (23) by noting that

$\mathbb{E}_{\xi}^{P}[\beta_{r}^{\sigma_{\mathsf{C}_{\eta}\times\\{0\\}}}]=\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}\mathbb{E}_{(x,n)}^{P}[\beta_{r}^{\sigma_{\mathsf{C}_{\eta}\times\\{0\\}}}]\,\xi(\mathrm{d}x\mathrm{d}n)$

and applying Lemma 4.

∎

Under almost minimal conditions (the minimal part will be discussed later in
this paragraph), the Metropolis-Hastings algorithm verifies the conditions of
Theorem 7 ensuring its geometric ergodicity. Indeed, consider a Metropolis-
Hastings algorithm with proposal kernel $A(x,\mathrm{d}y)$ and acceptance rate
$\alpha(x,y)$. Following Remark 3, it can be seen as an instance of the IMC
algorithm with the particular choice of the instrumental kernel
$Q(x,\mathrm{d}y)\propto\alpha(x,y)A(x,\mathrm{d}y)$ and
$\tilde{R}(x,\cdot)\sim\mathrm{Geom}(p(x))$ with
$p(x)=\int_{\mathsf{X}}\alpha(x,y)A(x,\mathrm{d}y)$ the integrated acceptance
rate, i.e. for all $x\in{\mathsf{X}}$:

$\tilde{R}(x,\mathrm{d}n)=p(x)(1-p(x))^{n}.$

Now assume that $p$ is lower bounded by a constant $a_{0}>0$. As shown by [22,
Theorem 3.1.], this condition is necessary for a Metropolis-Hastings algorithm
to be geometrically ergodic. Here, it is also sufficient, since if
$\alpha_{0}:=1-a_{0}<1$, for all $x\in{\mathsf{X}}$,
$\tilde{R}(x,\mathrm{d}n)\leqslant\alpha_{0}^{n}$, and for any choice of
$1<\beta_{0}<\frac{1}{\alpha_{0}}$,

$\int_{X}\beta_{0}^{n}\tilde{R}(x,dn)<\infty.$

This proves (H3) for this choice of $\beta_{0}$. Moreover, Remark 3 also shows
that (H1) and (H2) are satisfied.

##  5 Pseudo-marginal IMC

We develop in this section two different frameworks for pseudo-marginal
Importance Markov chain [23]. The first, simplest one is valid if we want to
replace $\pi(x)$ by an unbiased estimate $\hat{\pi}(x)$. This can be written
as a specific kernel $\tilde{R}$ using the same space as classic IMC.

The second framework tackles the issue of having the intrumental chain
$(\tilde{X}_{i})$ being itself a pseudo-marginal chain, i.e. in this case,
both $\pi$ and $\tilde{\pi}$ are computed through two unbiased estimates.

###  5.1 Pseudo-marginal within IMC

####  5.1.1 Adaptation of the kernel $\tilde{R}$ in the pseudo-marginal
setting

The first pseudo-marginal approach of the Importance Markov chain can be
directly implemented and fits within the framework we develop in this article.
Indeed, knowledge of the density of $\pi$ is never assumed, only the
unbiasedness assumption (H2) is needed (and the geometric control of
hypothesis (H3) for geometric ergodicity).

Assume that for $x\in{\mathsf{X}}$, the density $\pi(x)$ (with respect to some
measure $\mu$) is not directly computable but a nonnegative estimate
$\hat{\pi}(x)$ is available, drawn from a kernel $T_{\pi}(x,\cdot)$ such that
$\int_{\mathbb{R}^{+}}T_{\pi}(x,\mathrm{d}w)w=\pi(x)$ Then, one can replace
$\varrho_{\kappa}(x)$ in (18) by
$\hat{\varrho}_{\kappa}(x)=\kappa\frac{\hat{\pi}(x)}{\tilde{\pi}(x)}$ to get a
plug-in kernel $\tilde{R}_{\scriptscriptstyle\mathrm{pm}}$ that satisfies
(H2).

This can be formalized as follows. First, define an extended kernel
$\tilde{R}_{\psi}$ on
${\mathsf{X}}\times\mathbb{R}^{+}\times\mathcal{P}({\mathbb{N}_{0}})$ by

$\displaystyle\tilde{R}_{\psi}(x,w,\mathrm{d}n)=(1-\left\langle\kappa
w/\tilde{\pi}(x)\right\rangle)\delta_{\left\lfloor\kappa
w/\tilde{\pi}(x)\right\rfloor}(\mathrm{d}n)+\left\langle\kappa
w/\tilde{\pi}(x)\right\rangle\delta_{\left\lfloor\kappa
w/\tilde{\pi}(x)\right\rfloor+1}(\mathrm{d}n).$

Therefore, $\tilde{R}_{\psi}(x,\hat{\pi}(x),\cdot)$ corresponds to the plug-in
random kernel of $\tilde{R}_{\scriptscriptstyle\mathrm{opt}}$ of (18) using
the estimate $\hat{\pi}(x)$. By construction,
$\int_{{\mathbb{N}_{0}}}n\tilde{R}_{\psi}(x,w,\mathrm{d}n)=\kappa
w/\tilde{\pi}(x)$. We can now define the integrated kernel
$\tilde{R}_{\scriptscriptstyle\mathrm{pm}}$ by

$\tilde{R}_{\scriptscriptstyle\mathrm{pm}}(x,\mathrm{d}n)=\int_{\mathbb{R}^{+}}\tilde{R}_{\psi}(x,w,\mathrm{d}n)T_{\pi}(x,\mathrm{d}w),$

and $R_{\scriptscriptstyle\mathrm{pm}}$ using (11).

######  Lemma 5.

$\tilde{R}_{\scriptscriptstyle\mathrm{pm}}$ satisfies (H2).

###### Proof.

Let $x\in{\mathsf{X}}$,

$\displaystyle\int_{\mathbb{N}_{0}}n\tilde{R}_{\scriptscriptstyle\mathrm{pm}}(x,\mathrm{d}n)=\int_{\mathbb{N}_{0}}\int_{\mathbb{R}^{+}}n\tilde{R}_{\psi}(x,w,\mathrm{d}x)T_{\pi}(x,\mathrm{d}x)$
$\displaystyle=\int_{\mathbb{R}^{+}}\left(\int_{\mathbb{N}_{0}}n\tilde{R}_{\psi}(x,w,\mathrm{d}n)\right)T_{\pi}(x,\mathrm{d}w)$
$\displaystyle=\int_{\mathbb{R}^{+}}\frac{\kappa
w}{\tilde{\pi}(x)}T_{\pi}(x,\mathrm{d}w)=\varrho_{\kappa}(x).$

∎

As $\tilde{R}_{\scriptscriptstyle\mathrm{pm}}$ satisfies (H2), the whole
methodology developed in the paper applies to the pseudo marginal case. In
particular, the geometric ergodicity result (see next section) still holds if
the estimator $\hat{\pi}$ is bounded and under (H${}_{\mbox{\tiny rej}}$), as
(H3) will be satisfied.

####  5.1.2 Variance of $\tilde{R}_{\scriptscriptstyle\mathrm{pm}}$

The variance of $\tilde{R}_{\scriptscriptstyle\mathrm{pm}}$ can be expressed
as a function of the variance of $\tilde{R}$. We have the following setup:
$W\sim T_{\pi}(X,\cdot)$ and $N\sim\tilde{R}_{\psi}(X,W,\cdot)$. Then

$\displaystyle\mathbb{V}\mathrm{ar}(N|X)$
$\displaystyle=\mathbb{V}\mathrm{ar}\left(\mathbb{E}[N|X,W]\middle|X\right)+\mathbb{E}\left[\mathbb{V}\mathrm{ar}(N|X,W)\middle|X\right]$
$\displaystyle=\mathbb{V}\mathrm{ar}\left(\frac{\kappa
W}{\tilde{\pi}(X)}\middle|X\right)+\mathbb{E}[\mathbb{V}\mathrm{ar}(N|X,W)|X]$
$\displaystyle=\frac{\kappa^{2}}{\tilde{\pi}^{2}(X)}\mathbb{V}\mathrm{ar}(W|X)+\mathbb{E}\left[\left\langle\frac{\kappa
W}{\tilde{\pi}(X)}\right\rangle\left(1-\left\langle\frac{\kappa
W}{\tilde{\pi}(X)}\right\rangle\right)\middle|X\right].$

The variance can be decomposed into two terms: while the second one is similar
to the variance obtained from $\tilde{R}$ and can also be upper-bounded by
$1/4$, the first one is the direct contribution of the variance of the kernel
$T$ that generates the estimate. In particular, if
$T_{\pi}(X,\cdot)=\delta_{\pi(X)}$, we recover the same expression as in the
previous case.

###  5.2 Fully pseudo-marginal IMC

In this section, we will write that $\pi(\mathrm{d}x)=\pi(x)\mu(\mathrm{d}x)$
and $\tilde{\pi}(\mathrm{d}x)=\tilde{\pi}(x)\mu(\mathrm{d}x)$ for a common
measure $\mu$.

We suppose here that $(\tilde{X}_{k})$ is itself a pseudo-marginal chain where
the estimates of $\tilde{\pi}$ are drawn from a kernel
$T_{\tilde{\pi}}(x,\cdot)$ such that for all $x\in{\mathsf{X}}$,
$\int_{\mathbb{R}^{+}}T_{\tilde{\pi}}(x,\mathrm{d}u)u=\tilde{\pi}(x)$. This
constructs a two-component Markov chain $(\tilde{X}_{k},\tilde{U}_{k})$ on
${\mathsf{X}}\times\mathbb{R}^{+}$ that targets
$\mu(\mathrm{d}x)\,T_{\tilde{\pi}}(x,\mathrm{d}u)u$ [24].

In order to extend the Importance Markov chain to this case, we need once
again to increase the dimension of the chain. The space
${\mathsf{X}}\times\mathbb{R}^{+}$ is replaced by
${\mathsf{X}}\times\mathbb{R}^{+}\times\mathbb{R}^{+}$. In that case, the
second (resp. third) marginal corresponds to the nonnegative estimates of
respectively $\tilde{\pi}$ (resp. $\pi$). The third component $\tilde{V}_{k}$
is drawn from a kernel $T_{\pi}(\tilde{X}_{k},\cdot)$ such that for all
$x\in{\mathsf{X}}$, $\int_{\mathbb{R}^{+}}T_{\pi}(x,\mathrm{d}v)v=\pi(x)$.

This constructs a Markov chain $(\tilde{X}_{k},\tilde{U}_{k},\tilde{V}_{k})$
that targets the probability measure $\tilde{\Pi}$ defined by

$\tilde{\Pi}(\mathrm{d}x\,\mathrm{d}u\,\mathrm{d}v)=\mu(\mathrm{d}x)\,T_{\tilde{\pi}}(x,\mathrm{d}u)u\,T_{\pi}(x,\mathrm{d}v).$

The distribution $\tilde{\Pi}$ is the instrumental distribution for the
extended space. Its first marginal is $\tilde{\pi}$ as
$\tilde{\Pi}(A\times\mathbb{R}^{+}\times\mathbb{R}^{+})=\tilde{\pi}(A)$ for
any $A\in\mathcal{X}$. We can define our target distribution $\Pi$ on the same
extended space by

$\Pi(\mathrm{d}x\,\mathrm{d}u\,\mathrm{d}v)=\mu(\mathrm{d}x)\,\tilde{T}_{\tilde{\pi}}(x,\mathrm{d}u)\,T_{\pi}(x,\mathrm{d}v)v.$

We can perform an Importance Markov chain using the instrumental chain
$(\tilde{X}_{k},\tilde{U}_{k},\tilde{V}_{k})$ targeting the instrumental
density $\tilde{\Pi}$ and the target distribution $\Pi$. In this setting,
$\varrho_{\kappa}$ becomes:

$\varrho_{\kappa}(x,u,v)=\kappa\frac{\mathrm{d}\Pi}{\mathrm{d}\tilde{\Pi}}(x,u,v)=\kappa\frac{v}{u}.$

It remains to draw some random integer $\tilde{N}_{k}$ with conditional
expectation:

$\varrho_{\kappa}(\tilde{X}_{k},\tilde{U}_{k},\tilde{V}_{k})=\kappa\frac{\tilde{V}_{k}}{\tilde{U}_{k}},$

using for instance $\tilde{R}_{\scriptscriptstyle\mathrm{opt}}$.

##  6 Numerical experiments

###  6.1 Toy example: mixture of Gaussians

####  6.1.1 Setting

For starters, we apply Algorithm 2 on an multidimensional Gaussian mixture.
The target distribution $\pi$ writes

$\pi(x)=\sum_{i=1}^{n}\phi_{d}(x;\mu_{i},I_{d}),$

where $\phi_{d}(x;\mu,\Sigma)$ is the density of a Gaussian distribution in
dimension $d$ with mean $\mu$ and covariance matrix $\Sigma$. The means of
these distributions are random, i.i.d.,
$\mu_{i}\sim\mathcal{N}(0,10^{2}I_{d})$. The instrumental distribution
$\tilde{\pi}$ is chosen to be $\tilde{\pi}(x)=\pi(x)^{\beta}$ for a fixed
$\beta\in(0,1)$. The kernel $Q$ targeting $\tilde{\pi}$ is a No-U-turn Sampler
(NUTS) [25].

The parameter $\beta$ flattens the instrumental distribution, thus easy out
the way for the instrumental chain to move from one mode to another, when
compared with the original targetted $\pi$. Conversely, extremely small values
of $\beta$ lead to a very flat instrumental distribution from which it is hard
to reconstruct the original target $\pi$. Therefore, we test different values
of $\beta$ towards finding an optimal tradeoff. The kernel $\tilde{N}$ used to
draw the number of replicas is set to be the same as in (18), i.e. a shifted
Bernoulli.

Note that with such a choice for $\tilde{\pi}$, the ratio
$\frac{\pi}{\tilde{\pi}}$ becomes $\pi^{1-\beta}$ so that the computation of
$\varrho_{\kappa}$ is just a pointwise evaluation of $\pi$ plus basic float
operations. As the choice of $R$ also requires basic float operations (floor,
fractional part and multiplication) and a Bernoulli draw, the complexity of
Algorithm 2 is close to the complexity of directly running $Q$ targeting
$\pi$.

To assess the influence of $\beta$ we estimate the mean squared error in
approximating the expectation of $\pi$ by running $200$ chains for each
$\beta\in\\{0.004,0.01,0.04,0.1,1\\}$. The results are presented in table 1.
As expected, for an untempered instrumental distribution, i.e. $\beta=1$, the
MSE is high due to a low exploration of the space, and is minimal for
$\beta=0.04$.

$\beta$ | 0.004 | 0.010 | 0.040 | 0.100 | 1.000
---|---|---|---|---|---
MSE | 16.150 | 6.123 | 0.544 | 17.863 | 33.982
Table 1: MSE for different values of $\beta$ for the Gaussian mixture

####  6.1.2 Choice of $\kappa$ and analysis of the Effective Sample Size

The parameter $\kappa$ that appears in (6) is somewhat more arbitrary as the
other parameters as its value is directly impacted by the normalizing
constants of $\pi$ and $\tilde{\pi}$. Indeed, assume $\pi_{U}=\pi Z$ (resp.
$\tilde{\pi}_{U}$) is an unnormalized density and write $Z$ (resp. $\tilde{Z}$
) the unknown normalization constant. We can then compute the ratio
$\varrho_{\kappa,U}:=\kappa\frac{\pi_{U}}{\tilde{\pi}_{U}}=\kappa\frac{Z}{\tilde{Z}}\frac{\pi}{\tilde{\pi}}=\varrho_{\kappa\frac{Z}{\tilde{Z}}}$.
Thus, ignoring the normalizing constants leads to a multiplier term $\kappa$
compared to the case where both densities are normalized.

However it proves possible to overcome this issue by noticing that
$\mathbb{E}[\sum_{i}^{n}\tilde{N}_{i}]=\kappa n$, i.e. that, for an
instrumental chain of length $n$, the (expected) length final chain $(X_{i})$
is proportional to $\kappa$. So, one way to deal with this issue is to tune
$\kappa$ such that the length of the final chain is approximately $\alpha n$,
where $\alpha$ is fixed. This is easily solved by taking $\kappa=\frac{\alpha
n}{\sum_{i=1}^{n}\varrho_{U}(\tilde{X}_{i})}$, for
$\varrho_{U}=\frac{\pi_{U}}{\tilde{\pi}_{U}}$ .

######  Remark 8.

The diagnosis of the choice of $\kappa$ can be easily done for a fixed
instrumental chain via the computation of the number of replicas for different
values of $\kappa$ in parallel is cheap once the vector of the values taken by
$\varrho$ is stored.

The problem of tuning $\alpha$ remains. We define a metric that may help to
this effect, namely, the effective sample size (denoted
$\mathrm{ESS}_{\kappa}$, as it depends on $\kappa$). The ESS for an importance
Markov chain is similar to the ESS defined for importance sampling, at the
difference that here the weights are discretized. Formally,

$\mathrm{ESS}_{\kappa}:=\frac{(\sum_{i=1}^{n}\tilde{N}_{i})^{2}}{\sum_{i=1}^{n}\tilde{N}_{i}^{2}}$

and we also define the usual importance sampling ESS:

$\mathrm{ESS}_{\mathrm{IS}}:=\frac{(\sum_{i=1}^{n}\varrho(\tilde{X}_{i}))^{2}}{\sum_{i=1}^{n}\varrho(\tilde{X}_{i})^{2}}.$

######  Remark 9.

Note that the new definition of $\mathrm{ESS}_{\kappa}$ only considers the
replication step of the IMC algorithm, and does not take into account the
convergence of the instrumental chain to the instrumental target. See [9] for
a recent work on the effective sample size for IS with dependent proposals.
The authors add a term to the $\mathrm{ESS}$ to take into account the
correlation of the sample.

Writing
$w_{i}=\frac{\varrho(\tilde{X}_{i})}{\sum_{i=1}^{n}\varrho(\tilde{X}_{i})}$
for the self-normalized $\rho$ , we have that
$\mathbb{E}[\tilde{N}_{i}|\tilde{X}_{i}]=\kappa w_{i}$, so conditionally on
$\tilde{X}_{1:n}$,

$\mathrm{ESS}_{\kappa}\underset{\kappa\rightarrow\infty}{\longrightarrow}\mathrm{ESS}_{\mathrm{IS}},\,\,\,\,\mathbb{P}-\text{a.s.}$
(24)

As expected by the definition of kernel $R$, as $\kappa$ increases, the
stochastic part of the number of replicas vanishes and estimates built with
the chain will behave as with importance sampling, in the sense that :

$M_{n}^{-1}\sum_{i=1}^{M_{n}}h(X_{i})=\sum_{i=1}^{n}\frac{\tilde{N}_{i}}{M_{n}}h(\tilde{X}_{i})\underset{\kappa\rightarrow\infty}{\longrightarrow}\sum_{i=1}^{n}w_{i}h(\tilde{X}_{i}).$

Figure 1: Analysis of the effect of $\kappa$ on the ESS and the length of the
chain

Following remark 8, we plot the effect of $\kappa$ on $\mathrm{ESS}_{\kappa}$
in fig. 1. For a fixed instrumental chain of length $n$ and a fixed
$\beta=0.04$, we computed the $\mathrm{ESS}_{\kappa}$ for $10^{3}$ values of
$\kappa$ varying on log scale from $10^{-1}$ to $10^{4}$. The plot on the
right confirms a linear dependence between the length of the final chain and
$\kappa$. Both other plots ($\mathrm{ESS}_{\kappa}/n$ as a function of
$M_{n}/n$ and $\mathrm{ESS}_{\kappa}$ as a function of $\kappa$) have the same
shape. Overall, the ESS is increasing with $\kappa$ and reaches a stationary
regime for $\kappa$ large enough: it converges to
$\mathrm{ESS}_{\mathrm{IS}}$. Therefore, taking $\kappa$ too large will not
increase the quality of the estimate. In that specific case, this diagnosis
leads us to choose $\alpha\simeq 1$.

###  6.2 Independent IMC with normalizing flows

####  6.2.1 Settings

In this section, we compare the Metropolis Hastings algorithm with independent
proposals with the Importance Markov chain algorithm in the specific case
where $Q(x,\cdot)=\tilde{\pi}(\cdot)$ for all $x\in{\mathsf{X}}$.

The target $\pi$ in dimension $d$ is defined, for
$x=(x_{1},\dots,x_{d})\in\mathbb{R}^{d}$, by

$\displaystyle\pi(x)\propto$
$\displaystyle\exp\bigg{\\{}-\frac{1}{2}\left(\frac{\left\|x\right\|-2}{0.1}\right)^{2}+\sum_{i=1}^{d}\log\left(e^{-\frac{1}{2}(\frac{x_{i}+3}{0.6})^{2}}+e^{-\frac{1}{2}(\frac{x_{i}-3}{0.6})^{2}}\right)\bigg{\\}}.$

This distribution suffers from multimodality as it has $2$ modes per marginal,
for a total of $2^{d}$ modes.

The instrumental distribution $\tilde{\pi}$ is obtained by training a
normalizing flow targeting $\pi$. We recall that a normalizing flow is an
invertible map $T$ from $\mathbb{R}^{d}$ to $\mathbb{R}^{d}$. Taking a base
distribution $\mu$ on $\mathbb{R}^{d}$ (chosen in that case to be the standard
Gaussian), the map $T$ should be chosen such that the pushforward measure of
$\mu$ through $T$ denoted $T_{\sharp}\mu:=\mu(T^{-1}(.))$ is close to $\pi$.
This can be done by optimizing $T$ in a family of maps, for instance rational
quadratic splines (RQSplines) using neural networks [26]. The training is done
by minimizing the forward Kullback-Leibler divergence between $T_{\sharp}\mu$
and $\pi$. See [27] for further details on the training. To get a sample from
the flow, one can generate $x\sim\mu$ and derive $T(x)$. The density $\rho$ of
$T_{\sharp}\mu$ is given by, for $x\in\mathbb{R}^{d}$:

$\rho(x)=\mu\big{(}T^{-1}(x)\big{)}\left\lvert\mathrm{det}J_{T^{-1}}(x)\right\rvert.$

The flows are designed such that $T^{-1}$ and $J_{T^{-1}}(x)$ are easily
computable.

We used the Python package FlowMC [28] with a RQSpline model to train the
flow. Every training of a flow yields a different $\tilde{\pi}=T_{\sharp}\mu$
as the training is stochastic. For details of implementation, see B.1.

The Self Regenerative Markov Chain Monte Carlo (with no adaptation) [14, 15]
is close to the Independent IMC, for the special case where the distribution
of the number of replicas $\tilde{N}_{i}$ is written as

$\mathbb{P}(\tilde{N}_{i}=n|\tilde{X}_{i}=x)=\mathbb{P}(VS=n),$

where $V$ is a Bernoulli random variable with parameter $\alpha(x)$ and $S$ is
geometric with parameter $q(x)$. In [15], the author suggests
$\alpha(x)=\min(1,1/(\varrho_{\kappa}(x)))$ and
$q(x)=\min(1,\varrho_{\kappa}(x))$. This method, called optimal self-
regenerative chain (OSR), is the one we use as a comparison benchmark, with
the same tuning of $\kappa$.

We compare three methods: the independent Metropolis-Hastings (see [2] for
details), the independent importance Markov chain with kernel
$\tilde{R}_{\scriptscriptstyle\mathrm{opt}}$ of (18) and the OSR chain defined
above. In this case, the IMC algorithm is close to the rejection sampling
chain of [29].

The parameter $\kappa$ is tuned such that the length of the final chain is
equal to the length of the instrumental chain. The computational cost (and
running time) of all three algorithms is similar.

####  6.2.2 Results

##### Comparison with Metropolis Hastings and OSR

For each dimension $d\in\\{5,10,15,20,25\\}$, we trained 10 different flows on
the same target. For each flow, 30 i.i.d. samples of size $n=3\cdot 10^{4}$
were generated. We display in fig. 2 the boxplot of the effective sample size
(ESS) of the first marginal for each algorithm, computed using the bulk
method.

While the ESS of the OSR is a bit higher than the one obtained with IMH, both
are outperformed by the importance Markov chain. Regardless of the dimension,
the ESS of the IMC is approximately twice that of MH. Performances across
dimensions are quite similar, even if it is worth noting that the ESS is
slightly lower for the dimension 5, probably due to the training of the flow
and the fact that the hyperparameters of the flow are not optimized for each
dimension.

Figure 2: Comparison of the effective sample size between MH, OSR and IMC

##### Comparison with importance sampling

We compare in Table 2 the mean squared error (MSE) of the first four odd
moments of the first marginal obtained by either the importance Markov chain
or with the importance sampling (IS) estimate defined by

$\hat{I}_{IS}(h)=\frac{\sum_{i=1}^{n}\varrho_{\kappa}(\tilde{X}_{i})h(\tilde{X}_{i})}{\sum_{i=1}^{n}\varrho_{\kappa}(\tilde{X}_{i})}$

against the independent IMC estimate defined by

$\hat{I}_{IMC}(h)=\frac{1}{k}\sum_{j=1}^{k}h(X_{j})=\frac{1}{k}\sum_{i=1}^{n}\tilde{N}_{i}h(\tilde{X}_{i})\hat{I}_{IMC}(h)=\frac{1}{k}\sum_{i=1}^{k}h(X_{i})=\frac{1}{k}\sum_{i=1}^{n}\tilde{N}_{i}h(\tilde{X}_{i}).$

where we denote $k=\sum_{i=1}^{n}\tilde{N}_{i}$ the length of the final chain.
The performances of the IMC are very close to the ones of IS, while the gap
increases with the dimension. The last column shows the mean number points
that are replicated once or more by the IMC, which is the length of the chain
in the case of importance sampling. This has a strong implication: by storing
the importance chain under the representation $(\tilde{X}_{i},\tilde{N}_{i})$,
only (around) 16400 points (and their associated number of replicas) are
needed to be stored for the dimension 25, instead of $30000$ for importance
sampling. If the dimension and the number of points are large, this can be
useful to reduce the memory usage. In some contexts, the slight loss in the
MSE can be compensated by the gain in memory usage. Moreover, the IMC outputs
an actual sample and not a weighted one, sample that approaches the target
distribution at mostly a geometric rate.

dimension | kind | mean | third | fifth | seventh |  $\\#$ positive copies
---|---|---|---|---|---|---
5 | imc | 9.955e-05 | 2.057e-04 | 8.372e-04 | 4.638e-03 | 1.388e+04
is | 9.822e-05 | 2.024e-04 | 8.220e-04 | 4.552e-03 | 3.000e+04
10 | imc | 5.176e-05 | 4.146e-05 | 6.799e-05 | 1.653e-04 | 1.512e+04
is | 4.890e-05 | 3.976e-05 | 6.571e-05 | 1.607e-04 | 3.000e+04
15 | imc | 3.449e-05 | 1.403e-05 | 1.279e-05 | 1.830e-05 | 1.595e+04
is | 3.229e-05 | 1.341e-05 | 1.254e-05 | 1.804e-05 | 3.000e+04
20 | imc | 2.560e-05 | 7.580e-06 | 5.402e-06 | 5.940e-06 | 1.609e+04
is | 2.486e-05 | 7.218e-06 | 5.068e-06 | 5.425e-06 | 3.000e+04
25 | imc | 2.356e-05 | 5.585e-06 | 2.962e-06 | 2.483e-06 | 1.646e+04
is | 2.273e-05 | 5.180e-06 | 2.631e-06 | 2.112e-06 | 3.000e+04
Table 2: Mean square error for the first four odd moments of the first
marginal for both the IMC chain and the importance sampling estimate

##  7 Conclusion

The Importance Markov Chain is a meta-algorithm, in the sense that the
produced Markov chain $\\{X_{k}:\,k\in{\mathbb{N}_{0}}\\}$ is built upon
another one, namely $\\{\tilde{X}_{k}:\,k\in{\mathbb{N}_{0}}\\}$. This allows
the practitioner to change the target of the MCMC kernel used for the
sampling: instead of targeting the distribution of interest directly, our
algorithm targets another distribution that may have better properties, and
then transforms the obtained sample into an output sample following the
distribution of interest. Possible future developments involve building an
adaptative IMC sampler, by allowing the instrumental distribution to evolve
with time, or using several instrumental chains simultaneously.

## Acknowledgments

Charly Andral is supported by a grant from Région Ile-de-France. Christian
Robert is funded by the European Union under the GA 101071601, through the
2023-2029 ERC Synergy grant OCEAN.

##  Appendix A Postponed proofs

The following appendix contains supplementary information that either does not
constitute an essential part of the paper, but is helpful in providing a more
comprehensive understanding of the research problem, or is too cumbersome to
be included in the body of the paper.

###  A.1 Uniqueness of the invariant probability measure

###### (Proof of Proposition 4).

Let $\pi_{0}$ be an invariant probability measure for $P$ and denote by
$\pi_{0}^{k}$ the measure defined by
$\pi_{0}^{k}(\mathsf{A})=\pi_{0}(\mathsf{A}\times\\{k\\})$ for any
$\mathsf{A}\in\mathcal{X}$. To obtain that $\pi_{0}=\bar{\pi}$, we will first
express $\pi_{0}$ using the measure $\pi_{0}^{0}$ only. Let
$f\in\mathsf{F}_{b+}({\mathsf{X}})$. Applying (10) with
$h_{k}(x,n)=f(x)\mathbf{1}_{\\{n\neq k\\}}$ and integrating with respect to
$\pi_{0}$ yields for any $k\geqslant 0$,

$\pi_{0}^{k}f=\pi_{0}h_{k}=\pi_{0}Ph_{k}=\pi_{0}^{k+1}f+\int_{\mathsf{X}}\pi_{0}^{0}S(\mathrm{d}x)R(x,k)f(x).$
(25)

To simplify this equation, we first show that $\pi_{0}^{0}=\pi_{0}^{0}S$.
Indeed, using again $\pi_{0}=\pi_{0}P$ and (10) with the function
$h(x,n)=f(x)$,

$\displaystyle\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}\pi_{0}(\mathrm{d}x\mathrm{d}n)$
$\displaystyle
f(x)=\mathbb{E}^{P}_{\pi_{0}}[f(X_{0})]=\mathbb{E}^{P}_{\pi_{0}}[f(X_{1})]=\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}\pi_{0}(\mathrm{d}x\mathrm{d}n)\mathbf{1}_{\\{n\geqslant
1\\}}f(x)+\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}\pi_{0}(\mathrm{d}x\mathrm{d}n)\mathbf{1}_{\\{n=0\\}}Sf(x),$

which can be equivalently written as $\pi_{0}^{0}f=\pi_{0}^{0}Sf$. Plugging
$\pi_{0}^{0}=\pi_{0}^{0}S$ into (25), we get
$\pi_{0}^{k}f=\pi_{0}^{k+1}f+\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)R(x,k)f(x)$
and by straightforward induction,

$\pi_{0}^{k}f=\pi_{0}^{0}f-\sum_{\ell=0}^{k-1}\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)R(x,\ell)f(x)=\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)f(x)\left[1-\sum_{\ell=0}^{k-1}R(x,\ell)\right]=\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)f(x)\sum_{\ell=k}^{\infty}R(x,\ell).$

Hence, for any $h\in\mathsf{F}_{b+}({\mathsf{X}}\times{\mathbb{N}_{0}})$,

$\displaystyle\pi_{0}h$
$\displaystyle=\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}\pi_{0}(\mathrm{d}x\mathrm{d}n)h(x,n)=\sum_{k=0}^{\infty}\int_{\mathsf{X}}\pi_{0}^{k}(\mathrm{d}x)h(x,k)$
$\displaystyle=\sum_{k=0}^{\infty}\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)h(x,k)\sum_{\ell=k}^{\infty}R(x,\ell)=\sum_{\ell=0}^{\infty}\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)R(x,\ell)\sum_{k=0}^{\ell}h(x,k).$
(26)

Combining with (4.2.1), we can conclude the proof of Proposition 4 (ie
$\pi_{0}=\bar{\pi}$) provided that

$\pi_{0}^{0}(\mathrm{d}x)=\kappa^{-1}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x).$
(27)

All that follows consists in proving this identity. Denote by $\pi_{1}$ the
measure on $({\mathsf{X}}\times[0,1],\mathcal{X}\otimes{\mathcal{B}}([0,1]))$
defined by: for any function $h\in\mathsf{F}_{b+}({\mathsf{X}}\times[0,1])$,

$\pi_{1}h=\int_{{\mathsf{X}}\times[0,1]}\pi_{0}^{0}(\mathrm{d}x)\mathrm{d}u\mathbf{1}_{[0,\rho_{\tilde{R}}(x)]}(u)\rho_{\tilde{R}}(x)^{-1}\mathbb{E}_{(x,u)}^{G}\left[\sum_{k=0}^{\sigma_{\mathsf{D}}-1}h(X_{k},U_{k})\right],$
(28)

where
$\mathsf{D}=\left\\{(x,u)\in{\mathsf{X}}\times[0,1]:\,u\leqslant\rho_{\tilde{R}}(x)\right\\}$
and $G$ is defined in (2). We first show that $\pi_{1}=\pi_{1}G$.

$\displaystyle\pi_{1}Gh$
$\displaystyle=\int_{{\mathsf{X}}\times[0,1]}\pi_{0}^{0}(\mathrm{d}x)\mathrm{d}u\
\mathbf{1}_{[0,\rho_{\tilde{R}}(x)]}(u)\rho_{\tilde{R}}(x)^{-1}\mathbb{E}_{(x,u)}^{G}\left[\sum_{k=0}^{\sigma_{\mathsf{D}}-1}Gh(X_{k},U_{k})\right]$
$\displaystyle=\sum_{k=0}^{\infty}\int_{{\mathsf{X}}\times[0,1]}\pi_{0}^{0}(\mathrm{d}x)\mathrm{d}u\
\mathbf{1}_{[0,\rho_{\tilde{R}}(x)]}(u)\rho_{\tilde{R}}(x)^{-1}\mathbb{E}_{(x,u)}^{G}\left[h(X_{k+1},U_{k+1})\mathbf{1}_{\\{k+1\leqslant\sigma_{\mathsf{D}}\\}}\right]$
$\displaystyle=\int_{{\mathsf{X}}\times[0,1]}\pi_{0}^{0}(\mathrm{d}x)\mathrm{d}u\
\mathbf{1}_{[0,\rho_{\tilde{R}}(x)]}(u)\rho_{\tilde{R}}(x)^{-1}\mathbb{E}_{(x,u)}^{G}\left[\sum_{\ell=1}^{\sigma_{\mathsf{D}}}h(X_{\ell},U_{\ell})\right].$

This implies

$\displaystyle\pi_{1}Gh=\pi_{1}(h)+\int_{{\mathsf{X}}\times[0,1]}\pi_{0}^{0}(\mathrm{d}x)\mathrm{d}u\
\mathbf{1}_{[0,\rho_{\tilde{R}}(x)]}(u)$
$\displaystyle\rho_{\tilde{R}}(x)^{-1}\mathbb{E}_{(x,u)}^{G}\left[h(X_{\sigma_{\mathsf{D}}},U_{\sigma_{\mathsf{D}}})\mathbf{1}_{\\{\sigma_{\mathsf{D}}<\infty\\}}\right]$
$\displaystyle-\int_{{\mathsf{X}}\times[0,1]}\pi_{0}^{0}(\mathrm{d}x)\left(\int_{0}^{\rho_{\tilde{R}}(x)}h(x,u)\mathrm{d}u\right)\rho_{\tilde{R}}(x)^{-1}.$
(29)

Now, write

$\displaystyle\mathbb{E}_{(x,u)}^{G}\left[h(X_{\sigma_{\mathsf{D}}},U_{\sigma_{\mathsf{D}}})\mathbf{1}_{\\{\sigma_{\mathsf{D}}<\infty\\}}\right]$
$\displaystyle=\sum_{\ell=1}^{\infty}\mathbb{E}_{(x,u)}^{G}\left[h(X_{\ell},U_{\ell})\mathbf{1}_{\\{U_{\ell}\leqslant\rho_{\tilde{R}}(X_{\ell})\\}}\mathbf{1}_{\\{\sigma_{\mathsf{D}}\geqslant\ell\\}}\right]$
$\displaystyle=\sum_{\ell=1}^{\infty}\mathbb{E}_{(x,u)}^{G}\left[\left(\int_{0}^{\rho_{\tilde{R}}(X_{\ell})}h(X_{\ell},u^{\prime})\mathrm{d}u^{\prime}\right)\mathbf{1}_{\\{\sigma_{\mathsf{D}}\geqslant\ell\\}}\right]$
$\displaystyle=\sum_{\ell=1}^{\infty}\mathbb{E}_{(x,u)}^{G}\left[\left(\int_{0}^{\rho_{\tilde{R}}(X_{\ell})}h(X_{\ell},u^{\prime})\mathrm{d}u^{\prime}\right)\rho_{\tilde{R}}(X_{\ell})^{-1}\mathbf{1}_{\\{U_{\ell}\leqslant\rho_{\tilde{R}}(X_{\ell})\\}}\mathbf{1}_{\\{\sigma_{\mathsf{D}}\geqslant\ell\\}}\right]$
$\displaystyle=\mathbb{E}_{(x,u)}^{G}\left[\left(\int_{0}^{\rho_{\tilde{R}}(X_{\sigma_{\mathsf{D}}})}h(X_{\sigma_{\mathsf{D}}},u^{\prime})\mathrm{d}u^{\prime}\right)\rho_{\tilde{R}}(X_{\sigma_{\mathsf{D}}})^{-1}\right]$
$\displaystyle=\int_{\mathsf{X}}S(x,\mathrm{d}x^{\prime})\left(\int_{0}^{\rho_{\tilde{R}}(x^{\prime})}h(x^{\prime},u^{\prime})\mathrm{d}u^{\prime}\right)\rho_{\tilde{R}}(x^{\prime})^{-1}.$

Plugging this expression into (29) yields:

$\pi_{1}Gh=\pi_{1}(h)+\int_{\mathsf{X}}\pi_{0}^{0}S(\mathrm{d}x^{\prime})\left(\int_{0}^{\rho_{\tilde{R}}(x^{\prime})}h(x^{\prime},u^{\prime})\mathrm{d}u^{\prime}\right)\rho_{\tilde{R}}(x^{\prime})^{-1}-\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)\left(\int_{0}^{\rho_{\tilde{R}}(x)}h(x,u)\mathrm{d}u\right)\rho_{\tilde{R}}(x)^{-1}=\pi_{1}h,$

where we have used that $\pi_{0}^{0}S=\pi_{0}^{0}$. Hence $\pi_{1}$ is an
invariant measure for $G$. Since any invariant measure for $Q$ is proportional
to $\tilde{\pi}$, we deduce from (2) that $G$ admits a unique invariant
measure (up to a multiplicative constant) proportional to
$\tilde{\pi}(\mathrm{d}x)\mathbf{1}_{[0,1]}(u)\mathrm{d}u$ and hence,
$\tilde{\pi}(\mathrm{d}x)\mathbf{1}_{[0,1]}(u)\mathrm{d}u\propto\pi_{1}$. Now,
taking $h(x,n)=\mathbf{1}_{\mathsf{D}}(x,u)f(x)$ for any arbitrary
$f\in\mathsf{F}_{b+}({\mathsf{X}})$, we get, using (28),

$\displaystyle\int_{\mathsf{X}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)f(x)$
$\displaystyle=\int_{{\mathsf{X}}\times[0,1]}\tilde{\pi}(\mathrm{d}x)\mathbf{1}_{[0,1]}(u)\mathrm{d}u\
h(x,u)$
$\displaystyle\propto\pi_{1}h=\int_{{\mathsf{X}}\times[0,1]}\pi_{0}^{0}(\mathrm{d}x)\mathrm{d}u\
\mathbf{1}_{[0,\rho_{\tilde{R}}(x)]}(u)\rho_{\tilde{R}}(x)^{-1}\mathbf{1}_{\mathsf{D}}(x,u)f(x)=\pi_{0}^{0}f.$

Finally, there exists a constant $\gamma$ such that for any
$f\in\mathsf{F}_{b+}({\mathsf{X}})$,

$\pi_{0}^{0}f=\gamma\int_{\mathsf{X}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)f(x).$
(30)

Applying (26) with $h=\mathbf{1}$ and using sucessively (11), (H2) and the
identity above, we get

$1=\sum_{\ell=0}^{\infty}\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)R(x,\ell)(\ell+1)=\sum_{k=1}^{\infty}\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)\frac{\tilde{R}(x,k)k}{\rho_{\tilde{R}}(x)}=\int_{\mathsf{X}}\pi_{0}^{0}(\mathrm{d}x)\frac{\varrho_{\kappa}(x)}{\rho_{\tilde{R}}(x)}=\gamma\int_{\mathsf{X}}\tilde{\pi}(\mathrm{d}x)\varrho_{\kappa}(x)=\gamma\kappa.$

Combined with (30) we finally obtain
$\pi_{0}^{0}f=\kappa^{-1}\int_{\mathsf{X}}\tilde{\pi}(\mathrm{d}x)\rho_{\tilde{R}}(x)f(x)$
which proves (27) and concludes the proof. ∎

###  A.2 A martingale weak convergence result with random indexes.

######  Theorem 8.

Let $(\Omega,\mathcal{F},\mathbb{P})$ be a probability space and let
$(\mathcal{F}_{n})$ be a filtration on $\Omega$ such that
$\mathcal{F}_{n}\subset\mathcal{F}$ for any $n\in{\mathbb{N}_{0}}$. Let
$(M_{n})$ be a square-integrable $(\mathcal{F}_{n})$-martingale such that

$\frac{M_{n}}{\sqrt{n}}\overset{\mathbb{P}-law}{\leadsto}G\,\quad\mbox{and}\quad\frac{\mathbb{E}\left[M_{n}^{2}\right]}{n}\rightarrow\sigma^{2}$

and let $(k_{n})$ be a sequence of random integers such that $\frac{k_{n}}{n}\
\overset{\mathbb{P}-\mathrm{prob}}{{\longrightarrow}}\ \lambda\in(0,\infty)$.
Then

$(n\lambda)^{-1/2}M_{k_{n}}\overset{\mathbb{P}-law}{\leadsto}G.$

###### Proof.

Let $\alpha,\lambda^{-},\lambda^{+}$ be positive constants such that
$\lambda^{-}<\lambda<\lambda^{+}$. Define

$\displaystyle B_{n}=(n\lambda)^{-1/2}M_{k_{n}},$ $\displaystyle
C_{n}=(n\lambda)^{-1/2}M_{\left\lfloor n\lambda^{-}\right\rfloor},$

and note that
$C_{n}\overset{\mathbb{P}}{\rightsquigarrow}(\lambda^{-}/\lambda)^{1/2}G$.
Then we have for any $u\in\mathbb{R}$,

$\displaystyle\lvert\mathbb{E}[\mathrm{e}^{iuB_{n}}]-\mathbb{E}[\mathrm{e}^{iuG}]\rvert$
$\displaystyle\leqslant\lvert\mathbb{E}[\mathrm{e}^{iuB_{n}}]-\mathbb{E}[\mathrm{e}^{iuC_{n}}]\rvert+\lvert\mathbb{E}[\mathrm{e}^{iuC_{n}}]-\mathbb{E}[\mathrm{e}^{iuG}]\rvert$
$\displaystyle\leqslant 2\mathbb{P}\left(\lvert
B_{n}-C_{n}\rvert>\alpha\right)+\lvert\mathbb{E}[\lvert\mathrm{e}^{iuB_{n}}-\mathrm{e}^{iuC_{n}}\rvert\mathbf{1}_{\\{\lvert
B_{n}-C_{n}\rvert\leqslant\alpha\\}}]\rvert+\lvert\mathbb{E}[\mathrm{e}^{iuC_{n}}]-\mathbb{E}{\mathrm{e}^{iuG}}\rvert$
$\displaystyle\leqslant 2\mathbb{P}\left(\lvert
B_{n}-C_{n}\rvert>\alpha\right)+\sup_{\lvert\beta\rvert\leqslant\alpha}\lvert\mathrm{e}^{iu\beta}-1\rvert+\lvert\mathbb{E}[\mathrm{e}^{iuC_{n}}]-\mathbb{E}[\mathrm{e}^{iuG}]\rvert.$
(31)

Since $(M_{i}-M_{\left\lfloor
n\lambda^{-}\right\rfloor})_{i\geqslant\left\lfloor
n\lambda^{-}\right\rfloor}$ is a martingale and $x\to x^{2}$ is convex,
$\left((M_{i}-M_{\left\lfloor
n\lambda^{-}\right\rfloor})^{2}\right)_{i\geqslant\left\lfloor
n\lambda^{-}\right\rfloor}$ is a non-negative submartingale. Then, the first
term of the right hand side in (31) may be bounded by applying Doob’s maximal
inequality to the non-negative submartingale $\left((M_{i}-M_{\left\lfloor
n\lambda^{-}\right\rfloor})^{2}\right)_{i\geqslant\left\lfloor
n\lambda^{-}\right\rfloor}$,

$\displaystyle\mathbb{P}(\left\lvert B_{n}-C_{n}\right\rvert>\alpha)$
$\displaystyle\leqslant\mathbb{P}\left(k_{n}\notin[n\lambda^{-},n\lambda^{+}]\right)+\mathbb{P}\left(\left\lvert
B_{n}-C_{n}\right\rvert>\alpha,k_{n}\in[n\lambda^{-},n\lambda^{+}]\right)$
$\displaystyle\leqslant\mathbb{P}(k_{n}\notin[n\lambda^{-},n\lambda^{+}])+\mathbb{P}\left(\sup_{i\in[n\lambda^{-}:n\lambda^{+}]}\lvert
M_{i}-M_{\left\lfloor
n\lambda^{-}\right\rfloor}\rvert>(n\lambda)^{1/2}\alpha\right)$
$\displaystyle\leqslant\mathbb{P}(k_{n}\notin[n\lambda^{-},n\lambda^{+}])+\frac{\mathbb{E}\left[\left(M_{\left\lfloor
n\lambda^{+}\right\rfloor}-M_{\left\lfloor
n\lambda^{-}\right\rfloor}\right)^{2}\right]}{n\lambda\alpha^{2}}$
$\displaystyle=\mathbb{P}(k_{n}\notin[n\lambda^{-},n\lambda^{+}])+\frac{\sum_{k=\left\lfloor
n\lambda^{-}\right\rfloor}^{\left\lfloor
n\lambda^{+}\right\rfloor-1}\mathbb{E}\left[\left(M_{k+1}-M_{k}\right)^{2}\right]}{n\lambda\alpha^{2}}.$

Note that since $(M_{n})_{n\in\mathbb{N}}$ is a square-integrable martingale,
we have

$D_{n}=\frac{\mathbb{E}[M_{n}^{2}]}{n}=\frac{\mathbb{E}[M_{0}^{2}]+\sum_{k=0}^{n-1}\mathbb{E}[(M_{k+1}-M_{k})^{2}]}{n},$

and the previous bound writes:

$\mathbb{P}\left(\left\lvert
B_{n}-C_{n}\right\rvert>\alpha\right)\leqslant\mathbb{P}(k_{n}\notin[n\lambda^{-},n\lambda^{+}])+\frac{\left\lfloor
n\lambda^{+}\right\rfloor D_{\left\lfloor
n\lambda^{+}\right\rfloor}-\left\lfloor n\lambda^{-}\right\rfloor
D_{\left\lfloor n\lambda^{-}\right\rfloor}}{n\lambda\alpha^{2}}.$

Finally letting $n$ go to infinity and using successively that
$\frac{k_{n}}{n}\ \overset{\mathbb{P}-\mathrm{prob}}{{\longrightarrow}}\
\lambda$, $D_{n}\underset{n\rightarrow\infty}{\longrightarrow}\sigma^{2}$ and
$C_{n}\overset{\mathbb{P}}{\rightsquigarrow}(\frac{\lambda^{-}}{\lambda})^{1/2}G$,
we obtain

$\limsup_{n\rightarrow\infty}\lvert\mathbb{E}[\mathrm{e}^{iuB_{n}}]-\mathbb{E}[\mathrm{e}^{iuG}]\rvert\leqslant
2\sigma^{2}\frac{\lambda^{+}-\lambda^{-}}{\lambda\alpha^{2}}+\sup_{\lvert\beta\rvert\leqslant\alpha}\lvert\mathrm{e}^{iu\beta}-1\rvert+\lvert\mathbb{E}[\mathrm{e}^{iu(\lambda^{-}/\lambda)^{1/2}G}]-\mathbb{E}[\mathrm{e}^{iuG}]\rvert.$

Letting $\lambda^{+}\searrow\lambda$ and $\lambda^{-}\nearrow\lambda$, we get

$\limsup_{n\to\infty}\lvert\mathbb{E}[\mathrm{e}^{iuB_{n}}]-\mathbb{E}[\mathrm{e}^{iuG}]\rvert\leqslant\sup_{\lvert\beta\rvert\leqslant\alpha}\lvert\mathrm{e}^{iu\beta}-1\rvert,$

and letting $\alpha\to 0$, we finally obtain
$\limsup_{n\to\infty}\lvert\mathbb{E}[\mathrm{e}^{iuB_{n}}]-\mathbb{E}[\mathrm{e}^{iuG}]\rvert=0$.
Therefore $B_{n}\overset{\mathbb{P}-law}{\leadsto}G$ which concludes the
proof. ∎

###  A.3 Central Limit Theorem

####  A.3.1 Preliminary results

Let $\bar{Q}$ be the Markov kernel on
$({\mathsf{X}}\times{\mathbb{N}_{0}})\times(\mathcal{X}\otimes\mathcal{P}({\mathbb{N}_{0}}))$
defined by

$\bar{Q}(x,n;\mathrm{d}x^{\prime}\mathrm{d}n^{\prime})=Q(x,\mathrm{d}x^{\prime})\tilde{R}(x^{\prime},\mathrm{d}n^{\prime}),$

and $\hat{\pi}$ be the probability measure on ${\mathsf{X}}\times\mathbb{N}$
defined by

$\hat{\pi}(\mathrm{d}x\mathrm{d}n)=\tilde{\pi}(\mathrm{d}x)\tilde{R}(x,\mathrm{d}n).$

Let $S_{n}=\sum_{i=1}^{n}\tilde{N}_{i}$ and define
$k_{n}=\max\left\\{k\in\mathbb{N}^{*}:S_{k}\leqslant n\right\\}$ ensuring that
$S_{k_{n}}\leqslant n<S_{k_{n}+1}$.

######  Lemma 6.

Assume (H1). Then $\bar{Q}$ admits $\hat{\pi}$ as invariant probability
measure.

###### Proof.

Let $\mathsf{A}\in\mathcal{X}\otimes\mathcal{P}({\mathbb{N}_{0}})$. Then

$\displaystyle\hat{\pi}\bar{Q}(\mathsf{A})$
$\displaystyle=\int_{({\mathsf{X}}\times{\mathbb{N}_{0}})^{2}}\tilde{\pi}(\mathrm{d}x)\tilde{R}(x,\mathrm{d}n)Q(x,\mathrm{d}x^{\prime})\tilde{R}(x^{\prime},\mathrm{d}n^{\prime})\mathbf{1}_{\mathsf{A}}(x^{\prime},n^{\prime})$
$\displaystyle=\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}\tilde{\pi}(\mathrm{d}x^{\prime})\tilde{R}(x^{\prime},\mathrm{d}n^{\prime})\mathbf{1}_{\mathsf{A}}(x^{\prime},n^{\prime})$
$\displaystyle=\hat{\pi}(\mathsf{A}).$

∎

######  Lemma 7.

Assume (H1) and (H${}_{\mbox{\tiny lln}}$). Then, for every
${\xi^{\prime}}\in\mathsf{M}_{1}({\mathsf{X}}\times{\mathbb{N}_{0}})$ and
measurable function
$g:{\mathsf{X}}\times{\mathbb{N}_{0}}\rightarrow\mathbb{R}$ such that
$\hat{\pi}(\left\lvert g\right\rvert)<\infty$,

$\lim_{n\rightarrow\infty}n^{-1}\sum_{k=0}^{n-1}g(\tilde{X}_{k},\tilde{N}_{k})=\hat{\pi}(g),\quad\mathbb{P}_{\xi^{\prime}}^{\bar{Q}}-a.s.$
(32)

###### Proof.

The proof follows that of Theorem 5 and also relies on [21, Proposition 3.5].
Let $\bar{h}$ be a harmonic function for the kernel $\bar{Q}$, i.e. for all
$(x,n)\in{\mathsf{X}}\times\mathbb{N}$,

$\bar{Q}\bar{h}(x,n)=\bar{h}(x,n),$

and let us prove that $\bar{h}$ is constant. For all
$(x,n)\in{\mathsf{X}}\times\mathbb{N}$,

$\bar{h}(x,n)=\bar{Q}\bar{h}(x,n)=\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}Q(x,\mathrm{d}x^{\prime})\tilde{R}(x^{\prime},\mathrm{d}n^{\prime})\bar{h}(x^{\prime},n^{\prime}).$
(33)

The integral on the right hand side of the equation above does not depend on
$n$, therefore $\bar{h}$ is also independant of $n$ and we can write for all
$(x,n)\in{\mathsf{X}}\times\mathbb{N}$,

$\bar{h}(x,n)=\bar{h}(x,0)=:h_{0}(x).$

Then from (33) we have for all $x\in{\mathsf{X}}$,

$h_{0}(x)=\int_{{\mathsf{X}}}Q(x,\mathrm{d}x^{\prime})h_{0}(x^{\prime})=Qh_{0}(x).$
(34)

which proves that $h_{0}$ is harmonic for $Q$. Using [21, Proposition 3.5], we
get that $h_{0}$ is constant. Therefore so is $\bar{h}$ and using [21,
Proposition 3.5] again, we have completed the proof of the lemma. ∎

######  Lemma 8.

For all $\zeta\in\mathsf{M}_{1}({\mathsf{X}}\times{\mathbb{N}_{0}})$,
$\frac{k_{n}}{n}\
\overset{\mathbb{P}^{\bar{Q}}_{\zeta}-\mathrm{prob}}{{\longrightarrow}}\
\kappa^{-1}$.

###### Proof.

Let $\beta>\kappa^{-1}$, we will prove that
$\mathbb{P}^{\bar{Q}}_{\zeta}\left(\frac{k_{n}}{n}\geqslant\beta\right)\underset{n\rightarrow\infty}{\longrightarrow}0$.
From the definition of $k_{n}$ we have that $\left\\{k_{n}\geqslant\beta
n\right\\}=\left\\{\sum^{\left\lfloor\beta
n\right\rfloor}_{i=1}\tilde{N}_{i}\leqslant n\right\\}$, hence

$\displaystyle\mathbb{P}^{\bar{Q}}_{\zeta}\left(\frac{k_{n}}{n}\geqslant\beta\right)=\mathbb{P}^{\bar{Q}}_{\zeta}\left(\frac{1}{\left\lfloor\beta
n\right\rfloor}\sum^{\left\lfloor\beta
n\right\rfloor}_{i=1}\tilde{N}_{i}\leqslant\frac{n}{\left\lfloor\beta
n\right\rfloor}\right),$

which converges to 0 as $n$ goes to infinity since by applying Lemma 7 with
$g(x,\ell)=\ell$ we have

$\frac{1}{\left\lfloor\beta n\right\rfloor}\sum_{i=1}^{\left\lfloor\beta
n\right\rfloor}\tilde{N}_{i}\
\overset{\mathbb{P}^{\bar{Q}}_{\zeta}-a.s.}{{\longrightarrow}}\
\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}\tilde{\pi}(\mathrm{d}x)\tilde{R}(x,\mathrm{d}\ell)\ell=\int_{{\mathsf{X}}}\tilde{\pi}(\mathrm{d}x)\varrho_{\kappa}(x)=\kappa>\beta^{-1}.$

Similarly we prove that for any $\beta<\kappa^{-1}$,

$\mathbb{P}^{\bar{Q}}_{\zeta}\left(\frac{k_{n}}{n}<\beta\right)\underset{n\rightarrow\infty}{\longrightarrow}0,$

by using that $\left\\{k_{n}<\beta n\right\\}=\left\\{\sum^{\left\lfloor\beta
n\right\rfloor}_{i=1}\tilde{N}_{i}>n\right\\}$, which completes the proof. ∎

######  Lemma 9.

Assume (H1) and (H${}_{\mbox{\tiny lln}}$). Let
$\zeta\in\mathsf{M}_{1}({\mathsf{X}}\times{\mathbb{N}_{0}})$,
$f:{\mathsf{X}}\times\mathbb{N}\longrightarrow\mathbb{R}$ be a measurable
function, and $(k_{n})_{n}\in\mathbb{N}^{\mathbb{N}}$ be a sequence of random
variables such that $\frac{k_{n}}{n}\
\overset{\mathbb{P}^{\bar{Q}}_{\zeta}-\mathrm{prob}}{{\longrightarrow}}\
\kappa^{-1}$ and $\hat{\pi}f^{2}<\infty$. Then,

$\frac{f(\tilde{X}_{k_{n}},\tilde{N}_{k_{n}})}{\sqrt{n}}\
\overset{\mathbb{P}^{\bar{Q}}_{\zeta}-\mathrm{prob}}{{\longrightarrow}}\ 0.$

###### Proof.

Let $\epsilon>0$, we will prove that

$\mathbb{P}^{\bar{Q}}_{\zeta}\left(g(\tilde{X}_{k_{n}},\tilde{N}_{k_{n}})>\epsilon^{2}n\right)\underset{n\rightarrow\infty}{\longrightarrow}0,$

where $g=f^{2}$. Let $\alpha,\beta\in\mathbb{R}^{+}$ be two nonnegative
numbers such that $\alpha<\kappa^{-1}<\beta$ and
$\beta-\alpha<\frac{\epsilon^{2}}{\hat{\pi}(g)}$.

$\displaystyle\mathbb{P}^{\bar{Q}}_{\zeta}\left(g(\tilde{X}_{k_{n}},\tilde{N}_{k_{n}})>\epsilon^{2}n\right)$
$\displaystyle\leqslant\mathbb{P}^{\bar{Q}}_{\zeta}\left(\frac{k_{n}}{n}\notin\left[\alpha,\beta\right]\right)+\mathbb{P}^{\bar{Q}}_{\zeta}\left(g(\tilde{X}_{k_{n}},\tilde{N}_{k_{n}})>\epsilon^{2}n,\frac{k_{n}}{n}\in\left[\alpha,\beta\right]\right).$

From $\frac{k_{n}}{n}\
\overset{\mathbb{P}^{\bar{Q}}_{\zeta}-\mathrm{prob}}{{\longrightarrow}}\
\kappa^{-1}$ we get
$\limsup_{n}\mathbb{P}^{\bar{Q}}_{\zeta}\left(\frac{k_{n}}{n}\notin\left[\alpha,\beta\right]\right)=0$
and therefore,

$\limsup_{n}\mathbb{P}^{\bar{Q}}_{\zeta}\left(g(\tilde{X}_{k_{n}},\tilde{N}_{k_{n}})>\epsilon^{2}n\right)\leqslant\limsup_{n}\mathbb{P}^{\bar{Q}}_{\zeta}\left(g(\tilde{X}_{k_{n}},\tilde{N}_{k_{n}})>\epsilon^{2}n,\frac{k_{n}}{n}\in\left[\alpha,\beta\right]\right).$

Defining $A_{n}=\frac{1}{n}\sum_{i=1}^{n}g(\tilde{X}_{i},\tilde{N}_{i})$, we
have

$\displaystyle\mathbb{P}^{\bar{Q}}_{\zeta}\left(g(\tilde{X}_{k_{n}},\tilde{N}_{k_{n}})>\epsilon^{2}n,\frac{k_{n}}{n}\in\left[\alpha,\beta\right]\right)$
$\displaystyle\leqslant\mathbb{P}^{\bar{Q}}_{\zeta}\left(\sum^{\left\lfloor
n\beta\right\rfloor}_{i=\left\lfloor
n\alpha\right\rfloor}g(\tilde{X}_{i},\tilde{N}_{i})>\epsilon^{2}n\right)$
$\displaystyle=\mathbb{P}^{\bar{Q}}_{\zeta}\left(\frac{\left\lfloor
n\beta\right\rfloor}{n}A_{\left\lfloor n\beta\right\rfloor}-\frac{\left\lfloor
n\alpha\right\rfloor}{n}A_{\left\lfloor
n\alpha\right\rfloor}>\epsilon^{2}\right)\underset{n\rightarrow\infty}{\longrightarrow}0,$

since $A_{n}\ \overset{\mathbb{P}^{\bar{Q}}_{\zeta}-a.s.}{{\longrightarrow}}\
\hat{\pi}g$ by Lemma 7 and $\hat{\pi}(g)(\beta-\alpha)<\epsilon^{2}$. This
concludes the proof. ∎

######  Lemma 10.

Let $(Y_{k})_{k}$ be a Markov chain on ${\mathsf{Y}}$ generated by a kernel
$T$ on ${\mathsf{Y}}\times\mathcal{Y}$ with invariant measure $\mu$. Assume
that there exists a measurable function $j:{\mathsf{Y}}\rightarrow\mathbb{R}$
such that the Poisson equation on ${\mathsf{Y}}$ for the kernel $T$ associated
to the function $j$ admits a solution $J$, i.e. for all $y\in{\mathsf{Y}}$

$J(y)-TJ(y)=j(y)-\mu(j).$

Then, considering the filtration $\mathcal{F}_{n}=\sigma(Y_{0:n})$,

$\sum_{i=0}^{n-1}j(Y_{i})-\mu(j)=M_{n}-J(Y_{n})+J(Y_{0}),$

where $M_{n}:=\sum_{i=1}^{n}J(Y_{i})-TJ(Y_{i-1})$ is a
$\mathcal{F}_{n}$-martingale.

###### Proof.

A simple index shift gives that

$\sum_{i=0}^{n-1}j(Y_{i})-\mu(j)=\sum_{i=1}^{n-1}J(Y_{i})-TJ(Y_{i})=\sum_{i=1}^{n}J(Y_{i})-TJ(Y_{i-1})-J(Y_{n})+J(Y_{1})=M_{n}-J(Y_{n})+J(Y_{1}).$

∎

Let us state as a reminder the following theorem from [30] around which the
proof of the following lemma will revolve.

######  Theorem 9 (Theorem 13 of Douc-Moulines).

Let $(\Omega,\mathcal{F},\mathbb{P})$ be a probability space and let
$(\mathcal{F}_{n,i})_{i\leqslant n}$ be a filtration on $\Omega$. Assume
$\mathbb{E}\left[U_{n,i}^{2}\middle|\mathcal{F}_{n,i-1}\right]<\infty$ for any
$n\in{\mathbb{N}_{0}}$ and any $i=1,\dots,n$, and

$\displaystyle\sum_{i=1}^{n}\left(\mathbb{E}\left[U_{n,i}^{2}\middle|\mathcal{F}_{n,i-1}\right]-\mathbb{E}\left[U_{n,i}\middle|\mathcal{F}_{n,i-1}\right]^{2}\right)\longrightarrow\sigma^{2}$
for some $\sigma>0$ (H3)
$\displaystyle\sum_{i=1}^{n}\mathbb{E}\left[U_{n,i}^{2}\mathbf{1}_{\\{\left\lvert
U_{n,i}\right\rvert\geqslant\varepsilon\\}}\middle|\mathcal{F}_{n,i-1}\right]\longrightarrow
0$ for any $\varepsilon>0$ (H4)

Then, for any real $u$,

$\mathbb{E}\left[\exp\left(iu\sum_{i=1}^{n}(U_{n,i}-\mathbb{E}\left[U_{n,i}\middle|\mathcal{F}_{n,i-1}\right])\right)\middle|\mathcal{F}_{n,0}\right]\longrightarrow\exp(-(u^{2}/2)\sigma^{2}).$

######  Lemma 11.

Let $(Y_{k})_{k}$ be a Markov chain on ${\mathsf{Y}}$ with kernel $T$ on
${\mathsf{Y}}\times\mathcal{Y}$ admitting an invariant probability measure
$\mu$. Assume that for every $\nu\in\mathsf{M}_{1}({\mathsf{Y}})$ and any
measurable function $g:{\mathsf{Y}}\rightarrow\mathbb{R}$ such that
$\mu(\left\lvert g\right\rvert)<\infty$,

$\lim_{n\rightarrow\infty}n^{-1}\sum_{k=0}^{n-1}g(Y_{k})=\mu(g),\quad\mathbb{P}_{\nu}^{T}-a.s.$
(35)

Let $J:{\mathsf{Y}}\rightarrow\mathbb{R}$ be a measurable function such that
$\mu(J^{2})<\infty$. Consider the filtration $\mathcal{F}_{n}=\sigma(Y_{0:n})$
and the $\mathcal{F}_{n}$-martingale
$M_{n}=\sum_{i=1}^{n}J(Y_{i})-TJ(Y_{i-1})=\sum_{i=1}^{n}\Delta M_{i}$. Then,

$n^{-1/2}M_{n}\overset{\mathbb{P}^{M}_{\nu}-law}{\leadsto}\mathcal{N}(0,\sigma_{J}^{2}(h)),$

with $\sigma_{J}^{2}(h)=\mathbb{E}^{T}_{\mu}\left[\Delta M_{1}^{2}\right]$.

###### Proof.

Define $U_{n,i}=\frac{\Delta M_{i}}{\sqrt{n}}$ for any $n\geqslant i\geqslant
1$ and $\mathcal{F}_{n,i}=\mathcal{F}_{i}=\sigma(Y_{0:i})$ for any
$i,n\in\mathbb{N}$. Let us verify the hypotheses of Theorem 9 :

  * •

For (H3), we have

$\displaystyle\sum_{i=1}^{n}\left(\mathbb{E}_{\nu}^{T}\left[U_{n,i}^{2}\middle|\mathcal{F}_{n,i-1}\right]-\mathbb{E}_{\nu}^{T}\left[U_{n,i}\middle|\mathcal{F}_{n,i-1}\right]^{2}\right)$
$\displaystyle=\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}_{\nu}^{T}\left[\Delta
M_{i}^{2}\middle|\mathcal{F}_{i-1}\right]$
$\displaystyle=\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}^{T}_{Y_{i-1}}\left[\Delta
M_{1}^{2}\right]\underset{n\rightarrow\infty}{\longrightarrow}\mathbb{E}^{T}_{\mu}\left[\Delta
M_{1}^{2}\right],$

where the limit is obtained from (35) with the function
$g:y\mapsto\mathbb{E}^{T}_{y}\left[\Delta M_{1}^{2}\right]$.

  * •

For (H4), let $A>0$ be a positive integer,

$\displaystyle\sum_{i=1}^{n}\mathbb{E}_{\nu}^{T}\left[U_{n,i}^{2}\mathbf{1}_{\left\\{\left\lvert
U_{n,i}\right\rvert\geqslant\varepsilon\right\\}}\middle|\mathcal{F}_{n,i-1}\right]$
$\displaystyle=\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}_{\nu}^{T}\left[\Delta
M_{i}^{2}\mathbf{1}_{\left\\{\left\lvert\Delta
M_{i}\right\rvert\geqslant\varepsilon\sqrt{n}\right\\}}\middle|\mathcal{F}_{i-1}\right]$
$\displaystyle\leqslant\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}_{\nu}^{T}\left[\Delta
M_{i}^{2}\mathbf{1}_{\left\\{\left\lvert\Delta M_{i}\right\rvert\geqslant
A\right\\}}\middle|\mathcal{F}_{i-1}\right],$

for large enough values of $n$. Then, the Markov property gives us that

$\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}_{\nu}^{T}\left[\Delta
M_{i}^{2}\mathbf{1}_{\left\\{\left\lvert\Delta M_{i}\right\rvert\geqslant
A\right\\}}\middle|\mathcal{F}_{i-1}\right]=\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}^{T}_{Y_{i-1}}[\Delta
M_{1}^{2}\mathbf{1}_{\left\\{\left\lvert\Delta M_{1}\right\rvert\geqslant
A\right\\}}].$

Applying (35) on the right hand side with $g:y\mapsto\mathbb{E}^{T}_{y}[\Delta
M_{1}^{2}\mathbf{1}_{\left\\{\left\lvert\Delta M_{1}\right\rvert\geqslant
A\right\\}}]$ gives that

$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}^{T}_{Y_{i-1}}[\Delta
M_{1}^{2}\mathbf{1}_{\left\\{\left\lvert\Delta M_{1}\right\rvert\geqslant
A\right\\}}]\
\overset{\mathbb{P}^{\bar{Q}}_{{\xi^{\prime}}}-a.s.}{{\longrightarrow}}\ \mu
g=\mathbb{E}^{T}_{\mu}[\Delta M_{1}^{2}\mathbf{1}_{\left\\{\left\lvert\Delta
M_{1}\right\rvert\geqslant A\right\\}}].$

Now let $A\longrightarrow\infty$, $(\mathbb{E}^{T}_{\mu}[\Delta
M_{1}^{2}\mathbf{1}_{\left\\{\left\lvert\Delta M_{1}\right\rvert\geqslant
A\right\\}}])_{A}$ converges to 0 by monotone convergence. Hence

$\frac{M_{n}}{\sqrt{n}}\overset{\mathbb{P}^{T}_{\nu}-law}{\leadsto}\mathcal{N}(0,\sigma_{J}^{2}(h))\text{
with }\sigma_{J}^{2}(h)=\mathbb{E}^{T}_{\mu}\left[\Delta M_{1}^{2}\right].$

∎

####  A.3.2 Proof of Theorem 6

Let $h$ be a bounded measurable function and denote $h_{0}:=h-\pi h$. We want
to prove that :

$\frac{1}{\sqrt{n}}\sum_{i=0}^{n-1}h_{0}(X_{i})\overset{\mathbb{P}^{P}_{\chi}-law}{\leadsto}\mathcal{N}(0,\sigma^{2}(h)).$

Let $U_{n}=\frac{1}{\sqrt{n}}\sum_{i=0}^{n-1}h_{0}(X_{i})$. Let
${\xi^{\prime}}=\xi\otimes\mu$ where $\mu$ is any probability measure on
$\mathbb{N}$. The probability measure associated to the trajectories
$(X_{i})_{i}$ obtained from the Markov chains
$(\tilde{X}_{i},\tilde{N}_{i})_{i}$ generated by $\bar{Q}$ starting from
${\xi^{\prime}}$ is the same as the one associated to the sequence
$(X_{i})_{i}$ produced as the first component of the Markov chain
$(X_{i},N_{i})_{i}$ with kernel $P$ starting from $\chi$ defined by
$\chi(f)=\int\xi(\mathrm{d}x)S(x,\mathrm{d}x^{\prime})R(x^{\prime},n^{\prime})f(x^{\prime},n^{\prime})$.
We will work with the former probability distribution (i.e. under
$\mathbb{P}_{\xi^{\prime}}^{\bar{Q}}$) as it is best suited for our proof.
Denote
$V_{n}=\frac{1}{\sqrt{n}}\sum_{i=1}^{k_{n}-1}\tilde{N}_{i}h_{0}(\tilde{X}_{i})$
and let us start by proving that $U_{n}-V_{n}\
\overset{\mathbb{P}_{\xi^{\prime}}^{\bar{Q}}-\mathrm{prob}}{{\longrightarrow}}\
0$. By definition of $k_{n}$, $h_{0}(\tilde{X}_{k_{n}})$ appears less than
$\tilde{N}_{k_{n}}$ times in $U_{n}$ and therefore,

$\left\lvert U_{n}-V_{n}\right\rvert=(n-S_{k_{n}})\frac{\left\lvert
h_{0}(\tilde{X}_{k_{n}})\right\rvert}{\sqrt{n}}\leqslant\tilde{N}_{k_{n}}\frac{\left\lvert
h_{0}(\tilde{X}_{k_{n}})\right\rvert}{\sqrt{n}}.$

From Lemma 8, $\frac{k_{n}}{n}\
\overset{\mathbb{P}_{\xi^{\prime}}^{\bar{Q}}-\mathrm{prob}}{{\longrightarrow}}\
\kappa^{-1}$ and we can apply Lemma 9 to the function $(x,n)\mapsto
n\left\lvert h_{0}(x)\right\rvert$ to obtain
$\tilde{N}_{k_{n}}\frac{\left\lvert
h_{0}(\tilde{X}_{k_{n}})\right\rvert}{\sqrt{n}}\
\overset{\mathbb{P}_{\xi^{\prime}}^{\bar{Q}}-\mathrm{prob}}{{\longrightarrow}}\
0$. Let us write
$V_{n}=\frac{1}{\sqrt{n}}\sum_{i=1}^{k_{n}-1}f(\tilde{X}_{i},\tilde{N}_{i})$
where $f:(x,n)\mapsto nh_{0}(x)$. It now suffices to prove that
$V_{n}\overset{\mathbb{P}^{\bar{Q}}_{{\xi^{\prime}}}-law}{\leadsto}\mathcal{N}(0,\sigma^{2}(h))$.
We will procede in two steps :

  1. (i)

rewrite $V_{n}=\frac{1}{\sqrt{n}}M_{k_{n}}+\delta_{n}$ using a solution to the
Poisson equation associated to $f$, where $(M_{n})_{n\in\mathbb{N}}$ is a
martingale and $\delta_{n}\
\overset{\mathbb{P}_{\xi^{\prime}}^{\bar{Q}}-\mathrm{prob}}{{\longrightarrow}}\
0$ ;

  2. (ii)

prove that
$\frac{1}{\sqrt{n}}M_{n}\overset{\mathbb{P}^{\bar{Q}}_{{\xi^{\prime}}}-law}{\leadsto}\mathcal{N}(0,\sigma_{M}^{2}(h))$
and apply Theorem 8 to obtain
$\frac{1}{\sqrt{n}}M_{k_{n}}\overset{\mathbb{P}^{\bar{Q}}_{{\xi^{\prime}}}-law}{\leadsto}\mathcal{N}(0,\kappa^{-1}\sigma_{F}^{2}(h))$.

Starting with (i), let $H$ be the solution to the Poisson equation associated
to $\varrho h_{0}$ for the kernel $Q$ on ${\mathsf{X}}$ given by
(H${}_{\mbox{\tiny Poiss}}$), i.e. for all $x\in{\mathsf{X}}$,

$H(x)-QH(x)=\varrho(x)h_{0}(x).$

Then, $H_{\kappa}:=\kappa H$ is a solution to the Poisson equation associated
to $\varrho_{\kappa}h_{0}$ for the kernel $Q$ on ${\mathsf{X}}$, and

$F(x,n):=H_{\kappa}(x)+nh_{0}(x)-\varrho_{\kappa}(x)h_{0}(x)$ (36)

is a solution to the Poisson equation associated to $f$ for the Markov kernel
$\bar{Q}$ such that $\hat{\pi}F^{2}<\infty$. Indeed, for
$(x,n)\in{\mathsf{X}}\times\mathbb{N}$ we have

$\displaystyle\bar{Q}F(x,n)$
$\displaystyle=\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}Q(x,\mathrm{d}x^{\prime})\tilde{R}(x^{\prime},\mathrm{d}n^{\prime})H(x^{\prime})+\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}Q(x,\mathrm{d}x^{\prime})\tilde{R}(x^{\prime},\mathrm{d}n^{\prime})n^{\prime}h(x^{\prime})$
$\displaystyle\quad\quad-\int_{{\mathsf{X}}\times{\mathbb{N}_{0}}}Q(x,\mathrm{d}x^{\prime})\tilde{R}(x^{\prime},\mathrm{d}n^{\prime})\varrho_{\kappa}(x^{\prime})h(x^{\prime})$
$\displaystyle=QH_{\kappa}(x)$

since
$\int_{{\mathbb{N}_{0}}}\tilde{R}(x^{\prime},\mathrm{d}n^{\prime})n^{\prime}=\varrho_{\kappa}(x^{\prime})$,
and therefore

$\displaystyle
F(x,n)-\bar{Q}F(x,n)=H_{\kappa}(x)+nh_{0}(x)-\varrho_{\kappa}(x)h_{0}(x)-QH_{\kappa}(x)=nh_{0}(x)=f(x,n).$

Moreover, $\hat{\pi}F^{2}\leqslant
4\left(\kappa^{2}\hat{\pi}H^{2}+\hat{\pi}f^{2}+\hat{\pi}(\varrho_{\kappa}h_{0})^{2}\right)<\infty$
since $\hat{\pi}H^{2}=\tilde{\pi}H^{2}<\infty$ and $\hat{\pi}f^{2}<\infty$
from (H${}_{\mbox{\tiny Poiss}}$). Then,

$\displaystyle\hat{\pi}f^{2}$
$\displaystyle=\int_{\mathsf{X}}\left(\int_{\mathbb{N}_{0}}{n^{2}\tilde{R}(x,\mathrm{d}n)}\right)h_{0}(x)^{2}\tilde{\pi}(\mathrm{d}x)\geqslant\int_{\mathsf{X}}\left(\int_{\mathbb{N}_{0}}n\tilde{R}(x,\mathrm{d}n)\right)^{2}h_{0}(x)^{2}\tilde{\pi}(\mathrm{d}x)=\int_{{\mathsf{X}}}\varrho_{\kappa}(x)^{2}h_{0}(x)^{2}\tilde{\pi}(\mathrm{d}x),$

proving that
$\hat{\pi}\left(\varrho_{\kappa}h_{0}\right)^{2}=\tilde{\pi}\left(\varrho_{\kappa}h_{0}\right)^{2}<\infty$.
Now let
$M_{n}=\sum_{i=2}^{n}F(\tilde{X}_{i},\tilde{N}_{i})-\bar{Q}F(\tilde{X}_{i-1},\tilde{N}_{i-1})$
and consider the filtration
$\mathcal{F}_{n}=\sigma(\tilde{X}_{1:n},\tilde{N}_{1:n})$. From Lemma 10,
$(M_{n})_{n}$ is a $\mathcal{F}_{n}$-martingale, and

$\sum_{i=1}^{n-1}f(\tilde{X}_{i},\tilde{N}_{i})=M_{n}-F(\tilde{X}_{n},\tilde{N}_{n})+F(\tilde{X}_{1},\tilde{N}_{1}).$
(37)

Therefore, $V_{n}=\frac{1}{\sqrt{n}}M_{k_{n}}+\delta_{n}$ with
$\delta_{n}=\frac{F(\tilde{X}_{k_{n}},\tilde{N}_{k_{n}})}{\sqrt{n}}+\frac{F(\tilde{X}_{1},\tilde{N}_{1})}{\sqrt{n}}$.
Note that $\frac{F(\tilde{X}_{1},\tilde{N}_{1})}{\sqrt{n}}\
\overset{\mathbb{P}_{\xi^{\prime}}^{\bar{Q}}-\mathrm{prob}}{{\longrightarrow}}\
0$ trivially, and $\frac{F(\tilde{X}_{k_{n}},\tilde{N}_{k_{n}})}{\sqrt{n}}\
\overset{\mathbb{P}_{\xi^{\prime}}^{\bar{Q}}-\mathrm{prob}}{{\longrightarrow}}\
0$ as a consequence of Lemma 9. Using Lemma 11 combined to Lemma 7 with
$Y_{i}=(\tilde{X}_{i},\tilde{N}_{i})$,

$\frac{M_{n}}{\sqrt{n}}\overset{\mathbb{P}^{\bar{Q}}_{{\xi^{\prime}}}-law}{\leadsto}\mathcal{N}(0,\sigma_{F}^{2}(h))\text{
with }\sigma_{F}^{2}(h)=\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[\Delta
M_{1}^{2}\right].$

which proves (ii) and concludes the first part of the proof. Let us now turn
our attention to the expression of the variance claimed by the theorem :

$\sigma^{2}(h)=\kappa\tilde{\sigma}^{2}(\varrho
h_{0})+\kappa^{-1}\hat{\sigma}^{2}(h_{0},\kappa),$

with $\tilde{\sigma}^{2}(\varrho h_{0})=2\tilde{\pi}\left(\varrho
h_{0}H\right)-\tilde{\pi}\left((\varrho h_{0})^{2}\right)$ and
$\hat{\sigma}^{2}(h_{0},\kappa)=\int_{{\mathsf{X}}}h_{0}(x)^{2}\mathbb{V}\mathrm{ar}^{\tilde{R}(x,\cdot)}[N]\tilde{\pi}(\mathrm{d}x)$.
From the expression of $F$ in (36) and denoting $\Delta
H_{1}=H(\tilde{X}_{1})-QH(\tilde{X}_{0})$ we have

$\displaystyle\sigma^{2}(h)$
$\displaystyle=\kappa^{-1}\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[\left(\kappa\Delta
H_{1}+(\tilde{N}_{1}-\varrho_{\kappa}(\tilde{X}_{1}))h(\tilde{X}_{1})\right)^{2}\right]$
$\displaystyle=\kappa\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[\Delta
H_{1}^{2}\right]+\kappa^{-1}\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[(\tilde{N}_{1}-\varrho_{\kappa}(\tilde{X}_{1}))^{2}h(\tilde{X}_{1})^{2}\right]+2\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[\Delta
H_{1}(\tilde{N}_{1}-\varrho_{\kappa}(\tilde{X}_{1}))h(\tilde{X}_{1})\right].$
(38)

Let us take a look at the first term of the rhs :

$\displaystyle\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[\Delta H_{1}^{2}\right]$
$\displaystyle=\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[H(\tilde{X}_{1})^{2}\right]+\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[QH(\tilde{X}_{0})^{2}\right]-2\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[H(\tilde{X}_{1})QH(\tilde{X}_{0})\right]$
$\displaystyle=\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[H(\tilde{X}_{0})^{2}\right]-\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[QH(\tilde{X}_{0})^{2}\right],$

where we used that $\hat{\pi}$ is a stationary probability measure. Noting
that $H^{2}-QH^{2}=(H-QH)(H+QH)=\varrho h_{0}(2H-\varrho h_{0})$, we finally
obtain using (H${}_{\mbox{\tiny Poiss}}$)

$\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[\Delta
H_{1}^{2}\right]=2\tilde{\pi}\left(\varrho
h_{0}H\right)-\tilde{\pi}\left((\varrho h_{0})^{2}\right).$ (39)

Let us now rewrite the second term of (38) :

$\displaystyle\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[(\tilde{N}_{1}-\varrho_{\kappa}(\tilde{X}_{1}))^{2}h_{0}(\tilde{X}_{1})^{2}\right]$
$\displaystyle=\int_{\mathsf{X}}\left(\int_{\mathbb{N}_{0}}(n-\varrho_{\kappa}(x))^{2}\tilde{R}(x,\mathrm{d}n)\right)h_{0}(x)^{2}\tilde{\pi}(\mathrm{d}x)$
$\displaystyle=\int_{{\mathsf{X}}}h_{0}(x)^{2}\mathbb{V}\mathrm{ar}^{\tilde{R}(x,\cdot)}[N]\tilde{\pi}(\mathrm{d}x),$

where
$\mathbb{V}\mathrm{ar}^{\tilde{R}(x,\cdot)}[N]=\int(n-\varrho_{\kappa}(x))^{2}\tilde{R}(x,\mathrm{d}n)$
due to (H2). For the last term of 38, the same hypothesis gives that

$\displaystyle\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[\Delta
H_{1}(\tilde{N}_{1}-\varrho_{\kappa}(\tilde{X}_{1}))h_{0}(\tilde{X}_{1})\right]=\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[{\mathbb{E}^{\bar{Q}}_{\hat{\pi}}\left[\tilde{N}_{1}-\varrho_{\kappa}(\tilde{X}_{1})|\tilde{X}_{1},\tilde{X}_{0}\right]}\Delta
H_{1}h_{0}(\tilde{X}_{1})\right]=0,$

which concludes the second part of the proof.

###  A.4 Geometric ergodicity

####  A.4.1 Proof of Lemma 2 and Lemma 3

###### Proof of Lemma 2.

We start with the first point (i). From the definition of $S$ in (9) and
noting that $\mathsf{C}_{\eta}$ is a $(1,\varepsilon\nu)$-small set for $Q$,
we have for all $x\in\mathsf{C}_{\eta}$ and $\mathsf{A}\in\mathcal{X}$,

$\displaystyle S(x,\mathsf{A})$
$\displaystyle=\sum_{k=1}^{\infty}\mathbb{E}_{x}^{Q}\left[\rho_{\tilde{R}}(X_{k})\mathbf{1}_{\mathsf{A}}(X_{k})\prod_{i=1}^{k-1}\left(1-\rho_{\tilde{R}}(X_{i})\right)\right]$
$\displaystyle\geqslant
Q(x,\rho_{\tilde{R}}\mathbf{1}_{\mathsf{A}})\geqslant\varepsilon\nu(\rho_{\tilde{R}}\mathbf{1}_{\mathsf{A}}).$

Applying (10) with $n=0$, we deduce that for all $x\in\mathsf{C}_{\eta}$ and
all $\mathsf{B}\in\mathcal{X}\otimes\mathcal{P}({\mathbb{N}_{0}})$,

$\displaystyle P(x,0;\mathsf{B})$
$\displaystyle=\int_{\mathsf{B}}S(x,\mathrm{d}x^{\prime})R(x^{\prime},\mathrm{d}u^{\prime})$
$\displaystyle\geqslant\varepsilon\int_{\mathsf{B}}\nu(\mathrm{d}x^{\prime})\rho_{\tilde{R}}(x^{\prime})R(x^{\prime},\mathrm{d}u^{\prime})$
$\displaystyle=:\varepsilon\tilde{\nu}(\mathsf{B}),$

which shows (i). Let us now turn to (ii). Noting that
$\mathsf{C}_{\eta}^{+}:=\mathsf{C}_{\eta}\cap\\{\tilde{R}(.,1)>0\\}=\mathsf{C}_{\eta}\cap\\{R(.,0)>0\\}$
we have,

$\displaystyle\tilde{\nu}(\mathsf{C}_{\eta}\times\\{0\\})$
$\displaystyle=\int_{\mathsf{C}_{\eta}}\nu(\mathrm{d}x)\rho_{\tilde{R}}(x)R(x,0)$
$\displaystyle\geqslant\eta\int_{\mathsf{C}_{\eta}}\nu(\mathrm{d}x)R(x,0)$
$\displaystyle=\eta\int_{\mathsf{C}_{\eta}^{+}}\nu(\mathrm{d}x)R(x,0)>0,$

where the last inequality stems from $\nu(\mathsf{C}_{\eta}^{+})>0$ and
$R(x,0)>0$ for all $x\in\mathsf{C}_{\eta}^{+}$. Hence (ii). ∎

###### Proof of Lemma 3.

Let $\mathsf{A}$ be accessible for $Q$ such that
$\epsilon_{\mathsf{A}}:=\inf_{x\in\mathsf{A}}\rho_{\tilde{R}}(x)>0$. We first
show that $\mathsf{A}$ is accessible for $S$. Let $x\in{\mathsf{X}}$ and
$n\in{\mathbb{N}_{0}}$ such that $Q^{n}(x,\mathsf{A})>0$. Let us consider the
representation of $S$ using the kernel $G$ defined in (2). Define
$\mathsf{D}:=\left\\{(x,u)\in{\mathsf{X}}\times[0,1]:\,u\leqslant\rho_{\tilde{R}}(x)\right\\}$
and $\sigma_{\mathsf{D}}^{(m)}$ the $m$-th return time to the set
$\mathsf{D}$. Then, for any probability measure $\mu$ on $[0,1]$,

$\displaystyle\int_{\mathsf{X}}Q^{n}(x,\mathrm{d}y)\rho_{\tilde{R}}(y)\mathbf{1}_{\mathsf{A}}(y)$
$\displaystyle=\mathbb{P}^{G}_{\delta_{x}\otimes\mu}((X_{n},U_{n})\in\mathsf{D},X_{n}\in\mathsf{A})$
$\displaystyle=\mathbb{P}^{G}_{\delta_{x}\otimes\mu}(\exists
m\in[1:n],n=\sigma_{\mathsf{D}}^{(m)},X_{\sigma_{\mathsf{D}}^{(m)}}\in\mathsf{A})$
$\displaystyle\leqslant\mathbb{P}^{G}_{\delta_{x}\otimes\mu}(\exists
m\in[1:n],X_{\sigma_{\mathsf{D}}^{(m)}}\in\mathsf{A})$
$\displaystyle=\mathbb{P}^{S}_{x}(\exists m\in[1:n],X_{m}\in\mathsf{A})$
$\displaystyle\leqslant\sum_{m=1}^{n}\mathbb{P}^{S}_{x}(X_{m}\in\mathsf{A}).$

Hence,
$\sum_{m=1}^{n}S^{m}(x,\mathsf{A})\geqslant\epsilon_{\mathsf{A}}Q^{n}(x,\mathsf{A})>0$
and there exists an integer $m\leqslant n$ such that $S^{m}(x,\mathsf{A})>0$.
Thus, $A$ is an accessible set for $S$.

We now show that $\mathsf{A}\times\\{0\\}$ is accessible for $P$. Let
$(x,k)\in{\mathsf{X}}\times{\mathbb{N}_{0}}$. From the first part of the
proof, there exists $m\in{\mathbb{N}_{0}}$ such that $S^{m}(x,\mathsf{A})>0$.
Indeed, let $\mathsf{B}:=\mathsf{A}\times\\{0\\}$,
$\mathsf{F}:={\mathsf{X}}\times\\{0\\}$ and
$(Y_{n}=(X_{n},N_{n}))_{n\in{\mathbb{N}_{0}}}$ be a Markov chain of kernel
$P$. Start by noting that if $x\in{\mathsf{X}}$ and $N\sim R(x,\cdot)$,

$\mathbb{P}^{P}_{(x,0)}\left(\sigma_{\mathsf{F}}<\infty\right)=\mathbb{P}\left(N<\infty\right)=1$

since (H3) ensures that $\mathbb{E}\left[N\right]<\infty$. Then, using the
strong Markov inequality, we obtain by induction that

$\mathbb{P}^{P}_{(x,0)}\left(\sigma^{(m)}_{\mathsf{F}}<\infty\right)=1.$

Hence,

$\displaystyle\mathbb{P}^{P}_{(x,0)}(Y_{\sigma_{\mathsf{F}}^{(m)}}\in\mathsf{B})$
$\displaystyle=\mathbb{P}^{P}_{(x,0)}(Y_{\sigma_{\mathsf{F}}^{(m)}}\in\mathsf{B},\sigma_{\mathsf{F}}^{(m)}<\infty)$
$\displaystyle=\sum_{\ell=1}^{\infty}\mathbb{P}^{P}_{(x,0)}(Y_{\sigma_{\mathsf{F}}^{(m)}}\in\mathsf{B},\sigma_{\mathsf{F}}^{(m)}=\ell)$
$\displaystyle=\sum_{\ell=1}^{\infty}\mathbb{P}^{P}_{(x,0)}(Y_{\ell}\in\mathsf{B},\sigma_{\mathsf{F}}^{(m)}=\ell)$
$\displaystyle\leqslant\sum_{\ell=1}^{\infty}\mathbb{P}^{P}_{(x,0)}(Y_{\ell}\in\mathsf{B})$
$\displaystyle=\sum_{\ell=1}^{\infty}P^{\ell}(x,0;\mathsf{B}).$

Since by definition of $P$ in (10) we have

$\mathbb{P}^{P}_{(x,0)}(Y_{\sigma_{\mathsf{F}}^{(m)}}\in\mathsf{B})=\mathbb{P}^{S}_{x}(X_{m}\in\mathsf{A})>0,$

at least one of the terms in the sum above is positive and therefore there
exists $\ell\in{\mathbb{N}_{0}}$ such that

$P^{\ell}(x,0\,;\mathsf{A}\times\\{0\\})>0.$

Using the definition of $P$ in (10) once more, we have that
$P^{k}(x,k;\\{(x,0)\\})=1$ and so

$P^{k+\ell}(x,k;\mathsf{A}\times\\{0\\})>0,$

which concludes the proof.

∎

####  A.4.2 Proof of Lemma 4

###### (Proof of Lemma 4).

Let $\mathsf{B}:={\mathsf{C}_{\eta}\times\\{0\\}}$ and
$\mathsf{F}:={\mathsf{X}}\times\\{0\\}$. Let $\beta\in(1,\infty)$ be an
arbitrary constant and let $D<\infty$ be any positive constant (assuming it
exists) such that

$\sup_{x\in{\mathsf{X}}}\int_{\mathbb{N}_{0}}\beta^{n+1}R(x,\mathrm{d}n)\leqslant
D.$

We first show that for all $(x,n)\in{\mathsf{X}}\times{\mathbb{N}_{0}}$ and
$\beta>1$ we have:

$\mathbb{E}^{P}_{(x,n)}\left[\beta^{\sigma_{\mathsf{B}}}\right]\leqslant\beta^{n}\mathbb{E}_{x}^{S}\left[D^{\sigma_{\mathsf{C}_{\eta}}}\right].$
(40)

Let us start with the case $n=0$ :

$\mathbb{E}_{(x,0)}^{P}\left[\beta^{\sigma_{\mathsf{B}}}\right]=\sum_{\ell=1}^{\infty}g^{(\ell)}(x),$
(41)

where
$g^{(\ell)}(x):=\mathbb{E}_{(x,0)}^{P}\left[\beta^{\sigma_{\mathsf{F}}^{(\ell)}}\mathbf{1}_{\\{\sigma_{\mathsf{F}}^{(\ell)}=\sigma_{\mathsf{B}}\\}}\right]$.

  * •

For $\ell=1$,

$g^{(1)}(x)=\mathbb{E}_{(x,0)}^{P}\left[\beta^{N_{1}+1}\mathbf{1}_{\mathsf{C}_{\eta}}(X_{N_{1}+1})\right]=\int_{\mathsf{X}}\left(\int_{\mathbb{N}_{0}}\beta^{n+1}R(x,\mathrm{d}n)\right)\mathbf{1}_{\mathsf{C}_{\eta}}(x^{\prime})S(x,\mathrm{d}x^{\prime})\leqslant
D\times\mathbb{P}_{x}^{S}\left(\sigma_{\mathsf{C}_{\eta}}=1\right),$ (42)

where the second equality holds from the definition of $P$ in (10) ensuring
that $X_{N_{1}+1}=X_{1}\quad$ $\mathbb{P}_{(x,0)}^{P}$-a.s.

  * •

For $\ell>1$, let us note that
$\left\\{\sigma_{\mathsf{F}}^{(\ell)}=\sigma_{B}\right\\}=\left\\{\sigma_{\mathsf{F}}^{(\ell-1)}\circ\theta^{\sigma_{\mathsf{F}}}=\sigma_{B}\circ\theta^{\sigma_{\mathsf{F}}},X_{\sigma_{\mathsf{F}}}\notin\mathsf{C}_{\eta}\right\\}$
as $\mathbb{P}_{(x,0)}^{P}$-a.s. we have
$\sigma_{\mathsf{F}}^{(\ell)}=\sigma_{\mathsf{F}}^{(\ell-1)}\circ\theta^{\sigma_{\mathsf{F}}}+\sigma_{\mathsf{F}}$,
and
$\sigma_{B}=\sigma_{B}\circ\theta^{\sigma_{\mathsf{F}}}+\sigma_{\mathsf{F}}$
under the event
$\left\\{\sigma_{B}>\sigma_{\mathsf{F}}\right\\}\supset\\{\sigma_{\mathsf{F}}^{(\ell)}=\sigma_{B}\\}$
for $\ell>1$. Hence,

$\displaystyle g^{(\ell)}(x)$
$\displaystyle=\mathbb{E}_{(x,0)}^{P}\left[\beta^{\sigma_{\mathsf{F}}^{(\ell)}}\mathbf{1}_{\\{\sigma_{\mathsf{F}}^{(\ell)}=\sigma_{\mathsf{B}}\\}}\right]$
$\displaystyle\leqslant\mathbb{E}_{(x,0)}^{P}\left[\beta^{\sigma_{\mathsf{F}}}\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{\sigma_{\mathsf{F}}})\mathbb{E}_{(x,0)}^{P}\left[\beta^{\sigma_{\mathsf{F}}^{(\ell-1)}\circ\theta^{\sigma_{\mathsf{F}}}}\mathbf{1}_{\\{\sigma_{\mathsf{F}}^{(\ell-1)}\circ\theta^{\sigma_{\mathsf{F}}}=\sigma_{\mathsf{B}}\circ\theta^{\sigma_{\mathsf{F}}}\\}}|\mathcal{F}_{\sigma_{\mathsf{F}}}\right]\right]$
$\displaystyle=\mathbb{E}_{(x,0)}^{P}\left[\beta^{\sigma_{\mathsf{F}}}\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{\sigma_{\mathsf{F}}})\mathbb{E}_{(X_{\sigma_{\mathsf{F}}},0)}^{P}\left[\beta^{\sigma_{\mathsf{F}}^{(\ell-1)}}\mathbf{1}_{\\{\sigma_{\mathsf{F}}^{(\ell-1)}=\sigma_{\mathsf{B}}\\}}\right]\right]$
(43)
$\displaystyle=\mathbb{E}_{(x,0)}^{P}\left[\beta^{N_{1}+1}\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{1})\mathbb{E}_{(X_{1},0)}^{P}\left[\beta^{\sigma_{\mathsf{F}}^{(\ell-1)}}\mathbf{1}_{\\{\sigma_{\mathsf{F}}^{(\ell-1)}=\sigma_{\mathsf{B}}\\}}\right]\right],$
(44)

where (43) comes from the strong Markov property applied to the Markov chain
$(X_{i},N_{i})_{i\in{\mathbb{N}_{0}}}$ with the stopping time
$\sigma_{\mathsf{F}}$ and (44) comes from the definition of the kernel $P$
since $X_{1}$ is repeated $N_{1}+1$ times while the second component decreases
by one at each iteration until reaching zero. The definition of $P$ in (10)
gives

$\displaystyle\mathbb{E}_{(x,0)}^{P}\left[\beta^{N_{1}+1}\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{1})\mathbb{E}_{(X_{1},0)}^{P}\left[\beta^{\sigma_{\mathsf{F}}^{(\ell-1)}}\mathbf{1}_{\\{\sigma_{\mathsf{F}}^{(\ell-1)}=\sigma_{\mathsf{B}}\\}}\right]\right]$
$\displaystyle=\int_{\mathsf{X}}\left(\int_{\mathbb{N}_{0}}\beta^{n+1}R(x,\mathrm{d}n)\right)\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(x^{\prime})g^{(\ell-1)}(x^{\prime})S(x,\mathrm{d}x^{\prime})$
$\displaystyle\leqslant
D\int_{\mathsf{X}}\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(x^{\prime})g^{(\ell-1)}(x^{\prime})S(x,\mathrm{d}x^{\prime})$
$\displaystyle=D\times\mathbb{E}_{(x,0)}^{P}\left[\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{1})g^{(\ell-1)}(X_{1})\right].$

Hence,

$g^{(\ell)}(x)=\mathbb{E}_{(x,0)}^{P}\left[\beta^{\sigma_{\mathsf{F}}^{(\ell)}}\mathbf{1}_{\\{\sigma_{\mathsf{F}}^{(\ell)}=\sigma_{\mathsf{B}}\\}}\right]\leqslant
D\times\mathbb{E}_{(x,0)}^{P}\left[\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{\sigma_{\mathsf{F}}})g^{(\ell-1)}\left(X_{\sigma_{\mathsf{F}}}\right)\right],$

which used in conjunction with (42) gives

$\displaystyle g^{(\ell)}(x)$ $\displaystyle\leqslant
D^{\ell}\times\mathbb{E}_{(x,0)}^{P}\left[\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{\sigma_{\mathsf{F}}})...\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{\sigma^{(\ell-1)}_{\mathsf{F}}})\mathbb{P}_{X_{{\sigma_{\mathsf{F}}}^{(\ell-1)}}}^{S}\left(\sigma_{\mathsf{C}_{\eta}}=1\right)\right]$
$\displaystyle=D^{\ell}\times\mathbb{E}_{x}^{S}\left[\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{1})...\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{\ell-1})\mathbb{E}_{X_{\ell-1}}^{S}\left[\mathbf{1}_{\mathsf{C}_{\eta}}(X_{1})\right]\right]$
(45)
$\displaystyle=D^{\ell}\times\mathbb{E}_{x}^{S}\left[\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{1})...\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{\ell-1})\mathbb{E}_{x}^{S}\left[\mathbf{1}_{\mathsf{C}_{\eta}}(X_{\ell})|\mathcal{F}_{\ell-1}\right]\right]$
$\displaystyle=D^{\ell}\times\mathbb{E}_{x}^{S}\left[\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{1})...\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{\ell-1})\mathbf{1}_{\mathsf{C}_{\eta}}(X_{\ell})\right]$
$\displaystyle=D^{\ell}\times\mathbb{P}_{x}^{S}\left(\sigma_{\mathsf{C}_{\eta}}=\ell\right),$

where (45) comes from the definition of $S$ and $P$ in (9) and (10). Plugging
the above into (41) leads to the following upper bound,

$\displaystyle\mathbb{E}_{(x,0)}^{P}\left[\beta^{\sigma_{\mathsf{B}}}\right]\leqslant\sum_{\ell=1}^{\infty}D^{\ell}\mathbb{P}_{x}^{S}\left(\sigma_{\mathsf{C}_{\eta}}=\ell\right)=\mathbb{E}_{x}^{S}\left[D^{\sigma_{\mathsf{C}_{\eta}}}\right].$

If $x\in\mathsf{C}_{\eta}$, then
$\mathbb{E}^{P}_{(x,n)}\left[\beta^{\sigma_{\mathsf{B}}}\right]=\beta^{n}$,
showing (40) for $x\in\mathsf{C}_{\eta}$. If $x\notin\mathsf{C}_{\eta}$, then
we have $\sigma_{\mathsf{B}}=\sigma_{\mathsf{B}}\circ\theta^{n}+n$,
$\mathbb{P}_{(x,n)}$-a.s. Thus,

$\displaystyle\mathbb{E}^{P}_{(x,n)}\left[\beta^{\sigma_{\mathsf{B}}}\right]$
$\displaystyle=\beta^{n}\mathbb{E}^{P}_{(x,n)}\left[\beta^{\sigma_{\mathsf{B}}\circ\theta^{n}}\right]=\beta^{n}\mathbb{E}^{P}_{(x,0)}\left[\beta^{\sigma_{\mathsf{B}}}\right]\leqslant\beta^{n}\mathbb{E}_{x}^{S}\left[D^{\sigma_{\mathsf{C}_{\eta}}}\right],$

which completes the proof of the inequality (40).

For any $\eta\in(0,\eta_{0})$,

$\mathbb{P}_{x}^{S}\left(\sigma_{\mathsf{C}_{\eta}}>k\right)=\sum_{\ell_{1},\dots,\ell_{k}=1}^{\infty}\mathbb{E}_{x}^{Q}\left[\prod_{j=1}^{k}\left[\prod_{i=s_{j-1}+1}^{s_{j}-1}\left(1-\rho_{\tilde{R}}(X_{i})\right)\right]\rho_{\tilde{R}}(X_{s_{j}})\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{s_{j}})\right],$

where we have set $s_{j}=\sum_{i=1}^{j}\ell_{i}$ for $j\geqslant 1$ and
$s_{0}=0$. Using that
$\rho_{\tilde{R}}\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}\leqslant\eta\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}$
and
$1-\rho_{\tilde{R}}\leqslant(1-\eta_{0})^{\mathbf{1}_{\mathsf{C}_{\eta_{0}}}}$,

$\displaystyle\mathbb{P}_{x}^{S}\left(\sigma_{\mathsf{C}_{\eta}}>k\right)$
$\displaystyle\leqslant\sum_{\ell_{1},\dots,\ell_{k}=1}^{\infty}\eta^{k}\mathbb{E}_{x}^{Q}\left[(1-\eta_{0})^{\sum_{j=1}^{k}\sum_{i=s_{j-1}+1}^{s_{j}-1}\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(X_{i})}\prod_{j=1}^{k}\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{s_{j}})\right]$
$\displaystyle\leqslant\sum_{\ell_{1},\dots,\ell_{k}=1}^{\infty}\eta^{k}\mathbb{E}_{x}^{Q}\left[(1-\eta_{0})^{M_{[1:s_{k}-1]\setminus\\{s_{1},\ldots,s_{k-1}\\}}}\prod_{j=1}^{k-1}\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(X_{s_{j}})\right],$
(46)

where $M_{I}:=\sum_{i\in I}\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(X_{i})$ for any
$I\subset{\mathbb{N}_{0}}$. Moreover, since $\eta<\eta_{0}$, we have
$\mathsf{C}_{\eta_{0}}\subset\mathsf{C}_{\eta}$, hence

$\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(x)=(1-\eta_{0})^{\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(x)}\mathbf{1}_{\bar{\mathsf{C}}_{\eta}}(x)\leqslant(1-\eta_{0})^{\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(x)}.$
(47)

Finally, setting for any arbitrary $\alpha\in(0,1)$,

$\displaystyle
N_{\ell}=\sum_{i=1}^{\ell-1}\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(X_{i}),$
$\displaystyle
R_{\ell}:=\mathbb{E}_{x}^{Q}\left[(1-\eta_{0})^{N_{\ell}}\right],$
$\displaystyle
R_{\ell,1}=\mathbb{E}_{x}^{Q}\left[(1-\eta_{0})^{N_{\ell}}\mathbf{1}_{\\{N_{\ell}>\alpha(\ell-k)\\}}\right],$
$\displaystyle
R_{\ell,2}=\mathbb{E}_{x}^{Q}\left[(1-\eta_{0})^{N_{\ell}}V(X_{\ell})\mathbf{1}_{\\{N_{\ell}\leqslant\alpha(\ell-k)\\}}\right],$

and plugging (47) into (46) combined with $V\geqslant 1$, we get

$\mathbb{P}_{x}^{S}\left(\sigma_{\mathsf{C}_{\eta}}>k\right)\leqslant\sum_{\ell_{1},\dots,\ell_{k}=1}^{\infty}\eta^{k}R_{s_{k}}\leqslant\sum_{\ell_{1},\dots,\ell_{k}=1}^{\infty}\eta^{k}\left\\{R_{s_{k},1}+R_{s_{k},2}\right\\}\leqslant\sum_{\ell_{1},\dots,\ell_{k}=1}^{\infty}\eta^{k}\left[(1-\eta_{0})^{\alpha(s_{k}-k)}+R_{s_{k},2}\right].$
(48)

We now give an explicit upper bound for $R_{\ell,2}$ for $\ell\geqslant k$.
Under (H${}_{\mbox{\tiny dft}}$), $QV(x)\leqslant\lambda V(x)$ for
$x\notin\mathsf{C}_{\eta_{0}}$ and $QV(x)\leqslant b_{\infty}V(x)$ for
$x\in\mathsf{C}_{\eta_{0}}$. Therefore, for any $x\in{\mathsf{X}}$,

$(1-\eta_{0})^{\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(x)}QV(x)\leqslant(b_{\infty}(1-\eta_{0})\lambda^{-1})^{\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(x)}\lambda
V(x)\leqslant A^{\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(x)}\lambda V(x),$

where $A:=1\vee b_{\infty}(1-\eta_{0})\lambda^{-1}$. This implies by applying
the tower rule conditionally on $X_{1:\ell-1}$, then $X_{1:\ell-2}$ and so on,

$\mathbb{E}_{x}^{Q}\left[V(X_{\ell})\frac{(1-\eta_{0})^{N_{\ell}}}{A^{N_{\ell}}\lambda^{\ell-1}}\right]=\mathbb{E}_{x}^{Q}\left[V(X_{1})\prod_{i=1}^{\ell-1}\frac{(1-\eta_{0})^{\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(X_{i})}V(X_{i+1})}{A^{\mathbf{1}_{\mathsf{C}_{\eta_{0}}}(X_{i})}\lambda
V(X_{i})}\right]\leqslant\mathbb{E}_{x}^{Q}\left[V(X_{1})\right]=QV(x)\leqslant
b_{\infty}V(x),$ (49)

where the last inequality follows from (21). Since $A\geqslant 1$, we have
$1\leqslant A^{\alpha(\ell-k)}/A^{N_{\ell}}$ on
$\\{N_{\ell}\leqslant\alpha(\ell-k)\\}$ and hence

$R_{\ell,2}\leqslant
A^{\alpha(\ell-k)}\lambda^{\ell-1}\mathbb{E}_{x}^{Q}\left[V(X_{\ell})\frac{(1-\eta_{0})^{N_{\ell}}}{A^{N_{\ell}}\lambda^{\ell-1}}\right]\leqslant
A^{\alpha(\ell-k)}\lambda^{\ell-1}b_{\infty}V(x),$

where the last inequality follows from (49). Pick $\alpha$ small enough so
that $\lambda A^{\alpha}<1$. Plugging the inequality above (with $\ell$
replaced by $s_{k}=\ell_{1}+\cdots+\ell_{k}$) into (48) yields

$\displaystyle\mathbb{P}_{x}^{S}\left(\sigma_{\mathsf{C}_{\eta}}>k\right)$
$\displaystyle\leqslant\sum_{\ell_{1},\dots,\ell_{k}=1}^{\infty}\eta^{k}\left[(1-\eta_{0})^{\alpha(\ell_{1}+\dots+\ell_{k}-k)}+\left[\lambda
A^{\alpha}\right]^{\ell_{1}+\dots+\ell_{k}}A^{-\alpha
k}\lambda^{-1}b_{\infty}V(x)\right]$
$\displaystyle=\eta^{k}\left[\left(\frac{(1-\eta_{0})^{\alpha}}{1-(1-\eta_{0})^{\alpha}}\right)^{k}\frac{1}{(1-\eta_{0})^{\alpha
k}}+\left(\frac{\lambda A^{\alpha}}{1-\lambda
A^{\alpha}}\right)^{k}\frac{1}{A^{\alpha k}}\lambda^{-1}b_{\infty}V(x)\right]$
$\displaystyle=\eta^{k}\left[\left(\frac{1}{1-(1-\eta_{0})^{\alpha}}\right)^{k}+\left(\frac{\lambda}{1-\lambda
A^{\alpha}}\right)^{k}\lambda^{-1}b_{\infty}V(x)\right].$

Now set
$\gamma:=\max\left(\frac{1}{1-(1-\eta_{0})^{\alpha}},\frac{\lambda}{1-\lambda
A^{\alpha}}\right)$ and choose $\eta<\eta_{0}$ sufficiently small so that
$\eta\gamma<1$. Then,

$\mathbb{P}_{x}^{S}\left(\sigma_{\mathsf{C}_{\eta}}>k\right)\leqslant\eta^{k}\gamma^{k}(1+\lambda^{-1}b_{\infty}V(x)),$

and if $D\in(1,\eta^{-1}\gamma^{-1})$,

$\displaystyle\mathbb{E}_{x}^{S}\left[\frac{D^{\sigma_{\mathsf{C}_{\eta}}}-1}{D-1}\right]$
$\displaystyle=\mathbb{E}_{x}^{S}\left[\sum_{k=0}^{\infty}D^{k}\mathbf{1}_{\\{k<\sigma_{\mathsf{C}_{\eta}}\\}}\right]=\sum_{k=0}^{\infty}D^{k}\mathbb{P}_{x}^{S}\left(\sigma_{\mathsf{C}_{\eta}}>k\right)$
$\displaystyle\leqslant\sum_{k=0}^{\infty}D^{k}\gamma^{k}\eta^{k}(1+\lambda^{-1}b_{\infty}V(x))=\frac{1+\lambda^{-1}b_{\infty}V(x)}{1-D\gamma\eta}.$

From (H3) there exists $\beta_{0}\in(1,\infty)$ and $D_{0}<\infty$ such that

$\sup_{x\in{\mathsf{X}}}\int_{\mathbb{N}_{0}}\beta_{0}^{n+1}R(x,\mathrm{d}n)\leqslant
D_{0}.$

Let $r\in(0,1)$ and consider $\beta_{r}=\beta_{0}^{r}$. From Hölder’s
inequality,

$\int_{\mathbb{N}_{0}}\beta_{r}^{n+1}R(x,\mathrm{d}n)\leqslant\left(\int_{\mathbb{N}_{0}}\beta_{0}^{n+1}R(x,\mathrm{d}n)\right)^{r}\leqslant
D_{0}^{r}.$ (50)

Choose $r$ such that $D_{r}:=D_{0}^{r}\in(1,\eta^{-1}\gamma^{-1})$. We can now
apply the above inequality in combination with (40) using the couple
$(\beta_{r},D_{r})$ instead of $(\beta,D)$ and obtain for all
$(x,n)\in{\mathsf{X}}\times{\mathbb{N}_{0}}$,

$\displaystyle\mathbb{E}^{P}_{(x,n)}\left[\beta_{r}^{\sigma_{\mathsf{B}}}\right]\leqslant\beta_{r}^{n}\mathbb{E}^{S}_{x}\left[D_{r}^{\sigma_{\mathsf{C}_{\eta}}}\right]$
$\displaystyle\leqslant\beta_{r}^{n}\left[1+(D_{r}-1)\frac{1+\lambda^{-1}b_{\infty}V(x)}{1-D_{r}\gamma\eta}\right]\leqslant\beta_{\star}\beta_{r}^{n}V(x)$

since $V\geqslant 1$, and where

$\beta_{\star}=1+(D_{r}-1)\frac{1+\lambda^{-1}b_{\infty}}{1-D_{r}\gamma\eta}.$

The proof is completed.

∎

##  Appendix B Numerical experiments

###  B.1 Details of the hyperparameters for the normalizing flows

The normalizing flow is a RQSpline with $10$ layers, $8$ bins, and a
$(128,128)$ hidden size. The local sampler is a MALA algorithm with step size
$0.1$. Training consists of a total of 30 loops with a unique epoch each time.
10 global steps are implemented with a further 10 local steps between each,
and the optimizer (Adam) has a learning rate of $8\cdot 10^{-4}$ and a
momentum of 0.9. The seed is 1250. The code generating the figures is
available at https://github.com/charlyandral/importance_markov_chain.

## References

  * [1] N. Metropolis, S. Ulam, The Monte Carlo Method, Journal of the American Statistical Association 44 (247) (1949) 335–341.
  * [2] C. P. Robert, G. Casella, Monte Carlo Statistical Methods, 2nd Edition, Springer Texts in Statistics, Springer New York, New York, NY, 2010.  doi:10.1007/978-1-4757-4145-2.
  * [3] H. Kahn, "Modification of the Monte Carlo Method," in Proceedings, seminar on scientific computation, November, 1949, IBM, New York, NY, 1950.
  * [4] L. Devroye, Non-Uniform Random Variate Generation, Springer, New York Heidelberg, 1986.
  * [5] N. Metropolis, A. W. Rosenbluth, M. N. Rosenbluth, A. H. Teller, E. Teller, Equation of State Calculations by Fast Computing Machines, The Journal of Chemical Physics 21 (6) (1953) 1087–1092.  doi:10.1063/1.1699114.
  * [6] W. K. Hastings, Monte Carlo sampling methods using Markov chains and their applications, Biometrika 57 (1) (1970) 97–109.  doi:10.1093/biomet/57.1.97.
  * [7] L. Fosdick, Monte Carlo computations on the Ising lattice, Methods in computational physics 1 (1963) 245–280.
  * [8] Z. I. Botev, P. L’Ecuyer, B. Tuffin, Markov chain importance sampling with applications to rare event probability estimation, Statistics and Computing 23 (2) (2013) 271–285.  doi:10.1007/s11222-011-9308-2.
  * [9] I. Raices Cruz, J. Lindström, M. C. M. Troffaes, U. Sahlin, Iterative importance sampling with Markov chain Monte Carlo sampling in robust Bayesian analysis, Computational Statistics & Data Analysis 176 (2022) 107558.  doi:10.1016/j.csda.2022.107558.
  * [10] I. Schuster, I. Klebanov, Markov Chain Importance Sampling—A Highly Efficient Estimator for MCMC, Journal of Computational and Graphical Statistics 30 (2) (2021) 260–268.  doi:10.1080/10618600.2020.1826953.
  * [11] S. N. MacEachern, L. M. Berliner, Subsampling the Gibbs Sampler, The American Statistician 48 (3) (1994) 188.  doi:10.2307/2684714.
  * [12] W. A. Link, M. J. Eaton, On thinning of chains in MCMC, Methods in Ecology and Evolution 3 (1) (2012) 112–115.  doi:10.1111/j.2041-210X.2011.00131.x.
  * [13] A. B. Owen, Statistically Efficient Thinning of a Markov Chain Sampler, Journal of Computational and Graphical Statistics 26 (3) (2017) 738–744.  doi:10.1080/10618600.2017.1336446.
  * [14] S. K. Sahu, A. A. Zhigljavsky, Self-regenerative Markov chain Monte Carlo with adaptation, Bernoulli 9 (3) (2003) 395–422.  doi:10.3150/bj/1065444811.
  * [15] J. Gåsemyr, Markov chain monte carlo algorithms with independent proposal distribution and their relationship to importance sampling and rejection sampling, Preprint series. Statistical Research Report (2002) 1–25.
  * [16] W. H. Wong, F. Liang, Dynamic weighting in Monte Carlo and optimization, Proceedings of the National Academy of Sciences 94 (26) (1997) 14220–14224.  doi:10.1073/pnas.94.26.14220.
  * [17] J. S. Liu, F. Liang, W. H. Wong, A Theory for Dynamic Weighting in Monte Carlo Computation, Journal of the American Statistical Association 96 (454) (2001) 561–573.  doi:10.1198/016214501753168253.
  * [18] S. Malefaki, G. Iliopoulos, On convergence of properly weighted samples to the target distribution, Journal of Statistical Planning and Inference 138 (4) (2008) 1210–1225.  doi:10.1016/j.jspi.2007.05.030.
  * [19] R. Douc, C. P. Robert, A vanilla Rao–Blackwellization of Metropolis–Hastings algorithms, The Annals of Statistics 39 (1) (2011) 261 – 277.  doi:10.1214/10-AOS838.
  * [20] R. Douc, E. Moulines, P. Priouret, P. Soulier, Markov Chains, Operation Research and Financial Engineering, Springer, 2018.  doi:10.1007/978-3-319-97704-1.