TP - Apprentissage bayésien et chaîne de Markov cachée

Auteur

Stéphane Derrode, Centrale Lyon, Dpt Mathématiques & Informatique

Cadre — MSO 3.7, Apprentissage bayésien, séance #5 (4h en présentiel, CR 1 noté, déposé en fin de séance).

Objectifs

L’objectif de ce TP est de programmer certains des algorithmes vus en cours sur le modèle des chaînes de Markov cachées (Hidden Markov Chain, HMC). On applique la décision bayésienne (critère MPM) après avoir appris automatiquement les paramètres du modèle par EM (Expectation-Maximization), puis on segmente une image grâce au parcours de Peano.

Ce TP étend directement le tuto EM dans le modèle de mélange fait en autonomie entre les séances #3 et #5 : on remplace l’hypothèse i.i.d. des pixels par une dépendance markovienne entre pixels voisins. Le M-step et la classification finale (MPM) sont quasi inchangés ; seul le E-step se complique (un forward-backward au lieu d’une simple normalisation des responsabilités).

Format et pré-requis

Projet multi-fichiers .py à exécuter dans VS Code ou un terminal Python. Dépendances : numpy, scipy, matplotlib, Pillow.

Table des matières

Introduction
Simulation d’une HMC
Restauration supervisee
Restauration non supervisee
Application a la segmentation d’une image (unique partie a reporter dans votre CR)

Introduction¶

Téléchargez et décompressez le fichier zip disponible en suivant ce lien. Vous disposez ainsi de 4 fichiers Python (et le sous-répertoire Peano), formant le squelette du programme :

func.py : contient tous les algorithmes pour la décision bayésienne et pour l’apprentissage EM. Trois fonctions sont à programmer : getBeta, getGamma, et UpdateParameters. Les autres (getAlpha, getCtilde, getMPMClassif, EM_Iter, DrawCurves*, …) sont déjà fournies.
SimulHMC.py : permet de simuler N observations et N labels d’une HMC dont les paramètres sont fixés au début du programme principal. Les trajectoires simulées sont stockées dans sources/XY.out.
SupervisedHMCRestoration.py : restauration MPM lorsque les paramètres du modèle sont supposés connus (cas supervisé). Utilise getAlpha, getBeta, getGamma.
UnsupervisedHMCRestoration.py : restauration MPM avec apprentissage des paramètres par EM (cas non supervisé). Utilise tout le pipeline func.py.

Bon à savoir. Le squelette tourne sans crasher même si vous n’avez encore rien codé : getBeta et getGamma renvoient des zéros, et UpdateParameters renvoie des valeurs « neutres ». Vous pouvez donc lancer SimulHMC.py puis SupervisedHMCRestoration.py dès le début pour vérifier que tout est en place.

Simulation d’une HMC¶

Fichier principal : SimulHMC.py

Vérifiez que vous comprenez bien les algorithmes servant à simuler les états de la chaîne de Markov et les observations correspondantes.
Exécutez le programme SimulHMC.py. Dans le répertoire results, vous pouvez observer des extraits graphiques des 2 signaux générés.
Ouvrez le fichier sources/XY.out avec un éditeur de texte et observez le format des données enregistrées. Ces données seront utilisées comme données d’entrée lors des 2 prochains programmes.

Note : np.random.seed(42) est fixé au début de SimulHMC.py pour que les résultats soient reproductibles d’un lancement à l’autre.

Restauration supervisee¶

Fichier principal : SupervisedHMCRestoration.py

Vérifiez que les paramètres définis dans le programme principal sont bien identiques à ceux utilisés pour la simulation (on est ici dans le cas de la restauration à paramètres connus). Vérifiez aussi que le programme lit bien le fichier que vous avez simulé par le programme précédent (i.e. sources/XY.out).
Programmez la méthode getBeta(…) du fichier func.py, grâce à l’algorithme backward présenté en cours. Cet algorithme est utilisé dans le cadre du critère MPM (vous n’avez pas à programmer la phase forward ; elle est déjà programmée dans getAlpha).
Programmez la méthode getGamma(…) du fichier func.py.

Lancez le programme et vérifiez que vous obtenez un résultat proche du suivant :

I= [0.5 0.5]
t= [[0.95 0.05]
 [0.05 0.95]]
Confusion matrix for MPM:
 [[523.   8.]
 [  8. 461.]]
Global error rate for MPM: 0.016
By class error rate for MPM: [0.01506591 0.01705757]

Restauration non supervisee¶

Fichier principal : UnsupervisedHMCRestoration.py

Il s’agit ici d’apprendre automatiquement les paramètres du modèle à partir de l’algorithme EM (vu en cours et déjà mis en pratique dans le tuto EM dans le modèle de mélange), avant d’appliquer la décision bayésienne (critère MPM).

Analysez le squelette de l’algorithme : il est très similaire à celui du modèle de mélange (séance #3 / autonomie), avec en plus un forward-backward à chaque itération et une mise à jour de la matrice de transition \(t\).
Programmez la méthode UpdateParameters(…), appelée par la fonction EM_Iter(…) du fichier func.py. C’est le M-step : ré-estimer \((\mu_k, \sigma_k^2, \pi_k, t_{kl})\) à partir de \(\gamma\) et \(\tilde c\).

Lancez le programme. Après 30 itérations, vous devriez observer un résultat comparable à

--->iteration= 0
Initial estimations:
  Confusion matrix for MPM =
 [[272. 294.]
 [  0. 434.]]
  Global Error rate for MPM:  0.294
  Class Error rate for MPM:  [0.51943463 0.        ]
--->iteration= 1
--->iteration= 2
...
--->iteration= 29
Final estimations:
  Confusion matrix for MPM =
 [[525.   6.]
 [  8. 461.]]
  Global Error rate for MPM:  0.014
  Class Error rate for MPM:  [0.01129944 0.01705757]

Observez et commentez les courbes de convergence sauvegardées dans le répertoire results (*_EvolParam.png, *_EvolError.png).

Application a la segmentation d’une image (unique partie a reporter dans votre CR)¶

Nous allons segmenter une image carrée à niveaux de gris (typiquement 64×64 ou 128×128), en exploitant l’algorithme d’estimation des paramètres de la section précédente et le parcours de Peano.

Le parcours de Peano (famille des courbes de Hilbert) permet de convertir une image (intrinsèquement 2D) en un vecteur 1D (de taille égale au nombre de pixels dans l’image), tout en conservant au mieux le voisinage des pixels (cf. slides de cours). Il est ainsi possible d’appliquer les algorithmes précédents sur ce vecteur. Nous obtenons alors un vecteur segmenté (par MPM, suite à l’estimation des paramètres par EM) qu’il convient de reconstruire sous le format d’une image par le biais du parcours de Peano inverse.

Marche à suivre¶

Sélectionnez une image à niveaux de gris :
- convertir une image couleur via gimp ou un outil en ligne si besoin ;
- carrée et de dimensions en puissance de 2 (typiquement 64×64 ou 128×128, pour limiter les temps d’exécution) ;
- le répertoire Peano/sources/ contient déjà cible_64_bruit.png qui convient pour démarrer — mais il est préférable de travailler sur votre propre image.
Examinez les 2 programmes principaux dans le répertoire Peano/ (PeanoImage.py, InvPeanoImage.py) et exécutez-les sur votre image, pour convertir l’image en vecteur et le vecteur en image. Vérifiez que la double transformation préserve l’image.
Créez un nouveau programme principal (suggéré : SI_Peano_HMC.py) qui segmente l’image en :
1. lisant l’image et la convertissant en vecteur 1D via Peano(...) ;
2. apprenant les paramètres HMC par EM sur ce vecteur ;
3. classant par MPM ;
4. reconstruisant l’image segmentée via PeanoInverse(...).
Astuce import :
```
import sys
sys.path.append('./Peano')
from PeanoImage    import Peano
from InvPeanoImage import PeanoInverse
```

Consignes pour le compte-rendu (CR 1, noté)¶

Le CR est individuel, déposé sur Pedagogie3 avant la fin de séance, sous forme d’un fichier compressé contenant :

les sources : votre programme SI_Peano_HMC.py (et toute modification de func.py ou autre) ainsi que l’image originale que vous avez sélectionnée ;
un mini-rapport (.pdf ou .md ; pas de .doc / .docx), avec page de garde (noms, prénoms, date, titre du TP, module, encadrant).

Le rapport doit contenir :

Présentation de l’image choisie : taille, type de contenu, raison du choix.
Trois images côte-à-côte : originale · vecteur reconstitué par Peano-Inverse (sanity check) · résultat de segmentation après EM-HMC.
Courbes de convergence (*_EvolParam.png, *_EvolError.png si vous disposez de la vérité-terrain) commentées : quels paramètres bougent ? quand l’algorithme se stabilise-t-il ?
Tableau des paramètres estimés finaux (moyennes \(\mu_k\), variances \(\sigma_k^2\), matrice de transition \(t\), prior \(\pi\)) et interprétation physique (qu’est-ce que ces nombres disent de l’image ?).
Commentaire qualitatif de la segmentation obtenue — est-elle bonne ? Pourquoi ? Cas d’échec ?
(Optionnel — si vous en avez le temps)
- histogramme normalisé de l’image + mélange estimé superposé ;
- comparaison avec un algorithme baseline (seuillage d’Otsu ou KMeans).