Avant propos

Après avoir manipulĂ© en dĂ©tail les mĂ©thodes d'analyse factorielle que sont l'ACP et l'AFC, nous allons maintenant nous intĂ©resser Ă  la classification automatique. LĂ  oĂ¹ les analyses graphiques nous permettent parfois de dĂ©tecter de probables structures de groupes, nous n'avons pour l'instant pas utilisĂ© d'outils pour dĂ©finir automatiquement ces groupes (appelĂ©s aussi clusters). Les mĂ©thodes destinĂ©es Ă  accomplir cette tĂ¢che sont appelĂ©s des algorithmes de clustering et nous allons en Ă©tudier quelques uns, Ă  commencer par la Classification Ascendante HiĂ©rarchique (CAH) dans l'exercice 4 ci-dessous.

Exercice 4

Dans cet exercice, nous utilisons la base de données très ((beaucoup ?) trop ?) classique iris, présente par défaut dans la version de base de R. Ce dataset regroupe un total de 150 iris de 3 espèces différentes (50 chacune) observées sur 4 variables relatives à leur pétales et sépales. La vraie classification étant connue d'avance, nous allons pouvoir tester nos méthodes de clustering, pour ensuite comparer les groupes obtenus automatiquement avec les vraies espèces d'iris. Les algorithmes de clustering sont des méthodes dites d'apprentissage non-supervisé, un sous-ensemble des méthodes d'apprentissage automatique (machine learning en anglais). Comme habituellement, importons packages et données, et visualisons les données brutes.

#rm(list = ls())
library(tidyverse)
library(gridExtra)
library(FactoMineR)

raw_db = iris %>% as_tibble

raw_db
## # A tibble: 150 x 5
##    Sepal.Length Sepal.Width Petal.Length Petal.Width Species
##           <dbl>       <dbl>        <dbl>       <dbl> <fct>  
##  1          5.1         3.5          1.4         0.2 setosa 
##  2          4.9         3            1.4         0.2 setosa 
##  3          4.7         3.2          1.3         0.2 setosa 
##  4          4.6         3.1          1.5         0.2 setosa 
##  5          5           3.6          1.4         0.2 setosa 
##  6          5.4         3.9          1.7         0.4 setosa 
##  7          4.6         3.4          1.4         0.3 setosa 
##  8          5           3.4          1.5         0.2 setosa 
##  9          4.4         2.9          1.4         0.2 setosa 
## 10          4.9         3.1          1.5         0.1 setosa 
## # ... with 140 more rows

Pour pouvoir faire l'analyse dans des conditions 'rĂ©elles', nous devons appliquer les mĂ©thodes de clustering sur les 4 premières variables seulement, qui constitueront notre db d'Ă©tude. La variable Species va Ăªtre stockĂ©e Ă  part et nous servira seulement au moment d'Ă©valuer la qualitĂ© des rĂ©sultats. Jetons un premier regard Ă©galement Ă  un summary de la base de donnĂ©es ainsi obtenu.

Species = raw_db %>% pull(Species)
db = raw_db %>% select(- Species)
summary(db)
##   Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
##  Min.   :4.300   Min.   :2.000   Min.   :1.000   Min.   :0.100  
##  1st Qu.:5.100   1st Qu.:2.800   1st Qu.:1.600   1st Qu.:0.300  
##  Median :5.800   Median :3.000   Median :4.350   Median :1.300  
##  Mean   :5.843   Mean   :3.057   Mean   :3.758   Mean   :1.199  
##  3rd Qu.:6.400   3rd Qu.:3.300   3rd Qu.:5.100   3rd Qu.:1.800  
##  Max.   :7.900   Max.   :4.400   Max.   :6.900   Max.   :2.500

Pour une analyse plus détaillée et détecter d'éventuelles structures de groupes, nous allons commencer par appliquer une première ACP, avec laquelle nous sommes maintenant familliers.

res_pca = PCA(db, ncp = 4, graph = F)
 
vp = res_pca$eig %>% data.frame %>% rownames_to_column(var = 'Composante') %>% as_tibble

ggplot(vp) + geom_col(aes(x = Composante, y = percentage.of.variance)) +
             geom_point(aes(x = Composante, y = cumulative.percentage.of.variance)) +
             ggtitle("Pourcentage de variance expliquée cumulée") + 
             geom_hline(aes(yintercept = 90), color = 'red')

Nous voyons que deux axes suffisent à expliquer plus de 95% de la variance des données, nous allons donc retenir 2 composantes principales et, après avoir vérifié la contribution et qualité de représentation des individus et variables, les projeter sur ce premier plan factoriel.

summary(res_pca)
## 
## Call:
## PCA(X = db, ncp = 4, graph = F) 
## 
## 
## Eigenvalues
##                        Dim.1   Dim.2   Dim.3   Dim.4
## Variance               2.918   0.914   0.147   0.021
## % of var.             72.962  22.851   3.669   0.518
## Cumulative % of var.  72.962  95.813  99.482 100.000
## 
## Individuals (the 10 first)
##                  Dist    Dim.1    ctr   cos2    Dim.2    ctr   cos2    Dim.3
## 1            |  2.319 | -2.265  1.172  0.954 |  0.480  0.168  0.043 | -0.128
## 2            |  2.202 | -2.081  0.989  0.893 | -0.674  0.331  0.094 | -0.235
## 3            |  2.389 | -2.364  1.277  0.979 | -0.342  0.085  0.020 |  0.044
## 4            |  2.378 | -2.299  1.208  0.935 | -0.597  0.260  0.063 |  0.091
## 5            |  2.476 | -2.390  1.305  0.932 |  0.647  0.305  0.068 |  0.016
## 6            |  2.555 | -2.076  0.984  0.660 |  1.489  1.617  0.340 |  0.027
## 7            |  2.468 | -2.444  1.364  0.981 |  0.048  0.002  0.000 |  0.335
## 8            |  2.246 | -2.233  1.139  0.988 |  0.223  0.036  0.010 | -0.089
## 9            |  2.592 | -2.335  1.245  0.812 | -1.115  0.907  0.185 |  0.145
## 10           |  2.249 | -2.184  1.090  0.943 | -0.469  0.160  0.043 | -0.254
##                 ctr   cos2  
## 1             0.074  0.003 |
## 2             0.250  0.011 |
## 3             0.009  0.000 |
## 4             0.038  0.001 |
## 5             0.001  0.000 |
## 6             0.003  0.000 |
## 7             0.511  0.018 |
## 8             0.036  0.002 |
## 9             0.096  0.003 |
## 10            0.293  0.013 |
## 
## Variables
##                 Dim.1    ctr   cos2    Dim.2    ctr   cos2    Dim.3    ctr
## Sepal.Length |  0.890 27.151  0.792 |  0.361 14.244  0.130 | -0.276 51.778
## Sepal.Width  | -0.460  7.255  0.212 |  0.883 85.247  0.779 |  0.094  5.972
## Petal.Length |  0.992 33.688  0.983 |  0.023  0.060  0.001 |  0.054  2.020
## Petal.Width  |  0.965 31.906  0.931 |  0.064  0.448  0.004 |  0.243 40.230
##                cos2  
## Sepal.Length  0.076 |
## Sepal.Width   0.009 |
## Petal.Length  0.003 |
## Petal.Width   0.059 |
plot(res_pca, axes = c(1,2), choix = 'var')

On trouve donc une structure oĂ¹ la variable Sepa.width contribue fortement au deuxième axe, et est globalement orthogonale aux trois autres variables, qui elles contribuent majoritairement au premier axe. Voyons maintenant si la projection des individus dans le plan factoriel fait Ă©merger une structure de groupe par rapport aux caractĂ©ristiques que nous venons de dĂ©tailler.

plot(res_pca, axes = c(1,2), choix = 'ind')

Si le nombre exact de groupes apparents n'est pas parfaitement Ă©tabli (probablement 2 ou 3 Ă  première vue), un sĂ©paration apparait tout de mĂªme Ă©vidente et il est lĂ©gitime d'ensuite chercher Ă  appliquer un algorithme de clustering. Gardons Ă  l'esprit que, de mĂªme que pour un prĂ©cĂ©dent exercice, nous avons Ă©cartĂ© deux axes factoriels, qui auraient pu nous aider Ă  discriminer les individus d'une façon diffĂ©rente et peut Ăªtre informative.

Mettons de côté mainenant ces résultats pour nous concentrer sur la mise en place de la CAH, à l'aide de la fonction hclust() du package FactoMineR. Il a été vu en cours que la CAH permet de construire un dendrogramme, un graphique qui résume la structure des groupes qui sont successivement construits à chaque étape de regroupement. Sur ce grapique, la taille des segments horizontaux sont proportionels à la valeur de dissimilarité inter-groupe entre les deux groupes rassemblés. Puisque l'on cherche à regrouper ensemble des individus qui se ressemblent et séparer les autres, il faut éviter justement de rassembler des sous groupes avec un fort score de dissimilarité inter-groupe. C'est pourquoi nous cherchons à stopper à 'couper' le dendrogramme avant que de grandes 'branches' horizontales n'apparaissent. Voyons tout d'abord un exemple d'utilisation de la fonction hclust(), après avoir calculé les distances euclidiennes entre les individus (centrés réduits avec scale()) à l'aide de la fonction dist():

db %>% scale() %>% 
       dist(method = 'euclidean') %>% 
       hclust() %>% 
       plot()

Tout ceci est très intéressant, mais n'oublions pas qu'il existe divers critères de regroupement avec leur caractéristiques propres, qui sont proposés comme argument de la fonction hclust. Affichons donc ci-dessous les résultats en fonction du critère choisi:

dist = db %>% scale() %>% dist(method = 'euclidean')

par(mfrow=c(2,2))
hclust(dist, method="ward.D2") %>% plot(main = "Dendrogramme avec la methode de Ward")
hclust(dist, method = "single") %>% plot(main = "Dendrogramme avec le lien simple")
hclust(dist, method = "complete") %>% plot(main = "Dendrogramme avec le lien complet")
hclust(dist, method = "average") %>% plot(main = "Dendrogramme avec le lien moyen")

par(mfrow=c(1,1))

On cherche plutôt à avoir des groupes bien homogènes et à éviter les petits groupes isolés. Le critère de Ward semble le plus adapté dans cette situation. Reprenons son graph en proposant différents choix de 'coupûres' valables pour notre classification. Dans cet exemple il parait pertinent de selectionner 2 ou 3 groupes (voir 5 si on veut beaucoup de groupes).

dendro_ward = hclust(dist, method="ward.D2")
plot(dendro_ward)
abline(h=10, col="red")
abline(h=15, col="blue")
abline(h=5, col="green")

Nous pouvons aussi utiliser des critères de \(R^2\) ou de PseudoF pour dĂ©cider du nombre de groupes optimal. Ceux ci peuvent Ăªtre calculĂ©s en utilisant la fonction ci-dessous (vous pouvez simplement copier-coller le code pour utiliser ensuite directement la fonction).

R2.PseudoF = function(donnees,dendro, graph=T, cut=20){
  n=nrow(donnees)
  R2 = numeric(n)

  Iinter = 0
  g = apply(donnees, 2, mean)
  I = (n-1)/n*sum(diag(var(donnees)))

  for(i in 2:n){
      class = cutree(dendro,k=i)
      ncl = unique(class)
      d = numeric(length(ncl)) 
      nb = integer(length(ncl))
      for(j in ncl){
          nb[j] = sum(class==j)
          if(nb[j]>1){
              m = apply(donnees[class==j,], 2, mean)
          } else {
              m = donnees[class==j,]
          }
          d[j] = sum((m-g)^2)   
      }
      Iinter = (1/n)*sum(nb*d)
      R2[i] = Iinter/I
  }

  ncl = 2:(n-1)
  PseudoF = R2[ncl]/(1-R2[ncl])*(n-ncl)/(ncl-1)

  if(graph==T){
  par(mfrow=c(1,2))
  plot(1:20,R2[1:cut], type = 'b', xlab = "Nombre de classes", ylab = "Rsquare")
  title("Indice du R2")
  plot(ncl[1:20],PseudoF[1:cut], type = 'b', xlab = "Nombre de classes", ylab = "PseudoF")
  title("Indice du PseudoF")
  }
  resultat = list(R2,PseudoF)
  names(resultat)=c("Rsquare","PseudoF")
  resultat
}
nb_cluster = R2.PseudoF(db,dendro_ward)

Au vu de ces critères, on a envie de selectionner 3 groupes, car il y a une vraie rupture de pente à partir de 4 groupes et donc une grosse dégradation de l'inertie inter-classe/inertie intra-classe. Sachant qu'il existe effectivement 3 espèces différentes d'iris, ces résultats sont plutôt rassurants. Maintenant que l'on a selectionné un nombre de groupes qui semble adéquat, nous pouvons comparer les partitions obtenues avec les vraies données Species que l'on avait stocké. La fonction cutree() permet de selectionner la partion obtenue par notre CAH en donnant comme argument le nombre \(k\) de groupes retenu (ici \(k=3\) donc).

res_cah = cutree(tree = dendro_ward, k = 3)
res_cah
##   [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##  [38] 1 1 1 1 2 1 1 1 1 1 1 1 1 3 3 3 2 3 2 3 2 3 2 2 3 2 3 2 3 2 2 2 2 3 3 3 3
##  [75] 3 3 3 3 3 2 2 2 2 3 2 3 3 2 2 2 2 3 2 2 2 2 2 3 2 2 3 3 3 3 3 3 2 3 3 3 3
## [112] 3 3 3 3 3 3 3 3 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
## [149] 3 3

On voit donc que le résultat est donné sous forme d'un indice de groupe (1,2 ou 3) auquel est attribué chaque individu. Gardez à l'esprit que ce numéro de groupe ne correspond à rien de précis, car dans le cas de la classification non-supervisée (clustering), il n'y a justement pas d'étiquette (label) définie et l'interprétation de ce qui définit le groupe reste à déterminer. Dans ce cas précis cependant, nous savons que les iris appartiennent en réalité à 3 espèces différentes, et nous allons donc pouvoir évaluer la qualité de notre procédure de clustering. Pour celà, nous croisons pour chaque individu la valeur de Species avec celui de son groupe obtenu via la CAH à l'aide de la fonction table(). Ce qui permet de donner un tableau d'adéquation global entre nos résultats et la vérité:

table(Species,res_cah)
##             res_cah
## Species       1  2  3
##   setosa     49  1  0
##   versicolor  0 27 23
##   virginica   0  2 48

On voit donc que la CAH semble avoir correctement identifiĂ© les iris appartenant aux espèces setosa et virginica. Cependant les rĂ©sultats sont plus mitigĂ©s concernant les versicolors, qui semblent pouvoir Ăªtre regrouper aussi avec les virginica. Voyons graphiquement les raisons de ces difficultĂ©s de sĂ©paration. Tout d'abord, nous affichons ci-dessous les individus dans le plan factoriel, coloriĂ©s selon les vrais espèces sous jacentes. Nous utilisons ici le package factoextra, que certains d'entre vous ont dĂ©jĂ  experimentĂ©, qui offre des graphiques plutĂ´t Ă©lĂ©gants.

library(factoextra)
fviz_pca_ind(res_pca,
             habillage = Species, # Color by groups
             palette = c("#00AFBB","#E7B800","#FC4E07"),# Show off with your favorite hexadecimal color
             addEllipses = TRUE # Add concentration ellipses around the groups
             )

Nous pouvons voir ici clairement la cause des résultats de la CAH, avec notamment toute une partie des iris versicolor qui sont très proches du centre du groupe des virginica. Il est donc logique de les retrouver confondus avec ceux-ci, comme nous le confirme le graphique affichant les groupes obtenus via la CAH:

fviz_pca_ind(res_pca, habillage = as.factor(res_cah), palette = c("#1E3AB1", "#BCF634", "#F64E34"), addEllipses = TRUE)

Notez bien que ces groupes ont Ă©tĂ© obtenus de manière automatique, simplement en calculant des scores de similaritĂ©s, ici Ă  partir des distances, entre les individus. Essayez de remarquer quels individus ont Ă©tĂ© mal classifiĂ©s et pourquoi c'est gĂ©nĂ©ralement logique. Il existe de nombreuses approches diffĂ©rentes du problème de la classification non-supervisĂ©e. La CAH, et plus gĂ©nĂ©ralement les approches hiĂ©rarchiques, ne sont qu'un exemple parmis d'autres et nous verrons par la suite une des mĂ©thodes les plus simples et les plus populaires pour cette tĂ¢che: l'algorithme des k-means.

Exercice 5

Pour continuer sur la thĂ©matique des mĂ©thodes de clustering, nous allons introduire, toujours sur l'exemple du jeu de donnĂ©es iris, l'algorithme des kmeans. Cette procĂ©dure itĂ©rative consiste Ă  dĂ©finir (souvent alĂ©atoirement) des centres de groupe, puis alternativement rassembler les individus les plus proches du centre, et recalculer les nouvelles coordonnĂ©es des centres. Après quelques itĂ©rations, la procĂ©dure converge en gĂ©nĂ©ral vers une situation stable, avec des clusters que l'on espère pertinents. Pour comparaison, nous reprenons le mĂªme dĂ©but d'analyse que dans l'exercice prĂ©cĂ©dent, en important packages et donnĂ©es, et effectuant une ACP qui nous servira Ă  projeter les individus dans le premier plan factoriel.

rm(list = ls())
library(factoextra)
library(tidyverse)
library(gridExtra)
library(FactoMineR)

raw_db = iris %>% as_tibble()

Species = raw_db %>% pull(Species)

db = raw_db %>% select(- Species) %>% scale()

res_pca = PCA(db, graph = F)

Avant de lancer l'algorithme, il faut encore une fois pouvoir définir le nombre de classes attendues dans notre jeu de données. Dans le cas des kmeans, qui est un algorithme généralement peu couteux en calcul, il est classique de lancer l'algorithme pour différentes valeurs de nombre de groupes \(k\), et de regarder pour lequel on obtient les 'meilleurs' clusters. Il existe évidemment plusieurs manières de définir 'meilleurs' ici, et donc plusieurs critères, comme la méthode de la silouette, ou par un calcul d'inertie. En utilisant encore une fois le package factoextra, il est possible d'utiliser la fonction fviz_nb_clust qui nous calcul ces critères pour différentes valeures de \(k\), afin de décider graphiquement de la valeur à retenir. Par exemple, la méthode de la silouette ci-dessous nous indique de retenir 2 groupes:

fviz_nbclust(db, kmeans, method = 'silhouette',  k.max = 10)

Si l'on utilise maintenant un critère d'inertie, il faut repérer le 'coude' et cette fois, le choix est moins clair. On utilisera également l'algorithme des kmeans avec 3 groupes dans un second temps pour comparaison.

fviz_nbclust(db, kmeans, method = 'wss',  k.max = 10)

Comme nous l'avons dit, l'initialisation des kmeans est en gĂ©nĂ©rale alĂ©atoire, et les groupes obtenus sont souvent très sensibles aux positions initiales des centres. Ainsi, il est usuel de lancer l'algorithme un grand nombre de fois (ici 50 fois), et de retenir le meilleur rĂ©sultat. De plus, pour Ăªtre sur d'avoir des rĂ©sultats reproductibles d'une fois sur l'autre (ou d'un individu Ă  l'autre), il est intĂ©ressant d'utiliser la fonction set.seed() qui 'fixe' l'alĂ©atoire du script avec laquelle elle est executĂ©e. En utilisant le mĂªme chiffre en entrĂ©e de cette fonction, nous devrions tous obtenir les mĂªmes rĂ©sultats.

set.seed(42)
res_kmeans = kmeans(x = db, centers = 2, nstart = 50)
res_kmeans
## K-means clustering with 2 clusters of sizes 50, 100
## 
## Cluster means:
##   Sepal.Length Sepal.Width Petal.Length Petal.Width
## 1   -1.0111914   0.8504137    -1.300630  -1.2507035
## 2    0.5055957  -0.4252069     0.650315   0.6253518
## 
## Clustering vector:
##   [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##  [38] 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [75] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [112] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [149] 2 2
## 
## Within cluster sum of squares by cluster:
## [1]  47.35062 173.52867
##  (between_SS / total_SS =  62.9 %)
## 
## Available components:
## 
## [1] "cluster"      "centers"      "totss"        "withinss"     "tot.withinss"
## [6] "betweenss"    "size"         "iter"         "ifault"

Le résumé des résultats de l'algorithme nous donne de nombreuses informations. En particulier, parmi les informations pertinentes pour nous, il est intéressant d'extraire le vecteur du cluster d'appartenance pour chaque individu avec cluster. Mais également la taille de chaque groupe avec size, ou encore les coordonnées des centres de chaque groupe avec centers.

clusters = res_kmeans$cluster %>% as.factor()
clusters 
##   [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##  [38] 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [75] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [112] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
## [149] 2 2
## Levels: 1 2
res_kmeans$size
## [1]  50 100
res_kmeans$centers
##   Sepal.Length Sepal.Width Petal.Length Petal.Width
## 1   -1.0111914   0.8504137    -1.300630  -1.2507035
## 2    0.5055957  -0.4252069     0.650315   0.6253518

Il est ainsi possible de comparer les groupes obtenus avec les vraies espèces des iris:

table(Species, clusters)
##             clusters
## Species       1  2
##   setosa     50  0
##   versicolor  0 50
##   virginica   0 50

Et de visualiser les groupes obtenus en projetant sur le premier plan factoriel et coloriant les individus par groupe. On voit que la distinction est assez naturelle entre deux groupes, bien sépararés, qui sont tout à fait satisfaisants et logiques.

fviz_pca_ind(res_pca,
             habillage = clusters,
             palette = c("#00AFBB","#E7B800","#FC4E07"),
             addEllipses = TRUE)

A prĂ©sent, nous pouvons effectuer la mĂªme analyse en choisissant un nombre de groupes \(k=3\).

set.seed(42)
res_kmeans_3 = kmeans(x = db, centers = 3, nstart = 50)
res_kmeans_3
## K-means clustering with 3 clusters of sizes 50, 53, 47
## 
## Cluster means:
##   Sepal.Length Sepal.Width Petal.Length Petal.Width
## 1  -1.01119138  0.85041372   -1.3006301  -1.2507035
## 2  -0.05005221 -0.88042696    0.3465767   0.2805873
## 3   1.13217737  0.08812645    0.9928284   1.0141287
## 
## Clustering vector:
##   [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##  [38] 1 1 1 1 1 1 1 1 1 1 1 1 1 3 3 3 2 2 2 3 2 2 2 2 2 2 2 2 3 2 2 2 2 3 2 2 2
##  [75] 2 3 3 3 2 2 2 2 2 2 2 3 3 2 2 2 2 2 2 2 2 2 2 2 2 2 3 2 3 3 3 3 2 3 3 3 3
## [112] 3 3 2 2 3 3 3 3 2 3 2 3 2 3 3 2 3 3 3 3 3 3 2 2 3 3 3 2 3 3 3 2 3 3 3 2 3
## [149] 3 2
## 
## Within cluster sum of squares by cluster:
## [1] 47.35062 44.08754 47.45019
##  (between_SS / total_SS =  76.7 %)
## 
## Available components:
## 
## [1] "cluster"      "centers"      "totss"        "withinss"     "tot.withinss"
## [6] "betweenss"    "size"         "iter"         "ifault"

Encore une fois, nous pouvons extraire les informations utiles de ces nouveaux résultats.

clusters_3 = res_kmeans_3$cluster %>% as.factor()
clusters_3
##   [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##  [38] 1 1 1 1 1 1 1 1 1 1 1 1 1 3 3 3 2 2 2 3 2 2 2 2 2 2 2 2 3 2 2 2 2 3 2 2 2
##  [75] 2 3 3 3 2 2 2 2 2 2 2 3 3 2 2 2 2 2 2 2 2 2 2 2 2 2 3 2 3 3 3 3 2 3 3 3 3
## [112] 3 3 2 2 3 3 3 3 2 3 2 3 2 3 3 2 3 3 3 3 3 3 2 2 3 3 3 2 3 3 3 2 3 3 3 2 3
## [149] 3 2
## Levels: 1 2 3
res_kmeans_3$size
## [1] 50 53 47
res_kmeans_3$centers
##   Sepal.Length Sepal.Width Petal.Length Petal.Width
## 1  -1.01119138  0.85041372   -1.3006301  -1.2507035
## 2  -0.05005221 -0.88042696    0.3465767   0.2805873
## 3   1.13217737  0.08812645    0.9928284   1.0141287

Et comparer les groupes obtenus avec les vraies espèces. On peut remarquer qu'il existe encore deux classes entremelées qui compliquent les résultats, mais contrairement à la CAH, nous n'avons plus un groupe qui absorbe un autre.

table(Species, clusters_3)
##             clusters_3
## Species       1  2  3
##   setosa     50  0  0
##   versicolor  0 39 11
##   virginica   0 14 36

Nous pouvons comparer graphiquement les partitions obtenues (en haut) avec les vraies espèces (en bas) pour comprendre comment l'algorithme procède. Notamment ici, les groupes sont clairement définis comme l'ensemble des individus les plus proches du centre des clusters en distance euclidienne. Comme vous pouvez le voir, les résultats diffèrent légèrement de ceux que l'on avait pu obtenir avec la CAH, et de manière générale, les différents algorithmes de clustering ont chacun leur propre 'philosophie', qui les rendront plus ou moins adaptés en fonction des situations, et qui exprime une tendance à regrouper les individus selon des critères qui leur sont propres.

gg1  = fviz_pca_ind(res_pca,
                    habillage = clusters_3,
                    palette = c("#00AFBB","#E7B800","#FC4E07"),
                    addEllipses = TRUE)
gg2 = fviz_pca_ind(res_pca, 
                   habillage = Species,
                   palette = c("#1E3AB1", "#BCF634", "#F64E34"),
                   addEllipses = TRUE)

grid.arrange(gg1, gg2, ncol = 1, nrow = 2)