DellaData

Transformez vos données en connaissances

News et découvertes de mars 2022

Voici la récap des news et découvertes qui ont attiré mon attention au cours du mois de mars. Je les ai regroupées en 5 thématiques : Biostatistiques, Recherche reproductible, Recherche clinique, Data et Utilisation de R en sociologie. Bonne lecture !

Biostatistique

De l’intérêt de visualiser avant d’interpréter

Sur twitter je me suis abonnée au compte de Dr. Guy Prochilo (@GuyProchilo), qui est maître de conférence en méthodes de recherche en psychologie. Il partage des réflexions et du contenu autour des biostatistiques, des investigations cliniques, des essais contrôlés randomisés et de R. C’est toujours intéressant et pertinent.

Il a, par exemple, partagé le graphique ci- dessous qui montre qu’un même coefficient de corrélation linéaire peut correspondre à des situations extrêmement différentes (dont certaines incompatibles avec la linéarité). Cela illustre parfaitement qu’il faut toujours visualiser les données avant d’interpréter ce paramètre.

Et si vous voulez briller en salle café de votre labo 😉, voici le code pour reproduire les graphs:

library(devtools)
install_github("janhove/cannonball")
library(cannonball)
plot_r(r = 0.5, n = 35)

Obtenir des tables de résultats avec le package parameters

J’ai redécouvert le package parameters (que j’avais un peu oublié). Ce package permet d’obtenir des tables de résultats pour des tests, et des modèles de régressions (linéaire, anova, glm, modèles mixtes). Par exemple avec un test t:

library(funModeling)
library(tidyverse)
library(parameters)

# test t 
t.test(max_heart_rate~has_heart_disease, data=heart_disease) %>% 
    parameters()
## Welch Two Sample t-test
## 
## Parameter      |             Group | has_heart_disease = no | has_heart_disease = yes | Difference |         95% CI | t(272.27) |      p
## ----------------------------------------------------------------------------------------------------------------------------------------
## max_heart_rate | has_heart_disease |                 158.38 |                  139.26 |      19.12 | [14.33, 23.91] |      7.86 | < .001
## 
## Alternative hypothesis: true difference in means between group no and group yes is not equal to 0

Ici, avec une ANOVA à 2 facteurs, en ajoutant l’effect size eta squared :

ToothGrowth$dose <- as.factor(ToothGrowth$dose)
aov(len~dose*supp, data=ToothGrowth) %>% 
    parameters(eta_squared="partial",
               type=2)
## Parameter | Sum_Squares | df | Mean_Square |     F |      p | Eta2 (partial)
## ----------------------------------------------------------------------------
## dose      |     2426.43 |  2 |     1213.22 | 92.00 | < .001 |           0.77
## supp      |      205.35 |  1 |      205.35 | 15.57 | < .001 |           0.22
## dose:supp |      108.32 |  2 |       54.16 |  4.11 | 0.022  |           0.13
## Residuals |      712.11 | 54 |       13.19 |       |        |               
## 
## Anova Table (Type 2 tests)

En cas de régression, lorsque l’hypothèse d’homogénéité des résidus est rejetée, il est possible d’utiliser des estimateurs sandwich.

Dans l’exemple ci-dessous, on pourra voir que les résidus augmentent lorsque la réponse (Y) augmente.

mydf<- structure(list(Y = c(0.13, 0.15, 0.05, 0.11, 0.1, 0.05, 
    0.1, 0.02, 0.14, 0.07, 0.07, 0.12, 0.12, 0.33, 0.14, 0.2, 0.16, 
    0.64, 0.33, 1.11, 0.26, 0.36, 0.39, 0.59, 0.94),X= c(4.36, 
    3.91, 3.98, 4.24, 3.84, 4.32, 4.65, 4.89, 4.36, 5.1, 4.91, 6.19, 
    8.05, 4.15, 6.1, 5.32, 5.39, 6.67, 5.42, 10.42, 7.95, 5.99, 5.92, 
    7.2, 8.77)), class = "data.frame", row.names = c(NA, -25L))


library(ggplot2)
ggplot(mydf, aes(x=X ,y=Y))+
    geom_point()+
    geom_smooth(method="lm")

C’est ce que nous montre aussi le graph “Homogeneity of variance” , généré par la fonction check_model().

fit <- lm(Y ~ X, data = mydf)
library(performance)
check_model(fit)

Et le test de Breush Pagan rejette l’hypothèse d’homogénéité :

check_heteroscedasticity(fit)
## Warning: Heteroscedasticity (non-constant error variance) detected (p = 0.001).

Dans cette situation, on peut alors utiliser un estimateur sandwich de la matrice de variance co-variance pour estimer l’erreur standard des paramètres :

# avec l'estimateur classique
model_parameters(fit)
## Parameter   | Coefficient |   SE |         95% CI | t(23) |      p
## ------------------------------------------------------------------
## (Intercept) |       -0.46 | 0.13 | [-0.73, -0.20] | -3.61 | 0.001 
## X           |        0.13 | 0.02 | [ 0.08,  0.17] |  5.93 | < .001


# avec un estimateur sandwich de la matrice de variance-covariance
model_parameters(fit, vcov="HC")
## Parameter   | Coefficient |   SE |         95% CI | t(23) |     p
## -----------------------------------------------------------------
## (Intercept) |       -0.46 | 0.19 | [-0.85, -0.08] | -2.48 | 0.021
## X           |        0.13 | 0.04 | [ 0.05,  0.21] |  3.43 | 0.002

On voit ici que les estimations des erreurs standard des paramètres est plus élevée par l’estimateur sandwich.

Vous trouverez plus d’infos sur cette page : https://easystats.github.io/parameters/articles/model_parameters_robust.html, notamment des estimateurs de type sandwich pour les modèles mixtes :

library(lme4)
fit <- lmer(Reaction~Days+(1|Subject), data=sleepstudy)

check_heteroscedasticity(fit)
## Warning: Heteroscedasticity (non-constant error variance) detected (p < .001)

# estimation classique
model_parameters(fit)
## # Fixed Effects
## 
## Parameter   | Coefficient |   SE |           95% CI | t(176) |      p
## ---------------------------------------------------------------------
## (Intercept) |      251.41 | 9.75 | [232.17, 270.64] |  25.79 | < .001
## Days        |       10.47 | 0.80 | [  8.88,  12.05] |  13.02 | < .001
## 
## # Random Effects
## 
## Parameter               | Coefficient
## -------------------------------------
## SD (Intercept: Subject) |       37.12
## SD (Residual)           |       30.99




# estimateur de type sandwich
model_parameters(
  fit,
  vcov = "vcovCR",
  vcov_args = list(type = "CR1", cluster = sleepstudy$Subject)
)
## # Fixed Effects
## 
## Parameter   | Coefficient |   SE |           95% CI | t(176) |      p
## ---------------------------------------------------------------------
## (Intercept) |      251.41 | 6.82 | [237.94, 264.87] |  36.84 | < .001
## Days        |       10.47 | 1.55 | [  7.42,  13.52] |   6.77 | < .001
## 
## # Random Effects
## 
## Parameter               | Coefficient
## -------------------------------------
## SD (Intercept: Subject) |       37.12
## SD (Residual)           |       30.99

Remarque : C’est juste un essai, je ne sais pas si c’est pertinent !

Données manquantes

J’ai découvert cette publication par un tweet de Guy Prochilo.

Je l’ai parcourue et je la trouve vraiment pertinente avec des explications simples des 3 types de données manquantes (MAR : Missing at Random, MCAR : Missing Completely at random et MNAR : missing not at random) et des conséquences en termes de biais et de précision, et en fonction de 5 approches possibles (Listewise deletion, pairwise deletion, imputation, maximum likelihood, multiple imputation).

Cliquez ici pour accéder au pdf de la publication.

Recherche reproductible

Toujour sur Twitter, grâce à Yannick Morvan (@yannick_morvan), j’ai découvert 3 ressources, qui vous seront certainement très utiles pour travailler de manière reproductible et collaborative.

Bonnes pratiques pour l'analyse de données scientifiques

Dans cette publication, vous trouverez de très bons conseils répartis sur 6 axes :

Gestion des données : sauvegarde des données originales et intermédiaires, documenter toutes les étapes, créer des données ordonnées susceptibles d’être analysées.
Logiciel : écriture, organisation et partage des scripts et des programmes utilisés dans l’analyse.
Collaboration : faciliter la compréhension et la contribution des collaborateurs existants et nouveaux à un projet
Organisation du projet : organiser les artefacts numériques d’un projet pour faciliter la découverte et la compréhension
Suivi des modifications : enregistrement de la façon dont les divers composants de votre projet changent au fil du temps
Manuscrits : la rédaction de manuscrits de manière à laisser une piste d’audit et à minimiser la fusion manuelle des conflits

Apprendre à utiliser git pour versionner ses scripts

Ce document, rédigé par Daniel Lakens, explique en pas à pas, avec de nombreuses images, comment coupler git à RStudio pour versionner vos documents sous R, en local et à distance sur GitHub ou GitLab.

Ça rejoint le point 5. de la publication précédente.

Le package WORCS

Pour aller encore plus loin, vous pouvez organiser votre travail de recherche reproductible sous R avec le package WORCS (The Workflow for Open Reproducible Code in Science (WORCS).

Pour vous faire une idée des possibilités, je vous conseille cette courte vidéo :

Et pour une démonstration, cette autre vidéo :

Vous trouverez encore des informations sur la page du package : https://cjvanlissa.github.io/worcs/index.html

Recherche clinique

Comment bien construire un modèle prédictif

Voici une série de publications concernant la réalisation de modèles prédictifs, partagés sur Twitter par Gary Collins (@GSCollins)

Data

Si vous recherchez des données environnementales pour vous entraîner ou pour enseigner, l’organisation LTER (Long Term Ecological Research programm) vient de proposer un package contenant 9 jeux de données environnementales, destinés à l’enseignement.

Vous trouverez plus d’informations sur cette page https://lter.github.io/lterdatasampler/

Dans l’onglet article, vous trouverez des informations sur chaque jeu de données, ainsi que des courtes analyses, avec le code R correspondant.

Le package s’installe depuis un répertoire GitHub :

# install.packages("remotes")
remotes::install_github("lter/lterdatasampler")
library(lterdatasampler)
View(hbr_maples)

Utilisation de R en Sociologie

J’ai découvert le blog de Joel Gires (@JoelGires) qui est chercheur en sociologie. Vous pourrez trouver des scripts R, développés dans le cadre d’analyses sociologiques, principalement quantitatives.