######################################################################
# TP 1 : Analyse en Composantes Principales (ACP) avec le logiciel R
# (code R avec commentaires)
######################################################################
#===================================
# Jeu de donnees : eaux petillantes
#===================================
# Chargement du jeu de donnees
eaux <- read.table("../data/eaux.txt", header=TRUE,row.names=1,sep=";")
eaux
## intensite.emission.de.bulles saveur.salee appreciation.globale
## St Yorre 3.9 6.4 2.9
## Vichy 1.4 6.0 2.8
## Quezac 5.1 4.7 3.5
## Salvetat 2.9 4.1 3.4
## Perrier 8.2 4.9 2.8
# Chargement du package FactoMineR
library(FactoMineR) # permet de charger le package "PCAmixdata"
# afin de pouvoir l'utiliser par la suite
help(PCA)
# Mise en oeuvre de l'ACP
res <-PCA(eaux) # tous les calculs de l'ACP sont stockes dans l'objet "res"


# NB : par défaut les graphiques des plans factoriels 1-2
# sont affichés a l'écran :
# - projection des individus sur le plan 1-2
# - cercle des corrélations du plan 1-2
res <- PCA(eaux,graph=FALSE) # idem sans les graphiques
res # permet de voir la liste de l'ensemble des sorties numeriques disponibles
## **Results for the Principal Component Analysis (PCA)**
## The analysis was performed on 5 individuals, described by 3 variables
## *The results are available in the following objects:
##
## name description
## 1 "$eig" "eigenvalues"
## 2 "$var" "results for the variables"
## 3 "$var$coord" "coord. for the variables"
## 4 "$var$cor" "correlations variables - dimensions"
## 5 "$var$cos2" "cos2 for the variables"
## 6 "$var$contrib" "contributions of the variables"
## 7 "$ind" "results for the individuals"
## 8 "$ind$coord" "coord. for the individuals"
## 9 "$ind$cos2" "cos2 for the individuals"
## 10 "$ind$contrib" "contributions of the individuals"
## 11 "$call" "summary statistics"
## 12 "$call$centre" "mean of the variables"
## 13 "$call$ecart.type" "standard error of the variables"
## 14 "$call$row.w" "weights for the individuals"
## 15 "$call$col.w" "weights for the variables"
# COMMENTAIRES : par défaut, c'est une ACP centrée réduite qui est faite avec
# le même poids (1/n) pour tous les individus de l'étude (ici les eux minérales
# pétillantes). L'espace des individus est R^3 car p=3 variables quantitatives
# ont été mesurées sur chaque individu.
# L'espace des variables est R^5 car chaque variable a été mesurée sur n=5
# individus.
# NB : il s'agit d'un tout petit jeu de données (n=5, p=3) qui contient donc
# seulement 15 valeurs numériques. Cependant, on va voir que l'ACP permet de
# faire ressortir très rapidement et visuellement les grandes informations
# contenues dans ce jeu de données.
# Choix du nombre d'axes à retenir
round(res$eig,digit=2)
## eigenvalue percentage of variance cumulative percentage of variance
## comp 1 1.73 57.78 57.78
## comp 2 1.06 35.45 93.23
## comp 3 0.20 6.77 100.00
# Permet d'afficher les valeurs propres,
# les pourcentages de variances expliquées par
# axe, puis cumulés.
# Graphique de l'ébouli des valeurs propres
barplot(res$eig[,1], main="Eigenvalues", names.arg=1:nrow(res$eig))
abline(h=1,col=2,lwd=2)

# COMMENTAIRES : comme on est en ACP normée, on peut utiliser le critère
# de Kaiser qui dit que seuls les axes factoriels associés à une valeur propre
# (eigenvalue) plus grande que 1 sont des axes intéressants à retenir.
# En effet, en ACP normée, les variables ont été réduites (divisées par leur
# écart-type), ainsi leur "nouvelle" variance est 1.
# Les variables synthétiques (composantes principales) construites par l'ACP
# ne seront intéressantes que si elles ont une dispersion plus importante
# que celle de n'importe quelle variables de l'étude. Or par construction,
# la variance d'une composante principale est égale à la valeur propre de l'axe
# correspondant. D'où le critère de Kaiser qui préconise de ne retenir que
# les axes ayant une valeur propre supérieure à 1.
#
# NB 1 : si on a une valeur propre inférieure à 1, mais proche de 1, il est
# recommandé de regarder tout de même quelle est l'information apportée par
# l'axe correspondant.
#
# NB 2 : si on une valeur propre plus grande que 1 (et proche de 1) qui
# correspond à un axe factoriel pour lequel on n'arrive pas à faire des
# commentaires intéressants, il est recommandé de ne pas retenir cet axe
# factoriel pour la suite des interprétations.
#
# NB 3 : il y a au maximum p valeurs propres car on peut fabriquer au maximum
# p composantes principale (=nouvelles variables synthétiques, deux à deux non
# corrélées, combinaisons linéaires des variables initiales) à partir de p
# variables initiales.
#
# Ici on obtient les sorties numériques suivantes :
# Eigenvalue Proportion Cumulative
# dim 1 1.73 57.78 57.78
# dim 2 1.06 35.45 93.23
# dim 3 0.20 6.77 100.00
#
# On a donc les deux premières valeurs propres (les valeurs propres sont
# toujours rangées par ordre décroissant) supérieures à 1. Selon le critère
# de Kaiser, il convient donc de retenir 2 axes factoriels. Le premier axe
# factoriel permet d'expliquer 57,78% de l'inertie (ou de la variance), le
# second axe factoriel permet d'expliquer 35,45% d'inertie supplémentaire.
# Ainsi, en considérant le plan 1-2, on récupère 93,23% de l'information
# (c'est à dire en quoi les individus diffèrent = inertie avec un regard
# d'algèbre linéaire = variance avec un regard de statisticien).
#
# NB 4 : la troisième valeur propre est très inférieure à 1. Ce troisième axe
# factoriel ne porte donc pas d'information pertinente, les 6.77%
# d'inertie expliquée correspondent plutôt à du bruit (fluctuation aléatoire).
# Graphiques des individus et des variables sur le plan factoriel 1-2
?plot.PCA # permet d'afficher l'aide de la commande "plot.PCA"
plot(res,axes=c(1,2),choix="var") # on retrouve ici le cercle des corrélations (plan 1-2)

# Commentaires du cercle des corrélations du plan 1-2 :
#
# Il s'agit de la projection des variables (qui sont de "taille"longeur" 1 dans l'espace
# R^5 auquel elles appartiennent) sur le sous-espace de dimension 2 de R^5 qui
# apporte le plus d'information.
#
# 1) Les 3 variables sont bien représentées dans ce plan 1-2 (ce qui n'est
# pas une surprise car il contient 93,23% de l'information) car leurs
# projections sont proches de la circonférence du cercle.
#
# 2) Sur l'axe 1 (horizontal), on observe une opposition entre la variable
# "saveur.salee" (projetée à gauche) et la variable "appreciation.globale"
# (projetée à droite). Cela signifie que ces deux variables sont fortement
# corrélées linéairement et négativement entre elles : plus la saveur salée
# est forte, moins l'eau pétillante est appréciée globalement ; a contrario,
# les eaux les plus appréciées seront celles qui ont une faible saveur salée.
#
# 3) L'axe 2 (vertical) est essentiellement caractérisé par la variable
# "intensite.emission.de.bulles" qui est projetée vers le haut.
# Cette variable (en projection) est orthogonale aux deux autres variables.
# Cela signifie qu'il n'y a pas a priori des liaisons linéaires entre
# la variable "intensite.emission.de.bulles" et les deux variables
# "saveur.salee" et "appreciation.globale". Ainsi, l'intensité d'émission
# des bulles n'a pas a priori d'impact sur l'appréciation globale de l'eau
# minérale (ni sur la saveur salée).
#
# Il convient maintenant de projeter les individus (eaux pétillantes) sur le
# plan factoriel 1-2 associé au sous-espace de dimension 2 de R^3
# (espace des individus) qui apporte le plus d'information (maximisation de
# l'inertie du nuage projeté). Ceci nous permettra alors de caractériser
# les individus en fonction des variables.
plot(res,axes=c(1,2),choix="ind") # on retrouve ici le graphique des individus (plan 1-2)

# Commentaires du plan factoriel 1-2 des individus :
#
# 1) Sur l'axe 1 (horizontal), on voit une opposition très claire entre les
# aux pétillantes, à gauche, "St Yorre" et "Vichy", et, à droite, "Quezac" et
# "Salvetat". Au vu cercle des corrélations du plan 1-2, on en déduite que
# - "St Yorre" et "Vichy" sont des eaux qui ont une saveur salée
# importantes et qui n'ont pas été globalement appréciées.
# - C'est l'inverse pour "Quezac" et "Salvetat" qui sont des eaux
# pétillantes qui ont été globalement bien appréciées et qui n'ont pas
# une forte saveur salée.
# - L'eau pétillante "Perrier" se projette au centre (sur l'axe
# horizontal), il s'agit donc d'une eau qui a une saveur salée moyenne
# et dont l'appréciation globale a été aussi été moyenne
# (moyenne = par rapport aux cinq eaux pétillantes présentes dans
# l'étude).
#
# 2) Sur l'axe 2 (vertical), on va pouvoir caractériser les eaux pétillantes
# par rapport à leur intensité d'émission de bulles.
# - "Perrier" (en haut) diffère très clairement des quatre autres eaux
# pétillantes (en bas) : cela signifie que la caractéristique majeure
# de "Perrier" est que c'est l'eau pétillante ayant une très forte
# intensité d'émission de bulles (par rapport aux 4 autres eaux
# pétillantes).
# - Les eaux "Vichy" et "Salvetat" (en bas) sont celles présentant
# la moindre intensité d'émission de bulles.
# - "St Yorre" et "Quezac" sont dans la moyenne concernant cette
# caractéristique d'émission de bulles.
# Au delà des sorties graphiques, des sorties numériques sont disponibles
# et peuvent être utiles pour quantifier certaines observations graphiques.
# Sorties numeriques pour les individus et pour les variables
# Remarque : par défaut les calculs sont faits sur les 5 premiers axes.
# (NB : dans cette étude, on peut fabriquer au maximum que 3 axes factoriels).
res$ind # permet d'afficher les sorties numériques associées aux individus : coordonnées, contributions, cosinus carrés.
## $coord
## Dim.1 Dim.2 Dim.3
## St Yorre -1.422259675 -0.06130521 0.5185196
## Vichy -1.520947213 -0.90545543 -0.3700892
## Quezac 1.422493513 -0.13598570 0.5801695
## Salvetat 1.518799529 -0.83922841 -0.4280543
## Perrier 0.001913846 1.94197475 -0.3005456
##
## $cos2
## Dim.1 Dim.2 Dim.3
## St Yorre 8.812339e-01 0.001637300 0.11712881
## Vichy 7.074044e-01 0.250711156 0.04188440
## Quezac 8.507138e-01 0.007774444 0.14151179
## Salvetat 7.221493e-01 0.220488754 0.05736194
## Perrier 9.485212e-07 0.976607789 0.02339126
##
## $contrib
## Dim.1 Dim.2 Dim.3
## St Yorre 2.334112e+01 0.07067623 26.463285
## Vichy 2.669268e+01 15.41745729 13.481126
## Quezac 2.334880e+01 0.34774834 33.130132
## Salvetat 2.661735e+01 13.24460325 18.034795
## Perrier 4.226472e-05 70.91951488 8.890661
##
## $dist
## St Yorre Vichy Quezac Salvetat Perrier
## 1.515072 1.808341 1.542263 1.787257 1.965095
# NB : les coordonnées peuvent être utiles dans le cas où deux individus sont projetés
# au même endroit sur le plan factoriel et qu'il est difficile de lire leurs noms graphiquement.
round(res$ind$cos2,digits=3) # uniquement les cosinus carres (arrondis à 3 chiffres après la virgule)
## Dim.1 Dim.2 Dim.3
## St Yorre 0.881 0.002 0.117
## Vichy 0.707 0.251 0.042
## Quezac 0.851 0.008 0.142
## Salvetat 0.722 0.220 0.057
## Perrier 0.000 0.977 0.023
# Commentaires : les cosinus carrés nous donnent la qualité de représentation (en %)
# des individus sur les axes factoriels.
# Par exemple : 70.7 % de l'information de "Vichy" est recupérée sur l'axe 1, 25,1% sur l'axe 2.
# Vu que les axes sont orhtogonaux, on a 70,7+25,1=95,8% de l'information de "Vichy"
# qui est fournie par le plan 1-2.
# "Perrier" est extrêmement bien représentée sur l'axe 2 (97,7%).
# On observe ici que les n=5 individus sont très bien représentés sur le plan 1-2.
# Ceci n'est pas surprenant car la qualité globale de plan 1-2 est de 93,23%
# (d'inertie ou variance expliquée).
#
# NB : Ne jamais faire des commentaires sur des individus qui seraient mal représentés
# sur le plan factoriel que l'on regarde car on ne ferait que des fausses déductions...
round(res$ind$contrib, digits=1) # uniquement les contributions en % (arrondies à un chiffre après la virgule)
## Dim.1 Dim.2 Dim.3
## St Yorre 23.3 0.1 26.5
## Vichy 26.7 15.4 13.5
## Quezac 23.3 0.3 33.1
## Salvetat 26.6 13.2 18.0
## Perrier 0.0 70.9 8.9
# Commentaires : On voit que 4 individus contribuent formement à l'axe 1
# "St Yorre", "Vichy", "Quezac" et "Salvetat" (avec des contributions de
# l'ordre de 25%). "Pierrier" n'a pas du tout contribué à la fabrication de
# ce premier axe factoriel.
# Pour l'axe 2, c'est essentiellement "Perrier" qui a fortement contribué
# à sa fabrication (70,9%) et dans une moindre mesure "Vichy" et "Salvetat"
# (avec des contributions de l'ordre de 15%).
#
# NB : je rappelle que l'on ne souhaite pas que seul un individu (ou peu
# d'invidus lorsque n est grand) contribue à la fabrication du (ou des)
# premier(s) axe(s) factoriels car cela rendrait alors les résultats de l'étude
# instables dans le sens où la présence ou pas de cet individu dans
# les données va entrainer une modification dans la construction des premiers axes factoriels.
# Ici, ce n'est pas le cas. On a 4 individus sur 5 (80% des individus) qui
# contribuent au premier axe factoriel, et 3 sur 5 (60%) qui contribuent au
# second axe. C'est normal que "Perrier" qui n'a pas contribué à la fabrication
# de l'axe 1 (qui contient 57,78% de l'information) ait tendance à contribuer
# fortement à la constitution de l'axe 2 (qui contient 35,45% de l'information).
res$var # permet d'afficher les sorties numériques associées aux variables : coordonnées, contributions, cosinus carrés.
## $coord
## Dim.1 Dim.2 Dim.3
## intensite.emission.de.bulles 0.2644680 0.95324075 0.1462489
## saveur.salee -0.9462657 -0.08864987 0.3110022
## appreciation.globale 0.8763031 -0.38341534 0.2916943
##
## $cor
## Dim.1 Dim.2 Dim.3
## intensite.emission.de.bulles 0.2644680 0.95324075 0.1462489
## saveur.salee -0.9462657 -0.08864987 0.3110022
## appreciation.globale 0.8763031 -0.38341534 0.2916943
##
## $cos2
## Dim.1 Dim.2 Dim.3
## intensite.emission.de.bulles 0.06994335 0.908667926 0.02138873
## saveur.salee 0.89541882 0.007858799 0.09672238
## appreciation.globale 0.76790713 0.147007322 0.08508555
##
## $contrib
## Dim.1 Dim.2 Dim.3
## intensite.emission.de.bulles 4.035342 85.4385366 10.52612
## saveur.salee 51.660687 0.7389326 47.60038
## appreciation.globale 44.303971 13.8225308 41.87350
# NB : les coordonnées peuvent être utiles dans le cas où deux
# variables sont projetées au même endroit sur le cercle des
# corrélations et qu'il est difficile de lire leurs noms graphiquement.
round(res$var$cos2,digits=3) # uniquement les cosinus carres (arrondis à 3 chiffres après la virgule)
## Dim.1 Dim.2 Dim.3
## intensite.emission.de.bulles 0.070 0.909 0.021
## saveur.salee 0.895 0.008 0.097
## appreciation.globale 0.768 0.147 0.085
# Comme pour les individus, les cosinus carrés nous donnent la qualité (en %)
# de réprésentation des variables sur les axes factoriels.
# Par exemple, la variable "intensite.emission.de.bulles" est très mal
# représentée sur l'axe 1 (7% seulement de l'information de cette variable est
# portée par cet axe), mais elle est très bien représentée sur l'axe 2 (90,9%
# de l'information de cette variable est fournie par l'axe vertical).
# Globalement sur le plan principal (1-2), cette variable
# "intensite.emission.de.bulles" a une qualité de représentation égale à
# 7% + 90,9% = 97,9%, elle est donc quasiment parfaitement représentée sur
# le plan 1-2.
#
# NB : on peut sommer les qualités de représentation de l'axe 1 et de l'axe 2
# pour avoir la qualité de représentations ur le plan 1-2, car les axes sont
# construits orthogonaux entre eux, et il n'y a donc pas de redondance
# d'information entre les axes factoriels.
#
# Les deux autres variables "saveur.salee" et "appreciation.globale" sont
# très bien représentée sur l'axe 1 (resp. 89.5% et 76,8%) et mal représentée
# sur l'axe 2 (resp. 0.8% et 14,7%).
round(res$var$contrib,digits=1) # uniquement les contributions en % (arrondies à un chiffre après la virgule)
## Dim.1 Dim.2 Dim.3
## intensite.emission.de.bulles 4.0 85.4 10.5
## saveur.salee 51.7 0.7 47.6
## appreciation.globale 44.3 13.8 41.9
# Les contributions (en %) des variables à la fabrication des axes factoriels
# nous indiquent naturellement que seules les variables "saveur.salee" et
# "appreciation.globale" ont beaucoup contribué à la fabrication de
# l'axe 1 (resp. 51.6% et 44,3%). La variable "intensite.emission.de.bulles"
# est la plus forte contributrice à la fabrication de l'axe 2
# (avec un pourcentage de contribution de l'ordre de 85,4%).