Pierre Del Moral
email: pierre.del-moral@inria.fr

Arnaud Doucet
email: arnaud@stat.ubc.ca

Keywords: Genetic algorithms, particle filters, mean field particle models, Feynman-Kac formula, Boltzmann-Gibbs measures, nonlinear Markov chains, interacting processes, genealogical tree based algorithms, simulated annealing, central limit theorem Algorithmes génétiques, filtres particulaires, modèles de champ moyen, formules de Feynman-Kac, mesures de Boltzmann-Gibbs, chaînes de Markov non linéaires, processus en interaction, modèles d

Interacting particle methods are increasingly used to sample from complex high-dimensional distributions. They have found a wide range of applications in applied probability, Bayesian statistics and information engineering. Understanding rigorously these new Monte Carlo simulation tools leads to fascinating mathematics related to Feynman-Kac path integral theory and their interacting particle interpretations. In these lecture notes, we provide a pedagogical introduction to the stochastic modeling and the theoretical analysis of these particle algorithms. We also illustrate these methods through several applications including random walk confinements, particle absorption models, nonlinear filtering, stochastic optimization, combinatorial counting and directed polymer models.

Particle methods: An introduction with applications

Piere Del Moral , Arnaud Doucet July 14, 2009

Introduction

Interacting particle methods are a class of Monte Carlo methods to sample from complex highdimensional probability distributions and to estimate their normalizing constants. This class of algorithms approximate the target probability distributions by a large cloud of random samples termed particles. Practically, the particles evolve randomly around the space independently and to each particle is associated a positive potential function. Periodically we duplicate particles with high potentials at the expense of particles with low potentials which die. This intuitive genetic mutation-selection type mechanism has appeared in numerous applications ranging from nonlinear filtering [START_REF] Carpenter | An improved particle filter for non-linear problems[END_REF][START_REF] Moral | Feynman-Kac Formulae: Genealogical and Interacting Particle Systems with Applications[END_REF][START_REF] Doucet | On sequential Monte Carlo sampling methods for Bayesian filtering[END_REF][START_REF]Sequential Monte Carlo Methods in Practice[END_REF][START_REF] Doucet | A tutorial on particle filtering and smoothing: fifteen years later[END_REF][START_REF] Gordon | A novel approach to state estimation to nonlinear non-Gaussian state estimation[END_REF][START_REF] Kitagawa | Monte Carlo filter and smoother for non-Gaussian nonlinear state space models[END_REF][START_REF] Künsch | State-space and hidden Markov models[END_REF], Bayesian statistics [START_REF] Chopin | A sequential particle filter method for static models[END_REF][START_REF] Del Moral | Sequential Monte Carlo samplers[END_REF][START_REF] Fearnhead | Computational methods for complex stochastic systems: A review of some alternatives to MCMC[END_REF][START_REF] Liu | Monte Carlo Strategies in Scientific Computing[END_REF], combinatorial counting [START_REF] Aldous | Go with the winners algorithms[END_REF], molecular and polymer simulation [START_REF] Grassberger | Pruned-enriched Rosenbluth method: Simulations of θ polymers of chain length up to 1 000 000[END_REF], rare events simulation [START_REF] Cérou | Genealogical models in entrance times rare event analysis[END_REF][START_REF] Cérou | A non asymptotic variance theorem for unnormalized Feynman-Kac particle models[END_REF][START_REF] Glasserman | Multilevel splitting for estimating rare event probabilities[END_REF], quantum Monte Carlo methods [START_REF] Assaraf | Diffusion Monte Carlo methods with a fixed number of walkers[END_REF][START_REF] Hetherington | Observations on the statistical iteration of matrices[END_REF][START_REF] Melik-Alaverdian | Quantum Monte Carlo methods in statistical mechanics[END_REF] and genetic algorithms [START_REF] Goldberg | Genetic Algorithms in Search, Optimization and Machine Learning[END_REF][START_REF] Holland | Adaptation in Natural and Artificial Systems[END_REF] among others.

From a mathematical point of view, these methods can be interpreted as stochastic numerical approximations of Feynman-Kac measures. Feynman-Kac measures represent the distribution of the paths of a reference Markov process, weighted by a collection of potential functions. These functional models are natural mathematical extensions of the traditional change of probability measures, commonly used in importance sampling. The particle interpretation consists in evolving a population of particles mimicking natural evolution mechanisms. During the mutation stage, the particles evolve independently of one another, according to the same probability transitions as the ones of the reference Markov chain. During the selection stage, each particle evaluates the potential value of its location. The ones with small relative values are killed, while the ones with high relative values are multiplied. The corresponding genealogical tree occupation measure converges, as the population size tends to infinity, to the complete Feynman-Kac distribution on path space.

The two authors and their collaborators started working on the design and the mathematical analysis of these algorithms in the mid 90's (see for instance [START_REF] Moral | Feynman-Kac Formulae: Genealogical and Interacting Particle Systems with Applications[END_REF][START_REF]Sequential Monte Carlo Methods in Practice[END_REF], and references therein). Over the past few years, the popularity of these computationally intensive methods has dramatically increased thanks to the availability of cheap powerful computers. In particular in signal processing and machine learning, these algorithms are now widely used to solve nonlinear filtering problems. In this context, they are known as particle filters. The mathematical analysis of these algorithms offers a rigorous and unifying framework to analyze the convergence of numerous heuristic-like algorithms currently used in physics, statistics and engineering. It applies to any problem which can be translated in terms of functional Feynman-Kac type measures.

In this set of lecture notes, we provide a pedagogical introduction to the stochastic modeling and the theoretical analysis of these interacting particle algorithms. In a first section, section 1, we present several application areas and provide a detailed description of the corresponding interacting particle algorithms. Section 2 provides a brief treatise on Feynman-Kac modeling techniques. The last section, section 3, gives an overview of some convergence results, including variance and L p -mean error estimates, fluctuations and concentration properties. We have tried to give a brief "exposé" of the modern mathematical theory that is useful for the analysis of the asymptotic behavior of Feynman-Kac and particle models. To simplify the presentation and to clarify the main ideas behind these stochastic models, we have chosen to restrict these lectures notes to finite or countable state space models, avoiding any measure theory irrelevancies. In this simplified framework, we develop a rigorous mathematical analysis only involving vector and matrix operations. We emphasize that all of these particle models and the associated convergence results can be extended to general state-space models, including path-space models and excursion spaces on measurable state spaces.

We undertook this project for two main reasons: First, we felt that there was no accessible and pedagogical treatment on stochastic particle models and their application areas. One objective in writing these lecture notes is to throw some new light on some interesting links between physical, engineering, statistical and mathematical domains which appear disconnected at first glance.

INRIA

Second, the mathematical concepts and models are now at a point where they provide a very natural and unifying mathematical basis for a large class of Monte Carlo algorithms. We hope that this unifying point of view will help to develop fruitfully this field further.

Examples

The list of applications discussed here is by no means exhaustive and it just reflects the scientific interests of the authors.

Random walks confined in a set

We consider a symmetric random walk X n on the integers Z starting at the origin X 0 = 0. More formally, we take independent random variables U n , where P (U n = 1) = P (U n = -1) = 1/2 and we set X n = X 0 + 1≤p≤n U p . We fix A = {-a + 1, -a + 2, ..., a -1}, with a ∈ N. We want to compute the conditional distributions Law ((X 0 , . . . , X n ) | ∀0 ≤ p ≤ n, X p ∈ A ) (1.1) as well as the quantities

Z n := P (∀0 ≤ p < n, X p ∈ A)
This problem can be solved by simulation using the following particle algorithm. We start with N particles at the origin denoted by ξ i 0 = 0, with i = 1, . . . , N . Each of them evolve ξ i 0 ξ i 1 according to one transition of the random walk; more formally, we sample N independent copies (U i 1 ) 1≤i≤N of the random variables U 1 , and we set ξ i 1 = ξ i 0 + U i 1 . We denote

η N 1 (1 A ) = 1 N 1≤i≤N 1 A (ξ i 1 ) = 1 N Card 1 ≤ i ≤ N : ξ i 1 ∈ A
the proportion of points ξ i 1 in the set A. We define from the sample population ξ i 1 1≤i≤N a new population of N individuals ξ i 1 1≤i≤N as follows. For each i = 1, . . . , N , we perform the following operation: If ξ i 1 ∈ A, we set ξ i 1 = ξ i 1 . If ξ i 1 ∈ A, we pick randomly an individual ξi 1 among those ξ j 1 in the set A and we set ξ i 1 = ξi 1 . In other words, individuals within A do not move, while the individuals outside A are replaced by a randomly chosen individual among those in the set A. It may happen that all individuals ξ i 1 are outside of the set A. In this case, the algorithm stops and we set τ N = 1 to report the time of this event. If the algorithm has not stopped, we have a new configuration ξ i

1 1≤i≤N
of N individuals in the set A. We evolve ξ i 1 ξ i 2 according to one transition of the random walk; that is we sample N independent copies (U i 2 ) 1≤i≤N of the random variables U 2 , we set ξ i 2 = ξ i 1 + U i 2 and we define

η N 2 (1 A ) = 1 N 1≤i≤N 1 A (ξ i 2 ) = 1 N Card 1 ≤ i ≤ N : ξ i 2 ∈ A .
As before, we define from the sample population

ξ i 2 1≤i≤N a new population of N individuals ξ i 2 1≤i≤N
: individuals within A do not move, while the individuals outside the desired set are replaced by a randomly chosen individual among those in the set A. If all individuals ξ i 2 fall are outside of the set A, we set τ N = 2. Iterating this stochastic process, for every time n (< τ N ), we define a sequence of genetic type populations

ξ n := ξ i n 1≤i≤N ∈ Z N selection --------→ ξ n := ξ i n 1≤i≤N ∈ Z N mutation -------→ ξ n+1 ∈ Z N (1.2)
This stochastic algorithm can be interpreted as a genetic type model with mutation transitions given by the one of a symmetric random walk and an acceptance-rejection selection type transition associated with the potential indicator type function 1 A . Several estimates can be extracted from this interacting sampling algorithm.

First, we mention that the stopping time τ N tends to infinity as the size of the population N → ∞. More precisely, the probability that the algorithm stops at a given time n tends to zero exponentially fast, as N tends to infinity. More interestingly, the product of the proportions of surviving particles at each time step

Z N n := 0≤p<n η N p (1 A )
is asymptotically a consistent estimate of the quantity P n (A) and it is unbiased; that is we have lim

N →∞ Z N n = Z n and E Z N n = Z n (1.
3)

The convergence on the l.h.s. is an almost sure asymptotic convergence. It can be made precise by non asymptotic estimates including non asymptotic variance estimates and more refined exponential type deviations. If we interpret the selection transition as a birth and death process, then the important notion of the ancestral line of a current individual arises. More precisely, when a particle ξ i n-1 -→ ξ i n evolves to a new location ξ i n , we can interpret ξ i n-1 as the parent of ξ i n . Looking backwards in time and recalling that the particle ξ i n-1 has selected a site ξ j n-1 in the configuration at time (n -1), we can interpret this site ξ j n-1 as the parent of ξ i n-1 and therefore as the ancestor

ξ i n-1,n at level (n -1) of ξ i n .
Running back in time we can construct the whole ancestral line

ξ i 0,n ←-ξ i 1,n ←-. . . ←-ξ i n-1,n ←-ξ i n,n = ξ i n (1.4)
of each current individual. The occupation measures of the corresponding N -genealogical tree model converge as N → ∞ to the conditional distribution (1.1). In a sense to be given, for any function f on the set Z n+1 , we have the convergence, as N → ∞,

lim N →∞ 1 N N i=1 f (ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ) 1 τ N >n = E (f (X 0 , . . . , X n ) | ∀0 ≤ p < n, X p ∈ A ) (1.5)
This convergence result can be refined in various directions. For instance, we can prove that the ancestral lines are "almost" independent with a common distribution given by the limiting conditional distribution. This property is often called the propagation of chaos property in applied probability. It refers to the fact that the initial population consists of independent and identically distributed random variables and that this property "propagates" approximately despite the introduction of interactions. Many other results can be derived including the fluctuations and the exponential concentration of the occupation measures of the genealogical tree around the limiting conditional distribution.

Besides the fact that the particle model approximate the (rare event) probabilities (1.3) and the conditional distributions (1.5) in path spaces, it also contains some information about the top of the spectrum of the matrix Q defined below

∀(x, y) ∈ {-a, -a + 1, ..., a -1, a} Q(x, y) := G(x) M (x, y) with G(x) := 1 A (x) and M (x, y) = 1 2 1 x-1 (y) + 1 2 1 x+1 (y)
Indeed, if we let λ the top eigenvalue of Q and we denote by h the corresponding eigenvector s.t.

x h(x) = 1, then we have lim

N →∞ lim n→∞ 1 n 0≤p≤n log η N p (1 A ) = log λ INRIA
In addition, the value h(x) coincides with the long time proportion of visits of the algorithm to the state x. In other words, h(x) can be interpreted as the limiting distribution of the individuals within the set A; that is lim

N,n→∞ 1 n 0≤p≤n 1 N 1≤i≤N 1 x ( ξ i n ) 1 τ N >n = h(x) = lim N,n→∞ 1 N 1≤i≤N 1 x ( ξ i n ) 1 τ N >n
The particle approximation model discussed above is far from unique. Many other interacting sampling strategies can be introduced by a simple change of probability measure. For instance, we can replace the mutation or the free evolution of the individuals in the previous algorithm by local moves restricted to the desired set A. These mutation type transitions ξ n-1

ξ n can also be seen as transitions of a simple random walk on Z reflected at the boundaries of the set A. By construction all the individuals ξ i n at any time horizon n and for any index i = 1, . . . , N are in the desired set A. The corresponding selection transition ξ n ξ n is now defined as follows: Each individual ξ i n = x on the boundary x ∈ ∂A = {-a + 1, (a -1)} of the set A has a probability G(x) := 1/2 to stay in A, while the other individuals ξ i n (which are in the set A) have a probability G(x) = 1 to stay in A. The population ξ n is now defined as follows. For every index i, with a probability G(ξ i n ), we set ξ i n = ξ i n , otherwise we replace ξ i n be a new individual ξ i n = ξ j n randomly chosen in the whole population with a probability proportional to G(ξ j n ). If we now write

η N n (G) = 1 N 1≤i≤N G(ξ i n
), all the previous particle approximation results (corresponding to G(x) = 1 A (x)) we have presented remain valid for this new particle algorithm.

Particle absorption models

The sampling techniques described in section 1.1 are far from being restricted to random walks models confined to a set. These strategies apply to a variety of application areas including computational physics, nonlinear filtering, biology, as well as in rare event analysis. From the pure mathematical point of view, they correspond to interacting particle approximation models of Feynman-Kac measures in path spaces.

To introduce these models, we recall that the conditional distributions discussed in (1.1) can be represented in terms of the distributions of the free path evolutions

P n (x 0 , . . . , x n ) = Proba ((X 0 , . . . , X n ) = (x 0 , . . . , x n )) = 1 0 (x 0 ) M 1 (x 0 , x 1 ) . . . M n (x n-1 , x n ) (1.6)
of the simple random walk starting at the origin with elementary transitions given by the matrix M n := (M n (x, y)) x,y∈Z with entries given by

M n (x, y) := 1 2 1 x-1 (y) + 1 2 1 x+1 (y)
More formally, if we set

Q n (x 0 , . . . , x n ) := Proba ((X 0 , . . . , X n ) = (x 0 , . . . , x n ) | ∀0 ≤ p < n, X p ∈ A )
then we have

Q n (x 0 , . . . , x n ) = 1 Z n    0≤p<n G p (x p )    P n (x 0 , . . . , x n ) (1.7)
with the indicator potential functions G n (x) = 1 A (x) and P n (x 0 , . . . , x n ) being the distribution of a free path of length n of the symmetric random walk. In (1.7), Z n is the normalizing constant given by

Z n = P (∀0 ≤ p < n, X p ∈ A) = E   0≤p<n G p (X p )   RR n°6991
These path integration type models are called Feynman-Kac measures in reference to Feynman path integral formulation of quantum mechanics where the classical notion of a single deterministic trajectory for a system is replaced by a sum over all possible trajectories weighted by the contributions of all the histories in configuration space. The Feynman-Kac measures presented in (1.7) can be regarded as the distribution of the paths of a Markov particle evolving using the Markov transitions M n in an environment with absorbing obstacles related to potential functions G n , and starting with some initial distribution Law(X 0 ) = η 0 with η 0 (x 0 ) = 1 0 (x 0 ) in (1.6). To be more precise, we consider an auxiliary coffin or cemetery state c and we set E c = E ∪ {c}. We define an E c -valued Markov chain X c n with two separate killing/exploration transitions:

X c n killing ---------→ X c n exploration ---------→ X c n+1 (1.8)
This killing/exploration mechanism are defined as follows:

• Killing: If X c n = c, we set X c n = c. Otherwise the particle X c n is still alive. In this case, with a probability G n (X c n ), it remains in the same site so that X c n = X c n , and with a probability 1 -G n (X c n ) it is killed and we set X c n = c.

• Exploration: Once a particle has been killed, it cannot being brought back to life so if X c n = c then we set X c p = X p = c for any p > n. Otherwise, the particle X c n ∈ E evolves to a new location X c n+1 = x in E randomly chosen according to the distribution M n+1 (X c n , x).

In this physical interpretation, the measure Q n represent the conditional distributions of the paths of a non absorbed Markov particle. To see this claim, we denote by T the time at which the particle has been killed

T = inf {n ≥ 0 ; X c n = c} By construction, we have Proba(T > n -1) = Proba( X c 0 ∈ E, . . . , X c n-1 ∈ E) = E (n+1) η 0 (dx 0 ) G 0 (x 0 ) M 1 (x 0 , dx 1 ) . . . M n-1 (x n-2 , dx n-1 )G n-1 (x n-1 ) = E n-1 p=0 G p (X p )
This also shows that the normalizing constants Z n represent respectively the probability for the particle to be alive at time n -1. In other words, we have that

Z n = Proba(T > n -1)
Similar arguments yield that is the distribution of a particle conditional upon being alive at time n -1

Q n (x 0 , . . . , x n ) = Proba ((X c 0 , . . . , X c n ) = (x 0 , . . . , x n ) | T > n -1 )
The particle sampling technique of any distribution Q n associated with some Markov transition M n and some sequence of [0, 1]-valued potential function G n on some (countable) state space E is defined as before in terms of a genetic type algorithm with M n -mutations and G n -selection type transitions. More precisely, at every time step n, we sample the mutation-selection transitions as follows: During the mutation step, every individual performs a local random move according to the Markov transition M n . During the selection step, every individual evaluates its potential value G n (ξ i n ), with 1 ≤ i ≤ N . For every index i, with a probability G n (ξ i n ), we set ξ i n = ξ i n , otherwise

INRIA

we replace ξ i n be a fresh new individual ξ i n = ξ j n random chosen in the whole population with a probability proportional to G n (ξ j n ). More precisely, for any time horizon n and any function f on the set E n+1 , we have lim

N →∞ 1 N N i=1 f (ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ) 1 τ N >n = x0,...,xn f (x 0 , . . . , x n ) Q n (x 0 , . . . , x n ) (1.9)
Furthermore, the N -particle approximation of the normalizing constants Z n are given by

Z N n := 0≤p<n η N p (G p ) with ∀n ∈ N η N n (G n ) := 1 N 1≤i≤N G n (ξ i n ) (1.10)
For time homogeneous models (G n , M n ) = (G, M ) associated with a reversible matrix M w.r.t. to some measure λ on E, i.e. λ (x) M (x, y) = λ (y) M (y, x), the corresponding particle model also contains information about the top of the spectrum of the matrix Q defined through

∀(x, y) ∈ E Q(x, y) := G(x) M (x, y)
More precisely, if we let λ the top eigenvalue of Q in L 2 (λ) and we denote by h the corresponding eigenvector s.t.

x λ(x)h(x) = 1, then we have lim

N →∞ lim n→∞ 1 n 0≤p≤n log η N p (G) = log λ as well as lim N,n→∞ 1 n 0≤p≤n 1 N 1≤i≤N 1 x ( ξ i n ) 1 τ N >n = λ(x)h(x) = lim N,n→∞ 1 N 1≤i≤N 1 x ( ξ i n ) 1 τ N >n
For further details on this subject, we refer the reader to [START_REF] Moral | Feynman-Kac Formulae: Genealogical and Interacting Particle Systems with Applications[END_REF][START_REF] Del Moral | Particle motions in absorbing medium with hard and soft obstacles[END_REF][START_REF] Del Moral | Particle approximations of Lyapunov exponents connected to Schrödinger operators and Feynman-Kac semigroups[END_REF] and references therein.

Nonlinear filtering problems

We discuss here the application of these particle model to filtering problems. Suppose that at every time step the state of the Markov chain X n is partially observed according to the following schematic picture

X 0 -→ X 1 -→ X 2 -→ . . . ↓ ↓ ↓ Y 0 Y 1 Y 2 . . .
with some random variables Y n whose values only depend on the current state of the chain

Proba (Y n = y n | X n = x n ) = G(x n , y n )
We consider the following pair of events

A n (x) := {(X 0 , . . . , X n ) = (x 0 , . . . , x n )} and B n-1 (y) := {(Y 0 , . . . , Y n-1 ) = (y 0 , . . . , y n-1 )}
The filtering problem consists of computing the conditional distributions of the state variables A n (x) given the observations B n (y). By construction, given A n (x), the random variables are independent and identically distributed with a distribution given by

Proba (B n-1 (y) |A n (x) ) = 0≤p<n G(x p , y p ) RR n°6991
By direct application of Bayes' rule we have the following formula

Proba (A n (x) ∩ B n-1 (y)) = Proba (B n-1 (y) |A n (x) ) × Proba (A n (x)) =    0≤p<n G(x p , y p )    P n (x 0 , . . . , x n ) (1.11)
from which we conclude that

Proba (A n (x) | B n-1 (y) ) = 1 Z n (y)    0≤p<n G(x p , y p )    P n (x 0 , . . . , x n )
with the normalizing constants

Z n (y) := Proba(B n-1 (y)) = x0,...,xn    0≤p<n G(x p , y p )    P n (x 0 , . . . , x n )
These Feynman-Kac formulae express the conditional distributions of the path sequence (X 0 , . . . , X n ) as the distribution of the signal paths (X 0 , . . . , X n ) = (x 0 , . . . , x n ) weighted by the product of the likelihood functions G(x p , y p ) from the origin p = 0 up to time p = n. If we fix the observation sequence Y n = y n and set G n (x n ) := G(x n , y n ), these measures have exactly the same form as the one presented in (1.7). The corresponding particle approximations are often referred as particle filters in signal processing and statistics (see for instance [START_REF]Sequential Monte Carlo Methods in Practice[END_REF], and references therein). These particle algorithms can also be used to approximate the log-likelihood functions using (1.10); that is the log-likelihood L n (y) := log Z n (y) is approximated using

L N n (y) := log Z N n (y) = 0≤p<n log η N p (G p ).

Stochastic optimization algorithms

Suppose we want to compute the global minima of a given non negative cost function V on some finite state space E equipped with the counting measure λ(x) := 1

Card(E) . From the probabilistic point of view, this problem amounts of sampling random states according to the Boltzmann-Gibbs distributions associated with a large inverse temperature parameter β and given below

µ β (x) := 1 Z β e -βV (x) λ(x) with Z β := x e -βV (x) λ(x)
There is no loss of generality to assume that inf x V (x) = 0 and for any state

x ∈ V 0 := V -1 ({0}), V (x) ≥ δ for some δ > 0. It follows that we have Card(V 0 ) ≤ Z β ≤ Card(V 0 ) + Card(V c 0 ) e -βδ → β↑∞ Card(V 0 )
and therefore lim

β→∞ µ β (x) = µ ∞ (x) := 1 V0 (x)/Card(V 0 )
This simple observation shows that sampling according to µ β is roughly equivalent to that of sampling randomly an unknown state variable with minimal cost. For very large state spaces, it is typically impossible to sample from µ β directly.

INRIA

The celebrated simulated annealing algorithm to sample from µ ∞ consists of sampling approximately from a sequence of distributions µ βn where β n is a non-decreasing sequence going to ∞. The rationale is that it is "easier" to sample from µ β when β is small; if β = 0 then µ 0 is the uniform counting measure on E from which it is trivial to sample. For β n > 0, we sample approximately from each intermediate distribution µ βn using Markov chain Monte Carlo (MCMC) sampling techniques; that is we select a transition matrix M βn = (M βn (x, y)) x,y∈E with left eigenvector µ βn associated with the eigenvalue 1, that is

x µ βn (x)M βn (x, y) = µ βn (y)
The probabilistic interpretation of the above equation is as follows: pick randomly a state x with distribution µ βn (x) and take a random transition x y from the distribution M βn (x, y), then the probability of being at state y is again µ βn (y). The literature on MCMC methods discusses numerous choices of transitions M βn satisfying this property. The most famous is the Metropolis-Hastings transition associated to a symmetric transition matrix K(x, y) = K(y, x) and defined by

M βn (x, y) = K(x, y) min 1, e -βn(V (y)-V (x)) + 1 -z K(x, z) min 1, e -βn(V (z)-V (x)) 1 x (y)
Using the fundamental ergodic theorem for regular Markov chains, starting from any initial state x 0 , the n-th step of a run of the Markov chain with transitions M βn has a probability very close to µ βn (y) of being at the site y, for a large n. Practically, we select β 1 and we run the chain starting at X 0 = x 0 for a large enough number of runs n 1 such that the law of the state X n1 is close to µ β1

X 0 = x 0 M β 1 -→ X 1 M β 1 -→ . . . M β 1 -→ X n1 with n 1 large enough s.t. Law(X n1 ) ≃ µ β1
Notice that the choice of n 1 depends on β 1 : the larger β 1 is, the "peakier" µ β1 is and the larger n 1 is. When the chain is stabilized, we choose a β 2 > β 1 and we run the chain starting at X n1 for a new large enough number of time steps n 2 such that the law of the state X n1+n2 is close to µ β2

X n1 M β 2 -→ X n1+1 M β 2 -→ . . . M β 2 -→ X n1+n2 with n 2 large enough s.t. Law(X n1+n2 ) ≃ µ β2
The theoretical "optimal" inverse temperature parameter ensuring convergence in some sense of the Markov chain to µ ∞ is logarithmic. This amounts to say that we change by one unit the parameter β on every time interval with exponential length. This is unrealistic from a practical point of view.

We present now an alternative particle strategy for sampling random states according to the sequence of measures µ βn associated with a given non decreasing sequence of inverse temperature parameters β n . We suppose that β 0 = 0 so that µ β0 coincides with the uniform counting measure on the set E. We start with N independent individuals ξ 0 := ξ i 0 1≤i≤N randomly chosen in E according to µ β0 . We perform a selection transition ξ 0 ξ 0 := ξ i 0 1≤i≤N

using the potential functions G 0 defined by

G 0 (x) = exp (-(β 1 -β 0 )V (x))
In other words, every individual evaluates its potential value G 0 (ξ i 0 ). For every index i, with a probability G 0 (ξ i 0 ), we set ξ i 0 = ξ i 0 , otherwise we replace ξ i 0 be a new individual ξ i 0 = ξ j 0 randomly chosen in the whole population with a probability proportional to G 0 (ξ j 0 ). During the mutation step ξ 0 ξ 1 := ξ i 1 1≤i≤N , every selected individual ξ i 0 performs a local random move ξ i 0 ξ i 1 (independently of one another) according to the Markov transition M β1 . Then, we perform another

selection transition ξ 1 ξ 1 := ξ i 1 1≤i≤N
using the fitness functions G 1 defined below:

G 1 (x) = exp (-(β 2 -β 1 )V (x)) RR n°6991
After this selection stage we mutate each selected individual using the Markov transition M β2 , and so on. Iterating these transitions, we define a simple genetic model with mutations transitions M βn and selection fitness functions G n :

ξ n := ξ i n 1≤i≤N ∈ E N selection --------→ ξ n := ξ i n 1≤i≤N ∈ E N mutation -------→ ξ n+1 ∈ E N (1.12)
This algorithm was first proposed in [START_REF] Del Moral | Sequential Monte Carlo samplers[END_REF]. A variety of convergence results can be established for this algorithm. For instance, for any function f on E and any time horizon, we have lim

N →∞ 1 N 1≤i≤N f (ξ i n ) = x µ βn (x) f (x)
In addition, if we set

η N n (G n ) := 1 N 1≤i≤N G n (ξ i n )
, the unbiased N -particle approximation Z N βn of the normalizing constants Z βn is given by

Z N βn := 0≤p<n η N p (G p ) -→ N →∞ Z βn 1.

Combinatorial counting and sampling

Suppose we want to compute the cardinality of a given subset A of some finite state space E equipped with the counting measure λ(x) := 1

Card(E) . Once again, from a probabilistic point of view, this problem is equivalent to compute the normalizing constant of the following Boltzmann-Gibbs distribution

µ A (x) := 1 Z A 1 A (x) λ(x) with Z A := x 1 A (x) λ(x)
To sample from µ A and compute Z A , the idea consists of selecting a judicious sequence of decreasing subsets A n in such a way that it is easy to sample states in A n starting from the set A n-1 . We suppose that A 0 = E so that µ A0 coincide with the uniform counting measure on the set E. The algorithm is thus very similar to the one described above for optimization.

For any set A n , we introduce an MCMC transition matrix M An = (M An (x, y)) x,y∈E with left eigenvector µ An associated with the eigenvalue 1, that is

x µ An (x)M An (x, y) = µ An (y)
A simple Metropolis-Hasting type transition associated with a symmetric transition matrix K(x, y) = K(y, x) is given by

M An (x, y) = K(x, y) 1 An (y) + 1 - z K(x, z) 1 An (z) 1 x (y)
The N -particle stochastic algorithm is defined as follows. We start with N independent random individuals ξ 0 := ξ i 0 1≤i≤N randomly chosen in E with µ A0 . We perform a selection transition ξ 0 ξ 0 := ξ i 0 1≤i≤N

using the fitness functions G 0 = 1 A1 . In other words, every individual in the set A 1 stays in the same place ξ i 0 = ξ i 0 , otherwise we replace ξ i 0 be a fresh new individual ξ i 0 = ξ j 0 random chosen among the individuals ξ j 0 ∈ A 1 . When no individuals ξ j 0 are in the set A 1 , the algorithm stops and we set τ N = 0. Assuming that τ N > 0, during the mutation step ξ 0 ξ 1 := ξ i 1 1≤i≤N , every selected individual ξ i 0 performs a local random move ξ i using the fitness functions G 1 = 1 A2 . When no individuals ξ j 1 are in the set A 2 , the algorithm stops and we set τ N = 1. After this selection stage we mutate each selected individual using the Markov transition M A2 , and so on. For any function f on E and any time horizon n, we have lim

N →∞ 1 N 1≤i≤N f (ξ i n )1 τ N >n = x µ An (x) f (x)
In addition, if we set

η N n (G n ) := 1 N 1≤i≤N G n (ξ i n )
, the proportion of individuals in A n+1 after the n-th mutation, the unbiased N -particle approximation Z N An of the normalizing constants Z An is given by

Z N An := 0≤p<n η N p (G p ) -→ N →∞ Z An = Card(A n )/Card(E)

Genetic search algorithms

We consider an energy function or a cost criteria V : x ∈ E → (x) on some finite state space E where we assume inf x V (x) = 0 without loss of generality. The objective is to find the global minima points

x ⋆ ∈ E s.t. V (x ⋆ ) = inf x∈E V (x)
. Let V ⋆ denote the set of these points. We describe in Section 1.4 an interacting particle algorithm to solve this problem which relies on interacting simulated annealing type chains. We present here the more standard genetic algorithm with mutation and proportional selection.

To construct this algorithm, we introduce a collection of Markov transitions M n (x, y) from E into itself. This collection of transition matrices represents the probability M n (x, y) that a individual at site x evolves to a new state x during the n-th mutation transition.

The genetic algorithm with N individuals is defined as follows. We start with N independent random individuals ξ 0 := ξ i 0 1≤i≤N randomly chosen in E with some distribution, say η 0 . We perform a proportional type selection transition ξ 0 ξ 0 := ξ i 0 1≤i≤N

using the potential functions

G 0 ξ i 0 = exp (-β 0 V ξ i 0 )
, where β 0 ≥ 0 is an inverse temperature parameter. In other words, with probability G 0 (ξ i 0 ) every individual stays in the same place ξ i 0 = ξ i 0 ; otherwise, we replace ξ i 0 by a new individual ξ i 0 = ξ j 0 randomly chosen among the individuals ξ j 0 with a probability proportional to its weight G 0 (ξ i 0 ). Formally, we set

ξ i 0 = ǫ i 0 ξ i 0 + 1 -ǫ i 0 ξ i 0
where ǫ i 0 stands for a sequence of independent {0, 1}-valued Bernoulli random variables with distributions

G 0 (ξ i 0 ) := Proba ǫ i 0 = 1 | ξ 0 = 1 -Proba ǫ i 0 = 0 | ξ 0 and ξ 0 := ξ i 0 1≤i≤N
are independent, identically distributed and ξ j 0 , 1 ≤ j ≤ N -valued random variables with common distributions given for any index 1 ≤ i ≤ N by

∀1 ≤ j ≤ N Proba ξ i 0 = ξ j 0 | ξ 0 = G 0 (ξ j 0 )/ 1≤j≤N G 0 (ξ j 0 )
During the mutation step ξ 0 ξ 1 := ξ i 1 1≤i≤N , every selected individual ξ i 0 performs a local random move ξ i 0 ξ i 1 (independently of one another) according to the Markov transition M 1 . Then, we perform another proportional type selection transition ξ 1 ξ 1 := ξ i 1 1≤i≤N

using the potential

functions G 1 ξ i 1 = exp (-β 1 V ξ i 1 )
, where β 1 ≥ 0 is another inverse temperature parameter, and so on. We define in this way a sequence of genetic type populations ξ n , ξ n , as in (1.2) and the corresponding genealogical tree model (1.4) associated with the ancestral lines ξ i p,n 0≤p≤n of every i-th individuals after the n-th mutation. In the same way, running back in time we have the whole ancestral line

ξ i 0,n ←-ξ i 1,n ←-. . . ←-ξ i n-1,n ←-ξ i n,n = ξ i n (1.13)
of every i-th individual after the n-th selection.

For any function f on E (n+1) and any time horizon n, we can prove that lim

N →∞ 1 N N i=1 f ( ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ) = E f n (X 0 , . . . , X n ) exp -0≤p≤n β p V (X p ) E exp -0≤p≤n β p V (X p )
In other words, the proportion of paths ( ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ) taking some value (x 0 , . . . , x n ) is given by lim

N →∞ 1 N N i=1 1 (x0,...,xn) ( ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ) = 1 Z n+1 exp   - 0≤p≤n β p V (x p )   P n (x 0 , . . . , x n )
with the probability of a free evolution path involving only mutation transitions

P n (x 0 , . . . , x n ) = η 0 (x 0 )M 1 (x 0 , x 1 ) . . . M n (x n-1 , x n )
and Z n+1 is a normalizing constant. Suppose that every free evolution path has the same chance to be sampled, in the sense that

P n (x 0 , . . . , x n ) = P n (y 0 , . . . , y n )
for any paths admissible pair of paths (x 0 , . . . , x n ) and (y 0 , . . . , y n ). This condition is satisfied if η 0 is the uniform counting measure on E and the mutation transitions M n (x, y) correspond to local random choices of the same number of neighbors, starting from any state x. In this case, for any admissible path (x 0 , . . . , x n ) we have that lim

N →∞ 1 N N i=1 1 (x0,...,xn) ( ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ) = 1 Z ′ n exp   - 0≤p≤n β p V (x p )  
for some normalizing constant Z ′ n . When the inverse temperature parameter β p increases the r.h.s. probability mass quantity only charges admissible paths (x 0 , . . . , x n ) that minimize the path potential function

V n (x 0 , . . . , x n ) = inf (y0,...,yn) 0≤p≤n V (y p )
In other words at low temperature, the ancestral lines of the simple genetic model described above converge to the uniform measure on all the paths (x 0 , . . . , x n ) of length n that minimize the energy function V n . For time homogenous mutation transitions associated with stochastic matrices M n (x, y) = M (x, y) satisfying the following condition for some integer m ≥ 1 and any pair (x, y)

∈ E 2 M (x, x) > 0 and M m (x, y) ≥ ǫM m (x, z)
we also have the convergence result

lim n→∞ lim N →∞ 1 N N i=1 1 V ⋆ ( ξ i n ) = 1 INRIA
as soon as β n = C log (n + 1) for some constant C that depends on m and on the oscillations of the function V . This convergence result is also true for β n = C (n + 1) α , with any α ∈]0, 1, as soon as the above condition is met for m = 1. Further details on these concentration properties can be found in [START_REF] Del Moral | Annealed Feynman-Kac models[END_REF]. Related convergence results for fixed population sizes can be found in [START_REF] Cerf | Asymptotic convergence of genetic algorithms[END_REF]. To give a flavor of these results, let us suppose that the mutation transitions M n (x, y) also depend on the inverse temperature parameter and

M n (x, y) → n→∞ 1 x (y) as β n ↑ ∞
Intuitively speaking, the genetic mutations become rare transitions at low temperature. In this situation, we can prove that there exists a "critical population size" N ⋆ that depends on the energy function as well as on the free evolution model such that

∀N ≥ N ⋆ lim n→∞ Proba ∀1 ≤ i ≤ N ξ i n ∈ V ⋆ = 1 1.

Directed polymers and self avoiding walks

In biology and chemistry, flexible polymers describe the kinetic structure of macromolecules in a given chemical solvent. The polymer chain at time n is regarded as a sequence of random variables

X n = (X ′ 0 , . . . , X ′ n ) ∈ E n = E × . . . × E (n+1) times
taking values in some finite or countable state space E. Of course, in concrete situations the state space is not a finite set but some Euclidian state space E = R d , for some d ≥ 1. This countable state space condition is only made to avoid unnecessary technicalities in the presentation of these stochastic models in path spaces. The following discussion is mainly taken from section 12.5 in the book [START_REF] Moral | Feynman-Kac Formulae: Genealogical and Interacting Particle Systems with Applications[END_REF].

The elementary states X ′ p represent the different directed monomers of the macromolecules X n . The length parameter n represents the degree of polymerization. The monomers are connected by chemical bonds and interact with one another as well as with the chemicals in the solvent. The energy of a polymerization sequence

X 0 = X ′ 0 -→ X 1 = (X ′ 0 , X ′ 1 ) -→ X 2 = (X ′ 0 , X ′ 1 , X ′ 2 ) -→ . . . -→ X n = (X ′ 0 , . . . , X ′ n ) is defined in terms of a Boltzmann potential exp -β n p=0 V p (X ′ 0 , . . . , X ′ p ) (1.14)
The parameter β ∈ R + represents the inverse temperature of the solvent, and each potential function

V n : (x 0 , . . . , x n ) ∈ E n → V n (x 0 , . . . , x n ) ∈ R +
reflects the local intermolecular energy between the monomer X ′ n = x n in the polymer chain X n-1 = (x 0 , . . . , x n-1 ) during the nth polymerization

X n-1 = (x 0 , . . . , x n-1 ) -→ X n = ((x 0 , . . . , x n-1 ), x n )
The potential functions V n depend on the nature of the solvent and the physico-chemical structure of the polymer. At low temperature, β ↑ ∞, the interaction between monomers may be strongly repulsive at short distances and attractive or repulsive at larger ones. For instance, the monomers may tend to avoid being closed on each other. These excluded volume effects and repulsive interactions can be modeled by choosing a potential function satisfying the following condition:

V n (x 0 , . . . , x n ) = 0 ⇐⇒ x n ∈ {x 0 , . . . , x n-1 } (1.15)
In this situation, every self-interaction is penalized by a factor e -βVn so that the energy of an elementary polymerization is minimal if and only if the new monomer differs from the previous ones. In this context, the inverse temperature parameter β is sometimes called the strength of repulsion.

In the opposite case, at high temperature, β → 0, the interaction forces disappear. In this situation, it is commonly assumed that X ′ n is an E-valued Markov chain with elementary transitions M ′ n and initial distribution η ′ 0 . By the definition of the chain X n = (X ′ 0 , . . . , X ′ n ), this Markovian hypothesis implies that the Markov transitions M n have the form M n ((x 0 , . . . , x n-1 ), (y 0 , . . . , y n-1 , y n )) = 1 (x0,...,xn-1) (y 0 , . . . , y n-1 ) M ′ n (y n-1 , dy n ) for any paths (x 0 , . . . , x n-1 ) ∈ E n-1 and (y 0 , . . . , y n ) ∈ E n and for some Markov transition M ′ n from E into E.

The Boltzmann-Gibbs measures associated with these directed polymer models are the measures in the product spaces E n defined below

Q n (x 0 , . . . , x n ) := 1 Z n exp -β n-1 p=0 V p (x 0 , . . . , x p ) P n (x 0 , . . . , x n )
with the distribution of the free paths of the monomer chain

P n (x 0 , . . . , x n ) = η ′ 0 (x 0 )M ′ 1 (x 0 , x 1 ) . . . M ′ n (x n-1 , x n )
These measures can alternatively be rewritten in the following form

Q n (x 0 , . . . , x n ) := 1 Z n n-1 p=0 G p (x 0 , . . . , x p ) P n (x 0 , . . . , x n )
with the exponential weight functions

∀n ≥ 0 G n (x 0 , . . . , x n ) := exp (-βV n (x 0 , . . . , x n ))
To illustrate these models, let us suppose that X ′ n is a simple random walk on the d-dimensional lattice Z d starting at the origin and set

G n (x 0 , . . . , x n ) := 1 Z d -{x0,...,xn-1} (x n )
In this situation, we have that

Q n = Law X ′ 0 , . . . , X ′ n ∀0 ≤ p < q < n X ′ p = X ′ q and Z n = Proba ∀0 ≤ p < q < n X ′ p = X ′ q = Card (S n-1
) /(2d) n-1 where S n-1 is the set of self-avoiding walks of length n -1 starting at the origin.

The N particle approximation of these quantities is nothing but a simple genetic type evolution model with N path-valued particles

ξ i n := ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ξ i n := ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ∈ E n
During the selection stage, with a probability G n (ξ i n ) every path-valued individual stays in the same place ξ i n = ξ i n ; otherwise, we replace ξ i n be a new individual ξ i n = ξ j n randomly chosen among the individuals ξ j 0 with a probability proportional to its weight G n (ξ i n ). This mechanism is intended to favor minimal energy polymerizations. For instance, in the case of repulsive interaction (1.15), INRIA a given polymer with degree n, say (ξ i 0,n , . . . , ξ i n,n ), has more chance of being selected if the last monomer ξ i n,n added during the nth sampled polymerization differs from the previous ones; that is, if ξ i n,n ∈ {ξ i 0,n , . . . , ξ i n-1,n }. Each selected polymer is a path-valued individual ξ i n . During the mutation transition, ξ i n evolves randomly according to the transition M n+1 of the path-valued Markov chain X n .

ξ i n+1 = ((ξ i 0,n+1 , . . . , ξ i n,n+1 ), ξ i n+1,n+1 ) = (( ξ i 0,n , . . . . . . , ξ i n,n ), ξ i n+1,n+1 ) ∈ E n+1 = (E n × E) (1.16)
where ξ i n+1,n+1 is a random variable with distribution M ′ n ( ξ i n , .). Various asymptotic estimates can be derived. For instance, for any function f on E n = E (n+1) and any time horizon n, we have lim

N →∞ 1 N N i=1 f (ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ) = x0,...,xn Q n (x 0 , . . . , x n )f (x 0 , . . . , x n )
and unbiased estimates of the normalizing constants

Z N n := 0≤p<n η N p (G p ) -→ N →∞ Z n
where η N p (G p ) is the empirical mean value of the potential functions after the mutation stage

∀n ≥ 0 η N n (G n ) := 1 N 1≤i≤N G n ξ i 0,n , ξ i 1,n , . . . , ξ i n,n
2 A brief treatise on Feynman-Kac modeling

Stochastic matrices and Markov semigroups

We let E be a finite set. A Markov chain X = (X n ) n≥0 is defined by a matrix M with positive entries M (x 1 , x 2 ) ≥ 0 s.t.

x2 M (x 1 , x 2 ) = 1 for each state x 1 . Such matrices are called transition matrices or stochastic matrices or in the literature on Markov chains. For instance, for finite state spaces E = {1, . . . , d}, with cardinality d ≥ 1, we can identify a Markov transition M (x, y) from E into itself with a d × d matrix

   M (1, 1) • • • M (1, d) . . . . . . . . . M (d, 1) • • • M 1 (d, d)    .
Starting from x 1 we sample the chain transition by choosing randomly a state x 2 with distribution M (x 1 , x 2 ). The outcomes X 0 = x 0 , X 1 = x 1 , X 2 = x 2 , ... are referred as a realization of the chain X starting from X 0 = x 0 . By construction, we have

P (X n = x n | X n-1 = x n-1 ) = M (x n-1 , x n )
and more generally

P (X 1 = x 1 , X 2 = x 2 , . . . , X n = x n | X 0 = x 0 ) = M (x 0 , x 1 )M (x 1 , x 2 ) . . . M (x n-1 , x n )
These simple observations already show that the law of the random state X n at time n starting from X 0 = x 0 is given by the n-th power of the matrix M ; that is we have

P (X n = x n | X 0 = x 0 ) = x1,...,xn-1 M (x 0 , x 1 )M (x 1 , x 2 ) . . . M (x n-1 , x n ) = M n (x 0 , x n ) RR n°6991
Identifying functions f : x ∈ E → f (x) ∈ R with column type vector f = (f (x)) x∈E , the conditional expectations of he random variable f (X n ) given X 0 = x 0 are also given by a matrix operation

E (f (X n ) |X 0 = x 0 ) = xn M n (x 0 , x n )f (x n ) = M n (f )(x 0 )
where M n (f ) stands for the function given by the column vector M n (f ) = (M n (f )(x 0 )) x0∈E . For instance, for finite state spaces E = {1, . . . , d}, functions

f (x) = d i=1 f (i) 1 i (x)
are identify with the column vector

f = (f (1), . . . , f (d)) ⋆ :=    f (1) . . . f (d)    . (2.1)
In this situation, the function M (f ) defined by

M (f )(x) := E (f (X n ) | X n-1 = x)
is given by the finite sum

∀1 ≤ i ≤ d M (f )(i) = d j=1 M (i, j) f (j).
Using the vector notation we find that the function M (f ) can alternatively be defined by the following matrix operation:

M (f ) :=    M (f )(1) . . . M (f )(d)    =    M (1, 1) • • • M (1, d) . . . . . . . . . M (d, 1) • • • M (d, d)       f (1) . . . f (d)    = M f
In the same way, identifying probability measures η 0 on E with row type vectors η 0 = (η 0 (x 0 )) x0∈E s.t. η 0 (x 0 ) ≥ 0 and x0 η 0 (x 0 ) = 1, the distribution η n of the state X n of the chain starting randomly from X 0 with distribution η 0 is given by

η n (x n ) = P (X n = x n ) = x0 η 0 (x 0 )M n (x 0 , x n ) = η 0 M n (x n )
where η 0 M n stands for the probability on E given by the row vector

η 0 M n = (η 0 M n (x n )) xn∈E . Furthermore, we have η n (f ) = xn η n (x n )f (x n ) = E (f (X n ))
In matrix notation, we simply have that

η n = η 0 M n and η n (f ) = (η 0 M n )(f ) = η 0 ([M n (f )])
For instance, for instance, for finite state spaces E = {1, . . . , d}, the sum η n (f ) can be rewritten as

η n (f ) = η n f = (η n (1), . . . , η n (d))    f (1) . . . f (d)    = d i=1 f (i) η n (i).
INRIA and we have

(η 0 M ) := [(η 0 M )(1), . . . , (η 0 M )(d)] = [η 0 (1), . . . , η 0 (d)]    M (1, 1) • • • M (1, d) . . . . . . . . . M (d, 1) • • • M (d, d)    = η 0 M.
In much the same way, we find that

η n = η n-1 M = (η n-2 M )M = η n-2 M 2 = . . . = η 0 M n
By construction, the law of the random states X p and X n , with 0 ≤ p ≤ n are given by the following matrix operation

η n = η p M (n-p)
The collection of matrices M n , with n ≥ 0, is called the Markov semigroup associated with the transition probabilities M . For n = 0, we use the convention M 0 = Id the identity matrix with zero off-diagonal elements and unit diagonal elements. For instance, if E = {1, 2} then the Markov transition M is given by the 2 × 2-matrix

M = M (1, 1) M (1, 2) M (2, 1) M (2, 2)
We further assume that the entries M (i, j) ∈ [0, 1] are chosen s.t. c := M (1, 2) + M (2, 1) > 0. Using an elementary induction on the time parameter n, we prove

M n = 1 c M (2, 1) M (1, 2) M (2, 1) M (1, 2) + (1 -c) n c M (1, 2) -M (1, 2) -M (2, 1) M (2, 1) .
Nevertheless, except in some very particular cases, the matrix semigroup M n has no explicit expression. Time nonhomogenous models can be studied in the same line of arguments. For non homogeneous Markov transitions M n , we have

P (X n = x n | X n-1 = x n-1 ) = M n (x n-1 , x n )
The law of the chain X n is now given by the following matrix composition

η n = η 0 M 1 M 2 . . . M n
The semigroup M p,n associated with this sequence of measures is now given by the following formulae

η n = η p M p,n with the matrix M p,n := M p+1 M p+2 . . . M n (2.2)

The Bayes-Boltzmann-Gibbs transformation

We consider a non negative potential function G and a probability measure η on a finite or countable state space E. We recall that η(G) is the quantity defined by the sum

η(G) := x η(x)G(x)
For indicator potential functions G = I A , sometimes we slightly abuse notation and we set

η(1 A ) = x η(x)1 A (x) = x∈A η(x) := η(A) Definition 2.1
We let Ψ G (η) be the probability measure on E defined by

Ψ G (η)(x) := 1 η(G) G(x) η(x)
Note that this probability measure Ψ G (η) is well-defined if and only if η(G) > 0. The transformation Ψ G is called the Bayes-Boltzmann-Gibbs transformation ( abbreviate BBG-transformation) associated to the potential function G.

When G = 1 A , we notice that Ψ G (η) is the restriction of the measure η to the set A;

Ψ 1 A (η)(x) := η A (x) := 1 η(A) 1 A (x) η(x).
In this case, we have

A ⊂ B =⇒ Ψ 1 A (Ψ 1 B (η)) = Ψ 1 B (Ψ 1 A (η)) = η A∩B = Ψ 1 A∩B (η)
In nonlinear filtering problems discussed in section 1.3, the measure η represents the conditional law of a Markov chain X n given a series of observations (Y 0 , . . . , Y n-1 ) and G is the likelihood function G(x n ) = g n (y n , x n ) associated with an observation Y n = y n delivered by some sensor. In this context the BBG transformation describes the way a new observation is incorporated into the filtering distribution

η = Law(X n | (Y 0 , . . . , Y n-1 )) Ψ gn(yn,.) (η) = Law(X n | (Y 0 , . . . , Y n-1 , Y n ))
We shall return more precisely to this BBG transformation in the further developments of the article.

In matrix notation, if the state space is given by E = {1, . . . , d} we can also identify the multiplicative potential function G with the following diagonal matrix

G =        G(1) 0 • • • 0 0 0 G(2) 0 • • • 0 0 . . . . . . . . . . . . 0 0 • • • G(d -1) 0 0 0 • • • 0 G(d)        .
Using the row and column vector identifications, we find that ηG = (η(1), . . . , η(d)) G = (η(1)G(1), . . . , η(d)G(d))

with the normalizing constant given by the L 1 -norm of the vector ηG

η(G) = ηG 1 = d i=1 (ηG)(i) = ηG(1) = ηG    1 . . . 1    = d i=1 µ(x i ) G(x i ).
In this matrix notation, the BBG transformation takes the following form

Ψ G (η) = 1 ηG 1 ηG.
We end this section with a nonlinear Markov interpretation of the BBG-transformation.

Proposition 2.2 For any collection of non negative parameters ǫ(η, G) that may depend on G and on η, and such that ǫ(η, G) G(x) ≤ 1 for every x ∈ E s.t. η(x) > 0

INRIA

we have the following decomposition

∀y ∈ E Ψ G (η)(y) = (ηS η )(y) := x η(x) S η (x, y) (2.3)
with the collection of selection type Markov transitions on E defined by

S η (x, y) = ǫ(η, G) G(x) 1 x (y) + (1 -ǫ(η, G) G(x)) Ψ G (η)(y).
Before getting into the details of the proof of this proposition, for finite state spaces E = {1, . . . , d}, the BBG-transformation takes the following matrix form

(ηS η ) := [(ηS η )(1), . . . , (ηS η )(d)] = [η(1), . . . , η(d)]    S η (1, 1) • • • S η (1, d) . . . . . . . . . S η (d, 1) • • • S η (d, d)    = ηS η = Ψ G (η).
Proof:

To check this assertion, we choose a function f on E and observe that

S η (f )(x) := y S η (x, y) f (y) = ǫ(η, G) G(x) f (x) + (1 -ǫ(η, G) G(x)) Ψ G (η)(f )
from which we find that

η (S η (f )) := x η(x)S η (f )(x) = ǫ(η, G) η(Gf ) + (1 -ǫ(η, G) η(G)) Ψ G (η)(f ) = ǫ(η, G) η(Gf ) + Ψ G (η)(f ) -ǫ(η, G) η(G) η(Gf ) η(G) = Ψ G (η)(f ).
This ends the proof of the proposition.

Positive matrices and Feynman-Kac measures

Let E be a finite set. We consider a collection of matrices Q n := (Q n (x, y)) x,y∈E with non negative entries Q n (x, y) ≥ 0. To simplify the presentation, we further assume that the entries are strictly positive. Given a probability measure η 0 on E, we denote by Γ n the measure on E (n+1) defined for any path (x p ) 0≤p≤n ∈ E (n+1) of length n by the following formula:

Γ n (x 0 , . . . , x n ) := η 0 (x 0 ) Q 1 (x 0 , x 1 ) Q 2 (x 1 , x 2 ) . . . Q n (x n-1 , x n )
Notice that Γ n can alternatively be defined by the following recursive formula

∀n ≥ 1 Γ n (x 0 , . . . , x n ) = Γ n-1 (x 0 , . . . , x n-1 ) × Q n (x n-1 , x n )
with Γ 0 (x 0 ) = η 0 (x 0 ). We choose a function f n on path space E (n+1) and we set

Γ n (f n ) = x0,...,xn Γ n (x 0 , . . . , x n ) f n (x 0 , . . . , x n ) (2.4)
Feynman-Kac models provide a precise meaning of these sums. To be more precise, we set

M n (x, y) := Q n (x, y) z Q n (x, z) and G n (x) := z Q n+1 (x, z) = Q n (1)(x)
The r.h.s. term in the above formula expresses the sum in terms of a matrix type operation with the unit function 1, defined by 1(x) = 1 for any x ∈ E. By construction, it is readily checked that

Γ n (x 0 , . . . , x n ) :=η 0 (x 0 )    1≤p≤n Q p (1)(x p-1 )    ×    1≤p≤n Q p (x p-1 , x p ) Q p (1)(x p-1 )    =    0≤p<n G p (x p )    ×    η 0 (x 0 ) 1≤p≤n M p (x p-1 , x p )    We conclude that Γ n (x 0 , . . . , x n ) =    0≤p≤n G p (x p )    × P n (x 0 , . . . , x n )
with the law of the path (x p ) 0≤p≤n of a non homogeneous Markov chain X n with initial distribution η 0 = Law(X 0 ) and local Markov transitions

M n P n (x 0 , . . . , x n ) := η 0 (x 0 ) 1≤p≤n M p (x p-1 , x p )
In this interpretation, the quantities (2.4) can alternatively be rewritten in terms of a Feynman-Kac formula

Γ n (f n ) = E   f n (X 0 , . . . , X n ) 0≤p<n G p (X p )   It is also important to observe that Γ n (1) = x0,...,xn Γ n (x 0 , . . . , x n ) = E   0≤p<n G p (X p )  
The above constructions clearly works for any decomposition of Q n in the following form

Q n (x, y) = G n-1 (x) M n (x, y)
for some Markov transition M n and some potential function G n-1 . Note that in this case we have

Q n (1)(x) = y Q n (x, y) = G n-1 (x) y M n (x, y) = G n-1 (x) Definition 2.3
We let X n be a Markov chain with Markov transitions M n on some finite or countable state space E, and we consider a sequence of [0, 1]-valued potential functions G n on the set E.

The Feynman-Kac path measure associated to the pairs (M n , G n ) is the measure Q n defined for any f n on E (n+1) by the following formulae

Q 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.5)

INRIA

Notice that Q n can alternatively be defined by the following formulation:

Q n (x 0 , . . . , x n ) = 1 Z n    0≤p<n G p (x p )    P n (x 0 , . . . , x n )
with the normalizing constant Z n := Γ n (1), and the distribution of the free path evolution of the reference Markov chain X n given by:

P n (x 0 , . . . , x n ) = η 0 (x 0 ) 1≤p≤n M p (x p-1 , x p )
We also consider the measure Q n on E (n+1) defined as Q n by taking a product of potential functions from the origin p = 0, up to the final time horizon p = n

Q n (x 0 , . . . , x n ) = 1 Z n    0≤p≤n G p (x p )    P n (x 0 , . . . , x n )
with the normalizing constant Z n = Z n+1 .

Another important quantities are the n-th marginals η n and η n of the measure Q n and Q n defined by

η n (x n ) = x0,...,xn-1 Q n (x 0 , . . . , x n-1 , x n ) and η n (x n ) = x0,...,xn-1 Q n (x 0 , . . . , x n-1 , x n ) (2.6)
It is important to observe that this pair of measures can be alternatively defined for any function f on E by the following formulae

η n (f ) := γ n (f )/γ n (1) with γ n (f ) = E   f (X n ) 0≤p<n G p (X p )   (2.7)
and

η n (f ) = γ n (f )/ γ n (1) with γ n (f ) = E   f (X n ) 0≤p≤n G p (X p )   (2.8)
For n = 0, we use the convention ∅ = 1, so that η 0 = Law(X 0 ). The above assertions can be proved directly from the definitions of the measures Q n and Q n . In reference to filtering, the measure η n is often called the prediction measure at time n, while η n is referred as the updated measure at time n.

We end this section with the following elementary observation connecting the pair of measures

Q n and Q n .    0≤p≤n G p (x p )    P n (x 0 , . . . , x n ) = η 0 (x 0 ) G 0 (x 0 ) 1≤p≤n M p (x p-1 , x p ) G p (x p ) = η 0 (G 0 )    0≤p<n G p (x p )    P n (x 0 , . . . , x n )
with the probability measure P n defined below

P n (x 0 , . . . , x n ) := η 0 (x 0 ) 1≤p≤n M p (x p-1 , x p ) RR n°6991
and the pair ( G p , M p ) of potential and transitions given by

M p (x p-1 , x p ) = M p (x p-1 , x p )G p (x p ) M p (G p )(x p-1 ) and G p (x p ) = M p+1 (G p+1 )(x p )
These formulae imply that

Q n (x 0 , . . . , x n ) = 1 Z ′ n    0≤p<n G p (x p )    P n (x 0 , . . . , x n )
with some normalizing constant Z ′ n . In other words Q n coincides with the measure defined as Q n by replacing the pair (G n , M n ) by the pair of updated type quantities ( G n , M n ). From the pure mathematical point of view, this rather elementary observation shows that the two models Q n and Q n are equivalent.

Examples as Feynman-Kac models

We illustrate the rather abstract models presented in section 2.3 in the examples discussed in section 1.

Random walks confined to a set

We let M n (x, y) be the transitions of a symmetric random walk X n on the integers Z starting at the origin X 0 = 0. we also fix a interval A = {-a + 1, -a + 2, ..., a -1}, with a > 0, and we set G n := 1 A . In this situation, we have

Q n = Law ((X 0 , . . . , X n ) | ∀0 ≤ p < n, X p ∈ A ) and Q n = Law ((X 0 , . . . , X n ) | ∀0 ≤ p ≤ n, X p ∈ A )

Particle absorption models

For any [0, 1]-valued potential functions G n and any Markov transitions M n , the measures Q n have the following interpretations

Q n = Law ((X c 0 , . . . , X c n ) | T ≥ n ) and Q n = Law ((X c 0 , . . . , X c n ) | T > n ) where X c
n is the particle absorption model defined in (1.8), and T the random absorption time of the particle with M n -free evolution and absorption rates G n .

Nonlinear filtering problems

A given Markov chain X n with transitions M n is partially observed at every time step n by some noisy sensor that delivers a random observation

Y n = y n with distribution Proba (Y n = y n | X n = x n ) = G(x n , y n )
For every fixed sequence of observations Y n = y n , we set G n (x) = G(x, y n ). In this situation, the measures Q n and their normalizing constants Z n have the following interpretations

Q n = Law ((X 0 , . . . , X n ) | ∀0 ≤ p < n, Y p = y p )
and Z n = P ((Y 0 , . . . , Y n-1 ) = (y 0 , . . . , y n-1 ))

In the same way, we have

Q n = Law ((X 0 , . . . , X n ) | ∀0 ≤ p ≤ n, Y p = y p ) INRIA

Stochastic optimization algorithms

We consider a collection of MCMC transitions M β (x, y), indexed by some inverse temperature parameter β ≥ 0, on some finite state space E equipped with the counting measure λ(x) := 1/Card(E). We assume that M β (x, y) is a Markov transition with an invariant Boltzmann-Gibbs distribution µ β associated with a non negative cost function

V on E µ β (x) := 1 Z β e -βV (x) λ(x) with Z β := x e -βV (x) λ(x) (2.9)
The genealogical tree model associated with the genetic search model defined in (1.12), section 1.4, is defined by the line of ancestors of each individual ξ i n in the current population:

∀1 ≤ i ≤ N ξ i 0,n ←-ξ i 1,n ←-. . . ←-ξ i n-1,n ←-ξ i n,n = ξ i n (2.10)
We recall that this stochastic algorithm is defined in terms of a genetic type model with mutations transitions M βn and selection fitness functions G n = e -(βn+1-βn)V . For a more precise description of these ancestral lines we refer the reader to (1.4). In this situation, we can prove that Q n coincide with the limiting distribution of the ancestral lines of the genetic search algorithm, in the sense that

∀1 ≤ i ≤ N Q n = lim N →∞ Law ξ i 0,n , ξ i 1,n , . . . , ξ i n,n (2.11)
In addition, Q n can also be interpreted as the limiting occupation measure of the genealogical tree, in the sense that:

Q n (f n ) = lim N →∞ 1 N N i=1 f n (ξ i 0,n , ξ i 1,n , . . . , ξ i n,n ) (2.12)
Much more is true, we can prove that the ancestral lines are "almost" independent and identically distributed random paths with common distribution Q n . In much the same way, the measures Q n are approximated using the ancestral lines of the selected individuals ξ i n at time n.

In the end of section 2.5, we shall prove that the n-th marginal measures η n defined in (2.7) also coincide with the Boltzmann-Gibbs measures (2.9) associated with an inverse temperature parameter β n ; that is, we have that

η n (x n ) = µ βn (x) := 1 Z βn e -βnV (x) λ(x) and Z βn = 0≤p<n η p (G p ) = γ n (1) (2.13)

Combinatorial counting and sampling

Our next objective is to provide an interpretation of the measures Q n in terms of the genetic particle models presented in section 1.5. These models are defined as above by replacing in (2.9) the exponential Boltzmann-Gibbs measures by the measure λ restricted to a given non increasing sequence of subsets A n ⊂ A n-1 ⊂ E; that is, by the measures µ An defined below

µ An (x) := 1 Z An 1 An (x) λ(x) with Z An := x 1 An (x) λ(x)
We recall that the corresponding stochastic algorithm is defined in terms of a genetic type model with mutations transitions M An and selection potential functions G n = 1 An+1 , where M An stands for a sequence of MCMC transitions with invariant measures µ An . The genealogical tree model associated with this genetic model is defined as in (2.10) and the pair of convergence results (2.11) and (2.12) are satisfied for the measures Q n associated with the Markov transitions M An and the potential functions G n = 1 An+1 . Furthermore, the n-th time marginals measures η n defined as in (2.7) are given by

η n (f ) := γ n (f )/γ n (1) with γ n (f ) = E   f (X n ) 0≤p<n 1 Ap+1 (X p )  
In the end of section 2.5, we shall also prove that

η n (x n ) = µ An (x) := 1 Z An 1 An (x) λ(x) and Z An = 0≤p<n η p (1 Ap+1 ) = γ n (1) (2.14)

Genetic search algorithms

The interpretation of the measures Q n in terms of simple genetic models has already been developed in section 1.6. In this example, we recall that Q n can be seen as the asymptotic distribution of the ancestral lines of a simple genetic model with mutation transitions M n and selection potential functions G n . The n-th time marginal η n of the measures Q n defined in (2.6) or in (2.7) can be regarded as the asymptotic distribution of an individual as the size of the population N tends to infinity. We can also interpret η n as the limiting occupation measure of the genetic population at time n. These measures are often called the infinite population model in the literature on theoretical aspects of genetic models.

Directed polymers and self avoiding walks

The interpretation of the measures Q n in terms of directed polymer models has already been developed in section 1.7. We simply notice that the probability measures η n defined in (2.12) can be extended to path-valued Markov chains

X n = (X ′ 0 , . . . , X ′ n ) ∈ E n := E × . . . × E (n+1) times
where X ′ n stands for an auxiliary E-valued Markov chain. In this situation, we observe that the potential functions G n are defined on the path space E n . Whenever X ′ n is the simple random walk on the d-dimensional lattice E = Z d , for the choice of potential functions

∀n ≥ 0 G n (X n ) = G n (X ′ 0 , . . . , X ′ n ) = 1 Z d -{X ′ 0 ,...,X ′ n-1 } (X ′ n ) we have that Q n = Law (X ′ 0 , . . . , X ′ n ) ∀0 ≤ p < q < n X ′ p = X ′ q and Q n = Law (X ′ 0 , . . . , X ′ n ) ∀0 ≤ p < q ≤ n X ′ p = X ′ q
In this situation, we recall that

γ n (1) = Proba ∀0 ≤ p < q < n X ′ p = X ′ q = Card (S n-1 ) /(2d) n-1
where S n-1 is the set of self-avoiding walks of length (n -1) starting at the origin.

INRIA

Nonlinear recursive equations

This section is concerned with the analysis of the flow of probability measures η n and η n defined respectively in (2.7) and (2.8). First, using the Markov property of the signal, we readily prove that

γ n (f ) = E   E   f (X n ) 0≤p<n G p (X p ) |X 0 , . . . , X n-1     = E   E (f (X n ) |X 0 , . . . , X n-1 ) 0≤p<n G p (X p )  
from which we conclude that

γ n (f ) = E   E (f (X n ) | X n-1 ) 0≤p<n G p (X p )   Recalling that E (f (X n ) | X n-1 ) = M n (f )(X n-1 ), we conclude that γ n (f ) = γ n-1 M n (f ).
In matrix notation, we conclude that

γ n = γ n-1 M n (2.15)
On the other hand, we have

γ n (f ) = E   [f (X n )G n (X n )] 0≤p<n G p (X p )   = γ n (f G n ) (2.16) with the multiplicative function (f G n )(x n ) = f (x n )G n (x n ).
Our next objective is to derive the nonlinear recursions of the sequence of measures η n and η n . From previous formulae, we find that

η n (f ) = γ n (G n f ) γ n (G n ) = γ n (G n f )/γ n (1) γ n (G n )/γ n (1) = η n (G n f ) η n (G n )
This shows that η n is deduced from η n using the following formula

η n := Ψ Gn (η n )
with the BBG transformation:

Ψ Gn (η n )(x n ) := 1 η n (G n ) G n (x n ) η n (x n )
On the other hand, we have

η n (f ) = γ n-1 M n (f ) γ n-1 (1) = η n-1 M n (f )
from which we conclude that

η n := η n-1 M n
From these expressions, it is not difficult to establish the following proposition.

Proposition 2.4 The measures η n and η n satisfy the following equations:

η n := Ψ Gn-1 (η n-1 )M n and η n := Ψ Gn ( η n-1 M n )
Another important observation is that the unnormalized distributions γ n can be expressed in terms of the sequence of normalized measures from the origin up to the current time horizon. To be more precise, we have

γ n+1 (1) = γ n (G n ) = γ n (G n ) γ n (1) γ n (1) = η n (G n ) γ n (1)
This multiplicative recursion yields the following integration exchange formula.

γ n+1 (1) = E   0≤p≤n G p (X p )   = 0≤p≤n η p (G p ) (2.17)
More generally, we have the following proposition which can be easily checked.

Proposition 2.5 For any function f , and any time n ≥ 1, we have

γ n (f ) = η n (f ) γ n (1) with γ n (f ) = η n (f ) 0≤p≤n η p (G p )
We end this section with a two-line proof of the pair of formulae (2.13) and (2.13). In the first situation, combining the fixed point equation with the BGG transformation given below

η n = η n M βn and η n = Ψ Gn-1 (η n-1 ) with G n-1 = exp [-(β n -β n-1 )V ] we find that η n = η n M βn = Ψ Gn-1 (η n-1 )M βn
Also observe that

Z βn = λ e -βnV = η n-1 (G n-1 ) Z βn-1 = Z 0 × 0≤p<n η p (G p )
In much the same way, in the second situation, combining the fixed point equation with the BGG transformation given below η n = η n M An and η n = Ψ Gn-1 (η n-1 ) with G n-1 = 1 An we find that η n = η n M An = Ψ Gn-1 (η n-1 )M An

Nonlinear Markov chain models

In this section we design a nonlinear Markov interpretation of the flow of Feynman-Kac probability measures η n associated with a pair of potential functions and Markov transitions (G n , M n ). To simplify the presentation, we further assume that the potential functions G n take values in [0, 1]. More general models can be handled using the developments of section 2.2, dedicated to the transformations Ψ Gn associated with more general potential functions G n . Combining proposition 2.4 with (2.3), we readily find that

∀n ≥ 0 η n+1 = Ψ Gn (η n )M n+1 and Ψ Gn (η n ) = η n S n,ηn
with the Markov transitions

S n,ηn (x, y) = G n (x) 1 x (y) + (1 -G n (x)) Ψ Gn (η n )(y).

INRIA

This implies that η n+1 (z) = We let (X n ) n≥0 be a Markov chain on E, with initial distribution η 0 = Law(X 0 ) and elementary Markov transitions given by

P X n+1 = y | X n = x = K n+1,ηn (x, dy) with Law(X n ) = η n
By a simple induction on the time parameter, we can easily check that η n coincides with the distribution of the random states X n .

Definition 2.6 The law of the paths of this non linear Markov chain

P (X 0 , X 1 , . . . , X n ) = (x 0 , x 1 , . . . , x n ) = η 0 (x 0 )K 1,η0 (x 0 , x 1 ) . . . K n,ηn-1 (x n-1 , x n )
is called the Mc-Kean measure associated with the Markov transitions K n+1,ηn .

Feynman-Kac semigroup models

In section 2.2 we have seen that the distribution of a Markov chain on some finite space can be computed by simple matrix operations. In this case, if we denote the law of the chain η p = Law(X p ) at time p then we obtain the law at any time n ≥ p by an elementary composition of (n-p) transition matrices η n = η p M p,n with the matrix M p,n := M p+1 M p+2 . . . M n

For the flow of Feynman-Kac measures η n defined in (2.7) the situation is more involved. By proposition 2.4, the measure η n is connected to η n-1 by some transformation given by

η n = Φ n (η n-1 ) := Ψ Gn-1 (η n-1 )M n Definition 2.7
The semigroup of the flow of Feynman-Kac measures η n is given by the composition of the one step mappings Φ n ; that is, we have

∀0 ≤ p ≤ n Φ p,n = Φ n • Φ n-1 • . . . • Φ p+1
with the convention Φ n,n = Id the identity mapping. By construction, for any pair of indexes 0 ≤ p ≤ n, the measure η p and η n are connected with the following formula

η n = Φ p,n (η p )
Our next objective is to make more precise the definition given above. First, using (2.15) and (2.16) we find that

∀n ≥ 0 γ n (y) = ( γ n-1 M n )(y) = x γ n-1 (x)M n (x, y) and γ n (x) = γ n (x)G n (x)
from which we conclude that

γ n (y) = (γ n-1 Q n )(y) = x γ n-1 (x)Q n (x, y) with Q n (x, y) = G n-1 (x)M n (x, y)
In other words, we have proved that

γ n = γ n-1 Q n = γ n-2 Q n-1 Q n = . . . = γ p Q p,n with the matrix Q p,n := Q p+1 Q p+2 . . . Q n
and therefore, for any function f on the set E, we have

η n (f ) := γ n (f ) γ n (1) = γ p Q p,n (f ) γ p Q p,n (1) = γ p Q p,n (f )/γ p (1) γ p Q p,n (1)/γ p (1) = η p Q p,n (f ) η p Q p,n (1)
In summary, we have proved the following proposition.

Proposition 2.8 For any function f on E we have

∀0 ≤ p ≤ n Φ p,n (η p ) (f ) := η p Q p,n (f ) η p Q p,n (1)
3 Interacting particle algorithms

Mean field particle interpretation models

In section 2.6 we provide a nonlinear Markov interpretation of the flow of Feynman-kac measures η n . In this interpretation, η n can be seen as the distribution of the random states X n of a Markov chain (X n ) n≥0 on E, with initial distribution η 0 = Law(X 0 ) and elementary Markov transitions given by

P X n+1 = y | X n = x = K n+1,ηn (x, dy) with Law(X n ) = η n
At a first level of analysis, let us suppose that the elementary transitions X n-1 X n are easy to sample. In this situation,we can approximate the flow of measures by sampling a series of N independent copies (X i n ) 1≤i≤N of the chain X n . By the well-known law of large numbers, in some sense we have for any time horizon n and any function f

η N n (f ) := 1 N N i=1 f (X i n ) -→ N ↑∞ η n (f ) = E f (X n )
In the above displayed formula we have used the notation

∀x ∈ E η N n (x) = 1 N N i=1 1 X i n (x) =⇒ η N n (f ) = x η N n (x)f (x) = 1 N N i=1 f (X i n )
More generally, for any function f n on the path space E n+1 , we have

1 N N i=1 f n (X i 0 , . . . , X i n ) -→ N ↑∞ E f n (X 0 , . . . , X n ) INRIA ξ n ∈ E N selection --------→ ξ n ∈ E N mutation -------→ ξ n+1 ∈ E N (3.4)
The genetic type evolution of the system is summarized by the following synthetic picture:

ξ 1 n . . . ξ i n . . . ξ N n         S n,η N n ----------→          ξ 1 n Mn+1
----------→ . . .

ξ i n ----------→ . . . ξ N n ----------→ ξ 1 n+1 . . . ξ i n+1 . . . ξ N n+1        
with the selection Markov transition :

S n,η N n (ξ i n , x) := G n (ξ i n ) 1 ξ i n (x) + 1 -G n (ξ i n ) 1≤j≤N G n (ξ j n ) 1≤k≤N G n (ξ k n ) 1 ξ j n (x)
Using the Feynman-Kac modeling techniques presented in section 2.4, it is now easy to check that the corresponding mean field particle interpretations discussed above coincide with the algorithms discussed in section 1.

Some particle estimates

In section 3.1 we have designed a mean field particle interpretation of the flow of Feynman-Kac measures

η n (f ) := γ n (f )/γ n (1) with γ n (f ) = E   f (X n ) 0≤p<n G p (X p )   INRIA
We also recall from (2.17) that the normalizing constants γ n (1) are expressed in terms of the flow of measures (η p ) 0≤p<n with the following multiplicative formula

γ n (1) = E   0≤p<n G p (X p )   = 0≤p<n η p (G p )
This also implies that the unnormalized Feynman-Kac measures can be rewritten in the following form:

γ n (f ) = η n (f ) γ n (1) with γ n (f ) = η n (f ) 0≤p<n η p (G p )
The mean field particle interpretation of these measures presented in section 3.1 can be seen as genetic type algorithms with M n -mutation transitions and G n -potential functions, where M n represents the Markov transition of the reference Markov chain X n . Inversely, any genetic type model of that forms converge, in some sense, as the size of the population N tends to infinity to the above Feynman-Kac measures.

More formally, for any function f n on the product space E (n+1) , we have that

lim N →∞ 1 N N i=1 f n ξ i 0 , ξ i 1 , . . . , ξ i n = x0,...,xn η 0 (x 0 )K 1,η0 (x 0 , x 1 ) . . . K n,ηn-1 (x n-1 , x n ) f n (x 0 , . . . , x n )
and for any function

f on E lim N →∞ η N n (f ) = η n (f ) with η N n (f ) = 1 N N i=1 f n ξ i n (3.5)
The normalizing constants γ n (1) = lim N →∞ γ N n (1) can be computed using the unbiased particle estimates γ N n (1) defined below

γ N n (1) := 0≤p<n η N p (G p ) with ∀n ≥ 0 η N n (G n ) = 1 N G n (ξ i n )
To give a flavor of some non asymptotic variance estimates that can be obtained for such particle estimates, we present a recent result for time homogeneous models (G n , M n ) = (G, M ).

Theorem 3.1 ( [START_REF] Cérou | A non asymptotic variance theorem for unnormalized Feynman-Kac particle models[END_REF]) Suppose that the pair of potential-transitions (G, M ) are chosen so that

∀(x, x ′ , y) ∈ E 3 G(x) ≤ δ G(x ′ ) and M m (x, y) ≤ β M m (x ′ , y) (3.6)
for some m ≥ 1 and some parameters (δ,

β) ∈ [1, ∞[ 2 .
In this situation, any n ≥ 0, and any N > (n + 1)βδ m we have

E γ N n (1) -γ n (1) 2 ≤ γ n (1) 2 4 N (n + 1) β δ m (3.7)
More general results can be found in [START_REF] Moral | Feynman-Kac Formulae: Genealogical and Interacting Particle Systems with Applications[END_REF], under weaker regularity conditions. In this connection, we mention that the N -particle approximation of the unnormalized measures γ n is given by the unbiased estimates

γ N n (f ) := η N n (f ) × γ N n (1) -→ N →∞ γ n (f ) = η n (f ) × γ n (1)
At the end of this section, we discuss some path space particle estimates associated to these genetic algorithms. First, it is important to observe that the above constructions can be extended to Feynman-Kac models with reference Markov chains and potential functions on path spaces

X n := (X ′ 0 , . . . , X ′ n ) ∈ E n = (E ′ ) (n+1) and G n (X n ) := G n (X ′ 0 , . . . , X ′ n ) (3.8)
This observation is essential for the analysis of the directed polymers and the self avoiding walks discussed in section 1.7. In this situation, every particle at time n is a path-valued particle of length n. The selection transition consists of selecting the path-particle with high potential value, while the mutation transition simply consists of a path extension with an elementary move according to the auxiliary process X ′ n , with Markov transitions M ′ n on the state space E ′ . When the potential functions only depend on the terminal value of the paths

G n (X n ) := G ′ n (X ′ n )
for some fitness function G ′ n , we can check that the path-valued particle model represents the time evolution of the genealogical tree model associated with the time evolution of the individuals ξ i n evolving with M ′ n -mutations and G ′ n -selections. In this situation, if

ξ i n := (ξ i 0,n , ξ i 1,n , . . . , ξ i n,n )
stands for the i-th ancestral line of the current individual ξ i n,n after the n-th mutation, then for any function f n on E n , we have that lim

N →∞ 1 N N i=1 f n ξ i 0,n , ξ i 1,n , . . . , ξ i n,n = E f n (X ′ 0 , . . . , X ′ n ) 0≤p<n G ′ p (X ′ p ) E 0≤p<n G ′ p (X ′ p ) (3.9)
Note that the n-th time marginals converge to the Feynman-Kac measures associated to the pair (G ′ n , M ′ n ); that is, for any function f on the state space E ′ we have lim

N →∞ 1 N N i=1 f ξ i n = η ′ n (f ) = γ ′ n (f ) γ ′ n (1) with γ ′ n (f ) := E   f (X ′ n ) 0≤p<n G ′ p (X ′ p )   (3.10)
Another quantity of interest is the complete ancestral tree defined by the whole population of individuals from the origin up to the current time horizon

∀1 ≤ i ≤ N ξ i 0,0 , ξ i 1,1 , . . . , ξ i n-1,n-1 , ξ i n,n
The occupation measures of this complete ancestral converge to the Mc-Kean distributions of the nonlinear Markov chain with initial distribution η 0 = η ′ 0 and elementary transitions K ′ n+1,η ′ n (x, y) on E ′ defined as the transitions K n+1,ηn by replacing the triplet (η n , G n , M n ) by the triplet (η ′ n , G ′ n , M ′ n ). In other words, we have

lim N →∞ 1 N N i=1 f n ξ i 0,0 , ξ i 1,1 , . . . , ξ i n,n = x0,...,xn η ′ 0 (x 0 )K ′ 1,η ′ 0 (x 0 , x 1 ) . . . K ′ n,η ′ n-1 (x n-1 , x n ) f n (x 0 , . . . , x n )

An elementary proof

In this section, we provide a simple proof of the convergence result (3.5) using a rather elementary induction w.r.t. the time parameter n. To simplify the analysis we assume that the potential functions G n are chosen so that δ(G n ) := sup x,y (G n (x)/G n (y)) < ∞ and E is a finite set. The main result of this section is the following elementary variance estimates.

INRIA

L q -mean error and bias estimates

We extend in this section the variance estimates stated in proposition 3.2 to any L q -mean errors bounds. The forthcoming analysis is developed under the regularity conditions stated in the beginning of section 3.3.

Proposition 3.4

For any function f on E with osc(f ) ≤ 1, any time horizon n ≥ 0, any population size N , we have

E η N n (f ) -η n (f ) ≤ c(n)/N (3.15)
for some finite constant c(n) < ∞. In addition, for any integer q ≥ 1, we have the variance estimate

√ N E η N n (f ) -η n (f ) q 1 q ≤ b(q) c(n) (3.16)
for some finite constant c(n) < ∞, and the collection of constants b(q) defined below for any p ≥ 0:

b(2p) 2p := (2p) p 2 -p b(2p + 1) 2p+1 := (2p + 1) (p+1)
p + 1/2 2 -(p+1/2) with (q + p) p := (q + p)!/p! (3.17)

In the further development of this section, we provide a rather elementary proof of the above proposition based on a simple induction on the time parameter combined with Kintchine's type inequalities. The resulting estimates of the upper bound constants c(n) are far from being "sharp". Next, we present without proofs some more precise extensions of these results. Detailed proofs of these results can be found in [START_REF] Moral | Feynman-Kac Formulae: Genealogical and Interacting Particle Systems with Applications[END_REF].

1. When the flow of measures is stable, in the sense that it forgets exponentially fast its initial conditions, we can prove uniform estimates w.r.t. the time horizon. For instance, under the regularity condition (3.6) for some parameters (m, δ, β), we can prove that the L q -mean error estimates (3.16) holds with some constants c(n) such that sup n≥0 c(n) ≤ 2mδ 2m-1 β 3

2. Similar L q -mean error estimates (3.16) can be derived for the convergence of the occupation measures (3.9) associated with genealogical tree models. Up to a state space enlargement the convergence analysis of models associated with path-valued Markov chains (3.8) follows exactly the same line of arguments. In addition, when the flow of marginal measures η ′ n introduced in (3.10) is stable, we can prove L q -mean errors bounds (3.16) for the occupation measures of genealogical tree models with some constant that grows linearly w.r.t. the time parameter; that is, such that c(n) = c n, form some finite constant c < ∞.

3. The extensions of the L q -mean errors bounds (3.16) to general state space models and Zolotarev's type seminorms

η N n -η n F := sup f ∈F η N n (f ) -η n (f )
associated with some collection F of functions can also be found in [?].

4. Note that the mean value of the occupation measures coincides with the law of a given particle; that is, we have that

E η N n (f ) = 1 N 1≤i≤N E f (ξ i n ) = E f (ξ 1 n ) INRIA

Kintchine's type inequalities

The L m -mean errors analysis developed in this section is based on the following Kintchine type inequality for martingales with symmetric and independent increments. Lemma 3.5 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 ] m ′ /2 n 1 m ′ (3.18)
with the smallest even integer m ′ ≥ m, the bracket process [M ] n := 0≤p≤n ∆ 2 p , and the collection of constants b(m) defined in (3.17).

Proof:

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

(M n-1 + ∆ n ) 2m = 2m p=0 2m p M 2m-p n-1 (∆ n ) p
Using the symmetry condition, all the odd moments of ∆ n are null. Consequently, we find that

E (M n-1 + ∆ n ) 2m = m p=0 2m 2p E M 2(m-p) n-1 E ∆ 2p n
Using the induction hypothesis, we prove that the above expression is upper bounded by the quantity For odd integers we use twice the Cauchy-Schwarz inequality 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

INRIA

This implies

E η N n+1 (f ) -η n+1 (f ) m 1 m ≤ 1 √ N E V N n+1 (f ) m 1 m + δ(G n ) × E η N n ( fn ) -η n ( fn ) m 1 m ≤ b(m) √ N + δ(G n ) × E η N n ( fn ) -η n ( fn ) m 1 m
The end of the proof can be conducted using a simple induction on the time parameter.

Proof of the bias estimates

The proof of the bias stated in proposition 3.4 needs a little bit more work. Using decomposition (3.14), we have

Φ n+1 η N n (f ) -Φ n+1 (η n ) (f ) = 1 1 -[1 -η N n (g n )]
η N n ( fn ) -η n ( fn )

Using the fact that 1 1 -u = 1 + u + u 2 1 -u ,

we obtain 1 1 -[1 -η N n (g n )] = 1 + 1 -η N n (g n ) + 1 -η N n (g n ) 2 η N n (g n )
. for some finite constant c(n) < ∞. This implies

Now we have

E η N n+1 (f ) -η n+1 (f ) = E Φ n+1 η N n (f ) -Φ n+1 (η n ) (f ) = E η N n (
E η N n+1 (f ) -η n+1 (f ) ≤ E η N n ( fn ) -η n ( fn ) + c(n)/N
The end of the proof of the bias estimates stated in proposition 3.4 is now clear.

INRIA

Theorem 3.8 For any fixed time horizon n ≥ 0, the sequence of random fields V N n converges in law, as the number of particles N tends to infinity, to a sequence of independent, Gaussian and centered random fields V n ; with, for any f , and n ≥ 0, E(V n (f ) 2 ) = η n-1 K n,ηn-1 ([f -K n,ηn-1 (f )] 2 ) .

(3.20)

This fluctuation theorem can be used to analyze the fluctuations of the occupation measures η N n around their limiting value η n (cf. for instance [START_REF] Moral | Feynman-Kac Formulae: Genealogical and Interacting Particle Systems with Applications[END_REF], and the references therein). The extension of these results to general mean field models can be found in [START_REF] Moral | Concentration inequalities for mean field particle models[END_REF].

Theorem 3.9 For any n ≥ 0 and any function f n on E s.t. η n (f n ) = 0, we have the following convergence in law

lim N →∞ √ N η N n (f n ) = n q=0 V q Q q,n (f n ) Proof:
We consider the following matrices:

∀0 ≤ p ≤ n Q p,n (x, y) = Q p,n (x, y) η p Q p,n (1)
with η p Q p,n (1) =

x,y η p (x)Q p,n (x, y)

By construction, we have that

η p Q p,n = η n (3.21)
We also observe that from which we prove that Q p,n is a well defined semigroup of matrices, in the sense that

η p Q p,n (1)
Q p,n (x, y) = Q p,n (x, y) p≤q<n η q (G q ) =   p≤q<r η q (G q )     r≤q<n η q (G q )   (Q p,r Q r,n )(x, y) = Q p,r Q r,n (x, y)
The following easily checked decomposition is pivotal

η N p Q p,n (f ) = η N 0 Q 0,n (f ) + p q=1 η N q Q q,n (f ) -η N q-1 Q q-1,n (f ) = η N 0 Q 0,n (f ) + p q=1
η N q Q q,n (f ) -Φ q η N q-1 Q q,n (f )

+ p q=1 Φ q η N q-1 Q q,n (f ) -η N q-1 Q q-1,n (f )
Choosing f = f n s.t. η n (f n ) = 0, we find that

√ N η N p Q p,n (f n ) = p q=0 V N q Q q,n (f n ) + R N p,n (f n )
with the remainder term

R N p,n (f n ) := √ N p q=1 Φ q η N q-1 Q q,n (f n ) -η N q-1 Q q-1,n (f n ) INRIA
Using the fact that η N q-1 Q q-1,n (f n ) = η N q-1 Q q-1,q Q q,n (f n ) and Φ q η N q-1 Q q,n (f n ) = η N q-1 Q q-1,q Q q,n (f n ) η N q-1 Q q-1,q (1) = η N q-1 Q q-1,q Q q,n (f n ) η N q-1 Q q-1,q (1) we prove that

Φ q η N q-1 Q q,n (f n ) -η N q-1 Q q-1,n (f n ) = Φ q η N q-1 Q q,n (f n ) 1 -η N q-1 Q q-1,q (1)
By (3.21), we conclude that

Φ q η N q-1 Q q,n (f n ) 1 -η N q-1 Q q-1,q (1)
= Φ q η N q-1 Q q,n (f n ) -Φ q (η q-1 ) Q q,n (f ) η q-1 Q q-1,q (1) -η N q-1 Q q-1,q (1)

and therefore R N p,n (f n ) → N →∞ 0, almost surely. The end of the proof is a now a direct consequence of theorem 3.8.

0 ξ i 1 (

 1 independently of one another) in INRIA the set A 1 according to the Markov transition M A1 . Then, we perform another selection transition ξ 1 ξ 1 := ξ i 1 1≤i≤N

η

  n (x)S n,ηn (x, y) M n+1 (y, z)= x η n (x) y S n,ηn (x, y)M n+1 (y, z)from which we conclude thatη n+1 (z) = x η n (x)K n+1,ηn (x, z)with the Markov transitions K n+1,ηn defined as the matrix composition of the selection type transition S n,ηn with the free exploration Markov transition M n+1 given by K n+1,ηn (x, z) := y S n,ηn (x, y)M n+1 (y, z)

n= 2

 2 m-p) (2(m -p)) (m-p) E [M ]To take the final step, we use the fact that2m 2p 2 -(m-p) (2(m -p)) (m-p) = 2 -m (2m) m 2 -p (2p) p m p and (2p) p ≥ 2 pto conclude thatE (M n-1 + ∆ n ) 2m ≤ 2 -m (2m) -m (2m) m E ([M ] m n )

  fn ) -η n ( fn ) + E η N n ( fn ) -η n ( fn ) × 1 -η N n (g n ) +E η N n ( fn ) -η n ( fn ) × [1-η N n (gn)]

2 1 / 2 × 1 / 2 ≤

 1212 Using proposition 3.4, we haveE η N n ( fn ) -η n ( fn ) × 1 -η N n (g n ) ≤ E η N n ( fn ) -η n ( fn ) E 1 -η N n (g n ) 2 c(n)/Nand in much the same wayE η N n ( fn ) -η n ( fn ) × [1-η N n (gn)]

  = γ p Q p,n (1) γ p (1) = γ n (1) γ p (1) = p≤q<n η q (G q )

RR n°6991

The schematic picture of this elementary Monte Carlo algorithm is given below

. . .

X N n+1

Unfortunately, for most realistic models we can rarely compute explicitly the law η n of random states of the chain X n . Nevertheless, we can use the occupation measure of the N samples at every time step n to approximate the law η n .

Nevertheless, under some appropriate regularity conditions on the Markov transitions K n+1,ηn(x,y) we have

x, y) ≃ N ↑∞ K n+1,ηn (x, y) We illustrate this rather abstract discussion with the Markov transitions K n+1,ηn (x, z) discussed in section 2.6. We recall that in this case, K n+1,ηn (x, z) is given by the following matrix composition K n+1,ηn (x, z) := y S n,ηn (x, y)M n+1 (y, z)

with the Markov transitions

and the BBG transformation

Replacing η n by the occupation measure η N n (x) = 1 N N i=1 1 X i n (x), we find that

and the weighted BBG occupation measures

We use this mean field approximation strategy inductively to define a Markov chain

on the product state spaces E N . The initial state (ξ 1 0 , ξ 2 0 , . . . , ξ N 0 ) consists of N independent and identically distributed random variables with common distribution η 0 termed particles. The elementary transitions are defined by the following synthetic diagram

RR n°6991

The mean field particle interpretation of this nonlinear measure valued model is an E N -valued Markov chain ξ n = ξ 1 n , ξ 2 n , . . . , ξ N n ∈ E N with elementary transitions defined as

In other words, given the population of N particles ξ n = ξ 1 n , ξ 2 n , . . . , ξ N n at time n, the population of particles ξ n+1 = ξ 1 n+1 , ξ 2 n+1 , . . . , ξ N n+1 at time (n + 1) are independent random variables with respective distributions

The rationale behind this approximation is that η N n+1 is the empirical measure associated with N independent variables with distributions K n+1,η N n (ξ i n , x), so as soon as η N n is a good approximation of η n then, in view of (3.2), η N n+1 should be a good approximation of η n+1 . In the situation where K n+1,ηn is given by the (3.1), the flow of Feynman-Kac measures evolves according to the two-step updating/prediction transitions,

In the N -mean field particle interpretation model, this pair of recursions has been replaced by a two-step selection/mutation transition in product spaces Proposition 3.2 For any function f on E, any time horizon n ≥ 0, and any population size N , we have the variance estimate

for some finite constant c(n) < ∞, whose values do not depend on the parameter N .

The initial particle approximation

We observe that the initial population ξ 0 = ξ 1 0 , ξ 2 0 , . . . , ξ N 0 consists of N independent and identically distributed random variables with common distribution η 0 . For any function f on E, we have

with the random quantities

and the random occupation measures

Using the decomposition

it is now a simple exercise to check that

where osc(f ) stands for the oscillation of the function f

This shows that the initial distribution η 0 is approximated using the initial occupation measure η N 0 , and we have the variance estimates

The local sampling errors

We recall that the population of particle at time (n + 1) is the sequence of conditionally independent random variables

We observe that for any y ∈ E and any n ≥ 0, we have

On the other hand, for any function f on E, we have the decomposition

In the above displayed formulae, V N n (f ) stands for the random variable defined below

with the random measures

Rewritten in a different way, we have

while the limiting measures satisfy the following equation

These decompositions show that the random measures V N n can be interpreted as the local sampling errors induced by the mean field particle transitions. The above formulae also show that the Nparticle model can also be interpreted as a stochastic perturbation of the limiting system.

Using the fact that

It is rather elementary to check that

and

From these calculations, we readily prove the following local error estimates.

INRIA

Proposition 3.3 For any function f on E, any time horizon n ≥ 0, and any population size N , we have the almost sure estimates

The proof by induction

The main object of this section is to complete the proof of proposition 3.2.

By definition of the one step mappings, for any function f on E, any time horizon n ≥ 0, and any population size N , we have

Notice that the function fn is such that

From the above observations, we prove the following decomposition

This implies that

and therefore

We end the proof of proposition 3.2 using the following decomposition

Using the conditional unbiased property (3.12) and the local estimates stated in proposition 3.3, we prove that

The end of the proof of the variance estimates stated in proposition 3.2 is now clear. In section 3.3.1 we have seen that (3.16) is true at rank n = 0. Assuming that (3.16) is true at rank n, the above inequality implies that the result is also satisfied at rank (n + 1). This ends the proof of the variance estimates stated in the proposition.

For finite state space models, the bias estimates (3.15) states that

One again, the constant c(n) can be uniformly bounded w.r.t. the time parameter under appropriate regularity conditions. Furthermore, similar estimates holds true for the ancestral line of the corresponding genealogical tree models. Notice that in this situation η n is replaced by the Mc-Kean measure introduced in definition 2.6.

All of these estimates can be extended to general state space models and to the law of the first q particles in the population

where X i n is a sequence of independent and identically distributed random variables with common distribution η n . This property is called the propagation of chaos property in applied probability.

We end this section with a pair of Bernstein's type exponential concentration estimates recently obtained in [START_REF] Moral | Concentration inequalities for mean field particle models[END_REF]. For any N ≥ 1, any n ≥ 0, and any λ > 0, we have

and

.

with some parameters (r, σ, σ) such that

with the uniform local variance parameters:

In the above displayed formula the supremum is taken over all integers n ≥ 1, any probability measures and all functions f on E with osc(f ) ≤ 1.

To see how these parameters enter into the local fluctuation variance, we notice that for any probability measure we have

We also have the following formula

and their "carré du champ" function Γ Ln,µ (f, f ) defined for any x ∈ E by

we conclude that

This ends the proof of the lemma.

Initial and local sampling estimates

By construction, for any function f on the set E and any n ≥ 0 we have

For n = 0, we use the convention K n,η N n-1 (f )(ξ i n-1 ) = η 0 (f ), so that

Given ξ n-1 , we let (ζ i n ) 1≤i≤N be an independent copy of (ξ i n ) 1≤i≤N . We readily check that

For n = 0, the above conditional expectation is a simple expectation w.r.t. the initial random variables ξ 0 and its independent copy ζ 0 . This yields the formula

with the terminal value M

By a direct application of Kintchine's inequality, we prove the following proposition Proposition 3.6 For any N ≥ 1, m ≥ 1, n ≥ 0 and any test function f we have the almost sure estimate

The proof of the L m -mean errors stated in proposition 3.4 follows the same line of arguments as those developed in section 3.3.3. Once again, we use the decomposition

and the fact that

with the function fn defined in (3.13), to prove that

Central limit theorems

Let us come back to the interpretation of the mean field particle model as a stochastic perturbation of the limiting system. This perturbation model is summarized by the pair of formulae

with the centered random measures V N n defined below

These quantities can be thought as a stochastic perturbation of the limiting system, while the quantity 1 √ N is the usual precision of Monte Carlo approximation models. Our next objective is to analyze the local sampling mean error variance associated with the particle model. To clarify the presentation, we slightly abuse the notation, and we set

the function on E defined for any x by

Lemma 3.7 For any N ≥ 1, any n ≥ 0 and any functions f, g on E we have

Proof:

The first assertion is immediate. Using the fact that

we find that

This ends the proof of the lemma.

The next theorem indicates that these random perturbations behave asymptotically as Gaussian random perturbations. The details of the proof of this functional central limit theorem can be found is [START_REF] Moral | Feynman-Kac Formulae: Genealogical and Interacting Particle Systems with Applications[END_REF].