La base de données est la même que dans la séance précédente:
TME6_lettres.pkl
Les mécanismes de discrétisation et de tracé sont donc les mêmes (cf semaine 6)
Le but de la séance est de mettre en oeuvre un modèle de Markov caché (MMC). Pour cela, nous allons développer les algorithmes vus en TD.
On fait l'hypothèse que les états sont connus... Alors que ce n'est pas le cas. Mais il existe des stratégies simples (et parfois efficaces) pour attribuer arbitrairement des états sur les chaines. La plus classique est l'hypothèse gauche-droite, qui est bien adaptée aux signaux courts et non périodiques:
Sur un exemple:
X0 = [ 1. 9. 8. 8. 8. 8. 8. 9. 3. 4. 5. 6. 6. 6. 7. 7. 8. 9. 0. 0. 0. 1. 1.] S0 = [ 0. 0. 0. 0. 0. 0. 1. 1. 1. 1. 1. 1. 2. 2. 2. 2. 2. 3. 3. 3. 3. 3. 3.]
Au niveau de la mise en oeuvre, vous définirez la méthode def initGD(X,N):
qui prend un ensemble de séquences d'observations et qui retourne l'ensemble des séquences d'états. Pour chaque séquence x
, vous pouvez utiliser:
Etant donné la structure d'un MMC:
La nature des données nous pousse à considérer des lois de probabilités discrètes quelconques pour les émissions. L'idée est donc de procéder par comptage en définissant la matrice {$B$} comme suit:
Ce qui donne l'algorithme:
Donner le code de la fonction def learnHMM(allX, allS, N, K):
qui apprend un modèle à partir d'un ensemble de couples (seq. d'observations, seq. d'états)
Variante stabilisée: afin d'éviter les transitions à probabilité nulle et de régulariser l'ensemble du système, il est intéressant d'initialiser les matrices de transition à une valeur non nulle. Cette initialisation spéciale peut être faite de manière optionnelle en utilisant les arguments par défaut de python:
Validation sur les séquences de la classe A:
K = 10 # discrétisation (=10 observations possibles) N = 5 # 5 états possibles (de 0 à 4 en python) # Xd = angles observés discrétisés Pi, A, B = learnHMM(Xd[Y=='a'],q[Y=='a'],N,K) Pi=[ 1. 0. 0. 0. 0.] A=[[ 0.79 0.21 0. 0. 0. ] [ 0. 0.76 0.24 0. 0. ] [ 0. 0. 0.77 0.23 0. ] [ 0. 0. 0. 0.76 0.24] [ 0. 0. 0. 0. 1. ]] B=[[ 0.06 0.02 0. 0. 0. 0. 0. 0.04 0.49 0.4 ] [ 0. 0.04 0. 0.13 0.09 0.13 0.02 0.09 0.41 0.09] [ 0. 0. 0. 0.02 0.12 0.5 0.31 0.04 0. 0. ] [ 0.07 0. 0. 0. 0. 0. 0.26 0.33 0.2 0.15] [ 0.73 0.12 0. 0. 0. 0. 0. 0.02 0.02 0.12]]
Rappels:
1. Initialisation (avec les indices à 0 en python): {$$\begin{array}{ccccccccc} \delta_{0} (i) &=& \log \pi_{i} +\log b_{i} (x_{1}) \\ \Psi_{0}(i) &=& -1 \mbox{ Note: -1 car non utilisé normalement} \end{array}$$} 2. Récursion: {$$ \begin{array}{ccccccccc} \delta_{t} (j) &=&\displaystyle \left[\max_{i} \delta_{t-1}(i) + \log a_{ij}\right] + \log b_{j}(x_{t}) \\ \Psi_{t}(j) &=&\displaystyle \arg\max_{i\in [1,\ N]} \delta_{t-1} (i) + \log a_{ij} \end{array}$$} 3. Terminaison (indices à {$T-1$} en python) {$$ S^{\star} = \max_{i} \delta_{T-1}(i)$$} 4. Chemin {$$\begin{array}{ccccccccc} s_{T-1}^{\star} & = &\displaystyle \arg\max_{i} \delta_{T-1}(i) \\ s_{t}^{\star} & = & \displaystyle \Psi_{t+1}(s_{t+1}^{\star}) \end{array}$$}
L'estimation de {$\log p(x_0^{T-1} | \lambda)$} est obtenue en cherchant la plus grande probabilité dans la dernière colonne de {$\delta$}.
Donner le code de la méthode def viterbi(x,Pi,A,B):
validation en utilisant le modèle précédent {$\lambda_A$} (mêmes valeurs de N et K):
s_est, p_est = viterbi(Xd[0], Pi, A, B) s_est = [ 0. 0. 0. 0. 0. 0. 0. 0. 1. 2. 2. 2. 2. 2. 3. 3. 3. 3. 4. 4. 4. 4. 4.] p_est = -38.0935655456
En fonction de votre vitesse d'avancement, coder la version {$\alpha$} de l'estimation (exacte) des probabilités d'une séquence d'observations. Comparer avec l'approximation estimée précédemment (le résultat attendu est-il plus ou moins élevé avec la procédure exacte?).
validation
p = calc_log_pobs_v2(Xd[0],Pi, A, B) p = -34.8950422805
En utilisant la procédure de Baum-Welch simplifiée vue en TD et rappelée ci-dessous, proposer un code pour apprendre les modèles correspondant aux 26 lettres de l'alphabet.
Baum-Welch simplifié:
Le critère de convergence sera la vraisemblance.
{$$\log\mathcal L^k = \sum_{lettre}\sum_i \log p(\mathbf x_i^{lettre} | \lambda_{lettre}^k)$$}
Comme dans la séance précédente, il faut trouver un moyen d'évaluer les performances sans biaiser les résultats... Nous utilisons la même procédure.
separeTrainTest(y, pc)
documentée lors de la séance précédente.
Le code de cette partie est grossièrement identique à la séance passée. Pour plus de détails : semaine 6 Comme lors de la semaine précédente, dessiner la matrice de confusion pour détecter les lettres faciles/difficiles à distinguer.
Comme dans le TME précédent, proposer une procédure de génération de lettres.
generateHMM
qui génère une séquence d'observations (et une séquence d'états) à partir des sommes cumulées des modèles de lettres (cf usage ci-dessous)
Durant ces deux semaines, on a étudié des modèles de complexités croissantes: architecture plus complexe, nombre de paramètres en hausse. Cela nous amène à discuter le compromis biais-variance. On considère en général qu'un modèle pauvre:
A l'inverse, un modèle plus riche:
Les solutions de régularisation (comme celle proposée dans l'algorithme de comptage) cherche à optimiser le compromis, c'est à dire à isoler une partie de l'espace des hypothèses {$\mathcal H$} particulièrement intéressante:
En faisant tourner plusieurs fois les algorithmes de chaines de Markov simples et plusieurs fois les algorithmes de MMC avec différentes valeurs de N et K: