14

我正在处理一个不平衡的设计/样本，最初是从aov(). 我现在知道，对于我的 ANOVA 测试，我需要使用类型 III 平方和，这涉及使用拟合lm()而不是使用aov().

问题是使用lm(). 我所做的所有研究都表明simint在multcomp包中使用会起作用，但现在它已更新，该命令似乎不可用。它似乎也依赖于经过aov()计算。

基本上我为 R 找到的所有 Tukey HSD 测试都假设您aov()用于比较而不是lm(). 要获得不平衡设计所需的 III 型平方和，我必须使用：

mod<-lm(Snavg~StudentEthnicity*StudentGender)

Anova(mod, type="III")

我如何使用 Tukey HSD 测试与我的 mod 使用lm()？或者相反，使用 III 型计算我的 ANOVA 并且仍然能够运行 Tukey HSD 测试？

谢谢！

4

4 回答 4

9

试穿HSD.test_agricolae

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")

输出

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
于 2011-10-12T00:51:35.293 回答
5

作为初始说明，除非已更改，否则要获得类型 iii 平方和的正确结果，您需要为因子变量设置对比度编码。这可以在lm调用内部或使用options. 下面的示例使用options.

我会谨慎使用HSD.test具有不平衡设计的类似功能，除非文档说明了它们在这些情况下的使用。文档TukeyHSD提到它针对“轻度不平衡”设计进行了调整。我不知道HSD.test处理方式是否不同。您必须检查包的附加文档或该函数引用的原始参考。

附带说明一下，将整个HSD.test函数括在括号中会导致它打印结果。请参见下面的示例。

一般来说，我建议您使用灵活的emmeans(née lsmeans) 或multcomp软件包来满足您的所有事后比较需求。对于对相互作用进行平均分离或检查治疗之间的对比emmeans特别有用。[编辑：请注意，我是这些页面的作者。]

对于不平衡设计，您可能希望报告 EM（或 LS）均值而不是算术均值。请参阅SAEPER：什么是最小二乘均值？. [编辑：警告我是本页的作者。] 请注意，在下面的示例中，报告的边际均值emmeans不同于HSD.test.

另请注意，“Tukey”中的“Tukey”glht与 Tukey HSD 或 Tukey-adjusted 比较无关；正如输出所示，它只是为所有成对测试设置了对比。

然而，adjust="tukey"inemmeans函数确实意味着使用 Tukey-adjusted 比较，如输出所示。

以下示例部分改编自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)
于 2017-11-27T14:09:42.450 回答
3

我还发现您对构建用于它的模型或模型HSD.test()的方式非常细致。lm()aov()

HSD.test()当我使用以下编码思想时，我的数据没有输出lm()：

    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()，必须严格使用 MYaseen208 答案中提供的逻辑构建：lm()aov()

    model <- lm(yield~virus, data=sweetpotato)

希望这可以帮助那些没有从HSD.test().

于 2016-10-07T11:09:19.593 回答
1

我遇到了同样的问题，即 HSD.test 什么也没打印出来。您需要放入console=TRUE函数内部，因此它会自动打印出来。

例如：

HSD.test(alturacrit.anova, "fator", console=TRUE).
Hope it helps!
于 2014-10-28T03:33:08.620 回答