Comment faire un test Tukey HSD avec la commande Anova (car package)
j'ai affaire à un design/échantillon déséquilibré et appris à l'origine aov(). Je sais maintenant que pour mes tests D'ANOVA je dois utiliser la somme de carrés de Type III qui implique l'utilisation ajustement en utilisant lm() plutôt que d'utiliser aov().
le problème est d'obtenir des tests post-hoc (en particulier le HSD de Tukey) en utilisant lm(). Toutes les recherches que j'ai faites ont dit que l'utilisation simint dans le multcomp package fonctionnerait, mais maintenant qu'il est mis à jour cette commande semble ne pas être disponible. Il semble également invoquer en passant par aov() calculer.
essentiellement tous les tests Tukey HSD que j'ai trouvés pour r supposent que vous utilisez aov() pour la comparaison plutôt que lm(). Pour obtenir la somme de type III de carrés dont j'ai besoin pour le design déséquilibré je dois utiliser:
mod<-lm(Snavg~StudentEthnicity*StudentGender)
Anova(mod, type="III")
comment utiliser un test Tukey HSD avec mon mod en utilisant lm()? Ou à l'inverse, calculer mon ANOVA en utilisant le Type III et être encore en mesure d'exécuter un Tukey HSD test?
Merci!
4 réponses
HSD.testagricolae
library(agricolae)
data(sweetpotato)
model<-lm(yield~virus, data=sweetpotato)
comparison <- HSD.test(model,"virus", group=TRUE,
main="Yield of sweetpotato\nDealt with different virus")
Sortie
Study: Yield of sweetpotato
Dealt with different virus
HSD Test for yield
Mean Square Error: 22.48917
virus, means
yield std.err replication
cc 24.40000 2.084067 3
fc 12.86667 1.246774 3
ff 36.33333 4.233727 3
oo 36.90000 2.482606 3
alpha: 0.05 ; Df Error: 8
Critical Value of Studentized Range: 4.52881
Honestly Significant Difference: 12.39967
Means with the same letter are not significantly different.
Groups, Treatments and means
a oo 36.9
ab ff 36.33333
bc cc 24.4
c fc 12.86667
j'étais coincé avec le même problème du HSD.test n'imprimant rien. Vous avez besoin de mettre console=TRUE à l'intérieur de la fonction, de sorte qu'il imprime automatiquement.
Par exemple:
HSD.test(alturacrit.anova, "fator", console=TRUE).
Hope it helps!
j'ai trouvé HSD.test() aussi être très méticuleux au sujet de la façon dont vous avez construit lm() ou aov() modèle que vous utilisez pour cela.
Il n'y a pas de sortie de HSD.test() avec mes données quand j'avais utilisé l'idée suivante de codage pour lm():
model<-lm(sweetpotato$yield ~ sweetpotato$virus)
out <- HSD.test(model,"virus", group=TRUE, console=TRUE)
résultat a été que:
Name: virus
sweetpotato$virus
le résultat était tout aussi mauvais quand on utilisait la même logique pour aov()
model<-aov(sweetpotato$yield ~ sweetpotato$virus)
Pour obtenir la sortie de HSD.test()lm()
(ou également si en utilisant aov() pour le modèle )
doit être construit strictement en utilisant la logique présentée dans la réponse de MYaseen208:
model <- lm(yield~virus, data=sweetpotato)
J'espère que cela aidera quelqu'un qui n'obtient pas une sortie correcte de HSD.test().
comme note initiale, sauf si elle a été changée, pour obtenir les bons résultats pour la somme des carrés de type iii, vous devez définir le codage de contraste pour les variables de facteur. Cela peut être fait à l'intérieur du lm ou options. L'exemple ci-dessous utilise options.
je serais prudent sur l'utilisation de HSD.test et des fonctions similaires avec des conceptions déséquilibrées, à moins que la documentation ne traite de leur utilisation dans ces situations. La documentation pour TukeyHSD mentionne qu'il s'adapte pour les dessins "légèrement déséquilibrés". Je ne sais pas si HSD.test gère les choses différemment. Vous devez vérifier la documentation supplémentaire pour le paquet ou la référence originale Citée pour la fonction.
comme note d'accompagnement, enfermant le tout HSD.test la fonction entre parenthèses les fera imprimer les résultats. Voir l'exemple ci-dessous.
en général, je recommande d'abandonner Tukey.HSD et des fonctions similaires, et en utilisant le flexible emmeans(née lsmeans) ou multcomp paquets pour tous vos besoins de comparaison post-hoc. emmeans est particulièrement utile pour faire séparations moyennes sur les interactions ou examen des contrastes entre les traitements.
avec un design déséquilibré, vous pouvez vouloir rapporter les moyennes E. M. (ou L. S.) au lieu des moyennes arithmétiques. Voir SAEPER: Quel est le moyen le moins carré?. Remarque: dans l'exemple ci-dessous que le marginal moyen déclaré par emmeans sont différentes de celles rapportées par HSD.test.
notez aussi que le "Tukey" dans glht n'a rien à voir avec les comparaisons de Tukey HSD ou de Tukey-adjusted; il n'établit que les contrastes pour tous les tests en paires, comme le dit le résultat.
cependant, le adjust="tukey"emmeans fonctions signifie utiliser des comparaisons ajustées par Tukey, comme le dit la sortie.
L'exemple suivant est partiellement adaptée de ARCHBS: One-way Anova.
if(!require(car)){install.packages("car")}
library(car)
data(mtcars)
mtcars$cyl.f = factor(mtcars$cyl)
mtcars$carb.f = factor(mtcars$carb)
options(contrasts = c("contr.sum", "contr.poly"))
model = lm(mpg ~ cyl.f + carb.f, data=mtcars)
library(car)
Anova(model, type="III")
if(!require(agricolae)){install.packages("agricolae")}
library(agricolae)
(HSD.test(model, "cyl")$groups)
if(!require(emmeans)){install.packages("emmeans")}
if(!require(multcompView)){install.packages("multcompView")}
library(emmeans)
library(multcompView)
marginal = emmeans(model,
~ cyl.f)
pairs(marginal, adjust="tukey")
cld(marginal, adjust="tukey", Letters=letters)
if(!require(multcomp)){install.packages("multcomp")}
library(multcomp)
mc = glht(model,
mcp(cyl.f = "Tukey"))
summary(mc, test=adjusted("single-step"))
multcomp::cld(mc)