R中的3因子嵌套方差分析

Question

我正在尝试在一篇论文中复制 3 因子嵌套 ANOVA 分析：Underwood, AJ (1993) The Mechanics of spatially replicated sampling programs to detect environment impact in a variable world.

示例的数据（来自表 3，Underwood 1993）可以通过以下方式生成：

dat <-
structure(list(B = structure(c(2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L
), .Label = c("A", "B"), class = "factor"), C = structure(c(2L,
2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L,
1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L,
1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L,
1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L,
2L, 1L, 1L, 1L, 1L, 1L, 1L), .Label = c("C", "I"), class = "factor"),
    Times = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L,
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L,
    3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L,
    1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L,
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L,
    4L, 4L), .Label = c("1", "2", "3", "4"), class = "factor"),
    Locations = c(1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L,
    1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L,
    3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L,
    2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L,
    1L, 2L, 2L, 2L, 3L, 3L, 3L, 1L, 1L, 1L, 2L, 2L, 2L, 3L, 3L,
    3L), X = c(59L, 51L, 45L, 46L, 40L, 32L, 39L, 32L, 25L, 51L,
    44L, 37L, 55L, 47L, 41L, 31L, 38L, 45L, 41L, 47L, 55L, 43L,
    36L, 29L, 23L, 30L, 37L, 57L, 50L, 43L, 36L, 44L, 51L, 39L,
    29L, 23L, 38L, 44L, 52L, 31L, 38L, 45L, 42L, 35L, 28L, 52L,
    44L, 37L, 51L, 43L, 37L, 38L, 31L, 24L, 60L, 52L, 46L, 30L,
    37L, 44L, 41L, 34L, 27L, 53L, 46L, 39L, 40L, 34L, 26L, 21L,
    27L, 35L), Times.unique = structure(c(5L, 5L, 5L, 5L, 5L,
    5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L,
    7L, 7L, 7L, 7L, 7L, 7L, 7L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L,
    8L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L,
    2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L,
    4L, 4L, 4L, 4L, 4L, 4L, 4L), .Label = c("A_1", "A_2", "A_3",
    "A_4", "B_1", "B_2", "B_3", "B_4"), class = "factor")), .Names = c("B",
"C", "Times", "Locations", "Y", "Times.unique"), row.names = c(NA,
-72L), class = "data.frame")

dat

数据框dat有4个因素：

B - 有两个级别“A”和“B”（在v之前）

Times - 8 个级别，4 个在“B”之前，4 个在“A”之后，每个级别编码为 1:4。请注意，变量 Times.unique 是相同的，但每次（之前和之后）都有一个唯一的代码

位置 - 具有三个级别，每次都在之前和之后测量

C - 有两个级别控制 (C) 和 (I)。注意：两个位置是控制，一个是影响

虽然我很清楚如何使用混合模型 (lmer) 分析这样的设计，但我想准确地复制他的示例，以便我可以运行一些模拟来比较他的方法。

特别是我试图复制表 4 中“a”列下的 SS 值。他适合具有以下项的 SS 和 df 值的设计：

B -> SS = 66.13，df = 1

次(B) -> SS = 280.64, df = 6

位置 -> SS = 283.86，df = 2

B x 位置 -> SS = 29.26，df = 2

Times(B) x Locations-> SS = 575.45, df = 12

残差 -> SS = 2420.00, df = 48

总计 -> SS = 6208.34，df = 71

我假设 Times(B) 术语表示嵌套在处理前后“B”中的时间。对于这个例子，他忽略了位置来自控制和影响处理，并完全忽略了因素 C。

我已经尝试了所有我能想到的组合来重现这个嵌套的方差分析，在 B 中使用独特的时间编码和时间编码为 1:4（之前和之后）。我尝试使用 %in%、/ 和 Error() 参数，以及来自汽车的 Anova 来更改计算的 SS 类型。 %in% 和 / 嵌套拟合的示例包括：

aov(Y~B+Locations+Times%in%B+B:Locations+Times%in%B:Locations, data=dat)
aov(Y~B+Locations+B/Times+B:Locations+B/Times:Locations, data=dat)

我似乎无法准确复制安德伍德的 SS 值，尤其是对于两个交互项。有朋友让我在statistix中拟合模型，SS值可以精确重现，所以可以得到这个模型的上述SS值。

谁能帮我在 R 中拟合这个模型？我希望将其嵌入到更大的模拟中，并且真的需要能够在 R 中运行模型，以便准确再现 Underwood 1993 SS 值？

Answer 1

您的问题是 dat$Locations 是一个整数，而它应该是一个因子（三个唯一位置）。一个提示是，您的 ANOVA 线认为 Locations 只占用 1 df，而 Underwood 给它 2。

只需添加一行：

dat$Locations = factor(dat$Locations)

然后您的代码行完美地再现了安德伍德的结果：

aov(Y~B+Locations+B/Times+B:Locations+B/Times:Locations, data=dat)
#Call:
#   aov(formula = Y ~ B + Locations + B/Times + B:Locations + B/Times:Locations, 
#    data = dat)
#
#Terms:
#                        B Locations   B:Times B:Locations B:Locations:Times
#Sum of Squares    66.1250 2836.8611  280.6389     29.2500          575.4444
#Deg. of Freedom         1         2         6           2                12
#                Residuals
#Sum of Squares  2420.0000
#Deg. of Freedom        48