我正在处理一个不平衡的设计/样本，最初学习了aov()。我现在知道，对于我的方差分析测试，我需要使用III型平方和，其中包括使用lm()而不是aov()进行拟合。
问题是使用lm()进行事后测试(特别是Tukey的HSD)。我所做的所有研究都表明，在multcomp包中使用simint可以工作，但现在它更新了，该命令似乎不可用。它似乎也依赖于通过aov()来计算。
基本上，我为R找到的所有Tukey HSD测试都假设您使用aov()而不是lm()进行比较。为了获得III类平方和，我需要使用不平衡设计，我必须使用:
mod<-lm(Snavg~StudentEthnicity*StudentGender)
Anova(mod, type="III")
我如何使用Tukey HSD测试与我的mod使用lm() ?或者相反，使用III型计算我的方差分析，仍然能够运行Tukey HSD测试?
谢谢!
在agricolae中尝试HSD.test
library(agricolae)
data(sweetpotato)
model<-lm(yield~virus, data=sweetpotato)
comparison <- HSD.test(model,"virus", group=TRUE,
main="Yield of sweetpotatonDealt with different virus")
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
作为初始提示，除非已更改，否则要获得类型iii平方和的正确结果，需要为因子变量设置对比度编码。这可以在lm调用或options中完成。下面的示例使用options。
我会谨慎使用HSD.test和类似的功能与不平衡的设计，除非文档说明它们在这些情况下的使用。TukeyHSD的文档提到它对"轻度不平衡"进行了调整。设计。我不知道HSD.test是否有不同的处理方式。您必须检查包的附加文档或函数引用的原始参考。
作为旁注，将整个HSD.test函数括在括号中将导致它打印结果。
emmeans (nlsmeans)或multcomp包来满足您所有的事后比较需求。emmeans在进行相互作用的平均分离或检查处理之间的对比时特别有用。[编辑:警告，我是这些页面的作者。]
对于不平衡设计，您可能希望报告E.M.(或L.S.)均值而不是算术均值。参见SAEPER:最小二乘是什么?[编辑:警告，我是这个页面的作者。]注意，在下面的例子中，emmeans报告的边际均值与HSD.test报告的边际均值是不同的。
还要注意"Tukey"在glht中与Tukey HSD或Tukey调整的比较无关;它只是为所有成对测试设置对比，如输出所示。
然而，emmeans函数中的adjust="tukey"确实意味着使用tukey调整的比较，如输出所示。
以下示例部分改编自ARCHBS: One-way Anova。
### EDIT: Some code changed to reflect changes to some functions
### in the emmeans package
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")}
library(emmeans)
marginal = emmeans(model,
~ cyl.f)
pairs(marginal, adjust="tukey")
if(!require(multcomp)){install.packages("multcomp")}
library(multcomp)
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"))
cld(mc)
我发现HSD.test()对于您构建lm()或aov()模型的方式也非常细致。
当我使用以下lm()编码思想时，HSD.test()没有输出我的数据:
model<-lm(sweetpotato$yield ~ sweetpotato$virus)
out <- HSD.test(model,"virus", group=TRUE, console=TRUE)
输出只有:
Name: virus
sweetpotato$virus
当对aov()
model<-aov(sweetpotato$yield ~ sweetpotato$virus)
要获得HSD.test()的输出，请输入lm()(或者如果使用aov()模型)必须严格使用MYaseen208答案中提供的逻辑构造:
model <- lm(yield~virus, data=sweetpotato)
希望这能帮助那些没有从HSD.test()得到正确输出的人。
我遇到了HSD的相同问题。测试打印不出任何东西。您需要将console=TRUE放入函数中，以便它自动打印出来。
HSD.test(alturacrit.anova, "fator", console=TRUE).
Hope it helps!