Dans ce TP, nous allons étudier des données issues de l'article suivant :
Ces données son disponibles ici. Récupérez les et placez les dans le même répertoire que votre notebook. Elles représentent des mesures relatives de densité minérale osseuse spinale sur des adolescents nord-américains. Chaque valeur est la différence de mesures prises sur deux visites consécutives, divisées par la moyenne.
L'idée est de savoir si la population peut être décrite par un mélange de deux gaussiennes. Le but est donc d'appliquer à ces données l'algorithme EM, vérifier que le $K$ optimal de la méthode de sélection de modèle vue en cours est bien $K=2$ et proposer un clustering de ces observations.
import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as sps
import scipy.linalg as spl
On commence par créer le tableau de donnees "data" à partir du fichier.
data=np.loadtxt("../datasets/densitesOs.txt")
n,d = data.shape
1) De quelle dimension sont les données ? Affichez les données sous la forme d'un nuage de points.
2) On va maintenant appliquer l'algorithme EM aux données précédentes.
EM(K,nbpas,data)
prenant en entrée les données 'data', un entier K (nombre de gaussiennes du mélange) et un nombre de pas nbpas.On fixe le nombre de pas de EM et un intervalle pour les valeurs de K que l'on va tester.
nbpas=int(5e1) # nombre de pas dans EM
rangeK=range(1,5) # valeurs de K testees pour le nombre de classes (dans la selection de modele)
Pour la sélection de modèle, il faut calculer la log-vraisemblance et la complexité du modèle.
Une manière classique de choisir le nombre de composantes $K$ dans un intervalle d'entiers est de minimiser en $K$ $$-\ell(x;\widehat{\theta}_K)+ \frac{\operatorname{dim}_K \times \log n}{2},$$ avec
Le K optimal pour la sélection de modèle est l'argmin de cette fonction.
3) Ecrire une fonction logLikelyhood(K,alpha,mu,Cov,data)
permettant de calculer la vraisemblance du GMM donné par les paramètres (K,alpha,mu,Cov) sur les données data
.
4) Implémenter la sélection de modèle et l'appliquer aux données précédentes pour déterminer le K optimal. Commentez le résultat.
5) Afficher les données en les colorant en fonction de leur classe telle qu'estimée par EM.
6) Afficher la fonction $-\ell(x;\widehat{\theta}_K)+ \frac{\operatorname{dim}_K \times \log n}{2},$ en fonction de $K$
Refaire la même chose en utilisant la fonction sklearn.mixture.GaussianMixture
de la bibliothèque scikit-learn
. Voir la page suivante
import sklearn.mixture