Partitionnement

Objectifs du chapitre

À la fin de ce chapitre, vous serez capable de:

• formuler le partitionnement d'un jeu de données sans étiquettes comme la minimisation de l'inertie, et dire ce que ce critère favorise et ce qu'il ignore;
• dérouler à la main l'algorithme de Lloyd (les $K$-moyennes) et démontrer que ses deux étapes ne font jamais augmenter l'inertie, donc qu'il s'arrête;
• expliquer pourquoi le résultat dépend des centres de départ, et y remédier par des départs multiples et l'initialisation $K$-means++;
• choisir un nombre de groupes $K$ par la méthode du coude et le coefficient de silhouette, en sachant ce que ni l'un ni l'autre ne garantit;
• construire et lire le dendrogramme d'un partitionnement hiérarchique à lien simple et à lien complet;
• écrire un mélange gaussien, effectuer à la main une itération de l'algorithme EM, démontrer la formule de mise à jour des moyennes, et voir les $K$-moyennes comme la limite d'EM quand les affectations deviennent certaines.

Apprendre sans étiquettes

Une coopérative de consommation veut mieux connaître ses membres. Elle dispose, pour chacun, de ses achats de l'année: combien de fois il est venu, ce qu'il a acheté, à quelle heure, dans quel magasin. Elle n'a en revanche aucune «bonne réponse» à apprendre: personne n'a jamais classé ses membres en catégories, et la direction ne sait même pas combien de catégories il serait raisonnable de distinguer. La question n'est plus «prédire $y$ à partir de $\mathbf{x}$», comme dans les onze chapitres précédents, mais «ces données ont-elles une structure, et laquelle?».

C'est la situation de l'apprentissage non supervisé (unsupervised learning) de la définition 1.4: le jeu de données se réduit à $\mathcal{D}=\{{\mathbf{x}}_i{\}}_{i=1}^n$, sans étiquettes. Ce chapitre traite de sa tâche la plus ancienne et la plus répandue, le partitionnement (clustering): regrouper les exemples en groupes d'exemples semblables. Le chapitre 13 traitera de l'autre grande tâche non supervisée, la réduction de dimension.

La définition est volontairement vague sur «proches» et «éloignés»: c'est là que se logent tous les choix. Une méthode de partitionnement, c'est une façon précise de mesurer la ressemblance entre exemples et la qualité d'une partition. Nous en verrons trois. Les $K$-moyennes mesurent la compacité des groupes autour de leur centre; le partitionnement hiérarchique fusionne pas à pas les groupes les plus proches; les mélanges gaussiens supposent que les données proviennent de $K$ lois normales et cherchent lesquelles. Les trois répondent à la même question par des critères différents, et il arrive qu'ils ne donnent pas la même réponse.

Pourquoi partitionner? Pour résumer: onze millions de clients se décrivent mal un par un, mais assez bien par cinq profils types et leur effectif. Pour explorer: découvrir qu'un groupe de patients répond mal à un traitement, ou que les pannes d'une flotte de véhicules se répartissent en deux familles distinctes, ouvre une question qu'on ne savait pas poser. Pour préparer un apprentissage supervisé: le numéro du groupe peut devenir une caractéristique, ou l'on peut faire étiqueter à la main quelques exemples par groupe au lieu de tout étiqueter. Et pour compresser: remplacer chaque couleur d'une image par la plus proche de $K$ couleurs représentatives, ou chaque échantillon d'un signal par le plus proche de $K$ niveaux, est un partitionnement — c'est d'ailleurs le problème de quantification d'un signal qui a donné naissance à l'algorithme de Lloyd.

Jeu de données C — «Clients d'une coopérative». Onze membres P1 à P11, décrits par deux scores de comportement d'achat $x_1$ et $x_2$, déjà mis à l'échelle (la section sur les limites des $K$-moyennes reviendra sur ce point). Ce sont des données fictives, construites pour le cours, assez petites pour que chaque calcul se fasse à la main et que le résultat se lise sur une figure. Ici $n=11$ et $d=2$.

Un coup d'œil au tableau, ou à la figure 12.1 plus loin, suggère trois groupes: P1 à P4 en bas à gauche, P5 à P8 à droite, P9 à P11 en haut à gauche. Cette intuition sera confirmée par le calcul — mais un algorithme ne «voit» pas de figure, et nous verrons qu'il peut parfaitement passer à côté. Pour les programmes du chapitre:

# Jeu de données C: (x1, x2), de P1 à P11 dans l'ordre. Données fictives.
clients = [(1, 1), (2, 1), (1, 2), (2, 3), (6, 5), (7, 5),
           (6, 6), (8, 7), (3, 7), (2, 8), (4, 

Un critère: l'inertie

Pour transformer «des groupes compacts» en un problème de calcul, on associe à chaque groupe un centre, et l'on mesure la qualité d'une partition par la somme des carrés des distances de chaque point au centre de son groupe. C'est la même idée que le risque empirique quadratique du chapitre 1: on remplace chaque exemple par un représentant, et l'on paie le carré de l'écart.

Plus l'inertie est petite, plus les points sont serrés autour de leurs centres. Le critère est entièrement géométrique: il ne connaît que les distances euclidiennes, et il mesure l'écart au carré, ce qui lui donne la même sensibilité aux points éloignés que la perte carrée et la moyenne (exemple 1.2). La lettre est $K$, majuscule, et non $k$: au chapitre 7, $k$ comptait des voisins; ici $K$ compte des groupes, et les deux nombres n'ont rien à voir.

Le résultat suivant justifie le choix du barycentre comme centre: c'est, parmi tous les points de l'espace, celui qui minimise la somme des carrés. C'est la version vectorielle du théorème 1.2.

Démonstration. Écrivons ${\mathbf{x}}_i-\mathbf{c}=({\mathbf{x}}_i-\bar{\mathbf{x}})+(\bar{\mathbf{x}}-\mathbf{c})$ et développons le carré de la norme avec le produit scalaire:

Sommons sur $i$. Le terme croisé vaut $2({\sum }_i({\mathbf{x}}_i-\bar{\mathbf{x}}){)}^{\top }(\bar{\mathbf{x}}-\mathbf{c})$, et par définition du barycentre: il est nul. Il reste (12.3). Le premier terme du membre de droite ne dépend pas de , le second est positif et ne s'annule que pour .

L'identité (12.3) a une conséquence que nous utiliserons constamment: on peut écrire l'inertie comme le minimum d'une fonction de deux familles de variables, les affectations $z_1,\dots ,z_n$ et des centres ${\mathbf{m}}_1,\dots ,{\mathbf{m}}_K$ quelconques,

le minimum étant atteint quand chaque ${\mathbf{m}}_k$ est le barycentre ${\mathbf{\mu }}_k$ de son groupe. Minimiser $W$ sur les partitions revient donc à minimiser $J$ à la fois sur les affectations et sur les centres.

Combien de partitions faut-il examiner? Le nombre de façons de partager $n$ objets en $K$ groupes non vides est le nombre de Stirling de seconde espèce $S(n,K)$, qui se comporte comme $K^n/K!$. Pour le jeu C, $S(11,3)=28501$: un ordinateur les passe toutes en revue en une seconde, et c'est ainsi que nous avons vérifié que la partition en trois groupes de la figure 12.1 est la meilleure. Mais pour $n=100$ et , dépasse , et l'énumération est hors de question. On sait d'ailleurs que trouver une partition d'inertie minimale est un problème , même dans le plan quand fait partie des données: aucun algorithme efficace n'est connu, et l'on n'en attend pas. On se contente donc d'un algorithme qui une partition tant qu'il le peut.

L'algorithme de Lloyd

L'écriture (12.4) suggère une stratégie naturelle: minimiser $J$ alternativement en chacune des deux familles de variables, l'autre étant fixée. Les deux minimisations partielles sont faciles. À centres fixés, chaque terme $\Vert {\mathbf{x}}_i-{\mathbf{m}}_{z_i}{\Vert }^2$ ne dépend que de $z_i$: il suffit d'envoyer chaque point vers le centre le plus proche. À affectations fixées, le théorème 12.1 dit que le meilleur centre de chaque groupe est son barycentre.

L'algorithme porte le nom de Stuart Lloyd, qui l'a décrit en 1957 dans une note interne des laboratoires Bell sur la quantification des signaux, publiée en 1982; le nom de $K$-moyennes vient d'un article de MacQueen de 1967. La règle d'égalité n'est pas un détail de présentation: sur des données à coordonnées entières, les égalités de distance sont fréquentes, et deux programmes qui les tranchent différemment peuvent aboutir à des partitions différentes. Fixer la règle rend l'algorithme déterministe: les affectations sont une fonction des seuls centres. La démonstration de l'arrêt s'en sert.

Chaque itération coûte $O(nKd)$ opérations: on calcule $nK$ distances en dimension $d$, puis on fait $K$ moyennes. C'est très peu, et c'est la raison de la popularité de la méthode: elle passe sans difficulté à des millions de points.

Démonstration. Étape d'affectation. Les centres ${\mathbf{\mu }}^{(t)}$ sont fixés. Pour chaque $i$, le choix $z_i^{(t+1)}$ minimise $\Vert {\mathbf{x}}_i-{\mathbf{\mu }}_k^{(t)}{\Vert }^2$ sur , donc . En sommant sur , on obtient la seconde inégalité de (12.5).

Étape de mise à jour. Les affectations ${\mathbf{z}}^{(t+1)}$ sont fixées, et $J$ se découpe en une somme sur les groupes: $J={\sum }_k{\sum }_{i:z_i=k}\Vert {\mathbf{x}}_i-{\mathbf{\mu }}_k{\Vert }^2$, où le -ième terme ne dépend que de . Pour un groupe non vide, le théorème 12.1 dit que le barycentre minimise ce terme, strictement: il le rend plus petit que n'importe quel autre centre, en particulier que si . Un groupe vide contribue et garde son centre. En sommant, , avec égalité si et seulement si aucun centre n'a bougé.

Démonstration. Supposons qu'à l'itération $t$ les affectations aient changé, si bien que l'algorithme effectue une mise à jour. Deux cas se présentent. Si ${\mathbf{\mu }}^{(t+1)}={\mathbf{\mu }}^{(t)}$, l'affectation suivante est calculée à partir des mêmes centres que la précédente; comme l'affectation est une fonction déterministe des centres (c'est le rôle de la règle d'égalité), elle redonne ${\mathbf{z}}^{(t+1)}$, et l'algorithme s'arrête à l'étape 2. Sinon, le théorème 12.2 donne $W^{(t+1)}<J({\mathbf{z}}^{(t+1)},{\mathbf{\mu }}^{(t)})\le W^{(t)}$: l'inertie décroît strictement. Tant que l'algorithme ne s'arrête pas, les partitions successives ont donc des inerties strictement décroissantes, et elles sont — une partition qui reviendrait aurait la même inertie qu'à son premier passage. Or il n'y a qu'un nombre fini de partitions de points en au plus groupes. La suite ne peut donc pas se prolonger indéfiniment.

La démonstration borne le nombre d'itérations par le nombre de partitions, qui est astronomique; en pratique, l'algorithme s'arrête presque toujours après quelques dizaines d'itérations au plus. Mais le théorème dit aussi ce que l'algorithme ne fait pas: il s'arrête sur une partition qu'aucune des deux étapes ne sait améliorer, ce qui n'en fait pas une partition d'inertie minimale. C'est un minimum local de $J$, au sens où l'on ne peut l'améliorer ni en déplaçant les centres, ni en réaffectant les points un à un vers leur centre le plus proche.

La figure 12.1 donne une lecture géométrique de l'étape d'affectation: envoyer chaque point vers son centre le plus proche, c'est découper le plan en cellules de Voronoi des centres, comme le faisait la règle du plus proche voisin au chapitre 7, et regarder dans quelle cellule tombe chaque point. Chaque cellule est l'intersection des demi-plans «plus proche de ${\mathbf{\mu }}_k$ que de ${\mathbf{\mu }}_j$», tous délimités par des droites; elle est donc convexe. Retenez ce fait: il explique la principale limite des $K$-moyennes, que nous retrouverons plus loin.

Deux départs, deux partitions

L'exemple 12.1 a eu de la chance, ou plutôt ses centres initiaux étaient bien choisis: un dans chacun des trois groupes que l'œil distingue. Recommençons avec trois centres pris dans le même coin.

La figure 12.2 met les deux exécutions côte à côte. À droite, deux centres se sont partagé les quatre clients du bas et n'en sont jamais sortis; le troisième a dû, seul, représenter les sept autres, et il s'est posé au milieu d'eux, dans une région où il n'y a aucun client. C'est le visage typique d'un mauvais minimum local: un groupe coupé en deux d'un côté, deux groupes fondus en un de l'autre.

L'explorateur suivant permet de rejouer les deux exécutions, et une troisième, étape par étape. Le troisième départ, P9, P10 et P11, place lui aussi ses trois centres dans un même groupe — et atteint pourtant la partition optimale. Sa première affectation présente une égalité: P5 est à la même distance carrée, $13$, de P9 et de P11, et la règle l'envoie vers le centre d'indice le plus petit, P9.

Ce que l'explorateur montre, et qu'il faut retenir: l'inertie affichée ne remonte jamais, quel que soit le départ — c'est le théorème 12.2 —, et l'algorithme s'arrête toujours — c'est le théorème 12.3. Mais la valeur à laquelle il s'arrête dépend du départ. Être «tous dans le même groupe» n'est pas en soi une condamnation: depuis P9, P10 et P11, la première affectation donne au centre P9 tous les points du bas, et les centres se dispersent dès la première mise à jour. Depuis P1, P2 et P4, en revanche, les trois centres sont si serrés que deux d'entre eux se partagent durablement le même groupe.

Programmer l'algorithme demande deux fonctions, une par étape; la boucle qui les alterne et le critère d'arrêt sont fournis. Un point est un tuple de longueur quelconque, ce qui permet d'utiliser le même code en dimension 1, 2 ou 50.

Plusieurs départs et initialisation $K$-means++

Puisque le résultat dépend du départ, le premier remède est évident: lancer l'algorithme plusieurs fois depuis des départs différents, et garder la partition d'inertie la plus faible. L'inertie sert ici de juge, et c'est légitime: contrairement à l'erreur d'entraînement d'un modèle supervisé, elle est exactement le critère que l'on cherche à minimiser.

Sur le jeu C, le calcul se fait exactement. Il y a $\left(\genfrac{}{}{0ex}{}{11}{3}\right)=165$ façons de choisir trois clients distincts comme centres de départ. L'algorithme de Lloyd, lancé depuis chacune, aboutit à seulement quatre partitions différentes:

Un départ tiré au hasard atteint donc l'optimum avec probabilité $134/165\approx 0,8121$, et le manque avec probabilité $31/165\approx 0,1879$. Avec $R$ départs indépendants, la probabilité que tous le manquent est $(31/165{)}^R$: $0,006632$ pour $R=3$, pour . Les relances sont un remède puissant, et elles se parallélisent sans effort. Mais sur des données réelles, avec beaucoup de groupes, la proportion de bons départs peut être bien plus faible que %, et il vaut la peine de choisir les départs plus intelligemment.

L'idée est simple: un mauvais départ place deux centres dans le même groupe, alors qu'un bon les disperse. On choisit donc les centres un par un, en favorisant les points éloignés des centres déjà choisis.

Un point déjà choisi a $D=0$ et ne peut pas être tiré deux fois; un point proche d'un centre a peu de chances de l'être; un point lointain en a beaucoup. Le carré n'est pas arbitraire: c'est la contribution du point à l'inertie si l'on s'arrêtait là, et l'on tire donc les points proportionnellement à ce qu'ils «coûtent».

L'initialisation $K$-means++ a été proposée par Arthur et Vassilvitskii en 2007, avec une garantie qui a fait son succès.

Nous admettons ce résultat: sa démonstration, une récurrence sur les groupes optimaux qui ont déjà reçu un centre, dépasse le cadre du cours. Sa portée est réelle: aucune garantie de ce type n'existe pour un départ uniforme, et le facteur ne croît que comme le logarithme du nombre de groupes. Il faut aussi en mesurer la modestie. Pour $K=3$, le facteur vaut $8(\ln 3+2)\approx 24,79$; sur le jeu C, l'inertie moyenne de l'initialisation seule vaut $35,0111$, soit $2,94$ fois l'optimum, très en dessous de la borne. La borne est une assurance contre les catastrophes, pas une prévision. Et comme l'algorithme de Lloyd ne fait ensuite que diminuer l'inertie, elle vaut a fortiori pour la partition finale.

Choisir le nombre de groupes

Jusqu'ici, $K=3$ était donné. Dans la plupart des applications, il ne l'est pas, et c'est la question la plus difficile du partitionnement. Le premier réflexe — prendre le $K$ d'inertie minimale — ne marche pas: l'inertie optimale décroît toujours quand $K$ augmente (exercice 12.4), jusqu'à $W^{\ast }(n)=0$, atteint en mettant chaque point seul dans son groupe. Minimiser l'inertie sur $K$ reviendrait à ne rien regrouper du tout. Il faut un autre principe.

La méthode du coude. On trace $W^{\ast }(K)$ en fonction de $K$ et l'on cherche le $K$ à partir duquel ajouter un groupe ne rapporte plus beaucoup: le «coude» de la courbe. Sur le jeu C, l'énumération exhaustive des partitions donne les inerties optimales suivantes.

Passer de un à deux groupes divise l'inertie par trois, de deux à trois encore par plus de trois; ensuite, chaque groupe supplémentaire ne gagne plus que quelques unités, en isolant un point (P8 pour $K=4$, puis P4 pour $K=5$). Le coude est net, en $K=3$. Sur des données réelles, il l'est rarement: la courbe descend souvent en pente douce, et deux personnes y liront deux coudes différents. La méthode est un outil de lecture, pas une règle de décision.

Le coefficient de silhouette. Un critère plus local demande, pour chaque point, s'il est mieux dans son groupe que dans le groupe voisin.

Un $s(i)$ proche de $1$ signifie que le point est bien plus près de son groupe que du groupe voisin; proche de $0$, qu'il est à la frontière; négatif, qu'il serait mieux dans le groupe voisin. On choisit alors le $K$ dont la partition a la plus grande silhouette moyenne. Le critère est dû à Rousseeuw (1987). Contrairement à l'inertie, il ne décroît pas mécaniquement avec $K$: découper un groupe compact en deux fait chuter $b(i)$ pour tous ses points.

La silhouette demande les distances entre toutes les paires de points, soit $O(n^2)$ calculs: c'est bien plus cher que l'algorithme de Lloyd lui-même, et l'on se contente souvent, sur de grands jeux de données, de la calculer sur un échantillon.

Les limites des $K$-moyennes

Le critère d'inertie encode des hypothèses précises sur la forme des groupes. Quand elles sont fausses, l'algorithme ne se plaint pas: il renvoie la partition d'inertie minimale, qui n'est plus la bonne.

Des groupes convexes. Nous l'avons vu sur la figure 12.1: une partition stable de l'algorithme de Lloyd affecte chaque point à son centre le plus proche, donc chaque groupe est l'ensemble des points d'une cellule de Voronoi, et ces cellules sont convexes. Deux groupes en forme de croissants emboîtés, ou un anneau qui entoure un disque, ne peuvent pas être séparés: aucune frontière rectiligne ne coupe un anneau de son centre. Les $K$-moyennes couperont alors l'anneau et le disque par des droites, en morceaux sans signification. Le lien simple, plus loin, n'a pas cette limite.

Des groupes de même étendue. L'inertie compte des distances au carré, sans tenir compte de la dispersion propre à chaque groupe. Un groupe très étalé à côté d'un petit groupe dense est mal traité: la médiatrice entre les deux centres passe à mi-chemin, et elle ampute le grand groupe de ses points les plus proches du petit. Les mélanges gaussiens, qui donnent à chaque groupe sa propre variance, corrigent ce défaut.

Des points aberrants. Un barycentre est une moyenne, et une moyenne suit les valeurs aberrantes (exemple 1.2). Un seul client aux achats extravagants déplace le centre de son groupe, ou, pire, accapare un centre à lui seul. La variante des $K$-médoïdes (k-medoids), qui impose que chaque centre soit l'un des points du groupe, et qui peut utiliser d'autres distances, est plus robuste, au prix d'un calcul plus lourd.

L'échelle des caractéristiques. C'est la limite la plus fréquente en pratique, et le chapitre 7 l'a annoncée: comme les $k$ plus proches voisins, les $K$-moyennes reposent sur la distance euclidienne, qui est dominée par la caractéristique dont les valeurs s'étalent le plus. Supposons que le premier score du jeu C soit exprimé dans une unité dix fois plus petite, si bien que $x_1$ vaut $10,20,10,20,60,\dots$ au lieu de $1,2,1,2,6,\dots$. Les clients sont les mêmes, et le second score n'a pas changé. Mais la partition d'inertie minimale en trois groupes, obtenue par énumération, devient

P10 $=(2,8)$, client du groupe du haut, rejoint les clients du bas parce que son premier score, $2$, est le leur: la différence de $5$ à $7$ unités sur $x_2$ ne pèse plus rien face aux écarts de $10$ à $60$ sur $x_1$. La partition est désormais une partition selon presque seul.

Le partitionnement hiérarchique

Les $K$-moyennes demandent $K$ à l'avance et répondent par une seule partition. Une autre famille de méthodes construit toute une hiérarchie de partitions emboîtées, de $n$ groupes d'un point jusqu'à un seul groupe, et laisse choisir ensuite le niveau auquel couper.

Couper le dendrogramme par une droite horizontale à une hauteur $h$ donne une partition: les groupes sont les sous-arbres qui pendent sous la droite. Plus la coupe est basse, plus il y a de groupes. Pour que cette lecture ait un sens, il faut que les hauteurs croissent le long de l'arbre, de sorte qu'un groupe soit toujours formé plus haut que ses parties; c'est le cas pour les liens simple et complet. Si $A$ et $B$ forment la paire la plus proche, alors pour tout autre groupe $C$, $D(A\cup B,C)$ est le minimum (lien simple) ou le maximum (lien complet) de $D(A,C)$ et $D(B,C)$, deux nombres au moins égaux à : aucune fusion ultérieure ne peut avoir lieu plus bas. D'autres liens existent: le (), moyenne des distances entre paires, et le , qui fusionne les deux groupes dont la réunion fait le moins augmenter l'inertie (12.2) — c'est la version hiérarchique des -moyennes.

Déroulons le lien simple sur le jeu C. Les plus petites distances valent $1$: P1–P2, P1–P3, P5–P6 et P5–P7. Les premières fusions forment, à la hauteur $1$, les groupes $\{P1,P2,P3\}$ et $\{P5,P6,P7\}$. Viennent ensuite, à la hauteur $\sqrt{2}\approx 1,41$, P4 (à de P3) et le trio ; puis P8, à de P7. Il reste alors trois groupes — exactement ceux de la partition optimale des -moyennes. Le lien simple fusionne ensuite le groupe de droite et celui du haut à la hauteur , distance entre P7 et P11 , puis tout le monde à , distance entre P4 et P9.

Le lien complet procède dans le même ordre au début, mais à des hauteurs plus grandes, puisqu'il mesure la plus grande distance: $\{P1,P2,P3\}$ se forme à $\sqrt{2}$ (distance P2–P3), $\{P9,P10,P11\}$ à $2$, P4 rejoint le groupe du bas à , P8 celui de droite à . La coupe en trois groupes donne la même partition. Mais les deux dernières fusions ont lieu à (P8–P10, les deux points les plus éloignés des groupes de droite et du haut) et (P1–P8). La figure 12.3 montre les deux arbres sur la même échelle.

La différence de hauteur des deux arbres n'est pas une question d'échelle; elle traduit deux conceptions opposées de ce qu'est un groupe. Pour le lien simple, deux groupes sont proches dès que deux de leurs points le sont: c'est une notion de connexité. Couper le dendrogramme du lien simple à la hauteur $h$ revient à relier par un segment tous les couples de points à distance au plus $h$, et à prendre les composantes connexes du graphe obtenu. C'est pourquoi le lien simple sait suivre un anneau ou un croissant, de proche en proche, là où les $K$-moyennes échouent. C'est aussi son défaut, l'effet de chaîne (chaining): une file de points intermédiaires, même clairsemée, suffit à fondre deux groupes bien distincts. Sur le jeu C, il n'y a pas de telle file, mais on en voit l'amorce: le seul couple P7–P11 rapproche le groupe de droite de celui du haut à $2,83$, alors que leurs points sont en moyenne bien plus éloignés. Pour le lien complet, deux groupes ne sont proches que si tous leurs points le sont: il produit des groupes compacts, de diamètre limité, et il est moins sensible aux chaînes, mais un seul point aberrant suffit à repousser une fusion.

Le partitionnement hiérarchique a deux avantages: il ne demande pas $K$ à l'avance, et il est déterministe — pas de départ à tirer, pas de minimum local au sens des $K$-moyennes. Il a un coût: il faut les $n(n-1)/2$ distances entre paires, donc une mémoire en $O(n^2)$, et l'algorithme naïf fait $n-1$ fusions en parcourant à chaque fois toutes les paires de groupes, soit $O(n^3)$ opérations. Des implémentations soignées descendent à pour le lien simple, qui revient à construire un arbre couvrant de poids minimal (algorithme de Kruskal; Algorithmes, chapitre 8 «Arbres couvrants minimaux et union-find»), mais la barrière quadratique reste: au-delà de quelques dizaines de milliers de points, on lui préfère les -moyennes. Enfin, chaque fusion est définitive: une erreur faite tôt, à cause d'une égalité ou d'un point aberrant, n'est jamais corrigée.

Mélanges gaussiens et algorithme EM

Les $K$-moyennes affectent chaque point à un groupe, sans nuance: P4, un peu à l'écart de son groupe, compte autant que P3, en plein milieu. Et elles traitent tous les groupes comme des boules de même taille. Un modèle probabiliste lève les deux limites: il suppose que chaque point a été produit par l'un de $K$ mécanismes aléatoires, et il demande, pour chaque point, avec quelle probabilité il provient de chacun.

Le modèle

La variable $z_i$ est l'affectation du point $i$, exactement comme dans les $K$-moyennes, mais elle est désormais aléatoire et cachée: on observe $x_i$, jamais $z_i$. Le partitionnement devient un problème d'estimation: trouver les paramètres qui rendent les données les plus vraisemblables, par la méthode du (Probabilités et statistique, chapitre 10). La log-vraisemblance des points indépendants vaut

Le logarithme d'une somme ne se simplifie pas: les dérivées partielles de $\mathcal{L}$ couplent tous les paramètres, et il n'existe pas de solution explicite des équations $\nabla \mathcal{L}=0$. Si, en revanche, on connaissait les $z_i$, le problème serait trivial: chaque point serait attribué à une composante, et l'on estimerait chaque loi normale sur ses propres points par une moyenne et une variance empiriques, chaque poids par une proportion. L'algorithme EM exploite ce contraste: à défaut de connaître les $z_i$, il les devine, en probabilité.

Les responsabilités sont une affectation souple (soft assignment): au lieu de dire «P4 est dans le groupe 1», on dit «P4 est dans le groupe 1 avec probabilité $0,9$». Le dénominateur de (12.12) est la densité du mélange en $x_i$, et c'est lui qui apparaît dans la log-vraisemblance (12.11).

L'algorithme EM

Les formules (12.13) sont celles qu'on écrirait si l'on connaissait les $z_i$ — proportion, moyenne et variance empiriques de chaque groupe —, où chaque point compte dans chaque groupe pour la fraction $r_{ik}$ de lui-même. D'où viennent-elles? Si les $z_i$ étaient connus, la log-vraisemblance «complète» serait ${\sum }_i{\sum }_{k}1[z_i=k](\ln {\pi }_k+\ln \mathcal{N}(x_i;{\mu }_k,{\sigma }_k^2))$. Ne connaissant pas , l'étape E remplace l'indicatrice par son espérance conditionnelle sachant , qui est précisément ; l'étape M maximise la fonction obtenue,

les $r_{ik}$ étant fixés. Le nom de l'algorithme vient de ces deux étapes: une espérance, puis une maximisation. Il a été formulé dans cette généralité par Dempster, Laird et Rubin en 1977, pour tous les modèles à variables cachées.

Démonstration. Dans (12.14), $\ln \mathcal{N}(x_i;{\mu }_k,{\sigma }_k^2)=-\frac{1}{2}\ln (2\pi {\sigma }_k^2)-(x_i-{\mu }_k{)}^2/(2{\sigma }_k^2)$, et n'apparaît que dans les termes de la composante . Les termes de qui dépendent de se réduisent donc à

et maximiser $Q$ en ${\mu }_k$ revient à minimiser la somme pondérée $S(\mu )={\sum }_i{r}_{ik}(x_i-\mu {)}^2$, puisque . Posons et écrivons , comme dans la démonstration du théorème 12.1:

Le terme du milieu est nul, car ${\sum }_i{r}_{ik}(x_i-{\bar{x}}_k)={\sum }_i{r}_{ik}{x}_i-N_k{\bar{x}}_k=0$. Il reste , minimal en et seulement là, puisque .

La démonstration ne dépend ni de ${\sigma }_k^2$ ni de ${\pi }_k$: la moyenne optimale est la même quelle que soit la variance, ce qui permet de mettre à jour la moyenne en premier, puis la variance autour d'elle. L'exercice 12.5 établit de même les formules de ${\pi }_k$ et de ${\sigma }_k^2$. En dimension , la même démonstration, avec des normes au lieu de valeurs absolues, donne , et la covariance devient .

Nous admettons ce résultat, mais la raison tient en quelques lignes, et elle mérite d'être vue. Pour des paramètres quelconques $\mathbf{\theta }$ et des responsabilités $r_{ik}$ calculées avec ${\mathbf{\theta }}^{(t)}$, on écrit chaque terme de (12.11) comme le logarithme d'une moyenne pondérée par les $r_{ik}$:

par l'inégalité de Jensen, le logarithme étant concave. En sommant sur $i$, on obtient $\mathcal{L}(\mathbf{\theta })\ge Q(\mathbf{\theta })+H$, où $H=-{\sum }_{i,k}{r}_{ik}\ln r_{ik}$ ne dépend pas de . Pour , le rapport sous le logarithme vaut, par (12.12), la densité du mélange en quel que soit , et l'inégalité est une égalité. L'étape M choisit qui maximise , donc . EM monte ainsi le long d'une suite de minorants de , chacun collé à au point courant. Ce que nous n'établissons pas, c'est il converge: en général vers un point stationnaire de , qui peut être un maximum local, comme les -moyennes s'arrêtent sur un minimum local de l'inertie. Les mêmes remèdes s'appliquent: plusieurs départs, ou un départ fourni par les -moyennes.

Les $K$-moyennes, limite d'EM

Les deux algorithmes se ressemblent trop pour que ce soit un hasard: une étape qui attribue les points aux composantes, une étape qui recalcule les centres comme des moyennes. Le lien est exact.

Démonstration. Avec ces paramètres, ${\pi }_k\mathcal{N}({\mathbf{x}}_i;{\mathbf{\mu }}_k,\epsilon I)=\frac{1}{K}(2\pi \epsilon {)}^{-d/2}{\mathrm{e}}^{-\Vert {\mathbf{x}}_i-{\mathbf{\mu }}_k{\Vert }^2/(2\epsilon )}$, et le facteur commun disparaît de (12.12):

Si ${\mathbf{\mu }}_k$ est strictement la plus proche, chaque différence au numérateur est strictement positive, chaque exponentielle tend vers $0$ quand $\epsilon \to 0$, et $r_{ik}\to 1$; comme les responsabilités somment à 1, les autres tendent vers $0$. La responsabilité devient l'indicatrice du centre le plus proche, et la moyenne pondérée du théorème 12.5 devient le barycentre des points affectés.

Sur l'exemple 12.5, avec ${\mu }_1=2$ et ${\mu }_2=6$, le point $x=3,9$, un peu plus proche de ${\mu }_1$, reçoit la responsabilité pour , puis pour , pour et pour : l'affectation souple durcit à mesure que la variance diminue. Et les -moyennes lancées sur les mêmes sept valeurs depuis les centres et s'arrêtent sur les groupes et — le point , à égale distance, va au centre d'indice le plus petit — de centres et , tout près des moyennes et du mélange ajusté.

Le mélange gaussien est donc une généralisation des $K$-moyennes, et il en corrige plusieurs limites: chaque groupe a sa propre variance (sa propre forme ellipsoïdale en dimension $d$, par sa matrice ${\Sigma }_k$), son propre poids, et chaque point reçoit une probabilité d'appartenance plutôt qu'un verdict. Il le paie: davantage de paramètres — $Kd$ moyennes, $Kd(d+1)/2$ coefficients de covariance et $K-1$ poids libres —, donc davantage de données nécessaires, un calcul plus lourd, des maxima locaux et la dégénérescence de l'encadré précédent. Comme il définit une vraisemblance, il offre en revanche un critère pour choisir qui manquait aux -moyennes: comparer les modèles par leur vraisemblance pénalisée par le nombre de paramètres, ce que fait par exemple le ().

Synthèse

• Le partitionnement regroupe des exemples sans étiquettes selon un critère de ressemblance; il renvoie toujours des groupes, et leur sens doit être établi par ailleurs. Les $K$-moyennes minimisent l'inertie $W={\sum }_k{\sum }_{i\in C_k}\Vert {\mathbf{x}}_i-{\mathbf{\mu }}_k{\Vert }^2$, et le barycentre est le meilleur centre d'un groupe (théorème 12.1).

Exercices

Vous pouvez afficher le corrigé directement sous chaque énoncé après avoir cherché la solution.

Références

• James, G., Witten, D., Hastie, T. et Tibshirani, R., An Introduction to Statistical Learning, Springer, 2ᵉ éd., 2021, chap. 12 (apprentissage non supervisé: $K$-moyennes et partitionnement hiérarchique).
• Hastie, T., Tibshirani, R. et Friedman, J., The Elements of Statistical Learning, Springer, 2ᵉ éd., 2009, chap. 14 (partitionnement) et section 8.5 (algorithme EM).
• Bishop, C. M., Pattern Recognition and Machine Learning, Springer, 2006, chap. 9 (mélanges et EM, et les $K$-moyennes comme cas limite).
• Lloyd, S. P., «Least squares quantization in PCM», IEEE Transactions on Information Theory, 1982.
• Arthur, D. et Vassilvitskii, S., «k-means++: The advantages of careful seeding», Proceedings of the Eighteenth Annual ACM-SIAM Symposium on Discrete Algorithms, 2007.
• Dempster, A. P., Laird, N. M. et Rubin, D. B., «Maximum likelihood from incomplete data via the EM algorithm», Journal of the Royal Statistical Society, Series B, 1977.