7.2 Analyse en composantes principales

L’Analyse en Composantes Principales ou ACP (Principal Component Analysis ou PCA en anglais) est une méthode d’ordination de base qu’il est indispensable de connaître et de comprendre. La plupart des autres techniques d’ordination plus sophistiquées peuvent être vues comme des variantes de l’ACP. L’ACP sera utilisable lorsque :

  • Des relations linéaires sont suspectées entre les variables (si elles ne sont pas linéaires, penser à transformer les données auparavant pour les linéariser).

  • Ces relations conduisent à une répartition des individus (le nuage de points) qui forme une structure que l’on cherchera à interpréter.

  • Pour visualiser cette structure, les données sont simplifiées (réduites) de N variables à n (n < N et n = 2 ou 3 généralement). La représentation sous forme d’un nuage de points s’appelle une carte.

  • La réduction des dimensions se fait avec une perte minimale d’information au sens de la variance des données.

À vous de jouer !
h5p

7.2.1 Indiens diabétiques

Partons d’un exemple concret pour voir comment cela se passe en pratique. L’ACP est facilitée dans SciViews::R, grâce à sa section "explore". Donc, au début de votre script ou de votre document Quarto ou R Markdown, vous utiliserez l’instruction (ajoutant ou non que vous voulez travailler en français avec lang = "fr") :

SciViews::R("explore", lang = "fr")

Nous allons aborder un problème de santé humaine. Les Indiens Pimas sont des Amérindiens originaires du nord du Mexique qui sont connus pour compter le plus haut pourcentage d’obèses et de diabétiques de toutes les ethnies. Ils ont fait l’objet de plusieurs études scientifiques d’autant plus que les Pimas en Arizona développent principalement cette obésité et ce diabète, alors que les Pimas mexicains les ont plus rarement. Il est supposé que leur mode de vie différent aux États-Unis pourrait en être la raison. Voici un jeu de données qui permet d’explorer un peu cette question :

pima <- read("PimaIndiansDiabetes2", package = "mlbench")
pima
#      pregnant glucose pressure triceps insulin mass pedigree age diabetes
#   1:        6     148       72      35      NA 33.6    0.627  50      pos
#   2:        1      85       66      29      NA 26.6    0.351  31      neg
#   3:        8     183       64      NA      NA 23.3    0.672  32      pos
#   4:        1      89       66      23      94 28.1    0.167  21      neg
#   5:        0     137       40      35     168 43.1    2.288  33      pos
#  ---                                                                     
# 764:       10     101       76      48     180 32.9    0.171  63      neg
# 765:        2     122       70      27      NA 36.8    0.340  27      neg
# 766:        5     121       72      23     112 26.2    0.245  30      neg
# 767:        1     126       60      NA      NA 30.1    0.349  47      pos
# 768:        1      93       70      31      NA 30.4    0.315  23      neg

Ce jeu de données contient des valeurs manquantes. Le graphique suivant permet de visualiser l’importance des “dégâts” à ce niveau :

naniar::vis_miss(pima)

Moins de 10% des données sont manquantes, et c’est principalement dans les variables insulin et triceps. Si nous souhaitons un tableau sans variables manquantes, nous pouvons décider d’éliminer des lignes et/ou des colonnes (variables), mais ici nous souhaitons garder toutes les variables et réduisons donc uniquement le nombre de lignes avec sdrop_na().

pima <- sdrop_na(pima)
pima
#      pregnant glucose pressure triceps insulin mass pedigree age diabetes
#   1:        1      89       66      23      94 28.1    0.167  21      neg
#   2:        0     137       40      35     168 43.1    2.288  33      pos
#   3:        3      78       50      32      88 31.0    0.248  26      pos
#   4:        2     197       70      45     543 30.5    0.158  53      pos
#   5:        1     189       60      23     846 30.1    0.398  59      pos
#  ---                                                                     
# 388:        0     181       88      44     510 43.3    0.222  26      pos
# 389:        1     128       88      39     110 36.5    1.057  37      pos
# 390:        2      88       58      26      16 28.4    0.766  22      neg
# 391:       10     101       76      48     180 32.9    0.171  63      neg
# 392:        5     121       72      23     112 26.2    0.245  30      neg
Notre tableau est presque amputé de la moitié de ses données, mais il nous reste tout de même encore 392 cas sur 768, soit assez pour notre analyse.

Avant de nous lancer dans une ACP, nous devons décrire les données, repérer les variables quantitatives d’intérêt, et synthétiser les corrélations linéaires (coefficients de corrélation de Pearson) entre ces variables.

skimr::skim(pima)
Tableau 1.1: Data summary
Name pima
Number of rows 392
Number of columns 9
Key NULL
_______________________
Column type frequency:
factor 1
numeric 8
________________________
Group variables None

Variable type: factor

skim_variable n_missing complete_rate ordered n_unique top_counts
diabetes 0 1 FALSE 2 neg: 262, pos: 130

Variable type: numeric

skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
pregnant 0 1 3.30 3.21 0.00 1.00 2.00 5.00 17.00 ▇▂▂▁▁
glucose 0 1 122.63 30.86 56.00 99.00 119.00 143.00 198.00 ▂▇▇▃▂
pressure 0 1 70.66 12.50 24.00 62.00 70.00 78.00 110.00 ▁▂▇▆▁
triceps 0 1 29.15 10.52 7.00 21.00 29.00 37.00 63.00 ▅▇▇▃▁
insulin 0 1 156.06 118.84 14.00 76.75 125.50 190.00 846.00 ▇▂▁▁▁
mass 0 1 33.09 7.03 18.20 28.40 33.20 37.10 67.10 ▃▇▃▁▁
pedigree 0 1 0.52 0.35 0.09 0.27 0.45 0.69 2.42 ▇▃▁▁▁
age 0 1 30.86 10.20 21.00 23.00 27.00 36.00 81.00 ▇▂▁▁▁

Nous avons une variable facteur diabetes à exclure de l’analyse. De plus, la variable pregnant est une variable numérique discrète (nombre d’enfants portés). Nous l’éliminerons aussi de notre analyse.

La fonction correlation() du package {SciViews} nous permet d’inspecter les corrélations entre les variables choisies (donc toutes à l’exception de pregnant et diabetes qui ne sont pas quantitatives continues) :

pima_cor <- correlation(pima[, 2:8])
tabularise(pima_cor, digits = 2)
Tableau 7.1:

Matrice de coefficients de corrélation de Pearson r

glucose

pressure

triceps

insulin

mass

pedigree

age

glucose

1.000

0.2100

0.199

0.5812

0.2095

0.140

0.3436

pressure

0.210

1.0000

0.233

0.0985

0.3044

-0.016

0.3000

triceps

0.199

0.2326

1.000

0.1822

0.6644

0.160

0.1678

insulin

0.581

0.0985

0.182

1.0000

0.2264

0.136

0.2171

mass

0.210

0.3044

0.664

0.2264

1.0000

0.159

0.0698

pedigree

0.140

-0.0160

0.160

0.1359

0.1588

1.000

0.0850

age

0.344

0.3000

0.168

0.2171

0.0698

0.085

1.0000

plot(pima_cor)

Quelques corrélations positives d’intensités moyennes se dégagent ici, notamment entre mass (masse corporelle) et triceps (épaisseur du pli cutané au niveau du triceps), ainsi qu’entre glucose (taux de glucose dans le sang) et insulin (taux d’insuline dans le sang) et glucose et age. Par contre, le pedigree (variable qui quantifie la susceptibilité au diabète en fonction de la parenté) semble peu corrélé avec les autres variables.

Nous utiliserons la fonction pca() qui prend un argument data = et une formule du type ~ var1 + var2 + .... + varn, ou plus simplement, directement un tableau contenant uniquement les variables à analyser comme argument unique. Comme les différentes variables sont mesurées dans des unités différentes, nous devons les standardiser (écart type ramené à un pour toutes). Ceci est réalisé par la fonction pca() en lui indiquant scale = TRUE. Donc :

pima_pca <- pca(data = pima, ~ glucose + pressure + triceps + insulin + mass +
  pedigree + age, scale = TRUE)

Ou alors, nous sélectionnons les variables d’intérêt avec sselect() et appliquons pca() directement sur ce tableau, ce qui donnera le même résultat.

pima %>.%
  sselect(., glucose:age) %>.%
  pca(., scale = TRUE) %->%
  pima_pca

Le nuage de points dans l’espace initial à sept dimensions a été centré (origine ramenée au centre de gravité du nuage de points = moyenne des variables) par l’ACP. Ensuite une rotation des axes a été réalisée pour orienter son plus grand axe selon un premier axe principal 1 ou PC1 . Ensuite PC2 est construit orthogonal au premier et dans la seconde direction de plus grande variabilité du nuage de points, et ainsi de suite pour les autres axes. Ainsi les axes PC1, PC2, PC3, … représentent une part de variance de plus en plus faible par rapport à la variance totale du jeu de données. Ceci est présenté dans le résumé :

summary(pima_pca)
# Importance of components (eigenvalues):
#                          PC1   PC2   PC3   PC4    PC5   PC6    PC7
# Variance               2.412 1.288 1.074 0.878 0.6389 0.399 0.3098
# Proportion of Variance 0.345 0.184 0.153 0.126 0.0913 0.057 0.0443
# Cumulative Proportion  0.345 0.529 0.682 0.807 0.8988 0.956 1.0000
# 
# Loadings (eigenvectors, rotation matrix):
#          PC1    PC2    PC3    PC4    PC5    PC6    PC7   
# glucose   0.441 -0.455        -0.198         0.736       
# pressure  0.329  0.101 -0.613  0.206  0.654        -0.171
# triceps   0.439  0.488               -0.367        -0.644
# insulin   0.402 -0.418  0.263 -0.388  0.123 -0.642 -0.129
# mass      0.446  0.506        -0.181                0.711
# pedigree  0.198         0.625  0.711  0.251              
# age       0.325 -0.337 -0.384  0.471 -0.592 -0.168  0.179

Le premier tableau Importance of components (eigenvalues): montre la part de variance présentée sur chacun des sept axes de l’ACP (PC1, PC2, …, PC7). Le fait qu’il s’agit de valeurs propres (eigenvalues en anglais) apparaîtra plus clair lorsque vous aurez lu les explications détaillées plus bas. Ces parts de variance s’additionnent pour donner la variance totale du nuage de points dans les sept dimensions (propriété d’additivité des variances). Pour faciliter la lecture, la Proportion de Variance en % est reprise également, ainsi que les proportions cumulées. Ainsi, les deux premiers axes de l’ACP capturent ici 53% de la variance totale. Et il faudrait considérer les cinq premiers axes pour capturer presque 90% de la variance totale. Cependant, les trois premiers axes cumulent tout de même plus des 2/3 de la variance. Nous pouvons restreindre notre analyse à ces trois axes-là.

Le second tableau Loadings (eigenvectors, rotation matrix): est la matrice de transformation des coordonnées initiales sur les lignes en coordonnées PC1 à PC7 en colonnes. C’est les coefficients multiplicateurs appliqués aux variables initiales pour réaliser la rotation du système d’axe recherché. Nous pouvons donc y lire l’importante des variables initiales sur les axes de l’ACP, car plus ces coefficients sont proches de un en valeur absolue, plus la variable initiale prendra de l’importance sur les axes de la projection. Cela permet de déterminer ce que représente chaque axe principal. Par exemple, l’axe PC3 contraste essentiellement pressure et pedigree.

À vous de jouer !
h5p

Le graphique des éboulis sert à visualiser la “chute” de la variance d’un axe principal à l’autre, et aide à choisir le nombre d’axes à conserver (espace à dimensions réduites avec perte minimale d’information). Deux variantes en diagramme en barres verticales chart$screeplot() ou chart$scree() ou sous forme d’une ligne brisée chart$altscree() sont disponibles :

chart$scree(pima_pca, fill = "cornsilk")

chart$altscree(pima_pca)

La diminution est importante entre le premier et le second axe, mais plus progressive ensuite. Ceci traduit une structure plus complexe dans les données qui ne se réduit pas facilement à un très petit nombre d’axes. Nous pouvons visualiser le premier plan principal constitué par PC1 et PC2, tout en gardant à l’esprit que seulement 53% de la variance totale y est capturée. Donc, nous pouvons nous attendre à des pertes d’information non négligeables dans ce plan. Cela signifie que certains aspects n’y sont pas (correctement) représentés. Nous verrons qu’il est porteur, toutefois, d’information utile.

À vous de jouer !
h5p

Deux types de représentations peuvent être réalisées à partir d’ici : la représentation dans l’espace des variables, et la représentation dans l’espace des individus. Ces deux représentations sont complémentaires et s’analysent conjointement. L’espace des variables représente les axes initiaux projetés (comme l’ombre portée d’un cadran solaire) dans le plan choisi de l’ACP. Il se réalise à l’aide de chart$loadings(). Par exemple pour PC1 et PC2 nous indiquons choices = c(1, 2) (ou rien du tout, puisque ce sont les valeurs par défaut) :

chart$loadings(pima_pca, choices = c(1, 2))

Ce graphique s’interprète comme suit :

  • Plus la norme (longueur) du vecteur qui représente une variable est grande et se rapproche d’une longueur de un (matérialisée par le cercle gris), plus la variable est bien représentée dans le plan choisi. On évitera d’interpréter ici les variables qui ont des normes petites, comme pedigree ou pressure.

  • Des vecteurs qui pointent dans la même direction représentent des variables directement corrélées entre elles. C’est le cas de glucose, insulin et age d’une part, et par ailleurs aussi de mass et triceps.

  • Des vecteurs qui pointent en directions opposées représentent des variables inversement proportionnelles. Il n’y en a pas ici.

  • Des vecteurs orthogonaux représentent des variables non corrélées entre elles. ainsi le groupe glucose/insulin/age n’est pas corrélé avec le groupe mass/triceps.

  • Les PCs sont orientés en fonction des vecteurs qui représentent les variables initiales sur ce graphique. Puisque le vecteur qui représente la variable mass pointe vers le haut et la droite, les plus gros Indiens seront représentés aussi en haut à droite dans l’espace des individus que nous représenterons un peu plus tard. À l’inverse, les individus moins lourds seront représentés en bas à gauche. D’autre part, la variable age pointe vers le bas à droite. Donc, nous aurons un gradient des individus en fonction de leur âge, avec les plus jeunes à l’opposé en haut à gauche et les plus vieux en bas à droite dans le graphique des individus (voir plus loin). Et cela est corrélé également avec les variables insulin et glucose dans le sang.

Cela donne déjà une vision synthétique des différentes corrélations entre les variables. Naturellement, on peut très bien choisir d’autres axes, pour peu qu’ils représentent une part de variance relativement importante. Par exemple, ici, nous pouvons représenter le plan constitué par PC1 et PC3, puisque nous avons décidé de retenir les trois premiers axes :

chart$loadings(pima_pca, choices = c(1, 3))

Nous voyons que pedigree et pressure (inversement proportionnels) sont bien mieux représentés le long de PC3. Ici l’axe PC3 est plus facile à orienter : en haut les pedigrees élevés et les pressions artérielles basses, et en bas le contraire. Nous avons déjà lu cette information dans le tableau des vecteurs propres obtenu avec summary().

Le graphique entre PC2 et PC3 complète l’analyse, mais n’apportant rien de plus, il peut être typiquement éliminé de votre rapport.

chart$loadings(pima_pca, choices = c(2, 3))

À vous de jouer !
h5p

La seconde représentation se fait dans l’espace des individus. Ici, nous allons projeter les points relatifs à chaque individu dans le plan de l’ACP choisi. Cela se réalise à l’aide de chart$scores() (l’aspect ratio est le rapport hauteur/largeur qui peut être adapté) :

chart$scores(pima_pca, choices = c(1, 2), aspect.ratio = 3/5)

Ce graphique est peu lisible, tel quel. Généralement, nous représentons d’autres informations utiles sous forme de labels et/ou de couleurs différentes. Nous pouvons ainsi contraster les individus qui ont le diabète de ceux qui ne l’ont pas sur ce graphique et ajouter des ellipses de confiance à 95% autour des deux groupes pour aider à mieux les cerner à l’aide de stat_ellipse() :

chart$scores(pima_pca, choices = c(1, 2),
  labels = pima$diabetes) +
  stat_ellipse()

Ce graphique est nettement plus intéressant. Il s’interprète comme suit :

  • Nous savons que les individus plus âgés et ayant plus de glucose et d’insuline dans le sang sont dans le bas à droite du graphique. Or le groupe des diabétiques, s’il ne se détache pas complètement, tend à s’étaler plus dans cette région.

  • À l’inverse, le groupe des non-diabétiques s’étale vers la gauche, c’est-à-dire dans une région reprenant les individus les plus jeunes et aussi les moins gros.

On comprend mieux maintenant comment nous avons d’abord analysé le graphique des variables pour observer les corrélations entre ces variables, mais aussi, pour l’utiliser comme une boussole qui va orienter notre carte représentée par le graphique dans l’espace des individus. Une fois ce premier travail réalisé, la projection des individus dans le second graphique révèle une quantité complémentaire d’information très utile. Voyons ce que donne le graphique entre PC1 et PC3 (analyse du troisième axe) :

chart$scores(pima_pca, choices = c(1, 3),
  labels = pima$diabetes) +
  stat_ellipse()

Ici, la séparation se fait essentiellement sur l’axe horizontal (PC1). Donc, les différentes de pedigree (élevé dans le haut du graphique) et de tension artérielle (élevée dans le bas du graphique) semblent être moins liées au diabète. Le graphique PC3 versus PC2 peut aussi être réalisé, mais il n’apporte rien de plus dans le cas présent.

chart$scores(pima_pca, choices = c(2, 3),
  labels = pima$diabetes) +
  stat_ellipse()

Étant donné que les deux graphiques (variables et individus) s’interprètent conjointement, nous pourrions être tentés de les superposer. C’est réalisable, et ce graphique particulier s’appelle un biplot. Mais se pose alors un problème : celui de mettre à l’échelle les deux représentations pour qu’elles soient cohérentes entre elles. Ceci n’est pas facile et différentes représentations coexistent. L’argument scale = de la fonction chart$biplot() permet d’utiliser différentes mises à l’échelle. Enfin, ce type de graphique tend à être souvent bien trop encombré. Il est donc plus difficile à lire que les deux graphiques des variables et individus séparés. Voici ce que cela donne pour notre jeu de données exemple :

chart$biplot(pima_pca)

Bien moins lisible, en effet !

7.2.2 Biométrie d’oursin

Analysons à présent un autre jeu de données qui nous montrera l’importance de la transformation (linéarisation), du choix de réduire ou non (argument scale = dans la fonction pca()). Ces données nous permettront aussi de découvrir ce qu’est un effet saturant et comment s’en débarrasser. Il s’agit de la biométrie effectuée sur deux populations de l’oursin violet Paracentrotus lividus, une en élevage et une autre provenant du milieu naturel. Nous avons abondamment utilisé ce jeu de données en SDD I dans la section visualisation. Nous le connaissons bien, mais reprenons certains éléments essentiels ici…

urchin <- read("urchin_bio", package = "data.io", lang = "FR")
urchin
#        origin diameter1 diameter2 height buoyant_weight weight solid_parts
#   1: Pêcherie       9.9      10.2    5.0             NA 0.5215      0.4777
#   2: Pêcherie      10.5      10.6    5.7             NA 0.6418      0.5891
#   3: Pêcherie      10.8      10.8    5.2             NA 0.7336      0.6770
#   4: Pêcherie       9.6       9.3    4.6             NA 0.3697      0.3438
#   5: Pêcherie      10.4      10.7    4.8             NA 0.6097      0.5587
#  ---                                                                      
# 417:  Culture      16.7      17.2    8.5         0.5674 2.4300      2.2900
# 418:  Culture      16.5      16.5    7.9         0.5472 2.3200      2.1800
# 419:  Culture      16.8      16.7    8.2         0.4864 2.2200      2.1300
# 420:  Culture      17.3      17.2    8.5         0.4864 2.5200      2.3400
# 421:  Culture      17.0      16.6    7.9         0.4357 2.0500      1.9800
#      integuments dry_integuments digestive_tract dry_digestive_tract gonads
#   1:      0.3658              NA          0.0525              0.0079 0.0000
#   2:      0.4447              NA          0.0482              0.0090 0.0000
#   3:      0.5326              NA          0.0758              0.0134 0.0000
#   4:      0.2661              NA          0.0442              0.0064 0.0000
#   5:      0.4058              NA          0.0743              0.0117 0.0000
#  ---                                                                       
# 417:      1.8400            1.02          0.1661              0.0229 0.0215
# 418:      1.8000            1.01          0.0977              0.0147 0.0253
# 419:      1.6300            0.88          0.1704              0.0208 0.0154
# 420:      1.7200            0.89          0.1444              0.0167 0.0237
# 421:      1.4300            0.83          0.1462              0.0212 0.0266
#      dry_gonads skeleton lantern   test spines maturity  sex
#   1:     0.0000   0.1793  0.0211 0.0587 0.0995        0 <NA>
#   2:     0.0000   0.1880  0.0205 0.0622 0.1053        0 <NA>
#   3:     0.0000   0.2354  0.0254 0.0836 0.1263        0 <NA>
#   4:     0.0000   0.0630  0.0167 0.0180 0.0283        0 <NA>
#   5:     0.0000       NA      NA     NA     NA        0 <NA>
#  ---                                                        
# 417:     0.0034   0.9046  0.0750 0.3399 0.4896        0 <NA>
# 418:     0.0051   0.8965  0.0908 0.3189 0.4868        0 <NA>
# 419:     0.0020   0.7714  0.0877 0.2961 0.3876        0 <NA>
# 420:     0.0032   0.7938  0.0772 0.3077 0.4090        0 <NA>
# 421:     0.0051   0.7421  0.0723 0.2689 0.4009        0 <NA>

Ici aussi nous avons des valeurs manquantes :

naniar::vis_miss(urchin)

Ces valeurs manquantes sont rassemblées essentiellement dans les variables buoyant_weight, dry_integuments, les mesures relatives au squelette (skeleton, lantern, test et spines), et surtout au niveau de sex (impossible de déterminer le sexe des individus les plus jeunes). Si nous éliminons purement et simplement les lignes qui ont au moins une valeur manquante, nous perdons tous les individus jeunes, et c’est dommage. Nous allons donc d’abord éliminer la variable sex, ainsi que les quatre variables liées au squelette. Dans un second temps, nous appliquerons sdrop_na() sur ce qui reste :

urchin %>.%
  sselect(., -(skeleton:spines), -sex) %>.%
  sdrop_na(.) ->
  urchin2
urchin2
#        origin diameter1 diameter2 height buoyant_weight weight solid_parts
#   1: Pêcherie      16.7      16.8    8.4         0.5881   2.58        2.04
#   2: Pêcherie      19.9      20.0    9.2         1.0952   4.26        3.66
#   3: Pêcherie      19.9      19.2    8.5         0.6287   2.93        2.43
#   4: Pêcherie      19.3      19.8   10.2         0.7808   3.71        3.09
#   5: Pêcherie      18.8      20.0    9.3         0.7605   3.59        2.99
#  ---                                                                      
# 315:  Culture      16.7      17.2    8.5         0.5674   2.43        2.29
# 316:  Culture      16.5      16.5    7.9         0.5472   2.32        2.18
# 317:  Culture      16.8      16.7    8.2         0.4864   2.22        2.13
# 318:  Culture      17.3      17.2    8.5         0.4864   2.52        2.34
# 319:  Culture      17.0      16.6    7.9         0.4357   2.05        1.98
#      integuments dry_integuments digestive_tract dry_digestive_tract gonads
#   1:        1.77            1.06          0.0644              0.0159 0.0054
#   2:        3.11            1.95          0.2168              0.0489 0.0449
#   3:        2.13            1.17          0.1153              0.0231 0.0000
#   4:        2.52            1.45          0.1788              0.0381 0.0183
#   5:        2.40            1.41          0.2820              0.0520 0.0223
#  ---                                                                       
# 315:        1.84            1.02          0.1661              0.0229 0.0215
# 316:        1.80            1.01          0.0977              0.0147 0.0253
# 317:        1.63            0.88          0.1704              0.0208 0.0154
# 318:        1.72            0.89          0.1444              0.0167 0.0237
# 319:        1.43            0.83          0.1462              0.0212 0.0266
#      dry_gonads maturity
#   1:     0.0009        0
#   2:     0.0086        0
#   3:     0.0000        0
#   4:     0.0017        0
#   5:     0.0038        0
#  ---                    
# 315:     0.0034        0
# 316:     0.0051        0
# 317:     0.0020        0
# 318:     0.0032        0
# 319:     0.0051        0

Il nous reste 319 lignes des 421 initiales. Nous n’avons perdu qu’un quart des données, tout en nous privant seulement de quatre variables quantitatives liées au squelette (sex étant une variable qualitative, elle ne peut de toute façon pas être introduite dans l’analyse, mais elle aurait pu servir pour colorer les individus).

skimr::skim(urchin2)
Tableau 7.2: Data summary
Name urchin2
Number of rows 319
Number of columns 14
Key NULL
_______________________
Column type frequency:
factor 1
numeric 13
________________________
Group variables None

Variable type: factor

skim_variable n_missing complete_rate ordered n_unique top_counts
origin 0 1 FALSE 2 Cul: 188, Pêc: 131

Variable type: numeric

skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
diameter1 0 1 32.78 11.71 14.60 23.25 31.50 39.65 65.60 ▇▇▅▃▁
diameter2 0 1 32.71 11.67 15.00 23.45 31.60 39.60 65.60 ▇▇▆▃▁
height 0 1 16.78 6.25 7.30 11.10 16.20 21.50 32.20 ▇▆▆▅▂
buoyant_weight 0 1 4.27 3.84 0.31 1.35 3.18 5.67 17.73 ▇▃▁▁▁
weight 0 1 21.80 21.37 1.61 6.08 15.25 28.14 100.51 ▇▂▁▁▁
solid_parts 0 1 16.52 15.27 1.46 4.96 11.73 21.69 73.14 ▇▃▁▁▁
integuments 0 1 12.32 10.64 1.09 4.00 9.40 16.08 47.22 ▇▃▂▁▁
dry_integuments 0 1 7.16 6.30 0.58 2.22 5.42 9.42 28.80 ▇▃▂▁▁
digestive_tract 0 1 1.90 2.03 0.03 0.45 1.21 2.54 10.37 ▇▂▁▁▁
dry_digestive_tract 0 1 0.23 0.21 0.01 0.07 0.17 0.31 1.02 ▇▃▂▁▁
gonads 0 1 1.72 2.65 0.00 0.10 0.63 2.20 15.93 ▇▁▁▁▁
dry_gonads 0 1 0.51 0.82 0.00 0.03 0.17 0.64 5.00 ▇▁▁▁▁
maturity 0 1 0.37 0.71 0.00 0.00 0.00 0.00 2.00 ▇▁▁▁▂

Nous avons 12 variables quantitatives continues. Notez la distribution très asymétrique et similaire (voir colonne hist) de toutes ces variables. Les variables origin et maturity ne pourront pas être utilisées, mais seront éventuellement utiles pour colorer les points dans nos graphiques. Qu’en est-il des corrélations entre les 12 variables ?

urchin2_cor <- correlation(urchin2[, 2:13])
tabularise(urchin2_cor)
Tableau 3.2:

Matrice de coefficients de corrélation de Pearson r

diameter1

diameter2

height

buoyant_weight

weight

solid_parts

integuments

dry_integuments

digestive_tract

dry_digestive_tract

gonads

dry_gonads

diameter1

1.000

0.998

0.976

0.952

0.956

0.956

0.968

0.956

0.910

0.932

0.798

0.788

diameter2

0.998

1.000

0.974

0.952

0.957

0.957

0.969

0.956

0.912

0.933

0.798

0.788

height

0.976

0.974

1.000

0.928

0.924

0.926

0.944

0.934

0.875

0.909

0.759

0.747

buoyant_weight

0.952

0.952

0.928

1.000

0.986

0.992

0.994

0.999

0.921

0.939

0.882

0.868

weight

0.956

0.957

0.924

0.986

1.000

0.994

0.991

0.986

0.952

0.955

0.882

0.875

solid_parts

0.956

0.957

0.926

0.992

0.994

1.000

0.994

0.991

0.947

0.955

0.907

0.897

integuments

0.968

0.969

0.944

0.994

0.991

0.994

1.000

0.996

0.934

0.946

0.867

0.854

dry_integuments

0.956

0.956

0.934

0.999

0.986

0.991

0.996

1.000

0.920

0.938

0.874

0.860

digestive_tract

0.910

0.912

0.875

0.921

0.952

0.947

0.934

0.920

1.000

0.980

0.806

0.807

dry_digestive_tract

0.932

0.933

0.909

0.939

0.955

0.955

0.946

0.938

0.980

1.000

0.816

0.818

gonads

0.798

0.798

0.759

0.882

0.882

0.907

0.867

0.874

0.806

0.816

1.000

0.991

dry_gonads

0.788

0.788

0.747

0.868

0.875

0.897

0.854

0.860

0.807

0.818

0.991

1.000

# ou knitr::kable(urchin2_cor, digits = 2)
plot(urchin2_cor)

Toutes les corrélations sont positives, et certaines sont très élevées. Cela indique que plusieurs variables sont (pratiquement complètement) redondantes, par exemple, diameter1 et diameter2. Un effet principal semble dominer.

Si nous refaisons quelques graphiques, nous nous rappelons que les relations ne sont pas linéaires, par exemple, entre diameter1 et weight :

chart(data = urchin2, weight ~ diameter1) +
  geom_point()

Ce type de relation, dite allométrique, se linéarise très bien en effectuant une transformation double-logarithme, comme nous pouvons le constater sur le graphique suivant :

chart(data = urchin2, log(weight) ~ log(diameter1)) +
  geom_point()

Il est crucial de bien nettoyer son jeu de données avant une ACP, et aussi, de vérifier que les relations sont linéaires. Sinon, il faut transformer les données de manière appropriée. Rappelez-vous que l’ACP s’intéresse aux corrélations linéaires entre vos variables.

Attention toutefois à la transformation logarithmique appliquée sur des données qui peuvent contenir des zéros (par exemple, gonads ou dry_gonads). Dans ce cas, la transformation logarithme(x + 1) réalisée avec la fonction log1p() est plus indiquée. Nous allons ici transformer toutes les variables en log(x + 1). C’est assez fastidieux à faire avec smutate(), mais nous pouvons l’utiliser directement sur le tableau entier réduit aux variables quantitatives continues seules lors de l’appel à pca() comme suit :

urchin2 %>.%
  sselect(., -origin, -maturity) %>.% # Élimine les variables non quantitatives
  log1p(.) %>.% # Transforme toutes les autres en log(x + 1)
  pca(., scale = TRUE) -> # Effectue l'ACP après standardisation
  urchin2_pca

Nous avons standardisé les données puisqu’elles sont mesurées dans des unités différentes (longueurs en mm, masses en g). Voici ce que donne notre ACP :

summary(urchin2_pca)
# Importance of components (eigenvalues):
#                           PC1    PC2    PC3     PC4     PC5     PC6     PC7
# Variance               11.219 0.5010 0.1813 0.03862 0.02601 0.01657 0.00931
# Proportion of Variance  0.935 0.0418 0.0151 0.00322 0.00217 0.00138 0.00078
# Cumulative Proportion   0.935 0.9767 0.9918 0.99503 0.99720 0.99858 0.99936
#                            PC8     PC9    PC10    PC11    PC12
# Variance               0.00336 0.00210 0.00108 0.00082 0.00034
# Proportion of Variance 0.00028 0.00017 0.00009 0.00007 0.00003
# Cumulative Proportion  0.99964 0.99981 0.99990 0.99997 1.00000
# 
# Loadings (eigenvectors, rotation matrix):
#                     PC1    PC2    PC3    PC4    PC5    PC6    PC7    PC8   
# diameter1            0.295 -0.162        -0.441  0.237 -0.174  0.177  0.242
# diameter2            0.295 -0.166        -0.449  0.249 -0.178  0.157  0.124
# height               0.291 -0.218  0.154 -0.106 -0.902                     
# buoyant_weight       0.296         0.120  0.509  0.124                0.530
# weight               0.296 -0.149                                    -0.145
# solid_parts          0.297 -0.106  0.127                0.234        -0.594
# integuments          0.296 -0.159  0.157  0.153  0.125               -0.396
# dry_integuments      0.296 -0.115  0.160  0.455  0.114                0.105
# digestive_tract      0.288        -0.571 -0.187         0.485 -0.519  0.201
# dry_digestive_tract  0.283        -0.702  0.217        -0.278  0.513 -0.161
# gonads               0.271  0.568  0.226                0.575  0.430  0.143
# dry_gonads           0.259  0.697        -0.104        -0.465 -0.453 -0.111
#                     PC9    PC10   PC11   PC12  
# diameter1           -0.706                     
# diameter2            0.688  0.263              
# height                                         
# buoyant_weight                     0.265 -0.504
# weight               0.116 -0.912              
# solid_parts                 0.216  0.638       
# integuments                 0.148 -0.702 -0.371
# dry_integuments             0.128 -0.133  0.774
# digestive_tract                                
# dry_digestive_tract                            
# gonads                                         
# dry_gonads

Nous observons plus de 93% de la variance représentée sur le premier axe. Cela parait parfait ! Voici le graphique des éboulis :

chart$scree(urchin2_pca)

Ne vous réjouissez pas trop vite. Nous avons ici un effet saturant lié au fait que toutes les variables sont positivement corrélées entre elles. Cet effet est évident. Ici, c’est la taille des oursins. Nous allons conclure que plus un oursin est gros, plus ses dimensions et ses masses sont importantes. C’est trivial et d’un intérêt très limité, avouons-le.

Puisque l’ACP optimise la variance sur le premier axe, un effet saturant aura tendance à occulter d’autres effets intéressants. Nous pouvons nous en débarrasser en identifiant une des variables représentant le mieux cet effet, et en calculant les ratios entre toutes les autres variables et celle-là. Ainsi, nous passons de quantification de la taille sur toutes les variables à des ratios qui quantifient beaucoup mieux des effets de forme plus subtils.

Notez aussi les valeurs relativement faibles, mais homogènes de toutes les variables sur l’axe PC1 dans le tableau des vecteurs propres, avec des valeurs comprises entre 0,26 et 0,30. Le graphique des variables n’est pas beau du tout dans le premier plan de l’ACP, même si un effet différent relatif aux gonades apparaît tout de même sur l’axe PC2, il ne compte que pour 4,2% de la variance totale :

chart$loadings(urchin2_pca)
# Warning: ggrepel: 8 unlabeled data points (too many overlaps). Consider
# increasing max.overlaps

Recommençons tout de suite l’analyse en éliminant l’effet saturant. Nous pourrons considérer comme référence de la taille, par exemple, la masse immergée (buoyant weight) connue comme étant une grandeur pouvant être mesurée très précisément. Elle fait partie des variables les mieux corrélées sur l’axe PC1, représentant ainsi très bien cet effet saturant que nous voulons éliminer. Voici notre calcul :

urchin2 %>.%
  sselect(., -origin, -maturity, -buoyant_weight) %>.% # Élimination des variables inutiles
  (. / urchin2$buoyant_weight) %>.% # Division par buoyant_weight
  log1p(.) -> # Transformation log(x + 1)
  urchin3
head(urchin3)
#    diameter1 diameter2   height   weight solid_parts integuments
# 1:  3.380877  3.386644 2.726760 1.683990    1.497119    1.388714
# 2:  2.953357  2.958109 2.240741 1.587131    1.468302    1.345385
# 3:  3.485925  3.451231 2.675524 1.733496    1.582091    1.478861
# 4:  3.247200  3.271795 2.643585 1.749467    1.600897    1.441601
# 5:  3.247291  3.306831 2.582396 1.744070    1.595668    1.424509
# 6:  3.000850  2.973973 2.254397 1.690963    1.594759    1.406850
#    dry_integuments digestive_tract dry_digestive_tract      gonads  dry_gonads
# 1:        1.030481       0.1039141          0.02667720 0.009140213 0.001529182
# 2:        1.022630       0.1806157          0.04368131 0.040178983 0.007821777
# 3:        1.051165       0.1683868          0.03608357 0.000000000 0.000000000
# 4:        1.049797       0.2061975          0.04764294 0.023167059 0.002174887
# 5:        1.048737       0.3154008          0.06613980 0.028901124 0.004984271
# 6:        1.034084       0.3464496          0.05538973 0.016565715 0.003104904

Refaisons notre ACP sur urchin3 ainsi calculé :

urchin3_pca <- pca(urchin3, scale = TRUE)
summary(urchin3_pca)
# Importance of components (eigenvalues):
#                          PC1   PC2   PC3    PC4    PC5    PC6     PC7     PC8
# Variance               4.687 3.353 1.273 0.9666 0.3668 0.1724 0.10547 0.04761
# Proportion of Variance 0.426 0.305 0.116 0.0879 0.0333 0.0157 0.00959 0.00433
# Cumulative Proportion  0.426 0.731 0.847 0.9345 0.9678 0.9835 0.99308 0.99741
#                            PC9    PC10    PC11
# Variance               0.01943 0.00834 0.00068
# Proportion of Variance 0.00177 0.00076 0.00006
# Cumulative Proportion  0.99918 0.99994 1.00000
# 
# Loadings (eigenvectors, rotation matrix):
#                     PC1    PC2    PC3    PC4    PC5    PC6    PC7    PC8   
# diameter1           -0.425         0.145 -0.297               -0.101  0.101
# diameter2           -0.425 -0.101  0.143 -0.296               -0.103       
# height              -0.427         0.131 -0.300                       0.141
# weight               0.189 -0.428         0.216  0.497 -0.672 -0.145       
# solid_parts                -0.495  0.254 -0.117         0.335 -0.245 -0.662
# integuments         -0.142 -0.463  0.152  0.283  0.165  0.345  0.667  0.266
# dry_integuments     -0.259 -0.242  0.173  0.533 -0.669 -0.161 -0.277       
# digestive_tract      0.214 -0.370 -0.448 -0.102         0.388 -0.462  0.468
# dry_digestive_tract        -0.360 -0.485 -0.401 -0.429 -0.338  0.382 -0.155
# gonads               0.374         0.438 -0.257 -0.223                     
# dry_gonads           0.371         0.440 -0.269 -0.162 -0.122         0.446
#                     PC9    PC10   PC11  
# diameter1                   0.400  0.714
# diameter2                   0.425 -0.700
# height               0.197 -0.790       
# weight               0.109              
# solid_parts         -0.223              
# integuments                             
# dry_integuments                         
# digestive_tract      0.148              
# dry_digestive_tract                     
# gonads               0.725  0.130       
# dry_gonads          -0.589
chart$scree(urchin3_pca)

Maintenant que l’effet saturant est éliminé, la répartition des variances sur les axes se fait mieux. À partir du résumé de l’ACP, nous pouvons voir que l’axe PC1 contraste les diamètres avec les gonades, l’axe PC2 représente les masses somatiques (dans l’ordre inverse), et l’axe PC3 contraste de manière intéressante les masses du tube digestif avec celles des gonades (le tout en ratios sur la masse immergée, ne l’oublions pas). Les deux premiers axes reprennent 73% de la variance, mais il semble qu’un effet intéressant se marque également sur PC3 avec 85% de la variance totale sur les trois premiers axes.

Tout ceci est également visible sur les graphiques dans l’espace des variables (plans PC1 - PC2 et PC2 - PC3 représentés ici).

chart$loadings(urchin3_pca, choices = c(1, 2))

chart$loadings(urchin3_pca, choices = c(2, 3))
# Warning: ggrepel: 1 unlabeled data points (too many overlaps). Consider
# increasing max.overlaps

Enfin, dans l’espace des individus, avec l’origine reprise en couleur, nous observons ceci dans le premier plan de l’ACP :

chart$scores(urchin3_pca, choices = c(1, 2),
  col = urchin2$origin, labels = urchin2$maturity, aspect.ratio = 3/5) +
  theme(legend.position = "right") +
  stat_ellipse()
# Warning: The following aesthetics were dropped during statistical transformation: label
# ℹ This can happen when ggplot fails to infer the correct grouping structure in
#   the data.
# ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
#   variable into a factor?

Et pour le plan PC2 - PC3 :

chart$scores(urchin3_pca, choices = c(2, 3),
  col = urchin2$origin, labels = urchin2$maturity, aspect.ratio = 3/5) +
  theme(legend.position = "right") +
  stat_ellipse()
# Warning: The following aesthetics were dropped during statistical transformation: label
# ℹ This can happen when ggplot fails to infer the correct grouping structure in
#   the data.
# ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
#   variable into a factor?

Vous devriez pouvoir interpréter ces résultats par vous-même maintenant.

7.2.3 Visualisation de données quantitatives

7.2.3.1 Deux dimensions

Le nuage de points est le graphe idéal pour visualiser la distribution des données bivariées pour deux variables quantitatives. Il permet de visualiser également une association entre deux variables. Il permet aussi de visualiser comment deux ou plusieurs groupes peuvent être séparés en fonction de ces deux variables.

chart(data = pima, glucose ~ insulin %col=% diabetes) +
  geom_point()

7.2.3.2 Trois dimensions

Le nuage de points en pseudo-3D est l’équivalent pour visualiser trois variables quantitatives simultanément. Il est nécessaire de rendre l’effet de la troisième dimension (perspective, variation de taille des objets …). La possibilité de faire tourner l’objet 3D virtuel est indispensable pour concrétiser l’effet 3D et pour le visionner sous différents angles.

Le package {rgl} permet de réaliser ce genre de graphique 3D interactif (que vous pouvez faire tourner dans l’orientation que vous voulez à la souris). Dans un document R Markdown, il faut absolument configurer {knitr} qui est responsable de l’inclusion des graphiques R dans le document comme ceci avant toute chose :

library(rgl)
# Warning in rgl.init(initValue, onlyNULL): RGL: unable to open X11 display
# Warning: 'rgl.init' failed, running with 'rgl.useNULL = TRUE'.
knitr::knit_hooks$set(webgl = hook_webgl)

Ensuite, nous pouvons réaliser le graphique 3D {rgl} :

rgl::plot3d(pima$insulin, pima$glucose, pima$triceps,
  col = as.integer(pima$diabetes))

Enfin, nous l’introduisons dans le document à l’aide de ceci :

rgl::rglwidget(width = 800, height = 800)

À noter que cela ne fonctionne que pour les documents Quarto ou R Markdown au format HTML.

7.2.3.3 Plus de trois dimensions

Déjà à trois dimensions la visualisation devient délicate, mais au-delà, cela devient pratiquement mission impossible. La matrice de nuages de points peut rendre service ici, mais dans certaines limites (tous les angles de vue ne sont pas accessibles).

GGally::ggscatmat(pima, 2:6, color = "diabetes")
# Registered S3 method overwritten by 'GGally':
#   method from   
#   +.gg   ggplot2
# Warning: The dot-dot notation (`..scaled..`) was deprecated in ggplot2 3.4.0.
# ℹ Please use `after_stat(scaled)` instead.
# ℹ The deprecated feature was likely used in the GGally package.
#   Please report the issue at <https://github.com/ggobi/ggally/issues>.
# This warning is displayed once every 8 hours.
# Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
# generated.

Nous voyons qu’ici nous atteignons les limites des possibilités. C’est pour cela que, pour des données multivariées comportant beaucoup de variables quantitatives, les techniques de réduction des dimensions comme l’ACP sont indispensables.

7.2.4 ACP : mécanisme

Nous allons partir d’un exemple presque trivial pour illustrer le principe de l’ACP. Comment réduire un tableau bivarié en une représentation des individus en une seule dimension (classement sur une droite) avec perte minimale d’information ? Par exemple, en partant de ces données fictives :

Voici une représentation graphique 2D de ces données :

Si nous réduisons à une seule dimension en laissant tomber une des deux variables, voici ce que cela donne (ici on ne garde que Var1, donc, on projette les points sur l’axe des abscisses).

Au final, nous avons ordonné nos individus en une dimension comme suit :

C’est une mauvaise solution, car il y a trop de perte d’information. Regardez l’écart entre 7 et 9 sur le graphique en deux dimensions et dans celui à une dimension : les points sont trop près. Comparez sur les deux graphiques les distances 7 - 9 avec 9 - 8 et 1 - 2 versus 1 - 3. Tout cela est très mal représenté en une seule dimension.

Une autre solution serait de projeter le long de la droite de “tendance générale”, c’est-à-dire le long de l’axe de plus grand allongement du nuage de points.

Cela donne ceci en une seule dimension :

C’est une bien meilleure solution, car la perte d’information est ici minimale. Regardez à nouveau la distance entre 7 et 9 sur le graphique initial à deux dimensions et sur le nouveau graphique réduit à une dimension : c’est mieux qu’avant. Comparez aussi les distances respectives entre les paires 7 - 9 et 9 - 8, ainsi que 1 - 2 par rapport à 1 - 3. Tout cela est bien mieux représenté à présent.

L’ACP effectue précisément la projection que nous venons d’imaginer.

  • La droite de projection est appelée composante principale 1. La composante principale 1 présente la plus grande variabilité possible sur un seul axe.

  • Ensuite on calcule la composante 2 comme étant orthogonale (c.-à-d., perpendiculaire) à PC1 et présentant la plus grande variabilité non encore capturée par la composante 1.

  • Le mécanisme revient à projeter les points sur des axes orientés différemment dans l’espace à N dimensions (pour N variables initiales). En effet, mathématiquement, ce mécanisme se généralise facilement à trois, puis à N dimensions.

7.2.5 Calcul matriciel ACP

La rotation optimale des axes vers les PC1 à PCN se résout par un calcul matriciel. Nous allons maintenant le détailler. Mais auparavant, nous devons nous rafraîchir l’esprit concernant quelques notions.

  • Multiplication matricielle : \(\begin{pmatrix} 2 & 3\\ 2 & 1 \end{pmatrix} \times \begin{pmatrix} 1\\ 3 \end{pmatrix} = \begin{pmatrix} 11\\ 5 \end{pmatrix}\)

  • Vecteurs propres et valeurs propres (il en existe autant qu’il y a de colonnes dans la matrice de départ) :

\[ \begin{pmatrix} 2 & 3\\ 2 & 1 \end{pmatrix} \times \begin{pmatrix} 3\\ 2 \end{pmatrix} = \begin{pmatrix} 12\\ 8 \end{pmatrix} = 4 \times \begin{pmatrix} 3\\ 2 \end{pmatrix} \]

  • La constante (4) est une valeur propre et la matrice multipliée (à droite) est la matrice des vecteurs propres.

  • La rotation d’un système d’axes à deux dimensions d’un angle \(\alpha\) peut se représenter sous forme d’un calcul matriciel :

\[ \begin{pmatrix} \cos \alpha & \sin \alpha\\ -\sin \alpha & \cos \alpha \end{pmatrix} \times \begin{pmatrix} x\\ y \end{pmatrix} = \begin{pmatrix} x'\\ y' \end{pmatrix} \]

Dans le cas particulier de l’ACP, la matrice de transformation qui effectue la rotation voulue pour obtenir les axes principaux est la matrice rassemblant tous les vecteurs propres calculés après diagonalisation de la matrice de corrélation ou de variance/covariance (réduction ou non, respectivement). Le schéma suivant visualise la rotation depuis les axes initiaux X et Y (variables de départ) en bleu foncé vers les PC1, PC2 en rouge. Un individu p est représenté par les coordonnées {x, y} dans le système d’axes initial XY. Les nouvelles coordonnées {x’, y’} sont recalculées par projection sur les nouveaux axes PC1-PC2. Les flèches bleues sont représentées dans l’espace des variables, tandis que les points reprojetés sur PC1-PC2 sont représentés dans l’espace des individus selon les coordonnées primes en rouge.

7.2.5.1 Résolution numérique simple

Effectuons une ACP sur matrice var/covar sans réduction des données (mais calcul très similaire lorsque les données sont réduites) sur un exemple numérique simple.

  • Étape 1 : centrage des données (moyenne nulle de chaque colonne)

\[ \mathop{\begin{pmatrix} 2 & 1 \\ 3 & 4 \\ 5 & 0 \\ 7 & 6 \\ 9 & 2 \end{pmatrix}}_{\text{Tableau brut}} \xrightarrow{\phantom{---}\text{centrage}\phantom{---}} \mathop{\begin{pmatrix} -3.2 & -1.8 \\ -2.2 & \phantom{-}1.4 \\ -0.2 & -2.6 \\ \phantom{-}1.8 & \phantom{-}3.4 \\ \phantom{-}3.8 & -0.6 \end{pmatrix}}_{\text{Tableau centré (X)}} \]

  • Étape 2 : calcul de la matrice de variance/covariance

\[ \mathop{\begin{pmatrix} -3.2 & -1.8 \\ -2.2 & \phantom{-}1.4 \\ -0.2 & -2.6 \\ \phantom{-}1.8 & \phantom{-}3.4 \\ \phantom{-}3.8 & -0.6 \end{pmatrix}}_{\text{Tableau centré (X)}} \xrightarrow{\phantom{---}\text{var/covar}\phantom{---}} \mathop{\begin{pmatrix} 8.2 & 1.6 \\ 1.6 & 5.8 \end{pmatrix}}_{\text{Matrice carrée (A)}} \]

  • Étape 3 : diagonalisation de la matrice var/covar

\[ \mathop{\begin{pmatrix} 8.2 & 1.6 \\ 1.6 & 5.8 \end{pmatrix}}_{\text{Matrice carrée (A)}} \xrightarrow{\phantom{---}\text{diagonalisation}\phantom{---}} \mathop{\begin{pmatrix} 9 & 0 \\ 0 & 5 \end{pmatrix}}_{\text{Matrice diagonalisée (B)}} \]

  • La trace des deux matrices A et B (somme des éléments sur la diagonale) est égale à : 8.2 + 5.8 = 14 = 9 + 5.

  • 8.2 est la part de variance exprimée sur le premier axe initial (X)

  • 5.8 est la part de variance exprimée sur le second axe initial (Y)

  • 14 est la variance totale du jeu de données

  • La matrice diagonale B est la solution exprimant la plus grande part de variance possible sur le premier axe de l’ACP : 9, soit 64,3% de la variance totale.

  • Les éléments sur la diagonale sont les valeurs propres \(\lambda_i\) ! Vous vous rappelez les fameuses “eigenvalues” dans la sortie de summary(pima_pca).

  • Étape 4 : calcul de la matrice de rotation des axes (en utilisant la propriété des valeurs propres \(\text{A}.\text{U} = \text{B}.\text{U}\)).

\[ \mathop{\begin{pmatrix} 8.2 & 1.6 \\ 1.6 & 5.8 \end{pmatrix}}_{\text{Matrice A}} \times \text{U} = \mathop{\begin{pmatrix} 9 & 0 \\ 0 & 5 \end{pmatrix}}_{\text{Matrice B}} \times \text{U} \rightarrow \text{U} = \mathop{\begin{pmatrix} \phantom{-}0.894 & -0.447 \\ \phantom{-}0.447 & \phantom{-}0.894 \end{pmatrix}}_{\text{Matrice des vecteur propres (U)}} \]

  • La matrice des vecteurs propres (U) (“eigenvectors” en anglais) effectue la transformation (rotation des axes) pour obtenir les composantes principales.
  • L’angle de rotation se déduit en considérant que cette matrice contient des sinus et cosinus d’angles de rotation des axes :

\[ \begin{pmatrix} \phantom{-}0.894 & -0.447 \\ \phantom{-}0.447 & \phantom{-}0.894 \end{pmatrix} = \begin{pmatrix} \phantom{-}\cos(-26.6°) & \phantom{-}\sin(-26.6°) \\ -\sin(-26.6°) & \phantom{-}\cos(-26.6°) \end{pmatrix} \]

  • Étape 5 : représentation dans l’espace des variables. C’est une représentation dans un cercle de la matrice des vecteurs propres U sous forme de vecteurs.

  • Étape 6 : représentation dans l’espace des individus. On recalcule les coordonnées des individus dans le système d’axe après rotation.

\[ \mathop{\begin{pmatrix} -3.2 & -1.8 \\ -2.2 & \phantom{-}1.4 \\ -0.2 & -2.6 \\ \phantom{-}1.8 & \phantom{-}3.4 \\ \phantom{-}3.8 & -0.6 \end{pmatrix}}_{\text{Tableau centré (X)}} \times \mathop{\begin{pmatrix} \phantom{-}0.894 & -0.447 \\ \phantom{-}0.447 & \phantom{-}0.894 \end{pmatrix}}_{\text{Matrice des vecteur propres (U)}} \xrightarrow{\phantom{---}\text{X}.\text{U} = \text{X'}\phantom{---}} \mathop{\begin{pmatrix} -3.58 & \phantom{-}0.00 \\ -1.34 & \phantom{-}2.24 \\ -1.34 & -2.24 \\ \phantom{-}3.13 & \phantom{-}2.24 \\ \phantom{-}3.13 & -2.24 \end{pmatrix}}_{\text{Tableau avec rotation (X')}} \]

  • Ensuite, on représente ces individus à l’aide d’un graphique en nuage de points.

Tous ces calculs se généralisent facilement à trois, puis à N dimensions.

À vous de jouer !

Effectuez maintenant les exercices du tutoriel B07La_pca (Analyse en composantes principales).

BioDataScience2::run("B07La_pca")

Réalisez le travail B07Ia_acp_afc, partie I.

Travail individuel pour les étudiants inscrits au cours de Science des Données Biologiques II à l’UMONS (Q2 : analyse) à terminer avant le 2024-02-26 23:59:59.

Initiez votre projet GitHub Classroom

Voyez les explications dans le fichier README.md, partie I.

Pour aller plus loin
  • N’hésitez pas à combiner plusieurs techniques. Par exemple, vous pouvez représenter les groupes créés par classification ascendante hiérarchique ou des k-moyennes sur un graphique de l’ACP dans l’espace des individus en faisant varier les couleurs ou les labels des individus en fonction des groupes.

  • Une autre explication de l’ACP en utilisant quelques autres fonctions de R pour visualiser le résultat

  • Une explication détaillée de la PCA en anglais.

  • Une page qui reprend des explications et une série de vidéos qui présentent les différentes facettes de l’ACP (en anglais).