Échantillonnage stochastique efficace par modèle Bernoulli mélange de Gaussiennes pour la résolution des problèmes inverses parcimonieux 

Méthodes générales d’échantillonnage stochastique efficace

 Les approches d’échantillonnage stochastique de type Monte Carlo par chaînes de Markov (MCMC) permettent de tirer des échantillons à partir d’une distribution de probabilité donnée en construisant une chaîne de Markov. Des exemples classiques des méthodes MCMC sont l’échantillonneur de Gibbs [GG84] et l’algorithme Metropolis-Hastings [MU49, Has70]. Bien que ces méthodes aient prouvé leur utilité, aussi bien dans le cadre bayésien que d’autres domaines du calcul scientifique (analyse numérique, physique informatique…), elles peuvent parfois présenter quelques problèmes de convergence, particulièrement dans le cas de lois en grandes dimensions et/ou multimodales. Plus précisément, ces méthodes risquent d’être piégées dans un mode local (ou bassin d’attraction) les empêchant d’explorer efficacement la loi de probabilité cible. Dans la Figure 1.2 est donnée une illustration de ce problème de « piégeage ». En effet, le passage d’un mode vers un autre nécessite de passer par des états intermédiaires de faibles probabilités. Donc la probabilité de sortir du bassin d’attraction du mode locale Ps est très faible8 et nécessite en moyenne Ts = P −1 s itérations. De ce fait, la question de l’amélioration des propriétés de convergence des méthodes MCMC a été étudiée dans la littérature, et différentes méthodes ont été proposées pour construire des échantillonneurs efficaces. L’objectif est donc de construire un échantillonneur capable d’explorer la loi a posteriori de façon efficace dans un laps de temps raisonnable, permettant ainsi d’éviter de rester bloqué dans des maxima locaux. Dans cette section, nous allons présenter quelques-une des pistes étudiées dans le cadre bayésien pour la résolution des problèmes inverses. Pour faire la distinction avec les notations introduites dans la section précédente, et pour simplifier les explications, on considère ici le cas d’une loi cible a posteriori, dépendant d’une seule variable d’intérêt u ∈ R K : p(u|z) = G(z; Au,Γ)π(u) où z représente les données observées et G(·; µ,Γ) est la densité de probabilité de la loi gaussienne de moyenne µ et matrice de covariance Γ. Ici, A ∈ RM×K est un simple opérateur linéaire, Γ ∈ R K×K est une matrice de covariance quelconque et π(u) représente la loi a priori sur 8Cette probabilité Ps peut être considérée comme inversement proportionnelle au travail (au sens physique) nécessaire pour se déplacer, sous la courbe de la densité, afin de sortir du bassin d’attraction. Méthodes générales d’échantillonnage stochastique efficace p(x) x Figure 1.2 : Illustration du problème de piégeage des approches MCMC. u. Par ailleurs, afin d’alléger les notations et expressions mathématiques, nous écrivons p(u) au lieu p(u|z) 9 . Finalement, on introduit la fonction ψ tel que ψ(u) = ln (p(u)). On peut distinguer trois types d’approches pour l’échantillonnage efficace : les méthodes MCMC inspirées de l’optimisation, les méthodes d’augmentation de données, et les échantillonneurs partiellement marginalisés, que nous présentons ci-dessous. 

Méthodes MCMC inspirées de l’optimisation

Les méthodes MCMC inspirées de l’optimisation cherchent à guider la chaîne de Markov vers les régions les plus probables en intégrant des outils mathématiques et algorithmiques développés dans le contexte de l’optimisation. Nous présentons ici trois approches de ce type. Langevin Monte Carlo (LMC). L’algorithme de Langevin non ajusté (Unadjusted Langevin Algorithm) (ULA) est l’un des échantillonneurs les plus connus utilisant le gradient de la logdensité pour guider la chaîne de Markov. Celui-ci est basé sur la diffusion de Langevin qui est un processus de Markov à temps continu {u(τ ) : τ ∈ R+} et qui est défini comme la solution de l’équation différentielle stochastique suivante [RS02] : du(τ ) = 1 2 ∇ψ (u(τ )) dτ + dw(τ ), u(0) = u0, τ ∈ R+ (1.15) où w(τ ) décrit un mouvement brownien dans R K. Sous condition que la densité p(u) soit différentiable et non-nulle partout sur R K le processus {u(τ ) : τ ∈ R+} converge en distribution vers p(u) quand τ → ∞ [RT96]. Cependant, la simulation d’échantillons à partir de (1.15) est souvent impossible. En pratique, la discrétisation d’Euler-Mayurama de la diffusion de Langevin est alors utilisée pour aboutir à l’algorithme ULA : u (t) = u (t−1) + γ∇ψ  u (t−1) + p 2γb (1.16) où b est un bruit blanc gaussien et γ est un paramètre qui contrôle la discrétisation temporelle ainsi que la variance de b. Évidemment, cette discrétisation introduit un biais et donc ULA produit des échantillons qui ne sont pas distribués exactement suivant p(u). Dans l’algorithme 9La dépendance de la loi a posteriori aux données et autres paramètres auxiliaires sera explicitée quand nécessaire pour éviter les confusions. – 18 – Chapitre 1. Introduction Metropolis Adjusted Langevin Algorithm (MALA) ce biais est corrigé en ajoutant une étape Metropolis-Hastings (MH) d’acceptation/réjection qui garantit la convergence vers la loi cible p(u). L’algorithme MALA revient donc à un algorithme Metropolis-Hastings à marche aléatoire (Random-Walk) avec une loi de proposition gaussienne : u (t) ∼ N  u (t−1) + γ∇ψ  u (t−1) , 2γI  Ceci peut être interprété comme une proposition à partir d’une approximation quadratique locale de la fonction ψ(u) autour de u (t−1) avec une matrice hessienne Hψ = (2γ) −1 I. Par ailleurs, les performances de MALA peuvent être améliorées en introduisant une matrice d’incrément temporel P définie positive qui permet de prendre en compte la courbure de la loi a posteriori de manière plus adéquate, on parle alors de MALA pré-conditionné. Les choix classiques pour la matrice P sont la hessienne de ψ(u) ou la matrice d’information de Fisher. Toutefois, pour des problèmes de grande dimension le calcul (inversion et/ou factorisation) de ces matrices est généralement gourmand en coût de calcul, par conséquent des approximations locales de ces matrices sont privilégiées en pratique. Majorize-Minimize Metropolis-Hastings (3MH). S’appuyant sur l’algorithme MALA et les approches Majorize-Minimize, plus précisément les approches semi-quadratiques multiplicatives [GR92], une extension des méthodes Langevin Monte Carlo a été proposée dans [MCBBP20]. En particulier, sous certaines conditions (notamment la différentiabilité de ψ(u)), le logarithme de la loi a posteriori ψ(u) peut être approximé localement par une fonction majorante quadratique Q(u) dont la matrice hessienne P(t) = P(u (t) ) peut-être calculée directement à partir de ψ(u (t) ). u (t) ∼ N  u (t−1) + γ(P (t−1)) −1∇ψ  u (t−1) , 2γ(P (t−1)) −1  Cette stratégie permet ainsi d’adapter localement la matrice P(t) en fonction de l’itéré u (t) afin d’améliorer le taux d’acceptation de l’étape MH et donc l’efficacité de l’échantillonneur. En effet, il a été montré dans un cadre applicatif de restauration d’image que l’algorithme 3MH est plus efficace que l’algorithme MALA et sa version pré-conditionnée par la matrice Hessienne (i.e., Newton MCMC). Proximal Monte Carlo (PMC). Pour relaxer la condition de différentiabilité des méthodes à base de gradient, comme les méthodes que nous venons de voir, des approches s’inspirant des méthodes d’optimisation convexe non différentiable ont été récemment proposées. Ainsi, une adaptation de l’algorithme MALA pour échantillonner des distributions non-lisses (nonsmooth) a été proposée dans [Per16] donnant lieu à l’algorithme P-MALA où, comparé à (1.16) la discrétisation d’Euler-Mayurama fait intervenir l’opérateur proximal [CP11] du logarithme de la loi p(u) au lieu du gradient : u (t) =  1 − γ λ  u (t−1) + γ λ proxλ ψ(u (t−1)) + p 2γb (1.17) Rappelons que l’évaluation de l’opérateur proximal revient à résoudre le problème d’optimisation convexe10 suivant : proxλ ψ(u) = arg max v∈RK ψ(v) − 1 2λ ku − vk 2 10Quand p(u) est une densité log-concave i.e., ψ est concave. 

 Méthodes générales d’échantillonnage stochastique efficace 

De la même façon que l’opérateur proximal peut être interprété comme une descente de gradient sur l’enveloppe de Moreau, P-MALA peut être interprété comme un algorithme MALA classique effectué sur l’approximation de Moreau de la loi a posteriori p(u) qui est donnée par : pλ(u) = sup v∈RK p(v) exp  − 1 2λ ku − vk 2  (1.18) Toutefois, dans le cadre bayésien, l’évaluation de cet opérateur pour le logarithme de la loi a posteriori p(u) peut parfois s’avérer compliquée. Afin de surmonter cette difficulté, l’algorithme Moreau-Yosida Unadjusted Langevin Algorithm MYULA a été proposé dans [DMP18]. Celui-ci se base sur le fait qu’en inférence bayésienne, la loi a posteriori p(u) se compose d’un terme de fidélité aux données qui est généralement lisse (e.g., loi gaussienne) et un a priori potentiellement non-lisse. En d’autres termes, le logarithme de cette loi s’écrit ψ(u) = ψ1(u) +ψ2(u) avec ψ1(u) une fonction lisse (e.g., quadratique) et ψ2(u) une fonction potentiellement non-lisse. Ainsi, l’opérateur proximal est calculé uniquement sur la fonction non-lise ψ2(u) qui est généralement plus simple à évaluer que l’opérateur proximal de ψ(u). Similairement à P-MALA, MYULA peut également être interprété comme un algorithme MALA classique effectué sur la loi a posteriori donnée par la densité suivante : pλ(u|z) = G(z; Au,Γ)πλ(u) où πλ(u) est l’approximation de Moreau de l’a priori π(u). Notons que l’algorithme MYULA implique deux niveaux d’approximation la première vient de l’approximation de Moreau, et la seconde de la discrétisation temporelle de diffusion de la Langevin. Cependant, il est possible d’envisager, comme pour P-MALA, de rectifier les approximations par une étape d’acceptation/rejet. 

 MCMC et augmentation de données 

Dans la littérature, différentes études ont été menées autour des approches d’augmentation de données (Data Augmentation) ou de séparation de variables (Variable splitting)), qui consistent à introduire des variables auxiliaires (ou latentes) afin d’améliorer l’efficacité des méthodes MCMC. Ces méthodes partagent toutes le même principe suivant. Étant donné qu’il est difficile de tirer des échantillons directement à partir de la loi cible p(u), ces méthodes construisent, en introduisant des variables aléatoires auxiliaires v, une loi jointe p(u, v). Ensuite une chaîne de Markov {u (t) , v (t)}t=1,…,T est construite en tirant des échantillons u et v de façon alternée à partir des lois conditionnelles, tel que présenté dans Algorithme 1.2. Algorithme 1.2 Échantillonneur de Gibbs générique des approches augmentation de données. A chaque itération (t), à partir des échantillons courants (u (t) , v (t) ) : 1. tirer (u (t+1) , v (t+1)) comme suit : (a) v (t+1) ∼ p(v|u (t) ), (b) u (t+1) ∼ p(u|v (t+1)). Évidemment, la loi jointe p(u, v) doit être choisie tel que l’échantillonnage à partir des lois conditionnelles soit possible et facile à implémenter. Dans la littérature on trouve des schémas d’augmentation de données exacts, tel que p(u) = R p(u, v) dv comme dans le BLASSO [PC08] à titre d’exemple. Mais aussi des schémas d’augmentation de données approchés c’est-à-dire p(u) ≈ R p(u, v) dv comme par exemple le modèle Asymptotically Exact Data Augmentation (AXDA) proposé dans [VDC20] qui permet de s’affranchir de la question de la construction du schéma d’augmentation de données au prix d’une approximation.  Un exemple typique des méthodes MCMC avec variables auxiliaires est l’algorithme Monte Carlo Hamiltonien [DKPR87, GC11] qui peut être vu comme l’algorithme Metropolis-Hastings combinant les méthodes d’optimisation variationnelle et d’échantillonnage stochastique. Dans HMC la loi cible p(u) est augmentée par un vecteur de variables auxiliaires v ∈ R K distribué suivant une gaussienne centrée de matrice de covariance Q telle que : p(u, v) = p(u)G(v; 0, Q) Des échantillons de u et v sont ensuite générés suivant la dynamique hamiltonienne donnée par les équations différentielles paramétriques suivantes : du dt = −∇v ln (p(u, v)) = Q−1v dv dt = −∇u ln (p(u, v)) = ∇uψ(u) (1.19) Comme la résolution analytique des équations (1.19) est généralement impossible on fait appel en pratique à des schémas de discrétisation ciblant une approximation de p(u, v) puis cette approximation est corrigée par une étape de Metropolis-Hastings. Variables auxiliaires Monte Carlo (auxMC). L’échantillonnage à partir des lois a posteriori de grande dimension est une tache compliquée, notamment si la vraisemblance introduit une certaine corrélation (e.g., bruit corrélé). Le problème est d’autant plus compliqué que les variables d’intérêts sont considérées a priori corrélées (i.e., a priori non indépendant) Un exemple typique de ce type de situation est rencontré en restauration d’image quand A représente un filtre de lissage (blur operator) associé à un a priori prenant en compte certaines informations structurelles sur l’image par exemple un opérateur de gradient discret. Dans le but de dissocier les deux sources de corrélation et traiter chacune séparément, une stratégie basée sur l’introduction de variables auxiliaires à été proposé dans [MCBBP18] en particulier pour les lois a posteriori impliquant au moins un terme gaussien (la vraisemblance ou l’a priori, voire les deux). En s’inspirant des stratégies de minimisation semi-quadratique (en particulier les travaux dans [BBFAC04]) la fonction ψ(u) est augmentée par un vecteur de variable auxiliaires v ∈ R J tel que : ψ(u) = sup v ψ(u, v) = sup v ψ(u) + Q(u, v) (1.20) où Q(u, v) est une fonction concave quadratique en u et v, telle que la loi jointe p(u, v) s’écrit : p(u, v) = exp (ψ(u, v)) = p(u)G(v; Pu, Q) =⇒ p(u) = Z p(u, v) dv avec Q ∈ R K×K une matrice symétrique définie positive et P ∈ R J×K. L’efficacité de l’échantillonnage dépend bien entendu du choix de ces deux matrices qui est lié au problème traité et particulièrement, si la source de corrélation provient : (i) de la matrice de covariance du bruit Γ uniquement, (ii) du couplage de cette matrice avec l’opérateur A en d’autres termes ATΓ −1A. Par ailleurs, afin de ne pas compromettre le gain apporté par cette stratégie de variables auxiliaires, ces matrices doivent avoir une structure simple afin de permettre l’échantillonnage efficace de v qui est distribué suivant une loi gaussienne multidimensionnelle : v (t+1)|u (t) ∼ N (Pu (t) , Q). 

Méthodes générales d’échantillonnage stochastique efficace

 L’échantillonnage des variables d’intérêt u quand à lui se fait via la loi conditionnelle p(u|v) qui admet une forme plus simple à traiter que la loi p(u) initiale : Cas (i) : p(u|v) = G(Au; f1(v), µI)π(u), avec 0 < µ < 1/kΓ −1 kS Cas (ii) : p(u|v) = G(u; f2(v), µI)π(u), avec 0 < µ < 1/kATΓ −1AkS où k·kS est la norme spectrale et fi sont des fonctions de v dépendant des données. Observons que le terme gaussien ne fait plus intervenir les matrices sources de corrélation, car celles-ci sont couplées aux variables auxiliaires. Le cas (i) revient à un problème impliquant une « vraisemblance » avec un bruit non corrélé, alors que le cas (ii) peut s’interpréter comme un simple problème de débruitage sous l’a priori π(u). De plus, si cet a priori est indépendant π(u) = Q k πk(uk) l’échantillonnage de u peut alors se faire de façon indépendante pour chacune des composantes uk. Par ailleurs, des algorithmes plus élaborés peuvent être également considérés pour échantillonnage de u comme l’algorithme MALA pré-conditionné si l’a priori est différentiable ou si l’a priori remplit les conditions nécessaires, comme par exemple l’algorithme 3MH [MCPBB16]. Il a été montré empiriquement que cette stratégie d’augmentation de données est plus efficace que l’algorithme MALA classique ainsi que l’algorithme HMC. Modèles bayésiens hiérarchiques. Une autre approche d’augmentation de données dans le cadre bayésien s’appuie sur des modèles bayésiens hiérarchiques. Un exemple typique de ces méthodes est le BLASSO introduit dans [PC08] dans le cadre de la de régression linéaire parcimonieuse. Plus précisément, le BLASSO repose sur la représentation de l’a priori, supposé indépendant, en mélange de gaussiennes qui fait intervenir des variables auxiliaires v ∈ R K + . En d’autres termes π(u) s’écrit : π(u) = Y k πk(uk) = Y k Z R+ G(uk; 0, vk)πv(vk) dvk (1.21) Un échantillonneur de Gibbs peut alors être construit à partir de la loi a posteriori jointe de u et v : p(u, v) = G(z; Au,Γ)G(u; 0, V) Y k πv(vk) (1.22) où V = diag{v}. Cette stratégie permet de simplifier l’échantillonnage des variables d’intérêt u dont la loi conditionnelle a posteriori est gaussienne. Par ailleurs, comme l’a priori est indépendant, les variables auxiliaires vk peuvent être échantillonnées indépendamment. Il est important de noter que ce schéma d’augmentation est exact, dans le sens où il n’introduit aucun biais et ne nécessite donc pas d’étape de Metropolis-Hastings. En effet (1.21) implique que p(u) = Z RK + p(u, v) dv, donc les échantillons {u (t)}t=1,…,T générés conjointement avec {v (t)}t=1,…,T sont bien distribués suivant la loi cible p(u) lorsque l’échantillonneur à convergé. Ceci est un avantage par rapport aux algorithmes LMC et HMC. De plus, contrairement à l’approche auxMC, ici aucune supposition n’est faite sur la gaussianité de l’a priori ni de la vraisemblance. Dans [PC08] le cas particulier d’un a priori de Laplace est considéré avec des variables auxiliaires exponentielles, toutefois la stratégie d’augmentation de données s’applique pour n’importe quel a priori admettant la représentation (1.21). Par ailleurs, il est tout à fait possible d’envisager des cas où le modèle bayésien hiérarchique s’applique à la vraisemblance quand cette dernière n’est pas gaussienne comme par exemple la loi hyperbolique. S’appuyant sur ce travail, d’autres modèles hiérarchiques ont été étudiés dans [CGGK10], notamment des cas où l’a priori présente une certaine corrélation entre composantes adjacentes uk comme le Fused LASSO par exemple. Cependant, certaines applications requièrent la prise en compte de lois a priori pour lesquelles il n’est pas toujours évident de trouver une représentation hiérarchique de type (1.21). Dans ce cas, le modèle Asymptotically Exact Data Augmentation (AXDA) proposé dans [VDC20] permet de construire de façon systématique un modèle bayésien hiérarchique moyennant une approximation. Plus précisément, des variables auxiliaires v ∈ R K sont introduites pour définir une approximation asymptotiquement exacte de la loi cible p(u) ou de l’a priori π(u) (voire de la vraisemblance11) : Approximation de la loi cible : pρ(u) = Z Kρ(u, v)p(v) dv, ρ > 0 (1.23) Approximation de l’a priori : πρ(u) = Z Kρ(u, v)π(v) dv, ρ > 0 (1.24) où Kρ(u, v) = K(ρ −1 (u−v)) est un noyau (e.g., noyau de lissage [DE12]) permettant de contrôler l’écart entre les variables u et v. En particulier, le noyau Kρ doit converger vers la distribution de Dirac quand ρ → 0, afin que l’approximation soit asymptotiquement exacte : lim ρ→0 pρ(u) = lim ρ→0 Z Kρ(u, v)pρ(v) dv = p(u), ou lim ρ→0 πρ(u) = lim ρ→0 Z Kρ(u, v)πρ(v) dv = π(u). Il est important de souligner que le paramètre ρ joue un rôle crucial car celui-ci permet de contrôler l’approximation. En d’autres termes de faibles valeurs de ρ permettent de rapprocher pρ(u) de la loi cible p(u). Toutefois, plus ρ est faible, plus grande est la corrélation entre u et v ce qui peut compromettre considérablement l’efficacité de l’échantillonnage, nous reviendrons sur cet aspect plus en détails lorsque nous aborderons les approximations asymptotiquement exactes dans le chapitre suivant. Il a été montré empiriquement que l’échantillonneur de Gibbs construit à partir du modèle AXDA est efficace, notamment parce qu’il permet de simplifier les corrélations structurelles entre l’a priori et la vraisemblance et que l’approximation engendre un biais relativement faible notamment dans le cadre de la restauration d’images ou encore d’échantillonnage de gaussiennes en grande dimension [VDC20]. Notons que dans le cas où Kρ(u, v) est un noyau de lissage, le modèle AXDA peut-être une alternative intéressante pour les approches Proximal Monte Carlo en particulier quand l’évaluation de l’opérateur proximal est coûteuse et/ou ne peut pas être faite de façon exacte. 

 Échantillonneurs partiellement marginalisés 

Dans un autre registre, les échantillonneurs partiellement marginalisés ont été étudiés. Le principe est de générer des échantillons en combinant des étapes durant lesquelles l’échantillonnage se fait pour certains paramètres à partir de la loi a posteriori marginalisée et pour d’autres paramètres à partir de la loi a posteriori jointe. Dans [DP08] un certain nombre d’outils mathématiques ont été introduits qui permettent de construire un Partially Collapsed Gibbs Sampler (PCGS) afin d’améliorer les propriétés de mélange et donc la convergence des chaînes de Markov : 1. Réduction : consiste à ignorer certaines variables lors de l’échantillonnage, 11Notons que des travaux similaires on été conduit dans [RJLW20] pour l’approximation de la vraisemblance.  2. Permutation : consiste à interchanger l’ordre dans lequel les variables sont échantillonnées, 3. Marginalisation : consiste à échantillonner certaines variables à partir de lois conditionnelles marginalisées. Bien que la réduction est une étape qui permet de réduire considérable le coût de calcul par itération des échantillonneurs, celle-ci doit être manipulée avec soin afin de ne pas modifier la distribution stationnaire de la chaîne de Markov. La marginalisation et la permutation maintiennent toutes les deux la distribution stationnaire de la chaîne de Markov. Toutefois, la marginalisation permet d’améliorer considérablement la convergence des échantillonneurs, tandis que l’effet de la permutation est généralement faible. Prenons l’exemple simple d’une loi cible jointe p(u, θ). Celle-ci peut représenter une loi a posteriori où u ∈ R K est la variable d’intérêt et θ ∈ R est un hyper-paramètre inconnu du modèle qui nécessite d’être estimé également (i.e., π(u) devient π(u, θ)). Au lieu de considérer un échantillonneur de Gibbs classique à partir des lois a posteriori conditionnelles p(u|θ) et p(θ|u), il est possible, moyennant une combinaison des outils listés ci-dessus, de transformer l’échantillonneur de Gibbs en un échantillonneur de type PCGS plus efficace exploitant la loi a posteriori de θ marginalisée par rapport à u : p(θ) = Z R p(u, θ) du = Z R G(z; Au,Γ)π(u, θ) du. Un tel échantillonneur PCGS est présenté dans Algorithme 1.3. Algorithme 1.3 Échantillonneur de PCGS de la loi jointe p(u, θ) exploitant le principe de marginalisation A chaque itération (t), à partir des échantillons courants (u (t) , θ(t) ) : 1. tirer (u (t+1), θ(t+1)) comme suit : (a) θ (t+1) ∼ p(θ), (b) u (t+1) ∼ p(u|θ (t+1)), Remarquons d’abord que dans l’exemple de l’Algorithme 1.3 le paramètre θ est échantillonné à partir de sa loi marginalisée. Ceci implique que si p(θ) est une loi usuelle nous obtenons des échantillons {θ (t)} indépendants. De plus, s’il est possible d’échantillonner conjointement la variable d’intérêt u à partir de p(u|θ) ce PCGS revient à un échantillonneur Gibbs par bloc, c’est-à-dire que les deux étapes du PCGS sont équivalentes à un échantillonnage conjoint à partir de p(u, θ). Ainsi les échantillons {u (t) , θ(t)} produits par le PCGS sont des échantillons indépendants de la loi jointe cible. Dans le cadre d’applications concrètes, le PCGS n’est pas toujours un échantillonneur de Gibbs par bloc, cependant celui-ci reste nettement plus efficace que l’échantillonneur de Gibbs classique de la loi jointe (e.g., [GIL11, DT10]). Dans un cadre général, il existe plusieurs manières de construire un PCGS, notamment dans le cas où il y a plusieurs variables mises en jeu, et le choix de la stratégie de construction peut évidemment avoir un impact considérable sur l’efficacité de l’échantillonneur PCGS résultant, tant en termes de propriété de mélange de la chaîne de Markov qu’en termes de coût de calcul. En théorie, plus l’ensemble des variables marginalisées est grand plus l’échantillonneur PCGS est efficace. Cependant, la loi marginalisée par rapport à cet ensemble des variables doit être simple à manipuler et ne doit pas nécessiter un temps de calcul important pour être échantillonnée au point d’annihiler le gain apporté par la marginalisation. Par ailleurs, comme nous le verrons par la suite, la loi marginale n’est pas toujours calculable analytiquement ce qui empêche la – 24 – Chapitre 1. Introduction construction de l’échantillonneur PCGS. En effet, dans notre exemple de la loi jointe p(u, θ), la marginalisation de u n’est possible que si la loi π(u|θ) est conjuguée à la vraisemblance. 

 Conclusion sur l’échantillonnage stochastique efficace

 Nous avons vu que plusieurs méthodes ont été proposées dans littérature dans le but d’améliorer l’efficacité des méthodes d’échantillonnage stochastique. Certaines de ces méthodes s’inspirent des outils développés en optimisation, en particulier celles qui utilisent l’information du premier et/ou du second ordre du logarithme de la loi cible afin de guider les chaîne de Markov vers des régions plus probables. Évidement, ces méthode s’appliquent uniquement si la loi cible est différentiable. Les approches Proximal-MC permettent de relaxer cette condition de différentiabilité, toutefois l’évaluation de l’opérateur proximal peut parfois s’avérer coûteuse en temps de calcul. Par ailleurs, ces approches ne peuvent pas être considérées dans le cadre des a priori de type Bernoulli-D car ces deniers prennent en compte des variables discrètes. D’autre part, les méthodes d’augmentation de données, introduisent des variables auxiliaires afin de simplifier l’échantillonnage à partir des lois a posteriori conditionnelles. Le schéma d’augmentation de données peut être soit exact, sous certaines conditions sur la loi cible restreignant ainsi l’utilisation de ces approches uniquement à certaines lois de probabilité, soit en considérant un schéma d’augmentation de données systématique facile à mettre en œuvre, mais au prix d’une approximation de la loi cible. Finalement, les échantillonneurs partiellement marginalisés cherchent au contraire à réduire l’ensemble des variables mises en jeu de façon partielle (c’est-à-dire uniquement sur certaines étapes de l’échantillonneur) afin de réduire les corrélations entre les variables dans le but d’améliorer les propriétés de mélange des échantillonneurs. Bien que cette approche soit particulièrement efficace, la loi a posteriori marginale n’admet pas toujours une forme analytique exploitable. De plus, l’échantillonnage de cette loi marginale doit être suffisamment simple à implémenter pour ne pas compromettre le gain apporté par la marginalisation.