Thème NUM -Systèmes numériques Équipes-Projets CQFD

Keywords: Sequential Monte Carlo methods, Feynman-Kac formulae, mean field particle models, effective sample size, adaptive functional criteria

Sequential Monte Carlo (SMC) methods are a general class of techniques to sample approximately from any sequence of probability distributions. These distributions are approximated by a cloud of weighted samples which are propagated over time using a combination of importance sampling and resampling steps. This article is concerned with the convergence analysis of a class of SMC methods where the times at which resampling occurs are computed online using criteria such as the effective sample size. This is a popular approach amongst practitioners but there are very few convergence results available for these methods. It is shown here that these SMC algorithms correspond to a particle approximation of a Feynman-Kac flow of measures on adaptive excursion spaces. By combining a non-linear distribution flow analysis to an original coupling technique, we obtain functional central limit theorems and uniform exponential concentration estimates for these algorithms. The original exponential concentration theorems presented in this study significantly improve previous concentration estimates obtained for SMC algorithms.

Méthodes de Monte Carlo séquentielles et adaptatives

Résumé : Cet article porte sur la modélisation et l'analyse de la convergence de méthodes de Monte Carlo séquentielles et adaptatives. Ces techniques de simulation stochastiques sont fondées sur l'évolution aléatoire de systèmes de particules en interaction et associées à des critères d'adaptation. Nous présentons des interprétations de ces algorithmes en termes d'approximation de type champ moyen de flots de mesures de Feynman-Kac dans des espaces d'excursions. Les instants d'interaction sont calculés en lignes selon des critères d'adaptation empiriques tels des critères de tailles de ré-échantillonnage effectives et des critères d'entropie relative de type Boltzmann. Ces stratégies d'adaptation empiriques sont d'un usage courant en pratique souvent présentées comme de simples heuristiques, sans aucun fondement mathématique rigoureux. Nous présentons de nouvelles techniques d'analyse de type semi-groupe fondées sur des couplage de processus à valeurs mesures. Cette analyse permet d'obtenir des théorèmes de la limite centrale fonctionnels et des résultats de concentration exponentielle uniformes par rapport au paramètre temporel. Notons enfin que ces travaux sont les premières études rigoureuses sur la concentration et les fluctuations de ces classes d'algorithmes de simulation adaptatifs. Par ailleurs, les théorèmes de concentration exponentielle présentés dans cette étude améliorent de récents travaux sur ce thèmes dans le cadre des interprétations particulaires de flots de mesures de Feynman-Kac.

Mots-clés :

Méthodes de Monte Carlo séquentielles et adaptatives, formules de Feynman-Kac, méthodes particulaires de type champ moyen, taille d'échantillonnage effective, critères adaptatifs fonctionnels.

Introduction

Sequential Monte Carlo (SMC) methods are a generic class of simulation-based algorithms to sample approximately from any sequence of probability distributions. These methods are now extensively used in statistics, computer science, engineering, and physics; see [START_REF] Moral | Genealogical and interacting particle systems[END_REF][START_REF] Doucet | Sequential Monte Carlo methods in pratice[END_REF][START_REF] Liu | Monte Carlo methods for scientific computing[END_REF] for many applications. SMC methods approximate the target probability distributions of interest by a large number of random samples termed particles which evolve over time according to a combination of importance sampling and resampling steps.

In the resampling steps, new particles are sampled with replacement from a weighted empirical measure associated to the current particles; see Section 2.1.2 for more details. These resampling steps are crucial and, without them, it is impossible to obtain uniform (in time) convergence results for SMC estimates. However, resampling too often has a negative effect: it decreases the number of distinct particles, increases the variance of SMC estimates and thus should only be applied when necessary. Consequently, in most practical implementations of SMC methods, the times at which resampling occurs are selected by monitoring a criterion which assesses the quality of the current particle approximation. Whenever this criterion is above or below a given, user-set, threshold then resampling occurs. This approach was originally proposed in [START_REF] Liu | Blind deconvolution via sequential imputation[END_REF] and, subsequently, has been widely adopted.

For this class of adaptive SMC methods, the times at which resampling occurs are computed on-line using our current SMC approximation and thus are random. However, most of the theoretical results on SMC algorithms assume resampling occurs at deterministic times; see [START_REF] Douc | Limit theorems for weighted samples with applications to sequential Monte Carlo[END_REF] for an exception which is discussed later on. The objective of this article is to provide a thorough treatment of this problem.

Sharp convergence results are provided for this type of SMC algorithms; this is achieved through a coupling argument. Under weak assumptions, the random times at which resampling occurs converge, almost surely as the number of particles goes to infinity, toward some deterministic resampling times. We consider a reference SMC algorithm which resamples at these deterministic times; the times are not explicitly known so this algorithm cannot be implemented but all the convergence results developed in the literature hold for it. It is then shown that the difference, in probability, between this reference SMC algorithm and the adaptive algorithm is exponentially small in the number of particles. This approach allows us to transfer rather straightforwardly the convergence results of the reference SMC algorithm to the adaptive SMC algorithm. In particular, we establish multivariate and functional central limit theorems, yielding what seems to be the first results of this type for this class of algorithms. Moreover, for standard SMC methods and the adaptive SMC schemes presented here, we present new exponential concentration estimates, including time uniform estimates which significantly improve over those presented in [START_REF] Moral | Genealogical and interacting particle systems[END_REF]. These results allow us to calibrate the number of samples so as to achieve a given precision with a given probability, uniformly with respect to the time parameter.

Structure of the Article

The rest of the paper is organized as follows. Note that the proofs of some technical results are given in the appendix.

In Section 2, we detail some generic models (Section 2.1.2) and present the class of adaptive SMC algorithms studied here (Section 2.2). Various resampling criteria that are used by practitioners are discussed in Section 2.3. It is shown that they correspond to empirical approximations of limiting criteria defined in Section 2.4.1. This allows us to present our main coupling result in Section 2.4: the difference in probability between the reference SMC algorithm using these limiting criteria and the adaptive SMC algorithm using approximations of these limiting criteria is exponentially small in the number of particles.

In Section 3, a precise description of the sequence of distributions approximated by the reference SMC algorithm is given. We review Feynman-Kac distribution flows in Sections 3.1 and 3.2 then describe such flows on excursion spaces (Section 3.3). In Section 3.4 a general class of functional criteria which are used to define the deterministic resampling times are presented; the construction of these times is detailed in Section 3.5. We further detail the construction of these times for the criteria previously introduced in Section 2.4.1.

In Section 4, a theoretical analysis of the reference SMC algorithm is presented; this is defined formally in Section 4.1. A concentration analysis is given in Section 4.2 which significantly improves the exponential estimates for SMC approximations presented in earlier studies; e.g. those in Section 7.4.3 of [START_REF] Moral | Genealogical and interacting particle systems[END_REF]. These results are used to obtain a concentration result for the empirical criteria around their limiting values.

The results above are used, in Section 5, to bound the differences between the deterministic resampling times and their empirical approximations, up to an event with an exponentially small probability. This is achieved by introducing a randomized criterion. Finally we analyze the fluctuations of adaptive SMC algorithms in Section 6.

Adaptive SMC Algorithms and Main Results

Target distributions on path spaces

Some Notation

Standard notation used in this article is as follows. M(E), P(E), and B b (E), respectively denote, the set of bounded and signed measures, the subset of all probability measures on some measurable space (E, E), and the Banach space of all bounded and measurable functions f on E when equipped with the uni-

form norm f = sup x∈E |f (x)|. Osc 1 (E) is the convex set of E-measurable functions f with oscillations osc(f ) = sup (x,y)∈E 2 {|f (x) -f (y)|} ≤ 1. µ(f ) = µ(dx) f (x), is the integral of a function f ∈ B b (E), w.r.t a measure µ ∈ M(E).
In a slight abuse the notation, let µ(A) = µ(1 A ) with A ∈ E and 1 A the indicator of A. Recall that a bounded integral operator M from a measurable space (E, E) into an auxiliary measurable space (F, F), is an operator

f → M (f ) from B b (F ) into B b (E) such that the functions M (f )(x) = F M (x, dy) f (y)
are measurable and bounded, for any f ∈ B b (F ). A bounded integral operator M from a measurable space (E, E) into an auxiliary measurable space (F, F) INRIA also generates a dual operator µ → µM from M(E) into M(F ) defined by (µM )(f ) := µ(M (f )).

Motivation and Models

Most of the research in Monte Carlo simulation is centered around the computation of integrals, for a wide range of f γ(f

) := f (x) G(x) µ(dx) = E (f (X) G(X)) η(f ) := γ(f )/γ(1) (2.1)
where X is a random variable distributed according to µ on some measurable space (E, E), G is a given non-negative potential type function on E and f ∈ B b (E). If it is possible to draw independent copies (X i ) i≥1 from µ, an approximation of the pair of measures (γ, η) can be obtained as

γ N (•) := 1 N N i=1 G X i δ X i (•) and η N (•) := N i=1 G(X i ) N j=1 G(X j ) δ X i (•)
where δ a is the Dirac measure at the point a ∈ E. In practice, these approximations can suffer from several problems. For instance, assume that G = 1 A and µ(A) := P(X ∈ A) 10 -p , for some large parameter p ≥ 1. In this situation, it may not be possible to find a sample X i within the set A, so that γ N (•) = 0 and the measure η N (•) is not well defined. To overcome this difficulty, importance sampling techniques are often used. This strategy consists of finding an auxiliary importance distribution, say µ , such that µ µ with Radon-Nikodym derivative dµ dµ , and to set

γ(f ) = E (f (X ) G (X )) = f (x) G (x) µ (dx) with G (x) := G(x) dµ dµ (x)
where X ∼ µ and µ a distribution on E. This new model has exactly the same form as the one discussed above, so that we can use the same Monte Carlo strategy, replacing the pair of parameters (µ, G) by (µ , G ). For importance sampling, even if we have found a good importance distribution, another significant problem often arises. Namely, the Radon-Nikodym derivative reflects the relative differences between the desired reference measure and the importance sampling one. These weight functionals are generally degenerate (also known as weight degeneracy). For instance, assume X represents the path of a sequence of random variables (X 0 , . . . , X n ), and the potential functions G are a product of a sequence of functions from the origin up to a given time horizon n

X = (X 0 , . . . , X n ) and G(X) = 0<p<n G p (X p )
where G p is a potential function of the random state X p . As n increases, so does the relative variance of G, and accurate estimation of (2.1) is very difficult. In this situation, the pair of measures (γ, η) given below depend on the time parameter, and they are given for every test function f n on the path space by the following Feynman-Kac integral representation formulae

η n (f n ) = γ n (f n )/γ n (1) with γ n (f n ) = E f n (X 0 , . . . , X n ) 0<p<n G p (X p ) .
(2.2) From a pure mathematical point of view, we emphasize that every importance sampling/change of measure, of the type discussed above, reduces to the same type of Feynman-Kac model.

These probabilistic models and their SMC approximations appear in numerous fields including: probability, Bayesian statistics, theoretical chemistry, quantum physics and engineering sciences. For instance, in signal processing, (X n ) n≥0 is often a Markov chain and represents a signal process and the potential functions are the likelihood functions associated with some noisy observations. The corresponding measures η n are the conditional distributions of the random paths of the signal, given a series of observations. In particle physics, the law of an unabsorbed particle X n evolving in an absorbing medium E n associated with an energy well V is given by a Feynman-Kac integral formulae. In this case, there is a Boltzmann-Gibbs potential function

G n (X n ) = exp (-β V (X n )).
For more details on these models and their connections to molecular simulation and particle filtering techniques, see, for example, chapter 12 of [START_REF] Moral | Genealogical and interacting particle systems[END_REF].

We remark that in most of the application areas discussed above, the free evolution process X n and the potential functions G n are given by the physical problem at hand. In other instances, the choice of these parameters is more flexible. Let us assume that the target measure is defined in terms of a Boltzmann-Gibbs mass transformation of a reference 'static' measure λ, with respect to some potential function g n on some measurable space E; that is

π n (dx) = Ψ gn (λ)(dx) := 1 Z n g n (x) λ(dx) with Z n := λ(g n ).
For instance, in combinatorial counting problems λ often represents the uniform measure on a given finite set E and g n = 1 An the indicator function a subset A n ⊂ E. In global optimization problems, λ may represent a reference measure on a given space E and g n = exp (-β n V ) a Boltzmann-Gibbs factor associated with an inverse temperature parameter β and an energy type non-negative function V . Numerous algorithms have been proposed in the Markov chain Monte Carlo community to sample from these Boltzmann-Gibbs measures; that is, we design Markov kernels M n with a prescribed invariant measure π n = π n M n . For instance, Metropolis-Hastings transitions consist of sampling a candidate, and to accept or reject it according to an appropriate acceptance ratio which ensures that the target measure π n is reversible w.r.t these local transitions. The wellknown simulated annealing algorithm consists of finding a judicious schedule β n increasing sufficiently slowly to ensure that successive Metropolis-Hasting moves between these changes are close to equilibrium. Another strategy, consists of finding a sequence of multiplicative potential functions, so that

(g n+1 = g n × G n and π n = π n M n ) ⇒ π n+1 = Ψ G n (π n )M n+1 . ( 2
G n = exp (-(β n+1 -β n )V )
, and respectively

G n = 1 An+1 (2.4)
with a non-decreasing inverse temperature parameter β n and a non-increasing sequence of subsets A n .

Using the Markov property and the multiplicative structure of the Feynman-Kac measures (2.2), it can be proved that the measures π n given in (2.3) coincide with the n-time marginals of the path measures defined in (2.2); that is

π n (f ) ∝ E f (X n ) 0<p<n G p (X p ) and λ(g n ) := λ(g 0 )× 0≤p<n π p G p+1 .
(2.5) We underline that these integral representation formulae are valid for any choice of the parameters (G n , M n ) satisfying the pair of admissibility conditions presented in the L.H.S of (2.3). This formula is known as Crooks identity in statistical physics [START_REF] Crooks | Nonequilibrium measurements of free energy differences for microscopically reversible Markovian systems[END_REF]; see [START_REF] Del Moral | Sequential Monte Carlo samplers[END_REF] for a generalization. It is a discrete-time version of the celebrated Jarzynski identity [START_REF] Jarzynski | Equilibrium free energy differences from nonequilibrium measurements: a master equation approach[END_REF]; see [START_REF] Rousset | Equilibrium sampling from nonequilibrium dynamics[END_REF] for a continuous-time SMC implementation.

Adaptive sequential Monte Carlo

SMC methods are a popular class of methods for approximately sampling random variables according to the Feynman-Kac path measures (2.2) and computing their normalizing constants. These algorithms are sometimes referred to as population Monte Carlo, branching and interacting particle systems or diffusion/quantum Monte Carlo methods. The central idea is now given.

First, N particles evolve according to the free path evolution of the reference Markov chain X n . In a second stage, particles with low G n -potential value are killed, and those with a relative larger potential are duplicated. In this sense, we define a branching particle model with a genealogical tree occupation measure that approximates the desired measures (2.2) as N tends to ∞. The theoretical foundations of these algorithms are now rather well understood. From the pure mathematical point of view, they can be regarded as a mean field particle approximation of a non-linear Feynman-Kac flow of distributions.

In recent years, researchers have built new branching/resampling strategies. The idea is to evolve the particles up to some times that depend on criteria computed using the current particle approximation. These adaptive algorithms have several advantages over standard methods:

1. In the context of the Feynman-Kac path measures (2.2), these criteria can be used to design online strategies controlling the 'quality' of the current particle approximation. Popular criteria include the variance of the weights, or entropy like distances between the free evolution distribution and the weighted one.

2. In the context of Boltzmann-Gibbs measures (2.3), they allow us to define adaptive temperature schedules or adaptive subset sequences. For instance, we can decide to change the temperature parameter at some appropriate time periods. These time parameters and the next cooling temperature may depend on the performance of the exploration of the individuals in the energy landscape.

In terms of the latter point, some related investigations can be found in [START_REF] Jasra | Inference for Lévy driven stochastic volatility models via adaptive sequential Monte Carlo[END_REF].

Adaptive Framework

Adaptive SMC algorithms proceed as follows. Let (E n , E n ) n≥0 be a sequence of measurable spaces. The idea is to propagate forward a collection of N sample paths copies (Y

(N,i) n ∈ E n ) ∀1 ≤ i ≤ N Y (N,i) 1 := (Y (N,i) 0 , Y (N,i) 1 , . . . , Y (N,i) t N 1 )
of the path valued random variables (X 0 , . . . , X t N

1

) up to the first time t N 1 the importance weights

W 0,t N 1 (Y (N,i) 1 ) := 0<p≤t N 1 G p (Y (N,i) p )
(G 0 ≡ 1, although this need not be the case), of the N path samples become, in some respect, degenerate; see Section 2.3 for more details. Loosely speaking, at time t N 1 the weighted occupation measure of the system

η N 1 (•) := N i=1 W 0,t N 1 (Y (N,i) 1 ) N j=1 W 0,t N 1 (Y (N,j) 1 ) δ Y (N,i) 1 (•)
is an N particle approximation measure of the updated target measure given by

η 1 (f ) = γ 1 (f )/ γ 1 (1) with γ 1 (f ) = E   f X 0 , . . . , X t N 1 0<p≤t N 1 G p (X p )   .
The precise analysis of these approximations is one of the objectives of the article will be discussed later.

To prevent the weighted measure η N 1 from becoming degenerate after the selected time period t N 1 , these measures are replaced by a bootstrap type empirical measure

η N 1 (•) := 1 N N i=1 δ b Y (N,i) 1 (•) associated with N path particles Y (N,i) 1 := ( Y (N,i) 0,t N 1 , Y (N,i) 1,t N 1 , . . . , Y (N,i) t N 1 ,t N 1
) that are resampled; for instance randomly from the current weighted measure η N 1 . In a second stage we propagate forward, from every end point Y

(N,i) t N 1 ,t N 1 of the selected paths Y (N,i) 1 , a collection of N sample path copies (Y (N,i) t N 1 +1 , Y (N,i) t N 1 +2 , . . . , Y (N,i) t N 2 ) of the path-valued random variables (X t N 1 +1 , X t N 1 +2 , . . . , X t N 2 ) starting at X t N 1 = Y (N,i) t N 1 ,t N 1
, up to the first time t N 2 , the importance weights of the N path samples given by

W t N 1 ,t N 2 (Y (N,i) 2 ) := t N 1 <p≤t N 2 G p (Y (N,i) p )
INRIA again become degenerate. In the above expression, Y (N,i) 2

stands for the resulting path sequences of length t N 2 + 1 given by the formula

Y (N,i) 2 := Y (N,i) 1 , (Y (N,i) t N 1 +1 , Y (N,i) t N 1 +2 , . . . , Y (N,i) t N 2
) .

The following weighted occupation measure of the system

η N 2 (•) := N i=1 W t N 1 ,t N 2 (Y (N,i) 2 ) N j=1 W t N 1 ,t N 2 (Y (N,j) 2 ) δ Y (N,i) 2 (•)
is an approximation of the updated target measure at time t N 2 :

η 2 (f ) = γ 2 (f )/ γ 2 (1) with γ 2 (f ) = E   f X 0 , . . . , X t N 2 0<p≤t N 2 G p (X p )   .
As before, to prevent the degeneracy of the weighted measures, after the selected time period t N 2 , the weighted measure η N 2 is replaced by a bootstrap type empirical measure

η N 2 (•) := 1 N N i=1 δ b Y (N,i) 2 (•) associated with N particle paths Y (N,i) 2 := ( Y (N,i) 0,t N 2 , Y (N,i) 1,t N 2 , . . . , Y (N,i) t N 2 ,t N 2
) resampled from the current weighted measure η N 2 , and so on up to the desired final time horizon.

Remark 2.1 This adaptive excursion model is quite flexible. For instance, in the context of the Boltzmann-Gibbs model (2.3), it is possible to choose the same potential function G t N n +1 := G p between the interaction times t N n < p ≤ t N n+1 . In other words, the potentials are only changed at the random interaction times

t N n ∀t N n < p ≤ t N n+1 G p := G t N n +1 ↓ ∀t N n+1 < p ≤ t N n+2 G p := G t N n+1 +1 ↓ • • •
As an example, a temperature schedule can be used with the same decay rate between the interaction times t N n < p ≤ t N n+1 , and then changed only at the interaction times. Using this strategy, we have used only n different potential functions at time t N n+1 :

G t N 0 +1 G t N 1 +1 . . . G t N n +1
.

Without further work, the forthcoming analysis can be extended to situations where the potential functions G p p>t N n

, at a given level n, depend on the series of parameters (t N 0 , . . . , t N n ). More formally, it is possible to choose

∀t N n < p ≤ t N n+1 G p (x) := G p,(t N 0 ,...,t N n ) (x) (2.6) RR n°6700
for some collection of mappings G p,(t0,...,tn) , with

0 = t 0 ≤ • • • ≤ t n < p.
In much the same way, in the context of the Boltzmann-Gibbs model (2.3), we can use a sequence of Markov transitions M p p>t N n as well as an adaptive resampling criterion at a given level n that depends on (t N 0 , . . . , t N n ).

Some empirical criteria

A few empirical criteria have been proposed in the SMC literature to measure the degeneracy of the weighted samples; these criteria can be used to design adaptive resampling strategies. In all of these situations, the resampling times t N n n≥0 are random variables that depend on the current SMC approximation. Four criteria are now presented which are termed 1)-4), in the order that they appear. Due to their practical importance, and to better connect these strategies with existing results in the literature, we present a rather detailed discussion on the resulting adaptive resampling strategies.

Adaptive control of the weights

In this case, the particles explore the state-space with a series of local transitions

Y (N,i) t N n ,t N n -------→ Y (N,i) t N n +1 --- ---→ Y (N,i) t N n +2 --- ---→ . . . --- ---→ Y (N,i) s (2.7)
up to the first time s = t N n+1 (≥ t N n ) the proportion of particles with a weight of at least size s n is below a prescribed a n -threshold

C N t N n ,s := 1 N Card    1 ≤ i ≤ N : t N n <k≤s G k (Y (N,i) k ) ≥ s n    ≤ a n (2.8)
where 0 < a n ≤ 1. To our knowledge, this type of criteria was first discussed as an heuristic type resampling schedule in [START_REF] Moral | Measure valued processes and interacting particle systems. Application to non linear filtering problems[END_REF]. Next, we provide some comments on the choice of the adaptive parameters (s n , a n ). The forthcoming discussion applies to all the criteria discussed in the article.

For [0, 1]-valued potential functions G k , it is sensible to choose s n ∈ [0, 1]. In this situation, the random mappings s ∈ [t N n , ∞) → C N t N n ,s ∈ [0, 1]
are non-increasing functions starting at

C N t N n ,t N n = 1 (the standard convention ∅ = 1 is used).
It should be noted that when none of the sampled particles at time s = t N n + 1 have a weight greater than the prescribed threshold s n , then t N n+1 = t N n and the above sampling strategy is restarted. As discussed in Remark 2.1, our framework is general enough that it allows us to analyze, without further work, fully adaptive algorithms. For example, we may use a sequence of threshold parameters that depend on the previous times steps

s n = s n (t N 0 , . . . , t N n ) and a n = a n (t N 0 , . . . , t N n )
for some functions s and a n on R n+1 + . So, in cases where t N n+1 = t N n , it is possible to either choose to restart the complete run with the same set of potentials and

INRIA

Markov transitions used before and/or to decide to decrease the weight size parameter s n as well as the threshold level a n .

The choice of these threshold parameters depends on the problem at hand. From the previous discussion, if s n is excessively large, the algorithm may get stuck at a given time step t N n during several runs. In the further developments of the article, we provide a set of regularity conditions under which these empirical criteria approximate a limiting functional integration criteria. It will be shown, that the random times t N n converge, as N tends to ∞, to a series of deterministic times t n expressed in terms of these functional criteria. In this situation, for a large population size N we ensure that t N n+1 > t N n with high probability as soon as t n+1 > t n .

To illustrate this idea, suppose that (X k ) k≥0 is a Markov chain on E k = Z starting at the origin X 0 = 0. For instance, X k can be taken as a simple random walk on the integer lattice. If we choose the indicator potential functions and the unit weight size parameters

G k = 1 [k,∞) and s n = 1 (2.9)
then we find

t N n <k≤s G k (Y (N,i) k ) ≥ 1 ⇐⇒ ∀t N n < k ≤ s Y (N,i) k ≥ k for every 1 ≤ i ≤ N . Thus, C N t N n ,s = 1 N Card 1 ≤ i ≤ N : ∀t N n < k ≤ s Y (N,i) k ≥ k . (2.

10)

When a n = 1/2, we find that t N n+1 is such that

t N n+1 = inf s ≥ t N n : Card 1 ≤ i ≤ N : Y (N,i) t N n +1 ≥ t N n + 1, Y (N,i) t N n +2 ≥ t N n + 2 , . . . , Y (N,i) s ≥ s ≤ N/2 .
Roughly speaking, for the simple random walk model, this ensures that about half of the N particles are able to reach the next upper level set at every elementary time step. In this situation, the random times t N n+1 are more likely to take the value t N n + 1, or to return the value t N n . If we take a lower threshold parameter, say a n = 1/100, these random times t N n+1 are more likely to take large values. Indeed, in this case, about N/2 particles will be able to reach the first level, then about N/4 the second, N/8 the third one, and so on up to N/128 the seventh. Therefore, for a n = 1/100 it is expected that t N n+1 = t N n + 6. The random walk model discussed above can easily be extended to multisplitting branching excursion type models, absorbed at some lower level set B ⊂ Z ∩ R -. Here, the random path processes (Y

(N,i) k ) k≥t N
n are independent copies of a random excursion starting from a given level set t N n and ending at the absorbing set B. The next adaptive step t N n+1 coincides with the first upper level set the N particles are able to reach with a proportion less than a n (that is, before being absorbed by the set B). The resulting algorithm for adaptive multi-level splitting coincides with the adaptive branching algorithm for rare event simulation discussed in [START_REF] Cerou | Adaptive multilevel splitting for rare event analysis[END_REF].

Normalizing constants

In this scenario, the particles explore the space with a series of local transitions (2.7) up to the first time s = t N n+1 the average of the importance weights -which is an estimate of a normalizing constant -goes below a prescribed threshold a n

C N t N n ,s := 1 N N i=1 t N n <k≤s G k (Y (N,i) k ) ≤ a n . (2.11)
For indicator potential functions, this criterion coincides with the one discussed above; see (2.10). Furthermore, for (0

, 1]-valued potential functions G k , it is convenient in (2.11) to set a n ∈ [0, 1]; the random mappings s ∈ [t N n , ∞) → C N t N n ,s ∈ [0, 1] are non-increasing functions starting at C N t N n ,t N n = 1. Notice that a n = 1 =⇒ t N n+1 = t N n .
For (0, 1)-valued threshold parameters, we have the rather crude estimates

g k := inf G k ≤ sup G k := g k ⇓ t N n <k≤s g k := C N t N n ,s ≤ C N t N n ,s ≤ C N t N n ,s := t N n <k≤s g k . (2.12)
This pair of estimates can be improved by considering the infimum and the supremum of the potential product weight functions along the M k -admissible paths sequences (2.7) starting at the selected states Y

(N,i) t N n ,t N n . Note that (2.12) implies t N n+1 ≤ t N n+1 ≤ t N n+1 with t N n+1 := inf s ≥ t N n : C N t N n ,s ≤ a n and t N n+1 := inf s ≥ t N n : C N t N n ,s ≤ a n .
For instance, for the exponential Boltzmann-Gibbs potential given in (2.4), we have

V := inf V ≤ sup V := V =⇒ g k = e -(β k+1 -β k )V and g k = e -(β k+1 -β k )V .
This clearly implies that

C N t N n ,s = e - " βs+1-β t N n +1 " V and C N t N n ,s = e - " βs+1-β t N n +1 " V .
In addition, further assume that the temperature schedule (β k ) k>t N n after time t N n is associated with a fixed decay rate δ n > 0, then

∀t N n < k ∆β k+1 = (β k+1 -β k ) = δ n =⇒ β s+1 -β t N n +1 = s -t N n δ n =⇒      C N t N n ,s = e -(s-t N n ) δn V C N t N n ,s = e -(s-t N n ) δn V .

INRIA

From the previous discussion, it can be concluded that

t N n + 1 δ n V log (1/a n ) := t N n+1 ≤ t N n+1 ≤ t N n+1 = t N n + 1 δ n V log (1/a n ).
For small values of the inverse temperature decay rates δ n , it is ensured that the criterion (2.11) is met in at least

≤ (m n V ) -1 log (1/a n ) is needed; that is δ n ≤ (m n V ) -1 log (1/a n ) =⇒ t N n+1 -t N n ≥ m n .

Variance of the weights

In this situation, the particles explore the state-space with a series of local transitions (2.7) up to the first time s = t N n+1 the empirical particle variance is larger than some prescribed threshold

a n N i=1   t N n <k≤s G k (Y (N,i) k ) N j=1 t N n <k≤s G k (Y (N,j) k )   2 ≥ (1 + a n )/N
or equivalently

C N t N n ,s := 1 N N i=1   t N n <k≤s G k (Y (N,i) k ) 1 N N j=1 t N n <k≤s G k (Y (N,j) k )   2 -1 ≥ a n .
This is equivalent to resampling when the so-called effective sample size is below a prescribed threshold. This popular strategy was first proposed in [START_REF] Liu | Blind deconvolution via sequential imputation[END_REF].

For indicator potential functions, this criterion coincides with those above; as

(G k ) 2 = G k .

Boltzmann's entropy

Assume that the potential functions G k take values in (0, ∞). In this case, the particles explore the state-space with a series of local transitions (2.7) up to the first time s = t N n+1 such that the relative entropy of the empirical particle measures w.r.t its weighted versions are larger than some threshold a n

C N t N n ,s := - 1 N N i=1 t N n <k≤s log G k (Y (N,i) k ) ≥ a n . (2.13)
For the exponential Boltzmann-Gibbs potential given in (2.4), we have

C N t N n ,s := 1 N N i=1 t N n <k≤s (β k+1 -β k ) V (Y (N,i) k
).

In this situation, the random mappings

s ∈ [t N n , ∞) → C N t N n ,s ∈ [0, 1] RR n°6700
are non-decreasing functions starting at C N t N n ,t N n = 0 (using the standard convention ∅ = 0). Arguing as above, we have the rather crude estimates

V := inf V ≤ sup V := V ⇓ C N t N n ,s = β s+1 -β t N n +1 V ≤ C N t N n ,s ≤ C N t N n ,s := β s+1 -β t N n +1 V . This implies that t N n+1 ≤ t N n+1 ≤ t N n+1 with t N n+1 := inf s ≥ t N n : C N t N n ,s ≥ a n and t N n+1 := inf s ≥ t N n : C N t N n ,s ≥ a n .
When the temperature schedule (β k ) k>t N n is associated with a prescribed and fixed decay rate δ n > 0, it follows that

C N t N n ,s = s -t N n δ n V and C N t N n ,s = s -t N n δ n V
and therefore

t N n+1 := t N n + a n /(δ n V ) ≤ t N n+1 ≤ t N n+1 := t N n + a n /(δ n V ).
In this case, to ensure that the criterion (2.13) is met in at least m n iterations, we need to choose an inverse temperature decay rate δ n ≤ a n /(m n V ); that is

δ n ≤ a n /(m n V ) =⇒ t N n+1 -t N n ≥ m n .

Statement of some results

The following Section provides a guide of the major definitions and results in this article; these will be repeated at the relevant points in the paper.

A limiting reference SMC algorithm

Let (t n ) n≥0 be the deterministic sequence of time steps obtained by replacing the empirical criteria C N t N n ,s by their limiting values C tn,s , i.e. as N ↑ ∞; e.g. in the examples 1) and 2)

t n+1 := inf {t n ≤ s : C tn,s ≤ a n }.
In all situations, set t N 0 = t 0 = 0. In the examples 3) and 4), the limiting times (t n ) n≥0 are defined as above by replacing condition C tn,s ≤ a n by the condition C tn,s ≥ a n . In the four situations discussed above, these limiting criteria C tn,s are respectively given by

1) P n,b ηn   tn<k≤s G k (X k ) ≥ a n   2) E n,b ηn   tn<k≤s G k (X k )   INRIA and
3)

E n,b ηn tn<k≤s G k (X k ) 2 E n,b ηn tn<k≤s G k (X k ) 2 -1 4) - tn<k≤s E n,b ηn (log G k (X k ))
with the updated Feynman-Kac path measures η n given below

η n (f ) = γ n (f )/ γ n (1) with γ n (f ) = E   f X 0 , . . . , X tn 0<p≤tn G p (X p )   .
In the above formulae, E n,b ηn is the expectation w.r.t the law of a free evolution X k starting at time t n at the end point X tn = X tn,tn of a random path X n := ( X 0,tn , X 1,tn , . . . , X tn,tn ) with distribution η n ; that is

X tn,tn ∼ η n M tn +1 -------→ X tn+1 M tn+2 --- ---→ X tn+2 M tn+3 --- ---→ . . . M s --- ---→ X s (2.14) Definition 2.2 Let Y (N ) n := Y (N,1) n , Y (N,2) n , . . . , Y (N,N ) n and Y (N ) n := Y (N,1) n , Y (N,2) n , . . . , Y (N,N )
n denote the N particles associated to the adaptive SMC algorithm resampling at times t N n n≥0 and let

X (N ) n := X (N,1) n , X (N,2) n , . . . , X (N,N ) n and X (N ) n := X (N,1) n , X (N,2) n , . . . , X (N,N )
n denote the N particles associated to the reference SMC algorithm resampling at times (t n ) n≥0 . We also suppose that

X (N ) n , X (N ) n and Y (N ) n , Y (N ) n
coincide on every time interval 0 ≤ n ≤ m, as soon as the random times t N n = t n , for every 0 ≤ n ≤ m. This condition simply corresponds to the coupling of the two processes on the events ∩ 0≤n≤m {t N n = t n }.

An exponential coupling theorem

We give our main results, which hold under the following regularity condition:

(G ) ∀n ≥ 0 q n := sup (x,y)∈(E n ) 2 (G n (x)/G n (y)) < ∞. (2.15)
The analysis of more general models is technically more involved and often requires some additional assumptions on the free evolution process. We refer the interested reader to [START_REF] Moral | Genealogical and interacting particle systems[END_REF] for a more thorough discussion on these technicalities. The first result is a non-asymptotic exponential concentration estimate. Below, if constants are written with an argument, then they depend only on this given argument.

Theorem 2.3 For any n ≥ 0, f n ∈ Osc 1 (E 0 × • • • × E tn ), any N ≥ 1 and any 0 ≤ ≤ 1/2, there exist c 1 < ∞, 0 < c 2 (n) < ∞ such that we have the exponential concentration estimate P η N n -η n (f n ) ≥ ≤ c 1 exp{-N 2 /c 2 (n)}
for the empirical measures

η N n (•) = 1 N N i=1 δ " b X (N,i) n-1 , " X (N,i) t n-1 +1 ,X (N,i) t n-1 +2 ,...,X (N,i) tn
"" (•).

In addition, under appropriate regularity conditions on the potential transitions G k and on the Markov chain X k given in Section (4.2.1) the above estimates are valid for the marginal measures associated with the time parameters t n-1 ≤ p ≤ t n , for some constant c 2 (n) = c 2 whose values does not depend on the time parameter.

The second result is the following exponential coupling theorem.

Theorem 2.4 Assume the threshold parameters (a n ) n≥0 are sampled realizations of a collection of absolutely continuous random variables A = (A n ) n≥0 .

Then, for almost every realization of the sequence (a n ) n≥0 , (X

(N ) n , X (N ) n ) n≥0 and (Y (N ) n , Y (N )
n ) n≥0 are such that, for every m ≥ 0 and any N ≥ 1, there exist 0 < c 1 (m), c 2 (m) < ∞ and almost surely (m, A) ≡ (m) > 0 such that

P ∃0 ≤ n ≤ m (Y (N ) n , Y (N ) n ) = (X (N ) n , X (N ) n ) | A ≤ c 1 (m) e -N 2 (m)/c2(m) .
The proofs of Theorem 2.3 and Theorem 2.4 are detailed, respectively, in Sections 4.2.3 and 5.2. Some direct consequences of the exponential coupling Theorem 2.4 are as follows. The reference SMC algorithm (X (N ) n , X n ) n≥0 can be regarded as a 'standard' SMC algorithm as it resamples at deterministic times. A variety of convergence results have already been developed for these algorithms (see for instance [START_REF] Moral | Genealogical and interacting particle systems[END_REF]), including non asymptotic exponential estimates and L p -mean error bounds. Up to an event having an exponentially small occurrence probability, Theorem 2.4 allows us to transfer all of these estimates to the adaptive SMC algorithm that approximates these deterministic times.

Description of the models

For the convenience of the reader, in the first two sections 3.1, 3.2, we briefly introduce an abstract Feynman-Kac flow of distributions. The Feynman-Kac flow of distributions in excursion space associated to a prescribed sequence of resampling times is presented in Section 3.3.

The definition of the functional criteria dicussed in this work is provided in Section 3.4, and the precise definition of the associated 'limiting' deterministic time parameters (t n ) n≥0 discussed in the introduction is provided in Section 3.5.

INRIA

Feynman-Kac distributions flow

We consider a sequence of measurable state spaces S n , a probability measure η 0 on S 0 , and a sequence of Markov transitions M n (x n-1 , dx n ) from S n-1 into S n , with n ≥ 1. Let X n n≥0 be a Markov chain with initial distribution Law(X 0 ) = η 0 on S 0 and elementary transitions

P(X n ∈ dy | X n-1 = x) = M n (x, dy).
Let G n n≥0 be a sequence of non-negative and bounded potential-likelihood functions on S n . To simplify the presentation, and to avoid unnecessary technicalities, it is supposed that G 0 = 1 and for n ≥ 1, G n ∈ (0, 1) with (G)

q n := sup

(x,y)∈S 2 n (G n (x)/G n (y)) < ∞. (3.1)
This regularity condition can be relaxed in various ways; see [START_REF] Moral | Genealogical and interacting particle systems[END_REF]. The Boltzmann-Gibbs transformation Ψ n associated with a given potential function G n is the mapping

Ψ n : µ ∈ P(S n ) → Ψ n (µ) ∈ P(S n ) with Ψ n (µ)(dx) := 1 µ(G n ) G n (x) µ(dx).
Notice that Ψ n (µ) can be rewritten as a non-linear Markov transport equation Let (η n , η n ) n≥0 be the flow of probability measures, both starting at η 0 = η 0 , and defined for any n ≥ 1 by the following recursion

∀n ≥ 0 η n+1 = η n M n+1 with η n := Ψ n (η n ). (3.2)
Notice that (η n ) n≥0 and ( η n ) n≥0 satisfy the nonlinear recursion

η n+1 = Φ n+1 (η n ) := Ψ n (η n )M n+1 and η n+1 = Φ n+1 ( η n ) := Ψ n+1 ( η n M n+1 ) .
It is worth mentioning that the one step mapping Φ n+1 can be rewritten as follows

Φ n+1 ( η n ) := Ψ n ( η n ) M n+1
with the Boltzmann-Gibbs transformation Ψ n associated with a given potential function G n = M n+1 (G n ) and the Markov transition M n+1 given for any bounded measurable test function f on S n+1 , and any x ∈ S n by the formula

M n+1 (f )(x) = M n+1 (G n f )(x) M n+1 (G n )(x) .
From the pure mathematical point of view, this shows that the flow ( η n ) n≥0 has the same dynamical structure as the flow (η n ) n≥0 . Thus, w.l.o.g we shall concentrate our analysis on the flow of measures (η n ) n≥0 .

To get one step further, it can be readily checked that the flow of measures (η n ) n≥0 satisfies the non-linear transport equation

η n+1 = η n K n+1,ηn with K n+1,ηn := S n,ηn M n+1 .
Also recall that the solution (η n , η n ) of these recursive updating-prediction equations have the following functional representations

η n (f n ) = γ n (f n )/γ n (1) and η n (f n ) = γ n (f n )/ γ n (1) (3.3)
with the unnormalized Feynman-Kac measures γ n and γ n defined by the formulae

γ n (f n ) = E f n (X n ) 0<k<n G k (X k ) and γ n (f n ) = γ n (f n G n ). (3.4)

Feynman-Kac semi-groups

To analyze SMC methods, it is necessary to introduce the Feynman-Kac semigroup associated to the flow of measures (γ n ) n≥0 and (η n ) n≥0 . Let us start by denoting by Q n+1 (x n , dx n+1 ) the bounded integral operator from S n into S n+1 ; this is defined by the elementary multiplicative formula

Q n+1 (x n , dx n+1 ) = G n (x n ) M n+1 (x n , dx n+1 ).
Let (Q p,n ) 0≤p≤n be the corresponding linear semi-group defined by

Q p,n = Q p+1 Q p+2 . . . Q n
with the convention Q n,n = I, the identity operator. Note that Q p,n is alternatively defined by:

Q p,n (f n )(x p ) = E   f n (X n ) p≤k<n G k (X k ) X p = x p   . (3.5)
Using the Markov property, it follows that

γ n (f n ) = E   E   f n (X n ) p≤k<n G k (X k ) X p   0<k<p G k (X k )   = γ p (Q p,n (f n )) .
The last assertion shows that (Q p,n ) 0≤p≤n is the semi-group associated with the unnormalized measures (γ n ) n≥0 . Denote its normalized version by

P p,n (f n ) := Q p,n (f n ) Q p,n (1) . (3.6)
Finally, denote by (Φ p,n ) 0≤p≤n the non-linear semi-group associated to the flow of normalized measures (η n ) n≥0 :

Φ p,n = Φ n • . . . • Φ p+2 • Φ p+1
with the convention Φ n,n = I, the identity operator. From the previous discussion, note that (Φ p,n ) 0≤p≤n is alternatively defined in terms of (Q p,n ) 0≤p≤n with the following formulae:

Φ p,n (η p )(f n ) = γ p Q p,n (f n ) γ p Q p,n (1) = η p Q p,n (f n ) η p Q p,n (1)
.

(3.7)

INRIA

Path space and excursion models

Let X n n≥0 be a Markov chain taking values in some measurable state-spaces E n with elementary transitions

P X n ∈ dx n | X n-1 = x n-1 := M n (x n-1 , dx n )
and the initial distribution η 0 = Law(X 0 ). In addition, introduce a sequence of potential functions G n on the state spaces E n . To simplify the presentation, it is assumed that the G n take values in (0, 1). Associate, with a prescribed non-decreasing sequence of time parameters (t n ) n≥0 , the excursion-valued random variables X 0 = X 0 and ∀n ≥ 1

X n := (X tn-1+1 , X tn-1+2 , . . . , X tn ) and the random path sequences

X n := (X 0 , . . . , X n ) = (X 0 , X 1 , . . . , X tn ) ∈ E tn := E 0 × . . . × E tn .
Note that X n n≥0 forms a Markov chain

X n+1 := (X n , X n+1 ) = X n , X tn+1 , . . . , X tn+1 (3.8)
taking values in the excursion spaces S n := E tn . Furthermore, choosing, in (3.4), the potential functions

∀n ≥ 1 G n (X n ) = G n (X n ) := tn-1<k≤tn G k (X k ) (3.9)
then we readily find that

γ n (f n ) = E f n (X n ) 0<k<n G k (X k ) = E f n (X 0 , . . . , X n ) 0<k<n G k (X k ) = E   f n (X 0 , . . . , X tn ) 0<k≤tn-1 G k (X k )   .
By definition of the potential functions G n of the excursion Feynman-Kac model (3.9), it is easily proved that the condition (G) (eq. (3.1)) is satisfied, as soon as (G ) introduced in (2.15) holds true. More precisely, it holds that

(G ) ⇓ (G) with q n ≤ M -sup    tn-1<k≤tn (G k (x k )/G k (y k ))      ≤ tn-1<k≤tn q k   RR n°6700
where the essential supremum M -sup {.} is taken over all M k -admissible paths (x tn-1+1 , x tn-1+2 , . . . , x tn ) and (y tn-1+1 , y tn-1+2 , . . . , y tn ) of the underlying Markov chain X k . Moreover, the marginal γ n of γ n on the state space S n := (E tn+1 × . . . × E tn+1 ) is again given, for any f n ∈ B b (S n ), by

γ n (f n ) := E f n (X n ) 0<k<n G k (X k ) .

Functional criteria

In Section 3.3, we have assumed that a non-decreasing sequence of time parameters (t n ) n≥0 was available. We now introduce the functional critera which will be used to build this sequence. To connect the empirical criteria with their limiting functional versions, the latter need to satisfy some weak regularity conditions which are given below. Definition 3.1 We consider a sequence of functional criteria

∀n ≥ 0 ∀p ≤ q H (n) p,q : µ ∈ P (E q ) → H (n) p,q (µ) ∈ R +
satisfying the following Lipschitz type regularity condition

H (n) p,q (µ 1 ) -H (n) p,q (µ 2 ) ≤ | [µ 1 -µ 2 ] (h) | H (n) p,q (dh) (3.10)
for some collection of bounded measures H (n) p,q on B b (E q ) such that δ(H (n) p,q ) := osc(h) H (n) p,q (dh) < ∞.

To illustrate this rather abstract construction, we provide four examples of functional criteria corresponding those discussed in Section 2.4.1. All of them are defined in terms of the importance weight function W p,q : (x 0 , . . . , x q ) ∈ E q → W p,q (x 0 , . . . , x q ) :=

p<k≤q G k (x k ) ∈ (0, 1)
In all the situations described below, the measures H (n) tn,s in (3.10) have a finite support with one or two elements.

Weights thresholds control:

The functional defined below

H (n) p,q (µ) = µ A (n) p,q
with A (n) p,q := {x q ∈ E q : W p,q (x q ) ≥ s n } measures the proportion of the path sequences with a weight product larger than some parameter s n .

Normalizing constants:

In some sense, the functional defined below

H (n) p,q (µ) = µ(W p,q )
measures the degeneracy of the normalizing constants for large values of the parameter (q -p).

INRIA

3.

Variance of the weights: The functional defined below

H (n) p,q (µ) = µ W p,q µ(W p,q ) -1 2
coincides with the variance of the weights w.r.t the reference measure µ.

Boltzmann's entropy:

The functional defined below

H (n) p,q (µ) = Ent (dµ | W p,q dµ) := -µ (log W p,q )
measures the relative entropy distance between the reference measure µ and the updated weighted measure.

Resampling time construction

We now explain, precisely, how to define the sequence of resampling times (t n ) n≥0 .

Definition 3.2 For any pair of integers 0 ≤ p ≤ n and any η ∈ P(E p ), we denote by P η,(p,n) the law of a random path sequence of length n, starting with an initial distribution η on the path space E p and performing (n -p) elementary moves with the local probability transitions M k , with p < k ≤ n. More precisely, P η,(p,n) is the measure defined by

P η,(p,n) (d(x p , x p+1 , . . . , x n )) := η(dx p ) M p+1 (x p , dx p+1 ) . . . M n (x n-1 , dx n ) ∈ P E p × E p+1 × . . . × E n = P(E n )
where dx p = d(x 0 , . . . , x p ) denotes an infinitesimal neighborhood of a path sequence

x p = (x 0 , . . . , x p ) ∈ E p = E 0 × . . . × E p .
Given H

p,q , with n ≥ 0 and 0 ≤ p ≤ q, we define a non-decreasing sequence of deterministic time steps (t n ) n≥0 and a flow of Feynman-Kac measures (η n , η n ) as follows:

• Set n = 0 and t n = 0 and η 0 = η 0 = η 0 .

• Suppose that the time step t n is defined as well as the pair of measures (η n , η n ) ∈ P(E tn ) 2 at a given rank n. To define the time step t n+1 and the pair of measures (η n+1 , η n+1 ) ∈ P(E tn+1 ) 2 , we consider

H (n)
tn,s , with s ≥ t n and a given interval and

I n ⊂ R + of the form I n = [a n , ∞) or I n = (-∞,
η n+1 = Ψ n+1 (η n+1 ) (3.11)
with the Boltzmann-Gibbs transformation Ψ n+1 associated with the potential function G n+1 = W tn,tn+1 .

We end this construction with an important observation. By definition of the Markov transition M n+1 of the excursion model X n defined in Section 3.3, it can be checked that

η n+1 = P b ηn,(tn,tn+1) = η n M n+1 .
(3.12)

This yields the recursion (3.11) and (3.12) =⇒

η n+1 = Ψ n (η n ) M n+1 .
Hence the flow of measures η n and η n coincide with the Feynman-Kac flow of distributions defined in (3.3) with the Markov chain and potential function (X n , G n ) on excursion spaces defined in (3.8) and (3.9). The only difference between these models and the ones discussed in Section 3.3 comes from the fact that the sequence of time parameters (t n ) n≥0 is inductively defined in terms of a prescribed sequence of functional criteria. The SMC approximation of these distributions is discussed in Section 4.

Some applications

In this Section, we examine, in detail, the inductive construction of the deterministic resampling times t n introduced in Section 3.5 for the criteria 1)-4) in Section 2.3.

Adaptive control of the weights

Consider [0, 1)-valued potential functions. The resampling times are built to control the degeneracy of the Feynman-Kac weights through

H (n) tn,s P b ηn,(tn,s) = E n,b ηn 1 A (n) tn ,s X tn+1 , . . . , X s = P n,b ηn   tn<k≤s G k (X k ) ≥ s n  
where 0 < s n < 1 and

A (n) tn,s :=    (x tn+1 , . . . , x s ) ∈ E tn+1 × • • • × E s : tn<k≤s G k (x k ) ≥ s n    .
Here E n,b ηn denotes the expectation w.r.t the law P n,b ηn of a free evolution X k starting at time t n at the end point X tn of a random path (X 0 , . . . , X n ) with distribution η n on the space E n .

INRIA

The mappings

t n ≤ s → P n,b ηn   tn<k≤s G k (X k ) ≥ s n  
are clearly non-increasing functions starting at 1 for s = t n . If we choose an interval I n := (-∞, a n ] for some a n < 1 then the resampling time t n+1 coincides with the first time s the sequence of probabilities

s → P n,b ηn   tn<k≤s G k (X k ) ≥ s n  
goes below the level a n ; that is

t n+1 := inf    t n ≤ s : P n,b ηn   tn<k≤s G k (X k ) ≥ s n   ≤ a n    .
This strategy allows us to decompose the trajectories of the reference free motion X k into a sequence of excursion-processes

X 0 = X 0 and ∀n ≥ 1 X n := (X tn-1+1 , X tn-1+2 , . . . , X tn ) (3.13) with ∀s < t n+1 P n,b ηn   tn<k≤s G k (X k ) ≥ s n   > a n and s = t n+1 =⇒ P n,b ηn   tn<k≤tn+1 G k (X k ) ≥ s n   ≤ a n .
We illustrate this construction with the random walk model discussed in (2.9) (Section 2.3.1). In this situation, recall that X k is a Markov chain on E k = Z starting at the origin X 0 = 0. If the indicator potential functions,

G k = 1 [k,∞) , are chosen then E 0,b η0   0<k≤t1 G k (X k )   = P(X 1 ≥ 1 . . . , X k ≥ k).
When a 0 = 1/2, we find that t 1 coincides with the first level k the chain X k can reach with a probability less than 1/2 with X 0 = 0; more formally, we have

t 1 = inf {k ≥ 0 : P(X 0 ≥ 0, X 1 ≥ 1 . . . , X k ≥ k) ≤ 1/2}.
Similarly, when a n = 1/2, t n+1 coincides with the first level k the chain X k can reach with a probability less than 1/2 with X tn = x tn ; more formally, we have

t n+1 = inf k ≥ t n : P n,b ηn X tn+1 ≥ t n + 1, . . . , X k ≥ k ≤ 1/2 . RR n°6700

Normalizing constants

The resampling times are built to control the degeneracy of the normalizing constants estimates. To be more precise, set

H (n) tn,s P b ηn,(tn,s) = E n,b ηn   tn<k≤s G k (X k )   := η n (dx tn ) E tn,x tn   tn<k≤s G k (X k )   = η n (dx tn ) M tn+1 (x tn , dx tn+1 ) . . . . . . × M s (x s-1 , dx s ) tn<k≤s G k (x k ).
We have s = t n =⇒ H is clearly non-increasing starting at 1 at time s = t n . If we choose an interval

I n := (-∞, a n ]
for some a n < 1 then the resampling time t n+1 coincides with the first time s the quantities

s → E n,b ηn   tn<k≤s G k (X k )  
go below the level a n ; that is, we have

t n+1 := inf    t n ≤ s : E n,b ηn   tn<k≤s G k (X k )   ≤ a n    .
This criterion allows us to decompose the trajectories of the reference free motion X k into a sequence of excursion like processes

X 0 = X 0 and ∀n ≥ 1 X n := (X tn-1+1 , X tn-1+2 , . . . , X tn ) (3.14) with ∀s < t n+1 E n,b ηn   tn<k≤s G k (X k )   > a n and E n,b ηn   tn<k≤tn+1 G k (X k )   ≤ a n . INRIA 3.6.

Variance of the weights

The resampling times are built to control the variance of the weights. To be more precise, set

H (n) tn,s P b ηn,(tn,s) = E n,b ηn      tn<k≤s G k (X k ) E n,b ηn tn<k≤s G k (X k ) -1   2    = E n,b ηn tn<k≤s G k (X k ) 2 E n,b ηn tn<k≤s G k (X k ) 2 -1.
In this case, s = t n =⇒ H

G k = 1 A k we have ∀t n ≤ s ≤ r H (n) tn,s P b ηn,(tn,s) = E n,b ηn   tn<k≤s G k (X k )   -1 -1 ≤ H (n)
tn,r P b ηn,(tn,r) .

For another example in filtering problems, we refer the reader to [START_REF] Moral | Interacting particle filtering with discrete-time observations: asymptotic behaviour in the Gaussian case[END_REF]. One natural way to control these variances is to choose an interval

I n := [a n , ∞)
for some a n > 0 then the resampling time t n+1 coincides with the first time s the function

s → E n,b ηn   tn<k≤s G k (X k ) 2   / E n,b ηn   tn<k≤s G k (X k )   2
goes above the level 1 + a n ; that is

t n+1 := inf    t n ≤ s : E n,b ηn   tn<k≤s G k (X k ) 2   ≥ [1 + a n ] E n,b ηn   tn<k≤s G k (X k )   2      .
This construction allows us to decompose the trajectories of the reference free motion X k into a sequence of excursion like processes

X 0 = X 0 and ∀n ≥ 1 X n := (X tn-1+1 , X tn-1+2 , . . . , X tn ) (3.15)
with controlled Feynman-Kac weights variances

∀t n ≤ s < t n+1 E n,b ηn tn<k≤s G k (X k ) 2 E n,b ηn tn<k≤s G k (X k ) 2 < 1 + a n RR n°6700
and for

s = t n+1 E n,b ηn tn<k≤tn+1 G k (X k ) 2 E n,b ηn tn<k≤tn+1 G k (X k ) 2 = E n,b ηn G n+1 (X n+1 ) 2 E n,b ηn (G n+1 (X n+1 )) 2 ≥ 1 + a n .

Boltzmann's entropy

This criterion allows us to control the entropy like distance between the free motion trajectories and the weighted Feynman-Kac measures. To be more precise, set

H (n) tn,s P b ηn,(tn,s) = Ent P b ηn,(tn,s) | Q b ηn,(tn,s) = - tn<k≤s E n,b ηn (log G k (X k ))
with the weighted measures Q η,(p,n) defined by

Q η,(p,n) (d(x p , x p+1 , . . . , x n )) = P η,(p,n) (d(x p , x p+1 , . . . , x n )) × p<k≤n G k (x k ).
For (0, 1)-valued potential functions, the mappings

t n ≤ s → Ent P b ηn,(tn,s) | Q b ηn,(tn,s)
are clearly non-decreasing functions starting at 0 for s = t n . If we choose an interval I n := [a n , ∞) for some a n > 0 then the resampling time t n+1 coincides with the first time s the entropy distance goes above the level a n ; that is,

t n+1 := inf t n ≤ s : Ent P b ηn,(tn,s) | Q b ηn,(tn,s) ≥ a n .
This construction allows us to decompose the trajectories of the reference free motion X k into a sequence of excursion like processes

X 0 = X 0 and ∀n ≥ 1 X n := (X tn-1+1 , X tn-1+2 , . . . , X tn ) (3.

Convergence analysis of the reference SMC algorithm

In Section 4.1, we present the reference SMC algorithm which resamples at the deterministic times (t n ) n≥0 . The concentration analysis of this algorithm is developed in Section 4.2. In Section 4.2.1, we provide a brief discussion on the set of constants arising in our estimates and introduce assumptions ensuring a uniform control of these constants w.r.t the time horizon. In Section 4.2.2,

INRIA

several estimates on the bias and L m -mean errors are provided, whereas Section 4.2.3 gives an original non-asymptotic exponential concentration estimate.

In practice the times (t n ) n≥0 are unknown and are approximated, on-line, using the current SMC approximation. In Section 4.3, we use the concentration results obtained in Section 4.2 to analyze the convergence of these empirical approximations towards their limiting functional values.

A reference SMC algorithm

The SMC interpretation of the evolution equation (3.2) is the Markov chain

X (N ) n = X (N,1) n , X (N,2) n , . . . , X (N,N ) n ∈ S N n with elementary transitions P X (N ) n+1 ∈ dx n+1 | X (N ) n = N i=1 K n+1,η N n (X (N,i) n , dx i n+1 ) (4.1)
with

K n+1,η N n = S n,η N n M n+1 and η N n (•) := 1 N N j=1 δ X (N, j) n (•). (4.2)
Loosely speaking this integral decomposition shows that this SMC algorithm has the same updating/prediction nature as the one of the 'limiting' Feynman-Kac model. More precisely, the deterministic two-step updating/prediction transitions in distribution spaces

η n Sn,η n --------→ η n = η n S n,ηn = Ψ n (η n ) Mn+1 -------→ η n+1 = η n M n+1 (4.3)
have been replaced by a two-step resampling/mutation transition in a product space

X (N ) n ∈ E N tn resampling --------→ X (N ) n ∈ E N tn mutation -------→ X (N ) n+1 ∈ E N tn+1 (4.4)
where it is recalled that the state-space at the (n + 1)-th iteration is the path space given by

E tn+1 := E 0 × . . . × E tn+1 = E tn × E tn+1 × . . . × E tn+1 .
In our context, the SMC algorithm keeps track of all the paths of the sampled particles and the corresponding ancestral lines are denoted as follows

X (N,i) n = X (N,i) 0,tn , X (N,i) 1,tn , . . . , X (N,i) tn,tn X (N,i) n = X (N,i) 0,tn , X (N,i) 1,tn , . . . , X (N,i) tn,tn ∈ E tn := E 0 × . . . × E tn .
By definition of the reference Markov model X n given in (3.8), every path particle X

(N,i)

n+1 ∈ E tn+1 keeps track of the selected excursion X (N,i) n

∈ E tn and it

RR n°6700

evolves its terminal state X (N,i) tn,tn with (t n+1 -t n ) elementary random moves with the Markov transitions M k . More formally, we have that

X (N,i) n+1 =    X (N,i) 0,tn+1 , X (N,i) 1,tn+1 , . . . , X (N,i) tn,tn+1 , X (N,i) tn+1,tn+1 , . . . , X (N,i) tn+1,tn+1    = X (N,i) n , X (N,i) tn+1 , . . . , X (N,i) tn+1
with the (t n+1 -t n ) M k -exploration transitions given by the following synthetic diagram:

X (N,i) tn,tn M tn +1 ----→ X (N,i) tn+1 M tn +2 ----→ X (N,i) tn+2 M tn +3 ----→ . . . M t n+1 ----→ X (N,i) tn+1 .
From this discussion, it is worth mentioning a further convention that the particle empirical measures

η N n+1 (•) = 1 N N i=1 δ " b X (N,i) n , " X (N,i) tn +1 ,X (N,i) tn +2 ,...,X (N,i) t n+1
"" (•) are the terminal values at time s = t n+1 of the flow of random measures

t n ≤ s ≤ t n+1 → P N b η N n ,(tn,s) (•) = 1 N N i=1 δ " b X (N,i) n , " X (N,i) tn +1 ,X (N,i)
tn +2 ,...,X

(N,i) s

"" (•).

(4.5)

Concentration analysis

Introduction

This Section is concerned with the concentration analysis of the empirical measures η N n associated with (4.2) around their limiting values η n defined in (3.3). Our concentration estimates are expressed in terms of the following constants

q p,n = sup (x,y)∈S 2 p Q p,n (1)(x) Q p,n (1)(y)
and

β(P p,n ) := sup f ∈Osc1(Sn) osc(P p,n (f ))
where Q p,n is the linear semi-group (3.5) and P p,n its normalized version (3.6).

Note that these parameters can be expressed in terms of the mixing properties of the Markov transitions M n . It is clearly out of the scope of the present article to derive these estimates; we refer the interested reader to chapter 4 of [START_REF] Moral | Genealogical and interacting particle systems[END_REF]. Nevertheless, under appropriate mixing type properties we can prove that the series n p=0 q α p,n β(P p,n ) is uniformly bounded with respect to the final time horizon n, for any parameter α ≥ 0. Most of the results presented in this Section are expressed in terms of these series, including bias estimates, L m -mean error bounds, as well as exponential concentration estimates. As a result, these non-asymptotic results can be converted into uniform in time convergence results. To get a flavor of these uniform estimates, let us assume that the Markov transitions M k satisfy the following regularity property INRIA (M) m There exists an m ∈ Z + and a sequence (δ p ) p≥0 ∈ (0, 1)

N such that

∀p ≥ 0 ∀(x, y) ∈ S 2 p M p,p+m (x, •) ≥ δ p M p,p+m (y, •).
We also introduce the following quantities:

∀k ≤ l r k,l := sup k≤p<l G p (x p ) G p (y p )   ≤ k≤p<l q p  
with the collection of constants (q n ) ≥1 introduced in (3.1). In the above displayed formula, the supremum is taken over all admissible pair of paths with elementary transitions M p . Under the condition (M) m we have for any n ≥ m ≥ 1, and p ≥ 1, q p,p+n ≤ δ -1 p r p,p+m and β(P p,p+n )

≤ n/m -1 k=0 1 -δ 2 p r -1 p+km+1,p+ (k+1)m .
(4.6) The proof of these estimates relies on technical semi-group techniques; see chapter 4 of [START_REF] Moral | Genealogical and interacting particle systems[END_REF] for details. Several contraction inequalities can be deduced from these results. To understand this more closely, assume that (M) m is satisfied with m = 1, δ = ∧ n δ n > 0 and q = ∨ n≥1 q n . In this case, q p,p+n ≤ δ -1 q and β(P p,p+n

) ≤ 1 -δ 2 n ⇓ ∀α ≥ 0 n p=0 q α p,n β(P p,n ) ≤ q α /δ (2+α) .
More generally, assume (M) m is satisfied for some m ≥ 1 and that the parameters δ p and r k,l are such that ∧δ p := δ > 0 ∨ p r p,p+m := r < ∞ and ∨ p r p+1,p+m := r < ∞. (

In this situation, q p,p+n ≤ δ -1 r and β(P p,p+n ) ≤ 1 -δ 2 r -1 n/m and therefore

∀α ≥ 0 n p=0 q α p,n β(P p,n ) ≤ m r r α /δ (2+α) . (4.8)
It is well known that the mixing type condition (M) m is satisfied for any aperiodic and irreducible Markov chains on finite spaces, as well as for bi-Laplace exponential transitions associated with a bounded drift function. We also mention that this mixing condition is of course never met for the complete historical process X n = (X p ) 0≤p≤tn on S n = 0≤p≤tn E p discussed in Section 3.3. Nevertheless, under appropriate conditions on the Markov transitions M k , it is satisfied for the time marginal model associated with the excursion valued Markov chain model on tn-1<p≤tn E p . For instance, if

∀k ≥ 1 ∀(x, y) ∈ E p 2 M k (x, .) ≥ δ M k (y, .)
for some δ > 0, then condition (M) m is met with m = 1 and δ p = δ .

RR n°6700

4.2.2 Some L m -mean error bounds At this point, it is convenient to observe that the local sampling errors induced by the mean field particle model are expressed in terms of the collection of local random field models defined below. Definition 4.1 For any n ≥ 0 and any N ≥ 1, let V N n be the collection of random fields defined by the following stochastic perturbation formulae

η N n = η N n-1 K n,η N n-1 + 1 √ N V N n ⇐⇒ V N n := √ N η N n -η N n-1 K n,η N n-1
.

(4.9) For n = 0, the conventions K 0,η N -1 (x, dy) = η 0 (dy) and η N -1 K 0,η N -1 = η 0 are adopted.

In order to quantify high order L m -mean errors we need the following Kintchine type inequality for martingales with symmetric and independent increments. This result is more or less well known. For the convenience of the reader, a proof is provided in the appendix. Lemma 4.2 (Kintchine's inequality) Let M n := 0≤p≤n ∆ p be a real-valued martingale with symmetric and independent increments (∆ n ) n≥0 . For any integer m ≥ 1, and any n ≥ 0, we have

E (|M n | m ) 1 m ≤ b(m) E [M ]
E V N n (f n ) m F (N ) n-1 1 m ≤ b(m) osc(f n ). ( 4

.12)

Proof:

By construction, we have

V N n (f n ) = N i=1 ∆ (N ) n,i (f n ) with ∆ (N ) n,i (f n ) := 1 √ N f n (X (N,i) n ) -K n,η N n-1 (f n )(X (N,i) n-1 ) . Given X (N ) n-1 , we let (Y (N,i) n
) 1≤i≤N be an independent copy of (X

(N,i) n ) 1≤i≤N . We readily check that ∆ (N ) n,i (f n ) := E 1 √ N f n (X (N,i) n ) -f n (Y (N,i) n ) | X (N ) n-1

INRIA

with the pair of random measures (I N n , J N n ) given for any f n ∈ Osc 1 (S n ) by

I N n (f n ) := 2 n p=0 q p,n β(P p,n ) V N p U p,n,Φp(η N p-1 ) (f n ) J N n (f n ) := √ N n p=0 R p,n η N p , Φ p (η N p-1 ) (f n ). Now, observe that E W N n (f n ) = E(J N n (f n )). Using Proposition 4.3, for any f n ∈ Osc 1 (S n ) it can be checked that E I N n (f n ) m 1 m ≤ b(m) σ 2,n .
In much the same way, we find that

√ N E J N n (f n ) m 1 m ≤ b(2m) 2 σ 1,n . (4.13)
The end of the proof is now clear.

A concentration theorem

The following concentration theorem is the main result of this Section. It improves the exponential estimates presented in Section 7.4.3 in [START_REF] Moral | Genealogical and interacting particle systems[END_REF], by removing a multiplicative correction factor of order √ N in front of the exponential decay.

Theorem 4.6 For any n ≥ 0, f n ∈ Osc 1 (S n ), N ≥ 1, and any 0 ≤ ≤ 1/2,

P η N n -η n (f n ) ≥ ≤ 6 exp - N 2 8 σ 1,n (4.14)
where the constant σ 1,n is defined in Theorem 4.5.

In addition, suppose (M) m is satisfied for some m ≥ 1 and condition (4.7) holds true for some δ > 0 and some finite constants (r, r). In this situation, for any value of the time parameter n, for any f n ∈ Osc 1 (S n ), N ≥ 1 and for any ρ ∈ (0, 1) the probability that

η N n -η n (f n ) ≤ 4r δ 2 2mrr N δ log 6 ρ
is greater than (1 -ρ).

Proof:

We use the same notation as in the proof of Theorem 4.5. Recall that b(2m) 2m = E(X 2m ) for every centered Gaussian random variable with E(X 2 ) = 1 and

∀s ∈ [0, 1/2) E(exp sX 2 ) = m≥0 s m m! b(2m) 2m = 1 √ 1 -2s .
INRIA Using (4.13), for any f n ∈ Osc 1 (S n ) and 0 ≤ s < 1/(2σ 1,n ) it follows that

E exp s √ N J N n (f n ) ≤ m≥0 (sσ 1,n ) m m! b(2m) 2m = 1 1 -2sσ 1,n . (4.15)
To simplify the presentation, set f N p,n := U p,n,Φp(η N p-1 ) (f n ) and α p,n := 2 q p,n β(P p,n )

where U p,n,η (•) was introduced in Lemma 4.4. In this notation, it follows

V N p f N p,n = 1 √ N N i=1 f N p,n (X (N,i) p ) -K p,η N p-1 (f N p,n )(X (N,i) p-1 ) .
Recalling that E(e tX ) ≤ e t 2 c 2 /2 , for every real-valued and centered random variable X with |X| ≤ c, we prove that

E exp t α p,n V N p f N p,n X (N ) p-1 = N i=1 Sp K p,η N p-1 X (N,i) p-1 , dx e tαp,n √ N " f N p,n (x)-K p,η N p-1 (f N p,n )(X (N,i) p-1 ) « ≤ exp t 2 α 2 p,n 2
.

Iterating the argument, we find that

E e t I N n (fn) = E exp t n p=0 α p,n V N p f N p,n ≤ exp t 2 σ 2 n 2 (4.16) with σ 2 n := 4 n p=0 q 2 p,n β(P p,n ) 2 .
From these Laplace upper bounds, the proof of the exponential estimates now follows standard arguments. Indeed, for any 0 ≤ s < 1/(2σ 1,n ) and any > 0, by (4.15) we have

P √ N J N n (f n ) ≥ ≤ 1 1 -2sσ 1,n exp {-s}.
Replacing by N

P J N n (f n )/ √ N ≥ ≤ 1 1 -2sσ 1,n exp {-s N } .
Now, if choosing s = 3/(8σ 1,n ), then

P J N n (f n )/ √ N ≥ ≤ 2 exp {-N/(3σ 1,n )}.
To estimate the probability tails of I N n (f n ), we use (4.16) and the fact that > 0 and t ≥ 0

P I N n (f n ) ≥ ≤ exp -t - t 2 2 σ 2 n .
Now choosing t = /σ 2 n and replacing by √ N then we obtain

∀ > 0 P I N n (f )/ √ N ≥ ≤ exp - N 2 2σ 2 n .
Using the decomposition

η N n -η n = I N n / √ N + J N n / √ N
we find that for any parameter α ∈ [0, 1]

P η N n -η n (f n ) ≥ ≤ P I N n (f n )/ √ N ≥ α +P J N n (f n )/ √ N ≥ (1 -α) .
From previous calculations,

P η N n -η n (f n ) ≥ ≤ exp - N 2 α 2 2σ 2 n + 2 exp - N (1 -α) 3σ 1,n . (4.17)
Now choose α = (1 -)(≥ 1/2), then α 2 ≥ 1/4 and

P η N n -η n (f n ) ≥ ≤ exp - N 2 8σ 2 n + 2 exp - N 2 3σ 1,n ≤ 3 exp - N 2 8 (σ 1,n ∨ σ 2 n )
.

It remains to observe that q p,n ≥ 1 and

β(P p,n ) ≤ 1 =⇒ σ 2 n ≤ σ 1,n and η N n -η n (f n ) ≥ ⇐⇒ η N n -η n (f n ) ≥ or η N n -η n (f n ) ≤ - ⇐⇒ η N n -η n (f n ) ≥ or η N n -η n (-f n ) ≥ so that P η N n -η n (f n ) ≥ ≤ P η N n -η n (f n ) ≥ +P η N n -η n (-f n ) ≥ .
The end of the proof of (4.14) is now easily completed. We further assume that the mixing condition (M) m is satisfied for some m ≥ 1 and condition (4.7) holds true for some δ > 0 and some finite constants (r, r). In this situation, by (4.8) the following uniform concentration estimate holds

sup n≥0 P η N n -η n (f n ) ≥ ≤ 6 exp - N 2 δ 5 32 m r r 3 .
The proof of the theorem is concluded by choosing :=

n ( ) = b n 2a n ( )   1 + 4a n ( ) b n -1   = σ 2 n 3σ 1,n 1 1 + 6σ 1,n σ 2 n -1 (-→ ↓0 1)
and therefore

∀ ≥ 0 P η N n -η n (f ) ≥ ≤ 3 exp -N 2 2σ 2 n α 2 n ( ) .

Approximating the criteria

E s := E tn × E tn+1 × . . . × E s .
For instance, the following proposition is more or less a direct consequence of Theorem 4.6. Proposition 4.8 For any N ≥ 1, n ≥ 0, t n ≤ s ≤ t n+1 and any > 0, the concentration inequality:

P H (n) tn,s P N b η N n ,(tn,s) -H (n) tn,s P b ηn,(tn,s) ≥ ≤ (1 + N/2) exp - N 2 c(n)
holds for some finite constant c(n) < ∞ whose values only depend on the time parameter. In addition, when the measures H

(n) tn,s have a finite support, the concentration inequality:

P H (n) tn,s P N b η N n ,(tn,s) -H (n) tn,s P b ηn,(tn,s) ≥ ≤ c 1 (n) exp - N 2 c 2 (n)
also holds, with a pair of finite constants c 1 (n), c 2 (n) < ∞ whose values only depend on the time parameter.

Proof: By theorem 7.4.4 in [START_REF] Moral | Genealogical and interacting particle systems[END_REF], for any N ≥ 1, p ≥ 1, n ≥ 0 and any test function

f n ∈ Osc 1 (E tn ), we have sup N ≥1 √ N E η N n (f n ) -η n (f n ) p 1 p ≤ b(p) c(n)
with some finite constant c(n) < ∞ and with the collection of constants b(p) defined in (4.11). These estimates clearly imply that for any t n ≤ s ≤ t n+1 , and any test function

h n ∈ Osc 1 (E s ) sup N ≥1 √ N E P N b η N n ,(tn,s) (h n ) -P b ηn,(tn,s) (h n ) p 1 p ≤ b(p) c(n).
Under our Lipschitz regularity condition (3.10) on the criteria type functionals H

(n)

tn,s and using the generalized integral Minkowski inequality, it can be concluded that sup

N ≥1 √ N E H (n) tn,s P N b η N n ,(tn,s) -H (n) tn,s P b ηn,(tn,s) p 1 p ≤ b(p) c(n) δ(H (n)
tn,s ).

The proof of the exponential estimate follows exactly the same lines of arguments as the ones used in the proof of corollary 7.4.3 in [START_REF] Moral | Genealogical and interacting particle systems[END_REF]; thus it is omitted. The last assertion is a direct consequence of Theorem 4.6. This ends the proof of the proposition

An online adaptive SMC algorithm

The above proposition shows that the functional criteria H n ,(tn,s) , up to an exponentially small probability. Therefore, as we cannot compute the deterministic resampling times (t n ), it is necessary to approximate the reference particle model: Definition 4.9 We use the letters Y and Y to define the particle model associated to the following approximation of the reference particle model. More precisely, the particle systems Y

(N ) = Y (N,i) , Y (N ) = Y (N,i) , Y (N,i) s,t and Y (N,i) s,t are defined as X (N ) = X (N,i) , X (N ) = X (N,i) , and X (N,i) s,t and X (N,i) s,t
by replacing in the inductive construction of the deterministic sequence (t n ) n≥0 the measures P b ηn,(tn,s) by their current N -particle approximation measures

P N b η N n ,(t N n ,s) (•) := 1 N N i=1 δ " b Y (N,i) n , " Y (N,i) t N n +1 ,Y (N,i) t N n +2 ,...,Y (N,i) s «« (•).
Here n . We also assume that both models are constructed in such a way that they coincide on every time interval 0 ≤ n ≤ m, as soon as the random times t N n = t n , for every 0 ≤ n ≤ m. are inductively constructed using random times t N n based on the current empirical criteria

η N n (•) = 1 N N i=1 δ b Y (N,i) n (
H (n-1) t N n-1 ,s P N b η N n-1 ,(t N n-1 ,s) .
By construction, for the four functional criteria discussed in Section 3.6, we have H

(n) t N n ,s P N b η N n ,(t N n ,s) = C N t N n ,s where C N t N
n ,s are the empirical criteria discussed in 2.3.

Asymptotic analysis 5.1 A key approximation lemma

To go one step further in our discussion, it is convenient to introduce the following collection of events. Definition 5.1 For any δ ∈ (0, 1), m ≥ 0, a n ∈ R and N ≥ 1, we denote by Ω N m (δ, (a n ) 0≤n≤m ), the collection of events defined by:

Ω N m (δ, (a n ) 0≤n≤m ) := {∀0 ≤ n ≤ m , ∀t n ≤ s ≤ t n+1 H (

RR n°6700

To avoid this technical problem, one natural strategy is to introduce randomized criteria thresholds. We further assume that the parameters (a n ) n≥0 are sampled realizations of a collection of absolutely continuous random variables (A n ) n≥0 ; in practice, this element can be incorporated into the simulation, with no real change in the performance. The main simplification of these randomized criteria comes from the fact that the parameters m := inf We conclude that for almost every realization (A n ) 0≤n≤m = (a n ) 0≤n≤m the pair particle models (X

Functional central limit theorems

To demonstrate the impact of this functional fluctuation result and to connect this work with related studies, we provide a brief discussion on the proof of this multivariate central limit theorem.

Let us begin by recalling a rather well-known functional fluctuation theorem of the local errors associated with the mean field particle approximation introduced in (4.9). This result was initially presented in [START_REF] Moral | Branching and Interacting Particle Systems Approximations of Feynman-Kac Formulae[END_REF] in the case of elementary genetic type models. It has been extended to general mean field particle models in Section 9.3 in [START_REF] Moral | Genealogical and interacting particle systems[END_REF]. Theorem 6.1 For any fixed time horizon n ≥ 0, the sequence (V N p ) 0≤p≤n converges in law, as N tends to infinity, to a sequence of n independent, Gaussian and centered random fields (V p ) 0≤p≤n with, for any f p , g p ∈ B b (E p ), and 1 ≤ p ≤ n, E(V p (f p )V p (g p )) = η p-1 K p,ηp-1 ([f p -K p,ηp-1 (f p )][g p -K p,ηp-1 (g p )]). ( 6 (1)

n+1 (f n+1 ) η N n+1 T (1)
n+1

N ↑∞ η n+1,s (f n+1 ) :=

η n+1 T (1)
n+1 (f n+1 )

η n+1 T (1)
n+1 [START_REF] Cerou | Adaptive multilevel splitting for rare event analysis[END_REF] where η n+1 is the flow of Feynman-Kac measures on path spaces introduced in Section 3.3. Since the adaptive interaction time is taken such that t N n+1 = t n+1 , it holds that 1

N N i=1 δ b Y (N,i) n+1
N ↑∞ η n+1,tn+1 = η n+1 .

In other words, the marginal type functions are chosen The fluctuation analysis of these random fields strongly relies on the functional central limit theorem stated in corollary 6.4. In particular, the fluctuations of the random fields (6.3) depend on the fluctuations of a pair or random fields. From elementary stochastic calculus, the corresponding fluctuations cannot be deduced directly from the separate fluctuations of each components. For odd integers, apply the Cauchy-Schwarz inequality twice to deduce that

E(|M n | 2m+1 ) 2 ≤ E(M 2m n ) E(M 2(m+1) n ) ≤ 2 -(2m+1) (2m) m (2(m + 1)) (m+1) E [M ] m+1 n 2m+1 m+1 .
Since we also have that it can be concluded that (m+1) .

E(|M n | 2m+1 ) ≤ 2 -(m+1/2) (2m + 1) (m+1) m + 1/2 E [M ] m+1 n 1- 1 2
This ends the proof of the lemma.

Proof of Lemma 4.4

. 3 )

 3 Returning to the examples discussed above, we can choose the exponential Boltzmann-Gibbs potential and respectively the indicator potential functions INRIA defined by

Ψ

  n (µ)(dy) = (µS n,µ ) (dy) = Sn µ(dx) S n,µ (x, dy) with the accept-reject type selection transition S n,µ (x, dy) := G n (x) δ x (dy) + (1 -G n (x)) Ψ n (µ)(dy).

  tn,tn) = 1. When the potential functions take values in [0, 1), one concludes that the mapping s ≥ t n → H

  tn,tn) = 0. The variance mappings s → H (n) tn,s P b ηn,(tn,s) are generally increasing. For instance, in the case of indicator potential functions

  tn,s) | Q b ηn,(tn,s) ≥ a n .

Proposition 4 . 3

 43 where m stands for the smallest even integer m ≥ m and (b(m)) m≥1 is the collection of constants given below: b(2m) 2m := (2m) m 2 -m and b(2m + 1) 2m+1 := (2m + 1) (m+1) For any N ≥ 1, m ≥ 1, n ≥ 0 and any test function f n ∈ B b (S n ) we have the almost sure estimate

  By construction, the particle occupation measures P N b ηn,(tn,s) approximate the measures P b ηn,(tn,s) introduced in Definition 3.2; that is, in some senseP N b η N n ,(tn,s) N ↑∞ P b ηn,(tn,s) . Conversely, observe that P N b η N n ,(tn,s) , and respectively P b ηn,(tn,s) are the marginals of the measures η N n+1 , and respectively η n+1 w.r.t the (s -t n ) + 1 first coordinates. In other words, the measures P N b η N n ,(tn,s) , and respectively P b ηn,(tn,s) are the projections of the measures η N n+1 , and respectively η n+1 on the state space

  tn,s) can be approximated byH (n) tn,s P N b η N

  •) denotes the updated occupation measure of the particle model Y (N )

INRIA

  It is emphasized that the measures P N b η N n ,(t N n ,s) generally differ from the reference empirical measures P N b η N n ,(tn,s) introduced in (4.5). Indeed, the reference measures P N b η N n ,(tn,s) are built using the deterministic times t n based on the functional criteria H (n-1) tn-1,s P b ηn-1,(tn-1,s) , whilst the empirical measures P N b η N n ,(t N n ,s)

Lemma 5 . 2

 52 tn,s) -a n . On the event Ω N m (δ, (a n ) 0≤n≤m ), for any n ≤ m and for any t n ≤ s ≤ t n+1 , we haveH (n) tn,s P b ηn,(tn,s) > a n =⇒ H (n) tn,s P N b η N n ,(tn,s) > a n H (n) tn,s P b ηn,(tn,s) < a n =⇒ H (n) tn,s P N b η N n ,(tn,s) < a n .tn,s P N b η N n ,(tn,s) ≤ H (n) tn,s P b ηn,(tn,s) + δ H (n) tn,s P b ηn,(tn,s) -a n .

-N δ 2 2 m c 2 ≤ c 1

 221 tn,s) -a n are strictly positive for almost every realizations A n = a n of the threshold parameters.Theorem 5.4 For almost every realization of the random threshold parameters, and for any δ ∈ (0, 1), we have the following exponential estimates :P ∃0 ≤ n ≤ m t N n = t n | (A n ) 0≤n≤m ≤ c 1 (m) 1 + δ m N 2 e (m)for some constants c 1 (m), c 2 (m) < ∞ whose values only depend on the time parameter m. In addition, when the measuresH (n)tn,s have a finite support, for any δ ∈]0, 1/(2 m )[P ∃0 ≤ n ≤ m t N n = t n | (A n ) 0≤n≤m ≤ c 1 (m) e -N δ 2 2 m /c 2 (m)holds, for a possibly different pair of finite constants c 1 (m), c 2 (m) < ∞.Proof:Using Proposition 4.8, we obtain the rather crude estimateP Ω -Ω N m (δ, (A n ) 0≤n≤m ) | (A n ) 0≤n≤m = (a n ) 0≤n≤m (m) (1 + δ m N/2) exp -N δ 2 2m /c 2 (m) for a pair of finite constants c 1 (m), c 2 (m) < ∞ whose values only depend on the time parameter. The final line is a direct consequence of Proposition 4.8 and an application of Proposition 5.3 completes the proof.

n

  ) 0≤n≤m and (Y(N ) , Y(N ) n ) 0≤n≤m only differ on events Ω -Ω N m (δ, (a n ) 0≤n≤m) with exponentially small probabilities:P ∃0 ≤ n ≤ m (Y (N ) , Y n ) | (A n ) 0≤n≤m = (a n ) 0≤n≤m ≤ c 1 (m) (1 + δ m N/2) exp -N δ 2 2m /c 2 (m) .

6 A functional central limit theorem 6 . 1 .

 61 A direct approachIn this Section some direct consequences of the exponential coupling estimates, presented in Section 5.2, are discussed. For almost every realization (A n ) 0≤n≤m = (a n ) 0≤n≤m and for any test function f n ∈ Osc 1 (E tn ) the following decomposition holds (writing η N n for the online adaptive approximation introduced in Definition4n -η N n (f n ) 1 Ω-Ω N m (δ,(an) 0≤n≤m ) ≤ √ N P Ω -Ω N m (δ, (a n ) 0≤n≤m ) → N ↑∞ 0Thus we can conclude directly that, for almost every realization (A n ) 0≤n≤m = (a n ) 0≤n≤m , the random fields n -η n = and W N n := √ N η N n -η n converge in law, as N ↑ ∞, to the same Gaussian and centered random field W n .

. 1 ) 1 )

 11 Using arguments fairly similar to those in the proof of Lemma 4.4, we obtain the first order decomposition formula:[Φ n (µ) -Φ n (η)](f ) = (µ -η)D n,η + R n (µ, η) with the signed measure R n (µ, η) given by R n (µ, η)(f ) := -1 µ(G n,η ) [µ-η] ⊗2 (G n,η ⊗D n,η (f )) with G n,η := G n-1 /η(G n-1 )and the bounded integral operatorD n,η (f )(x) := G n,η (x) × M n (f -Φ n (η)(f )) (x).INRIA Indeed, if we choose the weighted functionsT (n+1 (f n+1 )(x 0 , . . . , x tn+1 ) := f n+1 (x 0 , . . . , x tn+1 ) tn<k≤s G k (x k )then in some senseη N n+1,s (f n+1 ) := η N n+1 T

T ( 0 )N

 0 n+1 (f n+1 )(x 0 , . . . , x tn+2 ) := f n+1 (x 0 , . . . , x tn+1 ) soη n+2 T (0) n+1 (f n+1 ) = η n+1 (f n+1 ) ∝ E   f n+1 (X 0 , . . . , X tn+1 ) 0<k≤tn+1 G k (X k ) ↑∞ η n+2 T (0) n+1 (f n+1 ) .From the previous discussion, for almost every realization (A n ) 0≤n≤m = (a n ) 0≤n≤m a central limit theorem is easily derived for the collection of random fields given below W (f n+1 ) -η n+1,s (f n+1 ) as well as for the mixture of random field sequences

INRIA

  to conclude thatE (M n-1 + ∆ n ) 2m ≤ 2 -m (2m)

Using ( 3 . 7 )

 37 , for any f ∈ B b (S n+1 ) we find that[Φ p,n (µ) -Φ p,n (η)](f ) = 1 µ(G p,n,η ) (µ -η)D p,n,η (f )with the bounded integral operatorD p,n,η (f )(x) := G p,n,η (x) × P p,n (f -Φ p,n (η)(f )) (x)where G p,n,η and P p,n stand for the potential function and the Markov operator is given byG p,n,η := Q p,n (1)/η(Q p,n (1)) and P p,n (f ) = Q p,n (f )/Q p,n (1). Now, since η (G p,n,η ) = 1, it is now immediate to check that [Φ p,n (µ) -Φ p,n (η)] = (µ -η)D p,n,η + R p,n (µ, η)with the measure R p,n (µ, η) given by R p,n (µ, η)(f ) := -1 µ(G p,n,η )[µ -η] ⊗2 (G p,n,η ⊗ D p,n,η (f )).

INRIA

  1 δnV log (1/a n ) time steps. Alternatively, if we want to ensure that (2.11) is met in at least m n time steps then δ n

  The resampling time t n+1 is defined as the first time s ≥ t n the quantity H

	Given t n+1 (> t n ), set
	η n+1 = P	b ηn,(tn,tn+1)
			a n ]
	such that
		H	(n) tn,tn P ηn,(tn,tn) ∈ I n . b
	(n) tn,s P
	RR n°6700

b ηn,(tn,s) hits the set I n ; that is

t n+1 := inf t n ≤ s : H (n) tn,s P b ηn,(tn,s) ∈ I n

  Returning to the end of the proof of Theorem 4.6, it can be observed that the exponential concentration estimates can be improved, marginally, by choosing in (4.17) the parameter α = α n ( ) ∈ [0, 1] such that a n ( ) α 2 = b n (1 -α), with a n ( ) := 2σ 2

	n	and b n = 1 3σ1,n . Elementary manipulations yield
	α
		1 √ N	4r δ 2	2mrr δ	log 6 ρ .

INRIA

Remark 4.7

RR n°6700

Acknowledgement

The third author was affiliated to the Institute of Statistical Mathematics (Tokyo, Japan), whilst this project was initiated.

This yields the formula

where M

(N ) n,N (f n ) is the terminal value of the martingale sequence defined below i ∈ {1, . . . , N } → M (N )

Using Kintchine's inequality the end of the proof of (4.12) is clear.

The proof of the following lemma is rather technical and is provided in the appendix. Lemma 4.4 For any 0 ≤ p ≤ n, any η, µ ∈ P(S p ) and any f n ∈ Osc 1 (S n ), we have the first order decomposition for the non-linear semi-group Φ p,n defined in (3.7):

and a collection of functions U p,n,η (f ), V p,n,η (f ), W p,n,η (f ) in Osc 1 (S p ), whose values only depends on the parameters (p, n, η).

We now present a bias estimate and some L m bounds of independent interest.

Theorem 4.5 For any n ≥ 0, f n ∈ Osc 1 (S n ) and any N ≥ 1,

In addition, for any m ≥ 1 we have

q p,n β(P p,n ).

Proof:

Using Lemma 4.4, we have the telescoping sum decomposition This ends the proof of the lemma.

Proposition 5.3 Assume that the threshold parameters a n are chosen so that H

(n) tn,s P b ηn,(tn,s) = a n , for any n ≥ 0. In this situation, for any δ ∈ (0, 1), m ≥ 0 and N ≥ 1, we have

Proof:

This result is proved by induction on the parameter m ≥ 0. Under our assumptions, for m = 0 we have t N 0 = t 0 = 0. Thus, by our coupling construction the pair of particle models coincide up to the time t N 1 ∧ t 1 . Therefore,

By Lemma 5.2, on the event Ω N m (δ, (a n ) 0≤n≤m ) we have t N 1 = t 1 . This proves the inclusion for m = 0 and m = 1. Suppose the result is true at rank m. Thus, on the event Ω N m (δ, (a n ) 0≤n≤m ) it is the case that t N n = t n , for any 0 ≤ n ≤ m. By our coupling construction, the pair of particle models coincides up to t N m+1 ∧ t m+1 ; that is,

Once again, by Lemma 5.2, on the event Ω N m+1 (δ, (a n ) 0≤n≤m+1 ) it also follows that t N m+1 = t m+1 . This ends the proof of the proposition.

Randomized criteria

The situation where the threshold parameters coincide with the adaptive criteria values H

(n) tn,s P b ηn,(tn,s) = a n cannot be handled using our analysis. Indeed, in this rather degenerate situation it is impossible, at least with our approach, to compare the random times t N n n≥0 with the reference deterministic ones (t n ) n≥0 . The analysis of this situation is clearly more involved since it requires us to control both the empirical approximating criteria and the mean field approximation model. For p = n, we use the convention D n,n = Id, the identity operator.

It should be noted that the semi-group D p,n introduced above can be explicitly described in terms of the semi-group Q p,n , with the following formula

The next technical lemma provides a key first order decomposition of the random fields W N n in terms of the local fluctuation errors. Its proof follows an elementary induction on the time parameter n, and is in the appendix. Note that R p is a remainder term that can be understood in the proof. Lemma 6.3 For any N ≥ 1 and any 0 ≤ p ≤ n, we have

Using the L m -mean error estimates presented in Section 4.2.2, it is easily proved that the sequence of remainder random fields R N n , introduced in (6.2), converge in law, in the sense of finite distributions, to the null random field as N tends to ∞. Therefore the fluctuations of W N n follows from Theorem 6.1.

Corollary 6.4 For any fixed time horizon n ≥ 0, the sequence of random fields (W N n ) n≥0 converges in law, as N tends to infinity, to the sequence of Gaussian and centered random fields (W n ) n≥0 given below

On the fluctuations of weighted occupation measures

We end this article with some comments on the fluctuations of weighted occupation measures on path spaces. Returning to the online adaptive particle model, given a pair of time steps (t N n , t N n+1 ) = (t n , t n+1 ) the N -particle measures

can be used to approximate the flow of updated Feynman-Kac path distributions ( η n+1,s ) tn≤s≤tn+1 given for any bounded test function

Related work

To our knowledge, the article [START_REF] Douc | Limit theorems for weighted samples with applications to sequential Monte Carlo[END_REF] is the only reference providing convergence results for a specific adaptive SMC scheme in a state-space models context. The authors establish a central limit theorem using an inductive proof w.r.t deterministic time periods; this does not account for the randomness of the latter quantities. They implicitly avoid the degenerate situation where the threshold parameter coincides with the limiting functional criterion. The fluctuations of the sum of the weighted and the resampled empirical approximations is deduced from the fluctuations of each separate components without the use of any multivariate type central limit theorem.

Alternative adaptive schemes can be found in [START_REF] Cornebise | Adaptive methods for sequential importance sampling with applications to state-space models[END_REF]. In that paper, the authors focus upon adapting the parameters of the instrumental kernel, so as to minimize the Kullback-Liebler or Chi-Square distance between this quantity and target distribution of interest. As such, this work is different from ours, other than the idea of using adaptive criteria within the simulation algorithm. An idea that combines this, in a different context, with the adaptation of target densities and dynamic resampling can be found in [START_REF] Jasra | Inference for Lévy driven stochastic volatility models via adaptive sequential Monte Carlo[END_REF].

Appendix

Proof of Lemma 4.2

We prove the lemma by induction on the parameter n. The result is clearly satisfied for n = 0. Suppose the estimate (4.10) is true at rank (n -1). To prove the result at rank n, use the binomial decomposition

Using the symmetry condition, all the odd moments of ∆ n are null. Consequently,

Then, via the induction hypothesis, the above expression is upper-bounded by the quantity

To take the final step, use the fact that

Using the fact that

Finally, for any f ∈ Osc 1 (S n ) observe that

with the pair of functions

This ends the proof of the lemma.

Proof of Lemma 6.3

The lemma is proved by induction on the time parameter n. For n = 0, we readily find that

) , with the convention Φ 0 (η N -1 ) = η 0 . Assuming that the decomposition is satisfied at rank n it may be proved that

Using the convention D n+1,n+1 = I it follows that

from which we conclude that (6.2) is satisfied at rank (n + 1). This ends the proof of the lemma.