Charger le package mixOmics ou l’installer si besoin :

library(mixOmics)
# install.packages("mixOmics")

Rappels

La variance expliquée

On note Variance expliquée,en \(\%\) de la variance totale de \(X\), dans le cas de la régression PLS (sparse et même DA), la grandeur, pour une matrice de covariables centrée \(X\), un poids associé \(u\)

\[ Var_{X,u}= \frac{||Xu||_2^2}{||X||_F^2}*100, \]

  • \(||Xu||_2^2=\{\text{Somme de tous les éléments de $Xu$ élevés au carré}\}\),

  • \(||X||_F^2=\{\text{Somme de tous les éléments de $X$ élevés au carré}\}\).

Petite simulation

Il a été demandé en cours de l’intérêt de construire une composante qui a une faible variance potentielle. L’exemple qui suit montre l’impérieuse nécessité de ne pas se concentrer seulement sur la variance expliquée lorsque l’on fait une analyse en PLS.

Génération des données

On regarde comment l’ajout d’une composante peut avoir un intérêt certain sur les résultats, en analyse discriminante, sur des données simulées.

n1 <- 50
n2 <- 30
n3 <- 20
n4 <- 10

mu <- 10

X1 <- cbind(rnorm(n1,mu),rnorm(n1,0),rnorm(n1,0),rnorm(n1,0))
X2 <- cbind(rnorm(n2,0),rnorm(n2,mu),rnorm(n2,0),rnorm(n2,0))
X3 <- cbind(rnorm(n3,0),rnorm(n3,0),rnorm(n3,mu),rnorm(n3,0))
X4 <- cbind(rnorm(n4,0),rnorm(n4,0),rnorm(n4,mu),rnorm(n4,mu))

X <- scale(rbind(X1,X2,X3,X4))
Y <- as.factor(c(rep("A",n1),rep("B",n2),rep("C",n3),rep("D",n4) ))

Construction du modèle

On construit un modèle à \(3\) composantes que l’on observe

r <- 3
mm <- plsda(X,Y,ncomp = r)
plot(as.data.frame(mm$variates$X),col=Y,pch=as.numeric(Y)+14)

On peut aussi représenter les densités marginales.

library(ggplot2)
toplot <- data.frame(cbind(mm$variates$X,as.factor(Y)))
alpha <- 0.4
cols = c("black","red","green","blue")
ggplot(toplot, aes(comp.1, fill = Y)) +
  geom_density(alpha = alpha)+
  scale_fill_manual(values=cols)

ggplot(toplot, aes(comp.2, fill = Y)) +
  geom_density(alpha = alpha)+
  scale_fill_manual(values=cols)

ggplot(toplot, aes(comp.3, fill = Y)) +
  geom_density(alpha = alpha)+
  scale_fill_manual(values=cols)

Calcul de la variance expliquée

On peut aussi s’attarder à calculer la variance expliquée par chaque axe, en \(\%\) de la variance totale de \(X\):

varTotX <- sum((X)^2)
var_comps <- diag(crossprod(mm$variates$X)/varTotX)*100
barplot(var_comps,pch=16,ylim=c(0,max(var_comps)),ylab="Variance explained (% of total X variance)");abline(h=0)

Ces valeurs sont aussi accessibles via l’argument explained_variance de tout modèle PLS de mixOmics. On y retrouve aussi les variances calculées pour Y.

Remarque :

Ces variances ne sont pas forcemment décroissantes par rapport à l’indice de la composante construite contrairement à l’ACP. Faire bien attention :)

Question On voit très nettement que les deux premières composantes apportent énormemment d’information. Mais qu’en est-il de la troisième ? Faut-il en conclure que cette composante n’est pas utile ?

Leave-One-Out validation croisée

errs <- rep(0,r)
for(h in 1:r){
  mm <- plsda(X,Y,ncomp = h)
  cvcv <- perf(mm,progressBar = F,validation = "loo")
  errs[h] <- cvcv$error.rate$overall[h,1]
}
names(errs) <- c("comp 1","comp 2","comp 3")
barplot(errs,pch=16,ylim=c(0,max(errs)),ylab="Error rate");abline(h=0)

Il est clair que ne construire que 2 composante amènerait à faire une grosse erreur potentielle, \(\sim 10\%\) d’erreur en leave-one-out.

La question est maintenant :

Est-ce intéressant de construire cette composante qui permet de réduire drastiquement l’erreur de prédiction mais qui apporte peu de variance ?

La réponse est bien entendu OUI. Se tromper dans \(\sim 10\%\) des cas c’est énorme. D’autant plus que les erreurs sont toutes commises pour un même groupe, le groupe bleu, si l’on ne construisait pas cette troisième composante, légère en terme de variance mais pas inutile du tout.

Ici on voit bien que la variance expliquée seule ne suffit pas pour faire le choix de modèle. La validation croisée s’avère être nécessaire. Ce n’est pas anecdotique, en grande dimension, il faut faire très attention, dans le cas de la PLS, aux conclusions que l’on peut faire en regardant seulement les variances expliquées par chacun des axes. La validation croisée permet seule de répondre à la question initiale de la prédiction du Y.