Matrice de corrélation, et R

Plus j’avance dans R, plus je trouve cet outil fabuleux.

En mai-juin dernier, une élève kiné vint me voir pour que je lui fasse des stats avec des petits p pour son mémoire final (je trouve ça parfaitement débile, mais c’est un autre débat).

Je n’avais pas touché un logiciel de statistiques (piraté, d’ailleurs) depuis ma thèse, et encore, j’avais mis à contribution fait faire tout le boulot à un interne du DIM de mon CHU.

Autrement dit, je partais de très très loin.

Je suis allé cherché le meilleur test adapté à ses données, puis je suis tombé sur un os, n’ayant plus de logiciel à disposition. Finalement, je m’en suis sorti avec je ne sais plus quel site qui proposait des feuilles Excel® permettant de faire le test voulu.

Je ne sais pas si le test était bon, en tout cas, le jury de distingués kinés, qui doit être aussi calé en statistiques que moi (huhuhuhuhuhu…), n’a pas élevé de protestations devant mon Wilcoxon.

Influencés par le Scrabble®, avec un W et un X dans un nom propre, ils ont dû être impressionnés, tout comme moi je le fus. La prochaine fois, je leur sortirai un test de Wald–Wolfowitz.

Au cours de mes recherches, j’ai découvert R, logiciel gratuit, communautaire, surpuissant et totalement imbitable.

Après un ou deux MOOC, le survol de dizaines de ressources disponibles sur la toile, et une rencontre, j’ai commencé à utiliser R, presque avec plaisir.

Je triche quand même un peu car, toujours à la suite de discussions sur Twitter, un pharma (@PotardDechaine ou @PierrikFaure ?) m’a fait découvrir RStudio, une interface qui rend R un peu moins aride.

Les choses simples sont difficiles sur R (ne serait-ce que pour visualiser ses données…), et les choses complexes sont très simples à faire. Pour corser le tout, certaines opérations ne sont possibles qu’après avoir téléchargé et installé des packages.

Bref, R se mérite, mais on peut faire des choses sympas.

Petit exemple.

Imaginez 65 individus chez lesquels ont été mesurées 7 variables quantitatives continues. Vous souhaitez savoir si certaines de ces variables sont corrélées entre elles.

Mon fichier texte (R ne lit que les fichiers textes, qui sont simples à obtenir à partir d’un fichier Excel®) s’appelle JMV.

Une seule ligne de commande, cor(JMV), permet de sortir une matrice de corrélation:

matcor8On ne peut pas faire beaucoup plus simple, non?

Par exemple, le coefficient de corrélation entre Var1 et Var2 est de 0.84 et des poussières.

On peut imaginer qu’un graphique sera plus parlant:

matcorr9Quelques lignes de code permettent d’obtenir un tel graphique:

>base <- JMV

>base.r <- abs(cor(base))

>cpairs(base, gap = .5)

Elles ne s’inventent pas, mais je n’ai pas eu trop de difficultés pour les trouver.

C’est joli, mais presque plutôt moins informatif que la matrice de corrélation, car on n’arrive toujours pas à discerner quelles sont les variables qui sont bien corrélées entre elles, et surtout si on pourrait en grouper certaines.

Point curieux, comme me l’a fait remarquer @Potarddechaine, s’agissant de coefficients de corrélation, ils devraient être compris entre -1 et +1, or l’échelle montre de toutes autres valeurs (67, 68, 69, 70…). Quelles sont-elles? Quels sont leurs réseaux?

Une autre série de lignes de code permet de classer ces corrélations et de les colorier en fonction de leur intensité (rouge si coefficient de corrélation>+/-0.5, bleu entre +/-0.5 et +/-0.18, blanc en dessous:

corrmtx3On commence à pouvoir grouper certaines variables en fonction de leurs corrélations (par exemple Var1, Var5 et Var2). Je vous avais déjà montré ce type de graphique qui m’avait beaucoup impressionné. Mais ça, c’était avant.

Aujourd’hui, j’ai découvert le package corrplot.

Et là, R devient sublime.

Voici quelques graphiques que j’ai tiré sans trop de difficultés de ma série de données:

corrmtx4corrplot1Par exemple, pour obtenir le graphique ci-dessus, le code à entrer dans R est « tout simplement »:

>corrplot(JMV_corr, method = « number »)

corrplot3

On peut ainsi regrouper les coefficients de corrélation (-1 à +1 😉 ), les illustrer de différents manières, même coupler les graphiques avec un test statistique, bref, rendre intelligibles des données en les rendant belles.

Un MOOC francophone: Virchow-Villerme/Fondamentaux en statistique

Je me suis inscrit au MOOC « Fondamentaux en statistique » de la plateforme FUN par curiosité, et je trouve l’expérience plutôt sympa. Les MOOC francophones se développent, j’avais déjà parlé de celui de la Lorraine University ici.

L’enseignant principal se nomme Avner Bar-Hen, et j’ai eu la très heureuse surprise d’échanger quelques tweets avec lui hier. Pour connaître la démarche ayant conduit à l’organisation de ce MOOC, je vous suggère ce lien.

La plateforme utilisée est une adaptation de celle de EdX, donc pas de mauvaise surprise. Certains devoirs doivent être déposés sur une plateforme externe, mais là-aussi, pas de souci particulier.

Les vidéos de cette première semaine sont claires nettes et précises.

J’ai été très surpris par le format très court des vidéos (autour de 8 minutes), alors que jusqu’à présent mes précédents MOOC m’avaient plutôt habitué à une durée de 15-30 minutes. Difficile de tirer une conclusion au bout de 6 vidéos, mais ce format ultra-court est pas mal et correspond bien à mon emploi du temps un peu trépidant.

Les quiz sont sympa, et l’exercice d’application m’a fait remettre le nez dans la syntaxe de R, que j’avais presque totalement oubliée après 2-3 semaines de non-utilisation.

Par contre, ce R, quel plaisir de sortir un graphique en boxplot avec une seule ligne de commande: boxplot().

Sur Excel, et bien, c’est tout simplement imbitable… L’exercice d’application m’a aussi fait découvrir la commande quantile().

Comme pour les MOOC que j’ai pu faire précédemment, il faut ne pas se laisser abuser par le « pas/peu de pré-requis ». ll faut avoir l’esprit curieux et savoir ouvrir un autre onglet pour aller chercher son bonheur sur Google. Différents fora de discussion et un wiki permettent et permettront à l’étudiant inscrit de se retrouver un peu si il est perdu.

Sinon, mon opinion des MOOC n’a pas tellement changé.

Les MOOC suivent clairement une courbe de type Technology Hype de Gartner.

Cette démarche est intéressante pour l’esprit curieux qui n’a pas besoin de diplôme, pour le reste (et l’avenir), Wait and See….

Pour conclure, j’aime beaucoup cette phrase de Avner Bar-Hen:

Les MOOCs : c’est une nouvelle ruée vers l’or, on se souviendra plus de la ruée que de l’or.

Ça, c’est chez nous…

Néanmoins, je ne suis pas certain que quelques universités américaines ne fassent pas quand même quelques dollars en vendant à plusieurs milliers d’étudiants avides de « diplômes » US des certificats « vérifiés ».

La preuve, le compte Instagram d’un statisticien américain de l’école bayésienne qui organise des MOOC, d’un cardiologue interventionnel du sud de la France, d’un baron de la drogue mexicain.

Un merveilleux exemple de Narcocorrido.

Mouhahahahahahahaha!

Le problème du char d’assaut allemand

Ce titre de note est tout à fait digne des Monty Python. C’est aussi un problème statistique (je sais, je suis un peu monomaniaque en ce moment…) assez rigolo.

Imaginons qu’une guerre moderne oppose, disons les bleus contre les rouges. Les rouges planifient une attaque massive qui va reposer essentiellement sur l’utilisation de blindés, mais ils savent que les bleus ont mis en production depuis quelques mois un char d’assaut qui surpasse tous les leurs.

Le point crucial, pour les rouges, est de savoir combien les bleus pourront aligner de ces chars de nouvelle génération lors de cette attaque.

Problème, de nombreux espions Bothan sont morts sans pouvoir ramener cette information.

Les espions rouges connaissent néanmoins la date de début de leur production, et aussi savent que chaque char d’assaut possède un numéro de série unique qui s’incrémente de 1 à chaque char produit: par exemple 1, 2, 3, ….

Les rouges ont réussi à capturer (« aléatoirement ») 5 chars plus ou moins intacts, avec 5 numéros de série. Par exemple 5, 48, 69, 110, et 16.

Pouvons-nous en déduire le nombre total de chars de nouvelle génération produits par les bleus, et ainsi leur production mensuelle?

Et bien, comme vous pouvez vous en douter, on peut estimer ce chiffre (avec une marge d’erreur potentiellement très raisonnable).

La formule à utiliser est même remarquablement simple:

n=m(1 + (1/k)) – 1

n est l’estimation du nombre de chars produits, k le nombre d’exemplaires capturés, et m le numéro de série le plus élevé observé sur ces exemplaires.

Dans notre cas, n=110(1+(1/5))-1 soit 131 blindés. Si la production a commencé il y a 2 mois, les bleus produisent donc environ 65 blindés par mois.

Bon, arrivés à ce point, vous devez vous dire que c’est n’importe quoi, que vous ayez ou pas des connaissances en statistiques, car ce chiffre est invérifiable.

(Ceux qui connaissent les statistiques peuvent me jeter des pierres, je ne fais pas bien mieux qu’un singe savant).

C’est là que cette histoire peut devenir drôle.

Imaginons que les bleus aient produit effectivement 342 chars. Imaginons que les rouges aient été un peu meilleurs, et qu’ils aient capturé 25 chars (plus l’échantillon capturé est important, meilleure est la précision).

Je vais demander à R de me donner 25 numéros de série aléatoires, parmi une population en comprenant 342.

En langage R, ça donne cela:

R GTP1Les lignes 1 et 2 permettent de créer une population pop de 342 numéros de série qui se suivent.

La ligne 3 permet d’extraire de pop un échantillon x de 25 numéros, que voici:

R GTP2Les lignes 4 et 5 définissent m et k.

Les deux dernières définissent n (est.n) qui est ici de…337.

Pas mal, non?

(pour ceux qui dorment, la bonne réponse est 342.)

Ok, vous ne me croyez toujours pas?

Je vais demander à R de créer 10.000 échantillons aléatoires x de 25 numéros tirés de pop (pop=342, k=25):

R GTP3

R GTP4La moyenne observée des populations prédites sur ces 10.000 échantillons est de… 342.07.

😉

Un petit coup de hist(sim.est.n) pour avoir un histogramme:

Rplot01(L’histogramme est pourri mais je suis incapable de configurer R pour le rendre acceptable.)

On peut aussi voir que la formule du problème du char d’assaut allemand (ce n’est pas son vrai nom, mais vous verrez à la fin pourquoi le char est allemand, et pas bleu) peut donner une approximation raisonnable de la population de pop en considérant la corrélation entre des populations pop et ces mêmes populations estimées par la formule, sur de nombreux essais.

Pour avoir une plus jolie corrélation, j’ai utilisé les paramètres suivants: k=25, pop autour de 100 et 1000 essais:

Rplot02Joli, non?

Pour en savoir plus sur le problème du char allemand (et pourquoi ce nom):

  • Le site d’ où j’ai tiré les codes de R.
  • Randall Munroe (xkcd) fait allusion succinctement à ce problème dans cette note.

Statistics One: le bilan

Je viens de terminer l’examen final de Statistics One, un MOOC de Princeton sur les statistiques. J’ai déjà parlé plusieurs fois de ce MOOC, mais une petite synthèse me paraît utile.

Comme beaucoup de MOOC anglo-saxons, Statistics One est gratuit. Par contre, aucun certificat n’est délivré à son issue. Donc ceux qui en font la chasse doivent passer leur chemin. Les quiz et les deux examens ont d’ailleurs un nombre illimité d’essais (en fait 100…), ce qui permet, même simplement par déduction, de répondre correctement aux QCM.

La note n’a donc aucune importance. Ce n’était pas le cas pour le MOOC de Stanford, où j’ai bien plus transpiré pour un joli certificat. De ce point de vue, j’ai pris bien plus de plaisir à obtenir ce bout de papier car j’ai retrouvé mes sensations (ça fait très sportif, pour un homme de salon comme moi) d’étudiant avant un examen/concours.

Les cours de Princeton sont intuitifs et pédagogiques. Les TD sur R, après un début brutal, sont en fait assez intéressants et représentent un bon début d’initiation (je reste très modeste) à l’utilisation de ce logiciel qui est quand même diabolique, même pour des pros des statistiques. Je suis toujours impressionné de lire dans des articles de Significance que des auteurs ont transpiré à grosses gouttes pour sortir un graphique.

Le seul point noir de ce MOOC est un certain laxisme de la part de l’équipe enseignante qui transparaît notamment dans les quiz. Pas mal de quiz comportaient une erreur ou une difficulté stupide qui empêchait d’avoir la note maximale.

Imaginez la scène se répétant presque chaque semaine: le quiz est un prolongement du TP avec quelques questions théoriques, mais surtout des questions sur l’utilisation pratique d’une fonction de R. Vous faites tourner R (en copiant/collant le TP, ce qui évite pas mal de migraines) et vous obtenez un résultat numérique. Vous êtes contents, mais le quiz répond que le résultat est faux. Perplexité, on recommence, c’est encore faux… Vous vous jetez sur les fora: en fait, il ne faut pas répondre 0.68 mais .68 (à l’américaine). Bon, ça, ça va encore. Mais d’autres fois, vous vous rendez compte qu’il faut répondre la valeur de la médiane, alors que c’est la moyenne qui est demandée, ou un p issu d’un test de Tukey pour la population A, alors que c’est la population B qui est demandée…

La plupart du temps, ce n’est pas le staff qui corrige le problème, mais d’autres étudiants qui trouvent le truc.

C’est énervant, mais rien de bien terrible.

Revoir des notions de bases et manipuler des statistiques, surtout via R, qui est peu commode, a été pour moi une bonne expérience.

Ces statistiques sont finalement tombées à un moment parfait pour moi. Un changement professionnel a fait exploser mon nombre de soucis à gérer et s’effondrer mon temps d’esprit libre. Je n’ai pour l’instant pas de regrets, et j’ai trouvé dans les statistiques un dérivatif sur lequel me fixer, totalement différent de ce que je fais toute la journée.

L’an prochain, je ferai un autre MOOC de statistiques, pour me rafraîchir l’esprit. Je vous raconterai!

(J’ai aussi trouvé pour le printemps 2014 un MOOC qui me parle bien: Unethical Decision Making in Organizations de l’Université de Lausanne.)