La méthode MCMC : Markoc Chain Monte Carlo (cours)

[@menadhara1352]

bonjours, question: comment utiliser q pour généré x ? merci pour votre réponse.

[@Achillionable]

Un moyen est en fait d'utiliser une sorte de Monte-Carlo pour ça. C'est brute force, mais ça marche Un petit code python pour illustrer : import numpy as np import matplotlib.pyplot as plt X=np.linspace(0,5,1000) def Y(X): ## ma distribution q pour obtenir x return(np.sin(X*np.pi/5)*np.abs(np.sin(1/3*X**2))) from scipy.integrate import quad I=quad(Y,0,5,args=(X)) ## j'en fait l'intégrale pour être sûr d'avoir une pdf avec intégrale=1 def Y(X): ## ma distribution q pour obtenir x, intégrale=1 return(np.sin(X*np.pi/5)/I*np.abs(np.sin(1/3*X**2))) L=[] for _ in range(10000000): x=np.random.rand()*5 if Y(x)>np.random.rand(): ## L'acceptation, pour avoir le nouveau x à partir de q L+=[x] H=np.histogram(L,bins=50) plt.hist(L,bins=50) plt.plot(X,Y(X)/np.max(Y(X))*np.max(H[0])) ## Je vérifie que mes tirages suivent en effet la distribution q. Et on observe que ça marche bien en effet. Du coup, à partir d'un q classique (chez moi Y), il suffit de faire ce petit tirage. ça fonctionnera forcément car la distribution q (ici Y) est une pdf, donc elle ne peut être supérieure à 1. Le plus efficace serait de vérifier que Y(x)> a ϵ [0,max(Y)], pour éviter les tirages inutiles. maxi=np.max(Y(X)) x=np.random.rand()*5 if Y(x)>np.random.rand()*maxi Et le prochain x' sera l'x qui respecte la condition if. Il y a probablement plus simple, mais ça fonctionne, j'ai directement pensé à ça.