Visualiser des catégories/classes avec R

Cette page n'a pas vocation a faire un cours sur les catégorisations. Vaste programme ! Pour voir quelques méthodes de classification/catégorisation : voir ce lien sur le site.

Elle tend à donner des solutions pour illustrer facilement différentes catégories de données sur un graphique ainsi que les relations existantes entre les différents paramètres associés et ce de la façon la plus automatique.

1/4 - Visualiser rapidement la variabilité par paramètre et pour chaque catégorie afin de se faire une idée

Prenons les données d'un exemple : des analyses sanguines, 3 catégories de malades A, B et C

Ce code est à copier-coller pour simuler un jeu de données à 5 paramètres et 3 catégories.

#Simulation d'un jeu de données à 3 catégories et 5 variables

# Glycémie

x_a <- rnorm(100,1,0.2);x_b <- rnorm(120,1.1,0.23);x_c <- rnorm(80,1,0.3);x <- c(x_a,x_b,x_c)

# Calcémie

y_a <- runif(100,2.2,2.6) ;y_b <- rnorm(120,2.3,0.5) ;y_c <- rpois(80,22)/10 ;y <- c(y_a,y_b,y_c)

# Protéine C-Réactive 6 mg/L

z_a <- rnorm(40,6,2) ;z_a <- c(z_a,rnorm(40,20,15));z_a <- c(z_a,rnorm(20,3,4));z_b <- rnorm(120,5.5,0.5) ;z_c <- rpois(80,50);z <- abs(c(z_a,z_b,z_c))

# TSH 0.4 à 4 mUI/l

w_a <- rnorm(50,2,1) ;w_a <- c(w_a,rnorm(50,2,1.1));w_a <- sort(w_a);w_b <- sort(abs(rnorm(120,1,2)));w_c <- rev(sort(abs(rnorm(80,2,1))));w <- abs(c(w_a,w_b,w_c))

# HDL 0.4 à 0.5

v_a <- rpois(100,50)/100;v_b <- sort(rpois(120,50)/100);v_c <- sort(rpois(80,50)/100)+0.2;v <- c(v_a,v_b,v_c)

# Catégories

categories<-c(rep("A",100),rep("B",120),rep("C",80))

# Compilation dans un tableau data.frame

compilation <- data.frame(categories,x,y,z,w,v);colnames(compilation) <- c("Catégories","Glycémie","Calcémie","Protéine C-Réactive","TSH","HDL");

head(compilation)

Affichons un histogramme pour les données de chaque catégorie pour chaque paramètre

Le mieux est de tracer des diagrammes en violons, à la fois boîtes à moustaches et histogramme.

On peut aussi faire des boîtes à moustaches (boxplot), et le mieux sera toujours des histogrammes (hist) si on a le temps.

ÉTAPE 1 - INSTALLER LA LIBRAIRIE VIOPLOT

install.packages("vioplot") # installation

library(vioplot) # exécution

ÉTAPE 2 - COPIER-COLLER LA FONCTION MULTIVIOPLOT

multivioplot <- function(data,cat,screen=c(16,10),col=c(1:100)) {

parameters <- colnames(data)

cat("Parameters : ",parameters,"\n\n")

nb_col <- length(parameters)

nb_graphics <- nb_col

nb_graphics_y <- floor(nb_graphics/sum(screen)*screen[2])

nb_graphics_x <- ceiling(nb_graphics/nb_graphics_y)

# Conversion des couleurs

# Fonction de conversion des couleurs en hexadécimal

col2hex <- function(col, alpha=1) rgb(t(col2rgb(col)), alpha=alpha*255, maxColorValue=255)

colors <- col2hex(col) # cf. col2rgb

# Affichage

layout(matrix(1:(nb_graphics_x*nb_graphics_y),nb_graphics_y,nb_graphics_x))

require(vioplot)

for (i in c(1:nb_col)) {

plot(c(0,(length(unique(cat))+1)),c(min(data[,i]),max(data[,i])),col="white",type="n",axes=F,xlab="",ylab="")

axis(2) # Ajoute l'axe des y à gauche

axis(1,at=c(1:length(unique(cat))),, labels=unique(cat),xlab="") # Ajoute l'axe des x en bas

mtext(parameters[i],side=1,col="black",line=2.5)

grid()

for (j in c(1:length(unique(cat)))) {vioplot(data[,i][cat==unique(cat)[j]],at=j,horizontal=F,add=T,col=colors[j]) }

}

}

ÉTAPE 3 - AFFICHER LES VIOPLOTS, A LA FOIS HISTOGRAMME ET BOITE A MOUSTACHES

multivioplot(compilation [,c(2:6)],compilation[,1],screen=c(3,3),col=c(1,2,3))

2/4 - Mettre en relations les différents jeux de données pour voir la réparition des données par catégorie

Voici deux exemples d'affichage avec multiplot() réalisé à partir des données ci-dessus (cf. commande multivioplot()).

La commande multiplot() est disponible sous R en téléchargeant le code dans la rubrique plus bas dans la page sous ces exemples.

• Voici le code de la commande multiplot() à copier-coller dans R. Cette commande nécessite la librairie/le package ellipse pour tracer les ellipses de confiance.

install.packages("ellipse")

multiplot <- function(data,cat,screen=c(16,10),col=c(1:100),pch=rep(16,100),grid=T,chull=F,chull.alpha = 0.5,ellipse=F,conf.ellipse=0.95) {

parameters <- colnames(data)

cat("Parameters",parameters,"\n\n")

nb_col <- length(parameters)

nb_graphics <- sum(c(1:(nb_col-1)))

nb_graphics_y <- floor(nb_graphics/sum(screen)*screen[2])

nb_graphics_x <- ceiling(nb_graphics/nb_graphics_y)

# Fonctions pour faire les nappes convexes et les ellipses

convexhull<-function(x, y, col="#00000000",border="black",lwd=lwd){

hpts <- chull(x = x, y = y)

hpts <- c(hpts, hpts[1])

polygon(x[hpts], y[hpts],col=col,border=border,lwd=lwd)

}

confellipse <-function(x=c(),y=c(),level=0.95,col,lwd=2,lty=1){

# Fonction pour éviter le problème d'interférence entre librairies : detach

#detach("package:car");detach("package:mixtools")

# Tracer une ellipse (intervalle de confiance à 50%)

xy_ell <- data.frame(x,y)

require("ellipse")

lines(ellipse(cov(xy_ell),centre=colMeans(xy_ell),level=level),type="l", lty=lty,lwd=lwd, col=col)

}

# Fonction de conversion des couleurs en hexadécimal

col2hex <- function(col, alpha=1) rgb(t(col2rgb(col)), alpha=alpha*255, maxColorValue=255)

# Mise en place d'une couleur par catégorie

#type_col <- is(col)[1]

type_pch <- is(pch)[1]

#if (type_col == "character") {colors <- rep("A",length(cat))

#} else {colors <-rep(1,length(cat))} #(type_col == "numeric" |"integer")

col <- col2hex(col) # cf. col2rgb

colors <- rep("#000000",length(cat))

if (type_pch == "character") {pchs <-rep("A",length(cat))

} else {pchs <-rep(1,length(cat))}

list_cat <- c(); list_col <- c() ; list_pch <- c()

cpt_temp = 0

# Répéter les lists col et pch si pas assez d'entrée // nombre de catégories

if (nb_col > length(col)) {col = rep(col,ceiling(nb_col/length(col)))}

if (nb_col > length(pch)) {pch = rep(pch,ceiling(nb_col/length(pch)))}

for (i in unique(cat)) {

cpt_temp <- cpt_temp+1

indices <- which(cat==i)

colors[indices] <- col[cpt_temp]

pchs[indices] <- pch[cpt_temp]

list_cat <- c(list_cat,i)

list_col <- c(list_col,col[cpt_temp])

list_pch <- c(list_pch,pch[cpt_temp])

}

compil_legend <- data.frame(list_cat,list_col,list_pch)

colnames(compil_legend) <- c("Catégories","Couleurs","Pch")

layout(matrix(1:(nb_graphics_x*nb_graphics_y),nb_graphics_y,nb_graphics_x))

max_i <- nb_col-1

for (i in c(1:max_i)) {

for (j in c((i+1):nb_col)) {

plot(data[,i],data[,j],xlab=parameters[i],ylab=parameters[j],col=colors,pch=pchs)

if (grid==T) {grid()}

if (chull ==T) {

for (k in c(1:length(compil_legend[,1]))) {

temp_col <- col2hex(compil_legend[k,2],alpha=chull.alpha)

convexhull(data[,i][cat==compil_legend[k,1]],data[,j][cat==compil_legend[k,1]],border=temp_col,lwd=2,col=temp_col)

}

}

if (ellipse ==T) {

for (k in c(1:length(compil_legend[,1]))) {

temp_col <- col2hex(compil_legend[k,2],alpha=chull.alpha)

confellipse(data[,i][cat==compil_legend[k,1]],data[,j][cat==compil_legend[k,1]],level=conf.ellipse,lwd=2,col=temp_col)

} } } }

cat("Légende:\n")

compil_legend}

3/4- Visualiser le jeu de données et les différentes catégories à partir d'une ACP ou d'une représentation de Fisher (LDA, Linear Discriminant Analysis)

ACP_visual() : pour visualiser en une ligne les 2 composantes principales, les catégories et le poids des facteurs.

• Il suffit de copier-coller la fonction ci-dessous et de l'utiliser comme dans l'exemple ci-dessous.

# Attribuer des couleurs aux catégories

color <- categories ;color[color=="A"] <- "black"

color[color=="B"] <- "red";color[color=="C"] <- "green"

# Visualiser

acp_visual(compilation[-1],col=color)

• Fonction à copier-coller :

install.packages("ade4");install.packages('plotrix') # ces 2 librairies sont nécessaires

acp_visual <- function(data.frame,col=1) {

require(ade4)

acp <- dudi.pca(data.frame, scannf= F,scale=T,center=T)

plot(acp$li[,1],acp$li[,2],pch=16,col=col,xlab="Composante 1",ylab="Composante 2")

# Une fois le graphique tracé

par(new = T) # pour ajouter par-dessus les vecteurs

plot(c(-1,1),c(-1,1),type="n",xlab="",ylab="",asp=1,axes=F) # pour tracer un graphique blanc orthonormé

# optionnel : tracer un cercle

require("plotrix")

draw.circle(0,0,1)

# Tracer les vecteurs ==> Remarque acp$c1 contient ces descriptions pour les acp faites avec ade4.

arrows(0,0,acp$c1[,1],acp$c1[,2],length=0.1,angle=17)

# Récupérer les noms des vecteurs

vecnames <- rownames(acp$c1)

# Tracer les noms des vecteurs

text(acp$c1[,1]+0.07,acp$c1[,2]+0.07,vecnames,cex=0.75)}

LDA_visual() : pour visualiser en une ligne les 2 composantes et les catégories dont la séparation aura été forcée par une analyse linéaire du discriminant.

• Il suffit de copier-coller la fonction ci-dessous et de l'utiliser comme dans l'exemple ci-dessous.

color <- categories ;color[color=="A"] <- "black"

color[color=="B"] <- "red";color[color=="C"] <- "green"

lda_visual(compilation[,-1],compilation[,1],col=color)

• Fonction à copier-coller :

lda_visual <- function(data.frame,categories,col=1) {

require("MASS")

my_lda <- lda(x=data.frame,data=data.frame,grouping=categories)

D1 <- my_lda$scaling[,1] # vecteur 1

D2 <- my_lda$scaling[,2] # vecteur 2

# Calcul de projeté des individus sur D1 et D2

xy <- as.matrix(data.frame)%*%as.matrix(cbind(D1,D2))# produit scalaire

plot(xy,col=col,pch=16,xlab="Composante 1 (LDA)",ylab="Composante 2 (LDA)")}

4/4- Appliquer un arbre de décision pour définir quels seront les critères qui permettront de distinguer les catégories

• Permet de créer un arbre séparant les catégories (compilation[,1] en fonction des critères 3 et 4 de compilation.

plot(compilation[,4],compilation[,3],col=compilation[,1],pch=16)

install.packages("tree") ; library("tree")

mytree <- tree(compilation[,1]~compilation[,3]+compilation[,4])

plot(mytree)

text(mytree)