R语言方差分析的注意事项

文章目录
R语言做方差分析很简单 , 就是一个函数 aov() , 包括但不限于单因素方差分析、多因素方差分析、协方差分析、重复测量方差分析等 , 都是这个函数 。
前面用一篇推文详细介绍了R语言中方差分析的各种实现方法:
R语言方差分析最全总结
R语言做方差分析和SPSS/SAS等传统统计软件不太一样 , 下面说一下需要注意的地方 , 主要是2个点:
均衡设计和非均衡设计
均衡设计是指不同组别之间的样本量相等 , 非均衡设计自然就是指不同组别之间样本量不相同 。
医学研究中大部分都是均衡设计 。
方差分析的3种类型
在计算方差分析中的平方和时 , 有3种类型(你可以简单理解为方差分析有3种类型) , SPSS/SAS在做方差分析的时候 , 默认是类型Ⅲ , 但是R语言中的aov()函数做方差分析时 , 默认是类型Ⅰ 。
R语言中做方差分析是公式表示的 , 比如:aov(y ~ A + B + A:B, data = http://www.kingceram.com/post/df) 。
表达式中效应的顺序在两种情况下会造成影响:(a)因子不止一个 , 并且是非平衡设计;(b)存在协变量 。出现任意一种情况时 , 等式右边的变量都与其他每个变量相关 。此时 , 我们无法清晰地划分它们对因变量的影响 。一般来说 , 越基础性的效应越需要放在表达式前面 。具体来讲 , 首先是协变量 , 然后是主效应 , 接着是双因素的交互项 , 再接着是三因素的交互项 , 以此类推 。对于主效应 , 越基础性的变量越应放在表达式前面 , 因此性别要放在处理方式之前 。有一个基本的准则:若研究设计不是正交的(也就是说 , 因子和/或协变量相关) , 一定要谨慎设置效应的顺序 。–《R语言实战》
也就是说:
3种类型的区别可以参考下面这张图:
R语言的aov()函数不能更改类型 , 但是我们通过其他R包实现更改类型 。比如car::Anova()或者包 。

R语言方差分析的注意事项

文章插图
示例
使用3个简单的小例子进行演示 。
one-way anova
首先是一个单因素(one-way anova)均衡设计的例子 , 来自课本例4-2 。
trt<-c(rep("group1",30),rep("group2",30),rep("group3",30),rep("group4",30))weight<-c(3.53,4.59,4.34,2.66,3.59,3.13,3.30,4.04,3.53,3.56,3.85,4.07,1.37,3.93,2.33,2.98,4.00,3.55,2.64,2.56,3.50,3.25,2.96,4.30,3.52,3.93,4.19,2.96,4.16,2.59,2.42,3.36,4.32,2.34,2.68,2.95,2.36,2.56,2.52,2.27,2.98,3.72,2.65,2.22,2.90,1.98,2.63,2.86,2.93,2.17,2.72,1.56,3.11,1.81,1.77,2.80,3.57,2.97,4.02,2.31,2.86,2.28,2.39,2.28,2.48,2.28,3.48,2.42,2.41,2.66,3.29,2.70,2.66,3.68,2.65,2.66,2.32,2.61,3.64,2.58,3.65,3.21,2.23,2.32,2.68,3.04,2.81,3.02,1.97,1.68,0.89,1.06,1.08,1.27,1.63,1.89,1.31,2.51,1.88,1.41,3.19,1.92,0.94,2.11,2.81,1.98,1.74,2.16,3.37,2.97,1.69,1.19,2.17,2.28,1.72,2.47,1.02,2.52,2.10,3.71)df <- data.frame(trt,weight)str(df)## 'data.frame': 120 obs. of2 variables:##$ trt: chr"group1" "group1" "group1" "group1" ...##$ weight: num3.53 4.59 4.34 2.66 3.59 3.13 3.3 4.04 3.53 3.56 ...
R语言中的方差分析结果比较:
【R语言方差分析的注意事项】# 1型fit <- aov(weight ~ trt, data = http://www.kingceram.com/post/df)summary(fit)##Df Sum Sq Mean Sq F valuePr(>F)## trt332.1610.71924.88 1.67e-12 ***## Residuals11649.970.431## ---## Signif. codes:0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# 2型car::Anova(fit, type=2)## Anova Table (Type II tests)## ## Response: weight##Sum SqDf F valuePr(>F)## trt32.156324.884 1.674e-12 ***## Residuals 49.967 116## ---## Signif. codes:0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# 3型car::Anova(fit, type=3)## Anova Table (Type III tests)## ## Response: weight##Sum SqDf F valuePr(>F)## (Intercept) 353.021 819.537