因此，假设有必要测试 anova 的正态性假设（参见1和2）
如何在 R 中进行测试？
我希望做类似的事情：
## From Venables and Ripley (2002) p.165.
utils::data(npk, package="MASS")
npk.aovE <- aov(yield ~ N*P*K + Error(block), npk)
residuals(npk.aovE)
qqnorm(residuals(npk.aov))
这不起作用，因为对于重复测量 anova 的情况，“残差”没有方法（也没有预测）。
那么在这种情况下应该怎么做呢？
可以从没有误差项的相同拟合模型中简单地提取残差吗？我对文献不够熟悉，不知道这是否有效，提前感谢您的任何建议。
您可能无法得到简单的响应，residuals(npk.aovE)但这并不意味着该对象中没有残差。执行str并查看在级别内仍然存在残差。我想您对“内部”级别最感兴趣
> residuals(npk.aovE$Within)
7 8 9 10 11 12
4.68058815 2.84725482 1.56432584 -5.46900749 -1.16900749 -3.90234083
13 14 15 16 17 18
5.08903669 1.28903669 0.35570336 -3.27762998 -4.19422371 1.80577629
19 20 21 22 23 24
-3.12755705 0.03910962 2.60396981 1.13730314 2.77063648 4.63730314
我自己的训练和实践不是使用正态性检验，而是使用 QQ 图和稳健方法的并行检验。
另一种选择是使用包的lme功能nlme（然后将获得的模型传递给anova）。您可以residuals在其输出上使用。
Venables 和 Ripley 在他们的书（第 284 页）后面关于随机和混合效应的部分中解释了如何对重复测量设计进行残差诊断。
为每个层的 aov 结果实现了残差函数（或 resid）：
从他们的例子：
oats.aov <- aov(Y ~ N + V + Error(B/V), data=oats, qr=T)
要获得拟合值或残差：
“因此fitted(oats.aov[[4]])和resid(oats.aov[[4]])是长度为 54 的向量，表示最后一层的拟合值和残差。”
重要的是，他们补充说：
“不可能将它们与原始实验的图唯一地联系起来。”
对于诊断，他们使用投影：
plot(fitted(oats.aov[[4]]), studres(oats.aov[[4]]))
abline(h=0, lty=2)
oats.pr <- proj(oats.aov)
qqnorm(oats.pr[[4]][, "Residuals"], ylab = "Stratum 4 residuals")
qqline(oats.pr[[4]][, "Residuals"])
正如另一位用户发布的那样，他们还表明该模型可以使用 lme 完成。
这里有两个选项，aov 和 lme （我认为第二个是首选）：
require(MASS) ## for oats data set
require(nlme) ## for lme()
require(multcomp) ## for multiple comparison stuff
Aov.mod <- aov(Y ~ N * V + Error(B/V), data = oats)
the_residuals <- aov.out.pr[[3]][, "Residuals"]
Lme.mod <- lme(Y ~ N * V, random = ~1 | B/V, data = oats)
the_residuals <- residuals(Lme.mod)
原始示例没有交互（Lme.mod <- lme(Y ~ N * V, random = ~1 | B/V, data = oats)），但它似乎正在使用它（并产生不同的结果，所以它正在做某事）。
就是这样……
但为了完整性：
1 - 模型摘要
summary(Aov.mod)
anova(Lme.mod)
2 - 重复测量 anova 的 Tukey 测试（3 小时寻找这个！！）。
summary(Lme.mod)
summary(glht(Lme.mod, linfct=mcp(V="Tukey")))
3 - 正态性和同方差图
par(mfrow=c(1,2)) #add room for the rotated labels
aov.out.pr <- proj(aov.mod)
#oats$resi <- aov.out.pr[[3]][, "Residuals"]
oats$resi <- residuals(Lme.mod)
qqnorm(oats$resi, main="Normal Q-Q") # A quantile normal plot - good for checking normality
qqline(oats$resi)
boxplot(resi ~ interaction(N,V), main="Homoscedasticity",
xlab = "Code Categories", ylab = "Residuals", border = "white",
data=oats)
points(resi ~ interaction(N,V), pch = 1,
main="Homoscedasticity", data=oats)
