Filtre particulaire

L'objectif d'un filtre à particule est d'estimer la densité postérieure des variables d'état compte tenu des variables d'observation. Le filtre de particules est conçu pour un modèle de Markov caché, où le système se compose de variables cachées et observables. Les variables observables (processus d'observation) sont liées aux variables cachées (state-processus) par une forme fonctionnelle connue. De même, le système dynamique décrivant l'évolution des variables d'état est également connu de façon probabiliste[1].

Un filtre de particules génériques estime la distribution postérieure des états cachés en utilisant le procédé de mesure d'observation. Considérez un espace d'état illustré dans le diagramme ci-dessous[2]

${\displaystyle {\begin{matrix}X_{0}&\longrightarrow &X_{1}&\longrightarrow &X_{2}&\longrightarrow &X_{3}&\longrightarrow &...&signal\\\downarrow &&\downarrow &&\downarrow &&...\\Y_{0}&&Y_{1}&&Y_{2}&&Y_{3}&&...&observation\end{matrix}}}$

Le problème de filtrage consiste à estimer séquentiellement les valeurs des états cachés ${\displaystyle X_{k}}$ , compte tenu des valeurs du processus d'observation ${\displaystyle Y_{0},...,Y_{k}}$, à tout moment étape ${\displaystyle k}$.

Toutes les estimations bayésiennes de ${\displaystyle X_{k}}$ suivent de la densité postérieure ${\displaystyle p(x_{k}\mid y_{0},y_{1},...,y_{k})}$. La méthodologie du filtre à particules fournit une approximation de ces probabilités conditionnelles en utilisant la mesure empirique associée à un algorithme de type génétique. En revanche, l'approche d'échantillonnage Méthode de Monte-Carlo par chaînes de Markov ou d'importance modéliserait la postérieure complète ${\displaystyle p(x_{0},x_{1},...,x_{k}\mid y_{0},y_{1},...,y_{k})}$.

Les méthodes de particules supposent souvent que ${\displaystyle X_{k}}$ et les observations ${\displaystyle Y_{k}}$ peuvent être modélisées sous cette forme:

• ${\displaystyle X_{0},X_{1},...}$ est un processus Markov sur ${\displaystyle \mathbb {R} ^{d_{x}}}$ (pour certains ${\displaystyle d_{x}\geqslant 1}$) Qui évolue en fonction de la densité de probabilité de transition ${\displaystyle p(x_{k}|x_{k-1})}$. Ce modèle est également souvent écrit de manière synthétique comme[3] - [2]

  ${\displaystyle X_{k}|X_{k-1}=x_{k}\sim p(x_{k}|x_{k-1})}$Avec une densité de probabilité initiale ${\displaystyle p(x_{0})}$.

• Les observations ${\displaystyle Y_{0},Y_{1},\cdots }$ prennent des valeurs dans un certain espace d'état sur ${\displaystyle \mathbb {R} ^{d_{y}}}$ (pour certains ${\displaystyle d_{y}\geqslant 1}$) sont conditionnellement indépendants à condition que ${\displaystyle X_{0},X_{1},\cdots }$ soient connus. En d'autres termes, chaque ${\displaystyle Y_{k}}$ dépend uniquement de ${\displaystyle X_{k}}$. En outre, nous supposons une distribution conditionnelle pour ${\displaystyle Y_{k}}$ donné ${\displaystyle X_{k}=x_{k}}$ sont absolument continus, et de manière synthétique, nous avons

  ${\displaystyle Y_{k}|X_{k}=y_{k}\sim p(y_{k}|x_{k})}$

Un exemple de système avec ces propriétés est[2]:

${\displaystyle X_{k}=g(X_{k-1})+W_{k}}$

${\displaystyle Y_{k}=h(X_{k})+V_{k}}$

Où les deux ${\displaystyle W_{k}}$ et ${\displaystyle V_{k}}$ sont des séquences mutuellement indépendantes avec la fonction de densité de probabilité s et ${\displaystyle g}$ et ${\displaystyle h}$ connues sont des fonctions connues. Ces deux équations peuvent être considérées comme des équations espace d'état et ressemblent aux équations d'espace d'état pour le filtre de Kalman. Si les fonctions g et h dans l'exemple ci-dessus sont linéaires, et si les deux ${\displaystyle W_{k}}$ et ${\displaystyle V_{k}}$ sont Gaussien, le Filtre de Kalman trouve la distribution de filtrage bayésienne exacte. Sinon, les méthodes basées sur le filtre de Kalman sont une approximation de premier ordre (EKF) ou une approximation de second ordre (UKF en général, mais si la distribution de probabilité est gaussienne, une approximation de troisième ordre est possible).

L'hypothèse selon laquelle la distribution initiale et les transitions de la chaîne de Markov sont absolument continues par rapport à la mesure de Lebesgue peuvent être assouplies. Pour concevoir un filtre à particules, nous devons simplement supposer que nous pouvons échantillonner les transitions ${\displaystyle X_{k-1}\rightarrow X_{k}}$ de la chaîne de Markov ${\displaystyle X_{k},}$ et calculer la probabilité Fonction ${\displaystyle x_{k}\mapsto p(y_{k}|x_{k})}$ (voir par exemple la description de la sélection de sélection génétique du filtre à particules donné ci-dessous). L'hypothèse absolument continue sur les transitions de Markov de ${\displaystyle X_{k}}$ ne sert qu'à dériver de manière informelle (et plutôt abusive) différentes formules entre les distributions postérieures en utilisant la règle de Bayes pour les densités conditionnelles.

Les filtres particulaires font l'hypothèse que les états ${\displaystyle x_{k}}$ et les observations ${\displaystyle y_{k}}$ peuvent être modélisées sous la forme suivante :

• La suite des paramètres ${\displaystyle x_{0},x_{1},\dots }$ forment une chaîne de Markov de premier ordre, telle que ${\displaystyle x_{k}|x_{k-1}\sim p_{x_{k}|x_{k-1}}(x|x_{k-1})}$ et avec une distribution initiale ${\displaystyle p(x_{0})}$.
• Les observations ${\displaystyle y_{0},y_{1},\dots }$ sont indépendantes conditionnellement sous réserve que les ${\displaystyle x_{0},x_{1},\dots }$ soient connus. En d'autres termes, chaque observation ${\displaystyle y_{k}}$ ne dépend que du paramètre ${\displaystyle x_{k}}$ : ${\displaystyle y_{k}|x_{k}\sim p_{y|x_{}}(y|x_{k})}$

Un exemple de ce scénario est ${\displaystyle \left\{{\begin{matrix}x_{k}=f(x_{k-1})+v_{k}\\y_{k}=h(x_{k})+w_{k}\end{matrix}}\right.}$

où à la fois ${\displaystyle v_{k}}$ et ${\displaystyle w_{k}}$ sont des séquences mutuellement indépendantes et distribuées à l'identique avec des fonctions de densité de probabilité connues et où ${\displaystyle f(\cdot )}$ et ${\displaystyle h(\cdot )}$ sont des fonctions connues. Ces deux équations peuvent être vues comme des équations de l'espace d'état et ressemblent à celles du filtre de Kalman.

Si les fonctions ${\displaystyle f(\cdot )}$ et ${\displaystyle h(\cdot )}$ étaient linéaires, et si à la fois ${\displaystyle v_{k}}$ et ${\displaystyle w_{k}}$ étaient des gaussiennes, alors le filtre de Kalman trouve la distribution de filtrage bayésien exacte. Dans le cas contraire, les méthodes à base de filtre de Kalman donnent une estimation de premier ordre. Les filtres particulaires donnent également des approximations, mais avec suffisamment de particules, les résultats peuvent être encore plus précis.

Les méthodes à particules, comme toutes les méthodes à base d'échantillonnages (telles que les MCMC), créent un ensemble d'échantillons qui approchent la distribution de filtrage ${\displaystyle p(x_{k}|y_{0},\dots ,y_{k})}$. Ainsi, avec ${\displaystyle P}$ échantillons, les valeurs espérées vis-à-vis de la distribution de filtrage sont approchées par : ${\displaystyle \int f(x_{k})p(x_{k}|y_{0},\dots ,y_{k})dx_{k}\approx {\frac {1}{P}}\sum _{L=1}^{P}f(x_{k}^{(L)})}$ où ${\displaystyle x_{k}^{(L)}}$ est la (L)-ième particule à l'instant ${\displaystyle k}$; et ${\displaystyle f(\cdot )}$, de la façon habituelle des méthodes Monte-Carlo, peut donner tous les données de la distribution (moments, etc.) jusqu'à un certain degré d'approximation.

En général, l'algorithme est répété itérativement pour un nombre donné de valeurs ${\displaystyle k}$ (que nous noterons ${\displaystyle N}$).

Initialiser ${\displaystyle x_{k}=0|_{k=0}}$ pour toutes les particules fournit une position de départ pour créer ${\displaystyle x_{1}}$, qui peut être utilisé pour créer ${\displaystyle x_{2}}$, qui peut être utilisé pour créer ${\displaystyle x_{3}}$, et ainsi de suite jusqu'à ${\displaystyle k=N}$.

Une fois ceci effectué, la moyenne des ${\displaystyle x_{k}}$ sur toutes les particules (ou ${\displaystyle {\frac {1}{P}}\sum _{L=1}^{P}x_{k}^{(L)}}$) est approximativement la véritable valeur de ${\displaystyle x_{k}}$.

L'échantillonnage avec rééchantillonnage par importance ou Sampling Importance Resampling (SIR) est un algorithme de filtrage utilisé très couramment. Il approche la distribution de filtrage ${\displaystyle p(x_{k}|y_{0},\ldots ,y_{k})}$ par un ensemble de particules pondérées : ${\displaystyle \{(w_{k}^{(L)},x_{k}^{(L)})~:~L=1,\ldots ,P\}}$.

Les poids d'importance ${\displaystyle w_{k}^{(L)}}$ sont des approximations des probabilités (ou des densités) a posteriori relatives des particules telles que ${\displaystyle \sum _{L=1}^{P}w_{k}^{(L)}=1}$.

L'algorithme SIR est une version récursive de l'échantillonnage d'importance. Comme en échantillonnage par importance, l'espérance de la fonction ${\displaystyle f(\cdot )}$ peut être approchée comme une moyenne pondérée : ${\displaystyle \int f(x_{k})p(x_{k}|y_{0},\dots ,y_{k})dx_{k}\approx \sum _{L=1}^{P}w^{(L)}f(x_{k}^{(L)}).}$

La performance de l'algorithme est dépendante du choix des distributions d'importances : ${\displaystyle \pi (x_{k}|x_{0:k-1},y_{0:k})}$.

La distribution d'importance optimale est donnée comme : ${\displaystyle \pi (x_{k}|x_{0:k-1},y_{0:k})=p(x_{k}|x_{k-1},y_{k}).}$

Cependant, la probabilité de transition est souvent utilisée comme fonction d'importance, comme elle est plus aisée de calculer, et cela simplifie également les calculs des poids d'importance subséquents : ${\displaystyle \pi (x_{k}|x_{0:k-1},y_{0:k})=p(x_{k}|x_{k-1}).}$

Les filtres à rééchantillonnage par importance (SIR) avec des probabilités de transitions comme fonction d'importance sont connues communément comme filtres à amorçage (bootstrap filters) ou algorithme de condensation.

Le rééchantillonnage permet d'éviter le problème de la dégénérescence de l'algorithme. On évite ainsi les situations où tous les poids d'importance sauf un sont proches de zéro. La performance de l'algorithme peut aussi être affectée par le choix de la méthode de rééchantillonnage appropriée. Le rééchantillonnage stratifié proposé par Kitagawa (1996) est optimal en termes de variance.

Un seul pas de rééchantillonnage d'importance séquentiel se déroule de la façon suivante :

1. Pour ${\displaystyle L=1,\ldots ,P}$, on tire les échantillons des distributions d'importances : ${\displaystyle x_{k}^{(L)}\sim \pi (x_{k}|x_{0:k-1}^{(L)},y_{0:k})}$
2. Pour ${\displaystyle L=1,\ldots ,P}$, on évalue les poids d'importance avec une constante de normalisation: ${\displaystyle {\hat {w}}_{k}^{(L)}=w_{k-1}^{(L)}{\frac {p(y_{k}|x_{k}^{(L)})p(x_{k}^{(L)}|x_{k-1}^{(L)})}{\pi (x_{k}^{(L)}|x_{0:k-1}^{(L)},y_{0:k})}}.}$
3. Pour ${\displaystyle L=1,\ldots ,P}$ on calcule les poids d'importance normalisés: ${\displaystyle w_{k}^{(L)}={\frac {{\hat {w}}_{k}^{(L)}}{\sum _{J=1}^{P}{\hat {w}}_{k}^{(J)}}}}$
4. On calcule une estimation du nombre effectif de particules comme ${\displaystyle {\hat {N}}_{\mathit {eff}}={\frac {1}{\sum _{L=1}^{P}\left(w_{k}^{(L)}\right)^{2}}}}$
5. Si le nombre effectif de particules est plus petit qu'un seuil donné ${\displaystyle {\hat {N}}_{\mathit {eff}}<N_{thr}}$, alors on effectue le rééchantillonnage:
   1. Tirer ${\displaystyle P}$ particules de l'ensemble de particules courant avec les probabilités proportionnelles à leur poids puis remplacer l'ensemble des particules courantes avec ce nouvel ensemble.
   2. Pour ${\displaystyle L=1,\ldots ,P}$ l'ensemble ${\displaystyle w_{k}^{(L)}=1/P}$.

Le terme Rééchantillonnage d'importance séquentiel (Sequential Importance Resampling) est aussi utilisé parfois pour se référer aux filtres SIR.

L'échantillonnage séquentiel par importance ou Sequential Importance Sampling (SIS) est similaire à l'échantillonnage avec rééchantillonnage par importance (SIR) mais sans l'étape de rééchantillonnage.

La version directe de l'algorithme est relativement simple en comparaison des autres algorithmes de filtrage particulaire et utilise la composition et le rejet. Pour produire un simple échantillon ${\displaystyle x}$ à ${\displaystyle k}$ de ${\displaystyle p_{x_{k}|y_{1:k}}(x|y_{1:k})}$:

(1) Fixer p=1

(2) Créer uniformément L depuis ${\displaystyle \{1,...,P\}}$

(3) Créer un test ${\displaystyle {\hat {x}}}$ depuis sa distribution ${\displaystyle p_{x_{k}|x_{k-1}}(x|x_{k-1|k-1}^{(L)})}$

(4) Créer les probabilités de ${\displaystyle {\hat {y}}}$ en utilisant ${\displaystyle {\hat {x}}}$ depuis ${\displaystyle p_{y|x}(y_{k}|{\hat {x}})}$ où ${\displaystyle y_{k}}$ est la valeur mesurée

(5) Créer une autre uniformément u depuis ${\displaystyle [0,m_{k}]}$

(6) Comparer u et ${\displaystyle {\hat {y}}}$

(a) Si u est plus grand alors répéter depuis l'étape (2)

(b) Si u est plus petite alors sauver ${\displaystyle {\hat {x}}}$ comme ${\displaystyle x{k|k}^{(p)}}$ et incrémenter p

(c) Si p > P alors arrêter

L'objectif est de créer P particules au pas ${\displaystyle k}$ en n'utilisant seulement que les particules du pas ${\displaystyle k-1}$. Cela requiert qu'une équation markovienne puisse être écrite (et calculée) pour créer un ${\displaystyle x_{k}}$ en se basant seulement sur ${\displaystyle x_{k-1}}$. Cet algorithme utilise la composition de P particules depuis ${\displaystyle k-1}$ pour créer à ${\displaystyle k}$.

Cela peut être plus facilement visualisé si ${\displaystyle x}$ est vu comme un tableau à deux dimensions. Une dimension est ${\displaystyle k}$ et l'autre dimension correspond au nombre de particules. Par exemple, ${\displaystyle x(k,L)}$ serait la L^ème particule à l'étape ${\displaystyle k}$ et peut être donc écrite ${\displaystyle x_{k}^{(L)}}$ (comme effectué plus haut dans l'algorithme).

L'étape (3) crée un potentiel ${\displaystyle x_{k}}$ basé sur une particule choisie aléatoirement (${\displaystyle x_{k-1}^{(L)}}$) a temps ${\displaystyle k-1}$ et rejette ou accepte cette particule à l'étape (6). En d'autres termes, les ${\displaystyle x_{k}}$ valeurs sont calculées en utilisant les ${\displaystyle x_{k-1}}$ calculées précédemment.

• filtre de Kalman, un estimateur analytique pour les distributions Gaussiennes
• estimation récursive bayésienne

• Sequential Monte Carlo Methods in Practice, par A Doucet, N de Freitas et N Gordon. Publié par Springer.
• On Sequential Monte Carlo Sampling Methods for Bayesian Filtering, par A Doucet, C Andrieu et S. Godsill, Statistics and Computing, vol. 10, no. 3, p. 197-208, 2000 CiteSeer link
• Tutorial on Particle Filters for On-line Nonlinear/Non-Gaussian Bayesian Tracking (2001); S. Arulampalam, S. Maskell, N. Gordon et T. Clapp; CiteSeer link

• Méthodes de Monte-Carlo Séquentielles (Filtrage Particulaire) à l'Université de Cambridge
• Animations MCL de Dieter Fox