Exemples d'erreurs dans les algorithmes MCMC

J'étudie une méthode de vérification automatique des méthodes Monte Carlo de la chaîne de Markov, et je voudrais quelques exemples d'erreurs qui peuvent se produire lors de la construction ou de la mise en œuvre de tels algorithmes. Points bonus si la méthode incorrecte a été utilisée dans un article publié.

Je suis particulièrement intéressé par les cas où l'erreur signifie que la chaîne a une distribution invariante incorrecte, bien que d'autres types d'erreurs (par exemple, chaîne non ergodique) seraient également intéressants.

Un exemple d'une telle erreur serait de ne pas produire de valeur lorsque Metropolis-Hastings rejette un mouvement proposé.

1. Probabilité marginale et estimateur de la moyenne harmonique

La vraisemblance marginale est définie comme la constante de normalisation de la distribution postérieure

L'importance de cette quantité vient du rôle qu'elle joue dans la comparaison des modèles via les facteurs de Bayes .

Plusieurs méthodes ont été proposées pour estimer cette quantité. Raftery et al. (2007) proposent l' estimateur de la moyenne harmonique , qui est rapidement devenu populaire en raison de sa simplicité. L'idée consiste à utiliser la relation

Par conséquent, si nous avons un échantillon de la partie postérieure, disons , cette quantité peut être approximée par$(\theta_1,...,\theta_N)$

Cette approximation est liée au concept d' échantillonnage d'importance .

Par la loi des grands nombres, comme discuté dans le blog de Neal , nous avons que cet estimateur est cohérent . Le problème est que le requis pour une bonne approximation peut être énorme. Voir le blog de Neal ou le blog de Robert 1 , 2 , 3 , 4 pour quelques exemples.$N$

Alternatives

Il existe de nombreuses alternatives pour approximer . Chopin et Robert (2008) présentent quelques méthodes basées sur l'échantillonnage de l'importance.$p({\bf x})$

2. Ne pas faire fonctionner votre échantillonneur MCMC assez longtemps (spécialement en présence de multimodalité)

Mendoza et Gutierrez-Peña (1999) déduisent la référence a priori / postérieure pour le rapport de deux moyennes normales et présentent un exemple des inférences obtenues avec ce modèle en utilisant un ensemble de données réelles. En utilisant les méthodes MCMC, ils obtiennent un échantillon de taille de la partie postérieure du rapport des moyennes qui est montré ci-dessous$2000$$\varphi$

Et obtenez l'intervalle HPD pour . Après une analyse de l'expression de la distribution postérieure, il est facile de voir qu'elle a une singularité à et que la postérieure devrait en fait ressembler à ceci (notez la singularité à )( 0,63 , 5,29 ) 0 $\varphi$ $(0.63,5.29)$$0$$0$

Ce qui ne peut être détecté que si vous exécutez votre échantillonneur MCMC assez longtemps ou si vous utilisez une méthode adaptative. Le HPD obtenu avec l'une de ces méthodes est de comme cela a déjà été rapporté . La longueur de l'intervalle HPD est considérablement augmentée, ce qui a des implications importantes lorsque sa longueur est comparée aux méthodes fréquentistes / classiques .$(0,7.25)$

3. D'autres questions telles que l'évaluation de la convergence, le choix des valeurs de départ, le mauvais comportement de la chaîne peuvent être trouvées dans cette discussion par Gelman, Carlin et Neal.

4. Échantillonnage d'importance

Une méthode d'approximation d'une intégrale consiste à multiplier l'intégrande par une densité , avec le même support, que l'on peut simuler à partir de$g$

Ensuite, si nous avons un échantillon de , , nous pouvons approximer comme suit( x 1 , . . . , x N ) $g$$(x_1,...,x_N)$$I$

Un problème possible est que devrait avoir des queues plus lourdes / similaires à / à ou le requis pour une bonne approximation pourrait être énorme. Voir l'exemple de jouet suivant dans R.f $g$$f$$N$

# Integrating a Student's t with 1 d.f. using a normal importance function   
x1 = rnorm(10000000)   # N=10,000,000
mean(dt(x1,df=1)/dnorm(x1))

# Now using a Student's t with 2 d.f. function
x2 = rt(1000,df=2)
mean(dt(x2,df=1)/dt(x2,df=2))

Darren Wilkinson sur son blog donne un exemple détaillé d'une erreur courante dans la marche aléatoire Metropolis-Hastings. Je recommande de le lire en entier, mais voici la version tl; dr.

Si la distribution cible est positive (comme les distributions gamma, etc. ) dans une dimension, il est tentant de rejeter immédiatement les propositions qui ont une valeur négative sur cette dimension. L'erreur est de jeter les propositions comme si elles ne s'étaient jamais produites et d'évaluer le ratio d'acceptation de Metropolis-Hastings (MH) des autres uniquement. C'est une erreur car cela revient à utiliser une densité de proposition non symétrique.

L'auteur suggère d'appliquer l'un des deux correctifs.

1. Comptez les «négatifs» comme échec d'acceptation (et perdez un peu d'efficacité).

2. Utilisez le bon rapport MH dans ce cas, qui est

$\pi$$\Phi$$\phi$ $\Phi(x) = \int_0^{\infty} \phi(y-x)dy$

Un cas très clair (lié à l'approximation de vraisemblance marginale mentionnée dans la première réponse) où la véritable convergence est l'exemple du problème du changement d'étiquette dans les modèles de mélange couplé à l'utilisation de l'estimateur de Chib (1995) . Comme l'a souligné Radford Neal (1999), si la chaîne MCMC ne converge pas correctement, dans le sens où elle explore certains modes de distribution cible, l'approximation Monte Carlo de Chib n'atteint pas la bonne valeur numérique.