From MAPSI

Cours: TME9

TME 9: Monte-Carlo par chaînes de Markov (MCMC)

Partie obligatoire

Nous allons tester expérimentalement ce qui a été vu en TD.

Estimation de {$\pi$} par Monte Carlo

Ecrivez une fonction tirage : float -> float x float prenant en un nombre paramètre réel {$m$} et retournant un point tiré aléatoirement selon la loi uniforme dans le carré {$ [-m, m] \times [-m,m] $}.

En utilisant la fonction précédente, écrivez une fonction monteCarlo : int -> float x float np.array x float np.array prenant en paramètre {$N$}, le nombre de tirages, et retournant un triplet composé d'une estimation de {$\pi$} par Monte Carlo ainsi que du tableau numpy des {$N$} abscisses et du tableau des {$N$} ordonnées tirées qui ont permis d'obtenir cette estimation.

Tracez le carré {$ [-1, 1] \times [-1,1] $}, le cercle de centre l'origine {$(0,0)$} et de rayon 1, et les points tirés avec une couleur différente selon qu'ils sont dans le cercle ou non. Vous pouvez utiliser le code suivant afin de faire l'affichage.

import matplotlib.pyplot as plt

plt.figure()

# trace le carré
plt.plot([-1, -1, 1, 1], [-1, 1, 1, -1], '-')

# trace le cercle
x = np.linspace(-1, 1, 100)
y = np.sqrt(1- x*x)
plt.plot(x, y, 'b')
plt.plot(x, -y, 'b')

# estimation par Monte Carlo
pi, x, y = monteCarlo(int(1e4))

# trace les points dans le cercle et hors du cercle
dist = x*x + y*y
plt.plot(x[dist <=1], y[dist <=1], "go")
plt.plot(x[dist>1], y[dist>1], "ro")
plt.show()

Vous devriez obtenir une figure similaire à celle ci-dessous:

Décodage par la méthode de Metropolis-Hastings

Nous travaillons sur des textes en anglais pour éviter les problèmes avec les accents.

1) Téléchargez le fichier countWar.pkl qui contient le modèle de bigrammes appris à partir du roman War and Peace de Tolstoï. Ce texte a été choisi car il contient plus de 3 millions de caractères, ce qui est assez long pour apprendre un modèle représentatif d'une langue. Téléchargez également le message secret.txt à décoder (ou celui-ci secret2.txt, même message codé différemment si vous avez des problèmes avec le premier fichier).

2) Exécutez le code suivant pour charger en python le modèle de bigrammes et le message secret.

import pickle as pkl

# si vos fichiers sont dans un repertoire "ressources"
with open("./ressources/countWar.pkl", 'rb') as f:
    (count, mu, A) = pkl.load(f, encoding='latin1')

with open("./ressources/secret.txt", 'r') as f:
    secret = f.read()[0:-1] # -1 pour supprimer le saut de ligne
 

3) Les fonctions de décodage seront représentées comme un dictionnaire où la clé et la valeur stockée sont toutes deux de type caractère (une lettre est codée/décodée en une autre lettre). Ecrivez une fonction swapF : (char,char) dict -> (char,char) dict qui prend en argument une fonction de décodage {$\tau_t$} et qui retourne une nouvelle fonction de décodage {$\tau$} construite en écheangeant deux lettres {$c_1$} et {$c_2$} comme indiqué à la question 3.3 du TD:

Vous pouvez tester votre fonction avec le dictionnaire suivant :

tau = {'a' : 'b', 'b' : 'c', 'c' : 'a', 'd' : 'd' }

4) Ecrivez une fonction decrypt : string x (char,char) dict -> string qui, étant donné une chaine de caractères mess et une fonction de décodage {$\tau$}, retourne la chaîne de caractères obtenue par décodage de mess par {$\tau$}.

Vous pouvez tester votre fonction avec le code suivant :

tau = {'a' : 'b', 'b' : 'c', 'c' : 'a', 'd' : 'd' }
decrypt ( "aabcd", tau )
decrypt ( "dcba", tau )

qui vous produira l'affichage :

 >>> decrypt ( "aabcd", tau )
 'bbcad'
 >>> decrypt ( "dcba", tau )
 'dacb'

5) Créer un dictionnaire (une table de hashage) associant à chaque caractère son indice dans mu/A. Le code est simplement le suivant:

chars2index = dict(zip(np.array(list(count.keys())), np.arange(len(count.keys()))))

chars2index['a'] donne simplement accès à l'indice correspondant à la lettre a dans mu ou A

Il est préférable d'utiliser le fichier des indices déjà généré. Lien : fichierHash.pkl

Attention aux tables de hash en python qui ne donnent plus accès à des listes de clés: il faut re-créer la liste quand on en a besoin.

6) Ecrivez une fonction logLikelihood : string x float np.array x float np.2D-array x (char,int) dict-> string qui, étant donné un message mess (chaîne de caractères), les tableaux mu et A créés à la question 2 à partir de pickle et du dictionnaire précédent chars2index, renvoie la log-vraisemblance du message mess par rapport au modèle bigramme (mu, A).

Vous pouvez tester votre fonction avec le code suivant :

logLikelihood( "abcd", mu, A, chars2index )
logLikelihood( "dcba", mu, A, chars2index )

qui vous produira l'affichage :

 >>> logLikelihood( "abcd", mu, A, chars2index )
 -24.600258560804818
 >>> logLikelihood( "dcba", mu, A, chars2index )
 -26.274828997400395

7) Codez la méthode de Métropolis-Hastings vue en TD sous la forme d'une fonction appelée MetropolisHastings(mess, mu, A, tau, N, chars2index) en utilisant les fonctions swapF, decrypt et logLikelihood:

La méthode est une simple boucle où on fait les étapes suivantes :

La fonction retourne le message décodé le plus vraisemblable. Vous afficherez la log-vraisemblance à chaque fois qu'elle est améliorée et le message décodé correspondant pour pouvoir observer l'évolution de l'algorithme.

Afin de tester votre fonction, vous pourrez exécuter le code suivant :

def identityTau (count):
    tau = {}
    for k in list(count.keys ()):
        tau[k] = k
    return tau
MetropolisHastings( secret2, mu, A, identityTau (count), 10000, chars2index)

Attention: ce code (un peu bête) ne fonctionne pas avec secret.txt, seulement avec secret2.txt!

7) Pour accélérer les calculs, nous allons partir d'une fonction de décodage prenant les fréquences d'occurences des lettres (c'est-à-dire la lettre la plus fréquente du message codé sera décodée vers la lettre la plus fréquente observée dans le roman de Tolstoï; ensuite la deuxième lettre la plus fréquente du message codé sera décodée vers la deuxième lettre la plus fréquente observée du roman; et ainsi de suite...) Vous pouvez utiliser le code suivant pour construire une telle fonction de décodage, nommée ici tau_init.

# ATTENTION: mu = proba des caractere init, pas la proba stationnaire
# => trouver les caractères fréquents = sort (count) !!
# distribution stationnaire des caracteres
freqKeys = np.array(list(count.keys()))
freqVal  = np.array(list(count.values()))
# indice des caracteres: +freq => - freq dans la references
rankFreq = (-freqVal).argsort()

# analyse mess. secret: indice les + freq => - freq
cles = np.array(list(set(secret2))) # tous les caracteres de secret2
rankSecret = np.argsort(-np.array([secret2.count(c) for c in cles]))
# ATTENTION: 37 cles dans secret, 77 en général... On ne code que les caractères les plus frequents de mu, tant pis pour les autres
# alignement des + freq dans mu VS + freq dans secret
tau_init = dict([(cles[rankSecret[i]], freqKeys[rankFreq[i]]) for i in range(len(rankSecret))])

MetropolisHastings(secret2, mu, A, tau_init, 50000, chars2index)

Le message sera normalement intelligible lorsque vous arrivez à une log-vraisemblance de plus de -3090. Toutefois, il reste en général des erreurs dans la traduction... A quel niveau?


Partie optionnelle


Estimation de {$\pi$} par MCMC

Codez la méthode MCMC int x float -> float vue en TD en utilisant la fonction tirage de l'exercice précédent. Votre méthode prendra en paramètre la longueur {$N$} de la chaîne de Markov générée par l'algorithme (le nombre d'itérations de l'algorithme) ainsi que le déplacement {$p$} maximal (en abscisse et en ordonnée) que l'on peut effectuer à chaque étape.

Retrieved from http://webia.lip6.fr/~mapsi/pmwiki.php?n=Cours.TME9
Page last modified on December 12, 2017, at 04:44 PM EST