R : Test de Kruskal-Wallis en boucle sur les colonnes spécifiées dans la base de données

Question

J'aimerais effectuer un test KW sur certaines variables numériques d'une base de données, en utilisant une variable de regroupement. Je préférerais faire cela dans une boucle, plutôt que de taper tous les tests, car il y a beaucoup de variables (plus que dans l'exemple ci-dessous).

Données simulées :

library(dplyr)
set.seed(123)
Data <- tbl_df(
data.frame(
muttype = as.factor(rep(c("missense", "frameshift", "nonsense"), each = 80)),
ados.tsc   = runif(240, 0, 10),
ados.sa    = runif(240, 0, 10),
ados.rrb   = runif(240, 0, 10))
) %>%
group_by(muttype)
ados.sim <- as.data.frame(Data)

Le code suivant fonctionne parfaitement en dehors de la boucle.

kruskal.test(formula (paste((colnames(ados.sim)[2]), "~ muttype")), data = 
ados.sim)

Mais ce n'est pas le cas à l'intérieur de la boucle :

for(i in names(ados.sim[,2:4])){  
ados.mtp <- kruskal.test(formula (paste((colnames(ados.sim)[i]), "~ muttype")), 
data = ados.sim)
}

J'obtiens l'erreur suivante :

Erreur dans terms.formula(formula, data = data) : terme non valide dans la formule du modèle

Quelqu'un sait-il comment résoudre ce problème ? Merci beaucoup !

Answer 1

Essayez :

results <- list()
for(i in names(ados.sim[,2:4])){  
  results[[i]] <- kruskal.test(formula(paste(i, "~ muttype")), data = ados.sim)
}

Cela permet également d'enregistrer les résultats dans une liste et d'éviter d'écraser les résultats au fur et à mesure qu'ils se présentent. ados.mtp à chaque itération, ce qui, je pense, n'est pas ce que vous vouliez faire.

Il convient de noter ce qui suit :

for(i in names(ados.sim[,2:4])){  
   print(i)
}
[1] "ados.tsc"
[1] "ados.sa"
[1] "ados.rrb"

C'est-à-dire, i vous donne déjà le nom de la colonne. Le problème dans votre code est que vous avez essayé de l'utiliser comme un entier pour le sous-ensemble, ce qui a transformé le résultat en NA .

for(i in names(ados.sim[,2:4])){  
   print(paste((colnames(ados.sim)[i]), "~ muttype"))
}
[1] "NA ~ muttype"
[1] "NA ~ muttype"
[1] "NA ~ muttype"

Et pour mémoire, tout cela pourrait également être fait de deux manières différentes, que je préfère souvent car elles facilitent légèrement l'analyse ultérieure :

Tout d'abord, stockez tous les objets de test dans un cadre de données :

library(tidyr)
df <- ados.sim %>% gather(key, value, -muttype) %>% 
      group_by(key) %>% 
      do(test = kruskal.test(x= .$value, g = .$muttype))

Vous pouvez ensuite subdiviser la base de données pour obtenir les résultats des tests :

df[df$key == "ados.rrb",]$test
[[1]]

    Kruskal-Wallis rank sum test

data:  .$value and .$muttype
Kruskal-Wallis chi-squared = 2.2205, df = 2, p-value = 0.3295

Il est également possible d'obtenir tous les résultats directement dans un cadre de données, sans stocker les objets de test :

library(broom)
df2 <- ados.sim %>% gather(key, value, -muttype) %>% 
       group_by(key) %>% 
       do(tidy(kruskal.test(x= .$value, g = .$muttype)))
df2
# A tibble: 3 x 5
# Groups:   key [3]
       key statistic   p.value parameter                       method
     <chr>     <dbl>     <dbl>     <int>                       <fctr>
1 ados.rrb 2.2205031 0.3294761         2 Kruskal-Wallis rank sum test
2  ados.sa 0.1319554 0.9361517         2 Kruskal-Wallis rank sum test
3 ados.tsc 0.3618102 0.8345146         2 Kruskal-Wallis rank sum test