Cours de R pour débutant pressé, deuxième épisode

Voi­ci donc la suite de mon pre­mier cours de R pour débu­tants pres­sés où nous avions vu les bases. Je vous invite à le lire si vous ne l’avez pas déjà fait.

Aujourd’hui le but est de s’approcher de la bio­lo­gie et de trai­ter d’un cas concret. Nous avons des don­nées de dif­fé­rents élé­ments fonc­tion­nels du génome de D. mela­no­gas­ter et leurs pro­por­tions de conser­va­tion et de contrainte évo­lu­tive. Nous sou­hai­tons repré­sen­ter la pro­por­tion dans laquelle chaque élé­ment est contraint ou pas. Pour cela, nous allons faire un gra­phique en barre ou bar­plot, avec barres d'erreurs. Ce que nous allons voir aujourd'hui peut très bien s'appliquer à d'autres types de don­nées : repré­sen­ter la quan­ti­té de gènes suite à une qPCR ou encore le nombre de fois qu'une sou­ris a choi­si le bon che­min dans le laby­rinthe sui­vant chaque condi­tion choi­sie. L'important n'est donc pas ici le fond, mais bien la forme.

À la fin de ce tuto­riel, vous devriez être capables de faire quelque chose comme ça :

R et les éditeurs de texte

Nous allons tout d'abord ouvrir R et un docu­ment (édi­teur de texte inté­gré à R), qui va vous ser­vir à sau­ve­gar­der le code de votre gra­phique, la console R n'étant pas for­cé­ment la plus pra­tique.

La console R :

Pour ouvrir un docu­ment, cli­quez sur "Fichier/​Nouveau Docu­ment" ou ctrl+N sur Win­dows et GNU/​Linux, cmd+N sur mac. Une bonne habi­tude à prendre est d'écrire son code dans un édi­teur de texte, nous per­met­tant de l'éditer, de le sau­ve­gar­der, d'y reve­nir, sans à chaque fois avoir à remon­ter toutes les lignes. Je vous deman­de­rai donc de suivre ce conseil. Sur­tout que si vous ne le faites pas, vous serez embê­tés plus bas sur cette page !

Je vou­drais que vous notiez avant de com­men­cer que j'utiliserai la nota­tion amé­ri­caine du sépa­ra­teur des déci­males, donc le point et non la vir­gule.

Chargement du tableau de données

Nous allons com­men­cer par créer une variable qui va conte­nir les don­nées. Pour cela, nous allons char­ger dans R un tableau que je vous ai pré­pa­ré : Table.txt

results = read.table ("Chemin_ou_vous_avez_sauve_Table/Table.txt", header = T);
attach (results);

Nous voyons ici plusieurs choses :

  • comment lire un tableau ("read.table") ;
  • comment sauver ces données dans une variable ("results =") ;
  • comment dire à "read.table" que le tableau qu'on lui donne a ce qu'on appelle un header (un entête), c'est-à-dire que les colonnes ont un nom (la première ligne du fichier est composée des noms des colonnes).

[important]Quand on utilise la fonction "read.table", R s'attend à ce que le fichier qu'on lui donne ait les colonnes séparées par un espace ; il est possible de lui spécifier un autre séparateur avec l'argument "sep".[/important]

Enfin, "attach" sert à attacher les données contenues dans "results", dans le but de pouvoir appeler directement les colonnes de Table.txt, simplement en les nommant.

[notice]Essayez par exemple de taper "Conserved", il vous répondra les valeurs contenues dans la colonne "Conserved"[/notice]

 

1er barplot : La partie haute

Maintenant que nos valeurs sont bien accrochées, dessinons ! "barplot" est une fonction qui permet de faire... roulement de tambours... un barplot ! Soit en bon français, un graphique en barres. C'est exactement ce qui nous intéresse ici. Commençons donc par la ligne suivante :

barplot (Conserved, ylim = c(-1,1), main = "Bases", col = "yellow", axes = F);

Ici, nous allons donc représenter les données de la colonne "Conserved" :

  • en définissant les limites sur l'axe y (ylim) entre -1 et 1 ("c(-1,1)") ;
  • en définissant un titre principal ("main") ;
  • en colorant les barres ("col") en jaune ;
  • le tout sans dessiner d'axes (axes = F) ;

Nous représentons ainsi la partie haute de notre premier graphe :

 

2ème barplot : La partie basse

Pour faire la partie basse, c'est aussi simple, il suffit de taper ceci :

barplot (-Unconserved, ylim = c(-1,1), col = "yellow4", add=T, axes = F);

Si vous avez bien suivi, vous devez avoir :

  • la partie du bas ("-Unconserved") ;
  • dans les mêmes limites ;
  • avec une couleur plus sombre ("yellow4") ;
  • toujours sans les axes

[notice]Vous devez vous demander ce que c'est que ce "add=T". Eh bien, c'est un paramètre très utile que l'on peut trouver dans beaucoup de fonctions, qui permet très simplement d'ajouter ("add") un nouveau graphique à un premier, sans l'effacer ! C'est l'équivalent de la fonction points ou lines pour la fonction plot (je vous invite à taper "?lines", "?points" ou "?plot" si vous vous posez des questions). Cet "add" présente l'avantage de ne pas avoir à apprendre ou retenir le nom d'une autre fonction pour faire la même chose ![/notice]

Bon, trêve de bavardages, la suite !

Les axes

axis (2, label = c(1, 0.5, "0.0", 0.5, 1), at = c(-1, -.5, 0, .5, 1))

La fonction axis permet de dessiner un axe, que l'on n'avait pas dessiné avec barplot. Pourquoi faire ceci ? Essayez les lignes précédentes en mettant "axes = T", vous verrez que barplot vous dessine un axe y dont les valeurs vont de -1 à 1. L'intérêt étant ici non pas d'avoir un axe de -1 à 1, mais de 1 à 0, et de 0 à 1, car on représente 2 éléments qui sont inverses, en valeur absolue.

Le 2 est l'axe que l'on veut dessiner, en retenant le code suivant : 1 = bas, 2 = gauche, 3 = haut et 4 = droite.

[notice]Je vous invite à vous amuser un peu avec ces axes.[/notice]

"Label" est la liste de noms que l'on veut donner aux points affichés sur l'axe (une liste se définit en mettant des éléments séparés par une virgule dans "c()"), et "at" est la liste des positions respectives de ces points.

axis (1, label = c("FB CDS","Cel. CDS","FB 5'UTR","Cel. 5'UTR","FB 3'UTR","Cel. 3'UTR","FB Exons",
      "Cel. Exons","FB Introns","FB Pseudogenes","PolII binding sites","TF Sites","Histones",
      "CNV","Origins"), las = 3, at = 0.7+1.2*(0:14), pos = -1.1)

Faisons de même pour l'axe des x, que l'on veut voir porter le nom des "Features" de Table.txt, mais de façon plus lisible :

  • pour cela, on n'utilise pas directement Features, mais on redonne une liste (vous vous souvenez de "c()") ;
  • on aurait aussi très bien pu définir d'abord le nom des Features et l'utiliser ;
  • ensuite, on indique les positions de ces noms (tapez "0.7+1.2*(0:14)" dans R pour comprendre). Vous pouvez trouver les positions empiriquement en tâtonnant, mais sachez qu'en toute logique, il y a un espace avant les barres (ici le 0.7), puis que les barres sont également espacées (de 1.2), et qu'il y a 15 barres (on compte de 0 à 14 car la première est en position 0.7).
  • enfin, "pos" détermine la position verticale de l'axe, que l'on veut ici en dessous du graphique du bas, donc placé à un peu plus de -1, par exemple -1.1.

 

Les barres d'erreurs

Plaçons maintenant les barres d'erreur, qui servent en général à représenter l'écart type, qui est en termes simples une indication de la dispersion des résultats. Vous pouvez cependant représenter ce que vous voulez, Variance ou autres. La fonction qui nous intéresse est "errbar", elle est contenue dans le package "Hmisc", qu'il vous faudra installer ("install.packages()") et appeler ("require()" ou "library()") :

install.packages("Hmisc")
library(Hmisc)
errbar (0.7+1.2*(0:14), Mean, yplus = Mean + 2*SD, yminus = Mean - 2*SD, add = T, pch = "")
errbar (0.7+1.2*(0:14), -Mean2, yplus = -(Mean2 + 2*SD2), yminus = -(Mean2 - 2*SD2), add = T, pch = "")

Nous donnons donc à errbar :

  • en premier la position sur l'axe x où nous voulons nos barres d'erreurs ;
  • ensuite les positions, qui sont les mêmes que celle des noms sur l'axe x ;
  • puis la position du milieu de la barre ("Mean"), la position haute de la barre d'erreur ("yplus"), la position basse ("yminus") ;
  • enfin on l'ajoute au graphique précédent ("add=T") ;
  • et on utilise "pch" qui est le symbole utilisé pour le milieu de la barre (nous n'en voulons pas, donc on le met égal à un ensemble vide).

[error]Vous devez avoir un erreur qui apparaît quand vous tapez le code de la première barre d'erreur :

Erreur dans y[Type == 1] : objet de type 'closure' non indiçable

Pour comprendre cette erreur, c'est très simple. Tapez dans votre console R :

Mean

Il vous renvoie :

function (object, ...)

UseMethod("Mean")

<environment: namespace:Hmisc>

Et pas la suite de moyennes attendu. Ceci est dû au fait que "Mean" est, en plus d'être notre liste de moyennes, une fonction interne à R et déjà utilisée. Donc R, quand on appelle "Mean", va tenter d'utiliser la fonction. Or on ne peut pas représenter graphiquement une fonction sans lui donner de données : on crée donc une erreur. Cette erreur nous permet de voir comment préciser que ce que l'on veut représenter appartient à un ensemble plus grand. Par exemple ici, notre ensemble est "results" duquel on veut représenter "Mean". Il faudra donc mettre dans notre commande non pas "Mean" mais "results$Mean", qui est en quelques sortes le chemin interne de notre liste.[/error]

[notice]Je vous invite à jouer avec "pch" en lui donnant des valeurs numériques, mais aussi des symboles ou des lettres ; pch = "+" est très utilisé.[/notice]

Les derniers détails

abline(h=0, lwd = 2)

Ajoutons maintenant une ligne droite pour délimiter l'espace des deux graphiques, avec la fonction "abline" (qui signifie ligne entre le point A et le point B). Nous voulons qu'elle soit horizontale, au point 0 de l'axe des y et qu'elle ait une épaisseur ("lwd = line width", épaisseur de ligne) de 2. Enfin, ajoutons un peu de texte, avec la fonction "mtext" :

mtext("Fraction", side = 2, adj = 0.5, padj = -6)
mtext("Conserved", side = 2, adj = 0.80, padj = -4)
mtext("Unconserved", side = 2, adj = 0.2, padj = -4)

On veut donc ajouter 3 mots, "Fraction", "Conserved" et "Unconserved", du côté gauche (2) du graphique :

  • "adj" sert à placer le texte parallèlement à l'axe auquel le texte se rapporte, donc ici en hauteur
  • "padj" place perpendiculairement à l'axe (d'où le p)

[important]Cette fois, cependant, on se place par rapport à l'extérieur du graphique et non plus par rapport aux valeurs de l'axe, car "mtext" place dans la marge (d'où le m). Alors là, vous allez me dire "Ouais mais moi ça sort de partout et ce n'est pas joli !"... et vous aurez raison. En effet, en écrivant des choses dans la marge, d'autant plus avec des valeurs inverses, sans avoir au préalable défini ces marges est dangereux, et vous fera bien souvent vous planter. [/important]

Nous définirons les marges ainsi :

par(oma = c(1,4,6,4), mar = c(9.1,4.1,1.1,4.1));

Ici, "oma" veut dire "outer margin", "mar" veut dire "margin", et les valeurs sont données toujours dans le même sens "bas, gauche, haut, droite". Maintenant, retapez votre code complet, ce qui devrait être facile si vous avez, au fur et à mesure de l'article, copié votre code dans votre éditeur de texte :

C'est-y pas beau tout ça ? Comment ? Vous voulez mettre ça dans un fichier .png ? Vous êtes bien exigeants... bon, eh bien c'est très simple, il suffit de taper ceci avant le graphe (personnellement, je le fais tout en haut) :

png ("chemin_ou_vous_voulez_le_sauver/nom.png", height = 1600, width = 1600, pointsize = 19);

[important]Pour les valeurs de hauteur ("height") et largeur ("width"), par défaut elles sont données en pixels que je trouve être une bonne unité pour se rendre compte de la qualité et de la lisibilité d'un graphique. Mais vous pouvez bien sur choisir d'autres unités. Enfin, la taille des points du texte est à faire varier en fonction des gouts, en vérifiant bien à chaque changement la position des divers éléments texte, et en sachant qu'elle est calculée en fonction de la résolution de l'image ("?png" pour plus d'informations).[/important]

Alors, ça fonctionne ? Non ? Le fichier est vide ? Ah oui, il faut le fermer, sinon R est toujours entrain d'écrire dedans. Pour cela, on tape :

dev.off();

R doit alors vous indiquer qu'il est revenu dans le cadre de dessin précédent (quartz, x11 ou autres).

Enfin, vu qu'on a attaché "results" au début et qu'on veut prendre de bonnes habitudes, il faut détacher "results" :

detach(result);

Si vous arrivez au même résultat, bravo ! Il ne vous reste plus qu'à vous attaquer à la dernière tâche ! Sinon, n'hésitez pas à poser des questions, je me ferai un plaisir d'y répondre !

À vous de jouer !

Votre défi, si vous l'acceptez, sera alors de faire comme sur la figure en haut de l'article, mettre 2 graphiques dans la même image, du texte en plus, et la tourner dans le bon sens. Pour cela 2 indices :

  • "mfrow" est un paramètre très utile de "par " ;
  • parfois, il vaut mieux se faciliter la vie et faire des choses simples avec un outil simple. Une rotation par exemple.

Le premier qui réussit aura gagné tout mon respect !

Allez, à la prochaine pour des choses encore plus folles !


Je tiens à remercier tous mes relecteurs, Hautbit, nallias, Guillaume Collet, Malicia, Yoann M. et nahoy !



Pour continuer la lecture :


Commentaires

6 réponses à “Cours de R pour débutant pressé, deuxième épisode”

  1. Vous êtes vrai­ment supers!!!! je suis sur qu'avec vous je serai un bio­in­for­ma­ti­cien avant d'avoir eu la licence ! Je par­le­rai de vous sur mon blog!!!!!!!!!!

    Mer­ci!!!!!!!!!!!

  2. Mer­ci à toi Nas­ser, ça me fait plai­sir de faire plai­sir ! Et mer­ci de par­ler de nous sur ton blog !

  3. Tous le plai­sir est pour moi !

  4. Super article 🙂

    Si vous aviez pu être à la place de mon prof il y a 8 ans, j'aimerai p‑e R au lieu de détes­ter ce lan­gage XD

  5. plus de tuto­rial depuis ?
    je suis bio­lo­giste et j'essai dans mes temps mort de me for­mé a R d'autre conseil ?
    site ?

    1. Bon­jour raf,

      mer­ci pour ton com­men­taire. Si tout va pour le mieux, le troi­sième épi­sode devrait arri­ver avant la fin de l'année. L'auteur habi­tuel a eu un très gros plan­ning ces der­niers mois 🙂
      Quant aux dif­fé­rents sites éven­tuels où tu pour­rais trou­ver des infor­ma­tions, je lui laisse le soin de te répondre mieux que moi je ne le ferai puisque je ne mani­pule pas ce "lan­gage" 🙂

Laisser un commentaire