Comprendre MCMC et l'algorithme Metropolis-Hastings

13

Au cours des derniers jours, j'ai essayé de comprendre comment fonctionne Markov Chain Monte Carlo (MCMC). En particulier, j'ai essayé de comprendre et de mettre en œuvre l'algorithme Metropolis-Hastings. Jusqu'à présent, je pense que j'ai une compréhension globale de l'algorithme, mais il y a quelques choses qui ne sont pas encore claires pour moi. Je souhaite utiliser MCMC pour adapter certains modèles aux données. Pour cette raison, je vais décrire ma compréhension de l'algorithme de Metropolis-Hastings pour ajuster une droite à certaines données observées :Df(x)=axD

1) Faites une première estimation pour . Définissez ce comme notre actuelle ( ). Ajoutez également à la fin de la chaîne de Markov ( ).a a a 0 a Caaaa0aC

2) Répétez les étapes ci-dessous plusieurs fois.

3) Évaluer la probabilité actuelle ( ) donné et . a 0 DL0a0D

4) Proposer un nouveau ( ) en échantillonnant à partir d'une distribution normale avec et . Pour l'instant, la est constante.a 1 μ = a 0 σ = s t e p s i z e s t e p s i z eaa1μ=a0σ=stepsizestepsize

5) Évaluer la probabilité nouvelle ( ) donné et . a 1 DL1a1D

6) Si est plus grand que , acceptez comme nouveau , ajoutez-le à la fin de et passez à l'étape 2.L 0 a 1 a 0 CL1L0a1a0C

7) Si est plus petit que générer un nombre ( ) dans la plage [0,1] à partir d'une distribution uniformeL 0 UL1L0U

8) Si est plus petit que la différence entre les deux probabilités ( - ), acceptez comme nouveau , ajoutez-le à la fin de et passez à l'étape 2.L 1 L 0 a 1 a 0 CUL1L0a1a0C

9) Si est plus grand que la différence entre les deux probabilités ( - ), ajoutez l' à la fin de , continuez à utiliser le même , passez à l'étape 2.L 1 L 0 a 0 C a 0UL1L0a0Ca0

10) Fin de la répétition.

11) Retirez certains éléments du début de (phase de rodage).C

12) prennent maintenant la moyenne des valeurs en . Cette moyenne est l'estimation .aCa

J'ai maintenant quelques questions concernant les étapes ci-dessus:

  • Comment puis-je construire la fonction de vraisemblance pour mais aussi pour toute fonction arbitraire?f(x)=ax
  • Est-ce une implémentation correcte de l'algorithme de Metropolis-Hastings?
  • Comment la sélection de la méthode de génération de nombres aléatoires à l'étape 7 peut-elle changer les résultats?
  • Comment cet algorithme va-t-il changer si j'ai plusieurs paramètres de modèle? Par exemple, si j'avais le modèle .f(x)=ax+b

Notes / Crédits: La structure principale de l'algorithme décrit ci-dessus est basée sur le code d'un atelier MPIA Python.

AstrOne
la source

Réponses:

11

Il semble y avoir des idées fausses sur ce que l'algorithme de Metropolis-Hastings (MH) est dans votre description de l'algorithme.

Tout d'abord, il faut comprendre que MH est un algorithme d'échantillonnage. Comme indiqué dans wikipedia

En statistique et en physique statistique, l'algorithme Metropolis – Hastings est une méthode Monte Carlo à chaîne de Markov (MCMC) permettant d'obtenir une séquence d'échantillons aléatoires à partir d'une distribution de probabilité pour laquelle l'échantillonnage direct est difficile.

Q(|)f()

  1. x0
  2. xQ(|x0)
  3. α=f(x)/f(x0)
  4. xfα
  5. x

Nk

Un exemple en R peut être trouvé dans le lien suivant:

http://www.mas.ncl.ac.uk/~ndjw1/teaching/sim/metrop/metrop.html

Cette méthode est largement utilisée dans les statistiques bayésiennes pour l'échantillonnage à partir de la distribution postérieure des paramètres du modèle.

f(x)=axx

Robert & Casella (2010), Introducing Monte Carlo Methods with R , Ch. 6, "Algorithmes de Metropolis – Hastings"

Il y a aussi beaucoup de questions, avec des pointeurs vers des références intéressantes, dans ce site sur la signification de la fonction de vraisemblance.

Un autre pointeur d'intérêt possible est le package R mcmc, qui implémente l'algorithme MH avec des propositions gaussiennes dans la commande metrop().

Habano
la source
Salut mon ami. Oui, je regarde MH dans le contexte de la régression linéaire. L'URL que vous m'avez donnée explique tout vraiment bien. Je vous remercie. Si je pose une autre question concernant MH, je posterai à nouveau une question. Merci encore.
AstrOne