R包rcompanion执行非参数双因素方差分析(Scheirer-Ray-Hare检验)
前文已经分别简介了单因素方差分析、单因素协方差分析的非参数方法在R语言中的执行过程,本篇将继续简介双因素方差分析(双因素ANOVA)的非参数方法,Scheirer-Ray-Hare检验。
本文使用的作图数据的网盘链接(提取码fjuq):
https://pan.baidu.com/s/1I2dX58-q5XJgexxfYTvpmg
示例数据说明
我们首先将示例数据读到R中,并从中挑选部分数据作为演示(测试数据同前文中的双因素ANOVA,以方便二者比较)。
soil <- read.table('soil.txt', sep = '\t', header = TRUE, stringsAsFactors = FALSE, check.names = FALSE)
group <- read.table('group.txt', sep = '\t', header = TRUE, stringsAsFactors = FALSE, check.names = FALSE)
soil <- merge(soil, group, by = 'sample')
#以 chao1 指数为例,同时将分组列转换为因子变量
chao1 <- soil[ ,c('sample', 'treat', 'times', 'chao1')]
chao1$treat <- factor(chao1$treat)
chao1$times <- factor(chao1$times)
str(chao1)
head(chao1)
假设存在这么一个关注的问题:
我们采集了来源于同一环境中的土壤(土壤类型一致),分为了3组,分别添加了三种类型的化学物质(a、b、c),并将土壤孵育了不同的时间(10、15、20、25 天)。在不同时间段收集样本后,通过16S测序,获得了土壤细菌群落的Alpha多样性指数,期望关注化学物质类型以及处理时间对土壤细菌群落的影响(关注交互作用)。
对应于上述挑选出的测试数据“chao1”:sample,试验样本名称,土壤类型完全相同;treat,在土壤中添加的三种化学物质(a、b、c),这列作为分组列,需要转换为因子变量类型;times,土壤孵育时间(天数),这列作为分组列,需要转换为因子变量类型(尽管它本来是数值类型的,但必需转化为因子类型后,才可作为分组变量用于方差分析);chao1,Alpha多样性指数中的Chao1指数,数值变量。
此处存在“化学物质类型”以及“土壤孵育时间”两组分组变量,对应于双因素,同时还需要关注二者的交互作用。在前文方差分析中,我们已经使用双因素ANOVA探究了化学物质类型以及处理时间是否对土壤细菌群落产生了显著影响。在这里,我们更换为Scheirer-Ray-Hare检验来实现。
Scheirer-Ray-Hare检验
双因素ANOVA建立在两组数据必须均服从正态分布,以及各组方差相等的基础上才能执行,如前文双因素ANOVA中的方法所述。若二者之一不符合,Scheirer-Ray-Hare检验将会是一种实用的方法。(尽管本示例数据符合方差分析的前提假设,但这里仅作为演示非参数检验,故请忽略这点)
Scheirer-Ray-Hare检验可使用rcompanion包scheirerRayHare()函数来实现,接下来我们使用scheirerRayHare()完成上述统计。式中表达方式类似于方差分析的aov()函数,写作scheirerRayHare(y~A*B)的样式,表示考虑所有可能的交互项:A、B以及A和B的交互(A:B),其中A、B分别为两组因子变量。表达式中各效应存在主次顺序,主因素(treat)在前,次因素(times)在后。
#rcompanion 包 scheirerRayHare(),详情使用 ?scheirerRayHare 查看帮助说明library(rcompanion)
scheirerRayHare(chao1~treat*times, data = chao1)
计结果屏幕输出如下,包含两组因素及其交互作用的检验结果。
通过对比前文的双因素ANOVA结果(本篇和前篇的所用数据是一致的),我们发现,前文检测到交互项显著,而这里却没有,侧面反映了这类非参数的方法较为保守。因此正如本系列开篇时所提到的,当数据满足ANOVA的条件时,尽可能使用ANOVA分析,能够有效地鉴别出非参数检验鉴别不到的差异;当无法适用ANOVA时,再考虑非参数的方法。
此外,对于双因素中交互项的可视化,同样可参考前文双因素ANOVA中的方法。
R语言执行非参数单因素方差分析(Kruskal-Wallis检验、Friedman检验)