R语言学习系列27-方差分析(2)

主效应、接着是双因素的交互项，再接着是三因素的交互项。若研究不是正交的，一定要谨慎设置效应的顺序。

有三种类型的方法可以分解y~A+B+A:B右边各效应对y所解释的方差：

类型I（序贯型）

效应根据表达式中先出现的效应做调整。A不做调整，B根据A调整，A:B交互项根据A和B调整。

类型II（分层型）

效应根据同水平或低水平的效应做调整。A根据B调整，B依据A调整，A:B交互项同时根据A和B调整。

类型III（边界型）

每个效应根据模型其他各效应做相应调整。A根据B和A:B做调整，A:B交互项根据A和B调整。

R默认调用类型I方法，其他软件（比如SAS和SPSS）默认调用类型III方法。car包中的Anova()函数（不要与标准anova()函数混淆）提供了使用类型II或类型III方法的选项，而aov()函数使用的是类型I方法。若想使结果与其他软件（如SAS和SPSS）提供的结果保持一致，可以使用Anova()函数。

三、单因素方差分析

1个因变量，1个影响因素：

总差异Yij = 平均差异μ + 因素差异αi + 随机差异εij

例1比较4种品牌的胶合板的耐磨性，各抽取5个样品，相同转速磨损相同时间测得磨损深度（mm），比较4个品牌胶合板的耐磨性有无差异？部分数据如下（ex27_ex1.Rdata）：

setwd(\办公资料/R语言/R语言学习系列/codes\ load(\head(datas) wear brand 1 2.30 A 2 2.32 A 3 2.40 A 4 2.45 A 5 2.58 A 6 2.35 B attach(datas)

table(brand)#各组的样本数 brand A B C D 5 5 5 5

aggregate(wear,by=list(brand),mean)#各组均值 Group.1 x 1 A 2.410 2 B 2.404 3 C 2.046 4 D 2.572

aggregate(wear,by=list(brand),sd)#各组标准差 Group.1 x 1 A 0.11269428 2 B 0.11760102 3 C 0.11216060

4 D 0.03271085 library(car)

qqPlot(lm(wear~brand,data=datas),simulate=TRUE)#用Q-Q图检验数据的正态性

leveneTest(wear~as.factor(brand),data=datas)#方差齐性检验

Levene's Test for Homogeneity of Variance (center = median) Df F value Pr(>F)

group 3 0.6987 0.5664 16

fit<-aov(wear~brand,data=datas)#单因素方差分析（检验组间差异） summary(fit)

Df Sum Sq Mean Sq F value Pr(>F)

brand 3 0.7398 0.24660 24.55 3.15e-06 *** Residuals 16 0.1607 0.01005 ---

Signif.codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

说明：方差齐性检验，原假设H0：方差齐，p值=0.5664>0.05, 故接受原假设，即方差齐。

单因素方差分析结果，brand是因素，Residuals是残差，各列依次为自由度、平方和、均方和、F统计量，p值=3.15e-06<0.05, 拒绝原假设，即不同品牌的磨损（均值）有显著差别。

library(gplots)

plotmeans(wear~brand,xlab=\品牌\磨损\#图形展示带95%置信区间的各组均值

通过前面的分析知道，不同品牌的磨损（均值）有显著差别，但并不知道哪个品牌与其它品牌有显著差别。TukeyHSD()函数提供了对各组均值差异的成对检验。

TukeyHSD(fit)

Tukey multiple comparisons of means 95% family-wise confidence level

Fit: aov(formula = wear ~ brand, data = datas) $brand

difflwrupr p adj

B-A -0.006 -0.18735345 0.1753535 0.9996826 C-A -0.364 -0.54535345 -0.1826465 0.0001610 D-A 0.162 -0.01935345 0.3433535 0.0886142 C-B -0.358 -0.53935345 -0.1766465 0.0001929 D-B 0.168 -0.01335345 0.3493535 0.0744337 D-C 0.526 0.34464655 0.7073535 0.0000019

说明：可以看出（H0：无差异），B与A的差异非常不显著，C与A、C与B、D与C的差异非常显著。

multcomp包中的glht()函数提供了更为全面的多重均值比较方法。

library(multcomp) attach(datas)

tuk<- glht(fit, linfct = mcp(brand = \#注意datas$brand必须是因子型 summary(tuk)

Simultaneous Tests for General Linear Hypotheses

Multiple Comparisons of Means: Tukey Contrasts

Fit: aov(formula = wear ~ brand, data = datas)

Linear Hypotheses:

Estimate Std. Error t value Pr(>|t|) B - A == 0 -0.00600 0.06339 -0.095 0.9997 C - A == 0 -0.36400 0.06339 -5.742 <0.001 *** D - A == 0 0.16200 0.06339 2.556 0.0886 . C - B == 0 -0.35800 0.06339 -5.648 <0.001 *** D - B == 0 0.16800 0.06339 2.650 0.0743 . D - C == 0 0.52600 0.06339 8.298 <0.001 *** ---

Signif.codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 (Adjusted p values reported -- single-step method)

plot(cld(tuk, level = 0.05), col = \

说明：标记相同字母（标记b）的品牌ABD认为是无显著差异，在同一亚组，而品牌C（标记a）与另外三个品牌有显著差异。

另外，也可以进行多重t检验，使用函数：

pairwise.t.test(x, g, p.adjust.method=,...)

其中，x为因变量，g为因子型的分组变量；

百度搜索“77cn”或“免费范文网”即可找到本站免费阅读全部范文。收藏本站方便下次阅读，免费范文网，提供经典小说综合文库R语言学习系列27-方差分析(2)在线全文阅读。