如何通过 R 中的交互进行 LSD 测试?

How to conduct LSD test with interactions in R?

我有一个字段数据

sowing_date<- rep(c("Early" ,"Normal"), each=12)
herbicide<- rep (rep(c("No" ,"Yes"), each=6),2)
nitrogen<- rep (rep(c("No" ,"Yes"), each=3),4)
Block<- rep(c("Block 1" ,"Block 2", "Block 3"), times=8)
Yield<- c(30,27,25,40,41,42,37,38,40,48,47,46,25,27,26,
          41,41,42,38,39,42,57,59,60)

DataA<- data.frame(sowing_date,herbicide,nitrogen,Block,Yield)

我进行了 3 向方差分析

anova3way <- aov (Yield ~ sowing_date + herbicide + nitrogen + 
            sowing_date:herbicide + sowing_date:nitrogen + 
            herbicide:nitrogen + sowing_date:herbicide:nitrogen + 
           factor(Block), data=DataA)
summary(anova3way)

3 个因素之间存在 3 向交互作用。所以,我想看看哪种组合的收益最高。

我知道如何比较均值差与单个因素,如下所示,但如果存在交互作用,我该怎么做?

library(agricolae)
LSD_Test<- LSD.test(anova3way,"sowing_date")
LSD_Test

例如,我想检查 3 种交互作用下的均值差,以及每两个因素之间的交互作用。

例如,我想在 R

中得到这个 LSD 结果

你能告诉我该怎么做吗?

非常感谢,

一种确实需要一些手动工作的方法是将实验参数编码为 -1 和 1,以便正确分离 2 和 3 参数相互作用。
对值进行编码后,您可以从方差分析模型中提取剩余自由度和误差平方和,并将其传递给 LSD.test 函数。

参见下面的示例:

sowing_date<- rep(c("Early" ,"Normal"), each=12)
herbicide<- rep (rep(c("No" ,"Yes"), each=6),2)
nitrogen<- rep (rep(c("No" ,"Yes"), each=3),4)
Block<- rep(c("Block 1" ,"Block 2", "Block 3"), times=8)
Yield<- c(30,27,25,40,41,42,37,38,40,48,47,46,25,27,26,
          41,41,42,38,39,42,57,59,60)
DataA<- data.frame(sowing_date,herbicide,nitrogen,Block,Yield)

anova3way <- aov (Yield ~ sowing_date * herbicide * nitrogen + 
                 factor(Block), data=DataA)
summary(anova3way)

#Encode the experiment's parameters as -1 and 1
DataA$codeSD <- ifelse(DataA$sowing_date == "Early", -1, 1)
DataA$codeherb <- ifelse(DataA$herbicide == "No", -1, 1)
DataA$codeN2 <- ifelse(DataA$nitrogen == "No", -1, 1)

library(agricolae)

LSD_Test<- LSD.test(anova3way, c("sowing_date"))
LSD_Test

#Manually defining the treatment and specifying the
# degrees of freedom and Sum of the squares (Frin the resduals from the ANOVA)
print(LSD.test(y=DataA$Yield, trt=DataA$sowing_date, DFerror=14, MSerror=34.3))

#Example for a two parameter value
print(LSD.test(y=DataA$Yield, trt=(DataA$codeSD*DataA$codeherb), DFerror=14, MSerror=34.3))

print(LSD.test(y=DataA$Yield, trt=(DataA$codeSD*DataA$codeherb*DataA$codeN2), DFerror=14, MSerror=34.3))


#calaculate the means and std (as a check)
#DataA %>% group_by(sowing_date) %>% summarize(mean=mean(Yield), sd=sd(Yield))
#DataA %>% group_by(codeSD*codeherb*codeN2) %>% summarize(mean=mean(Yield), sd=sd(Yield))

您需要手动跟踪最终报告中哪个 run/condition 与 -1 和 1 一致。

编辑: 所以我上面的回答显示了基于交互的整体效果。例如除草剂和氮肥相互作用如何影响产量。

根据您想要确定哪种组合提供最大收益的评论,您再次使用 LSD.test() 函数,但传递参数名称向量。

LSD_Test<- LSD.test(anova3way, c("sowing_date", "herbicide", "nitrogen"))
LSD_Test

从输出的组部分可以看到Normal、Yes 和Yes 是最佳产量。 “组”列将标识唯一的结果集群。例如,最后两行提供了类似的收益。

...
$groups
                  Yield groups
Normal:Yes:Yes 58.66667      a
Early:Yes:Yes  47.00000      b
Normal:No:Yes  41.33333      c
Early:No:Yes   41.00000     cd
Normal:Yes:No  39.66667     cd
Early:Yes:No   38.33333      d
Early:No:No    27.33333      e
Normal:No:No   26.00000      e
...