经管之家App
让优质教育人人可得
立即打开
什么是方差分析？
方差分析是基于平方和分解的一种统计方法，其目的在与推断两组或多组资料的均值是否相等，检验两个或多个均值的差异是否在统计意义上显著。
（1）单因素方差分析（One-way ANOVA）
平方和分解：
总平方和(SST)=组内平方和(SSE)+组间平方和(SSB) （均是离均差平方和）
自由度 n-1 = n-g + g-1
离均差平方和自能反应变异的绝对大小，变异程度除与离均差平方和的大小有关外，还与其自由度有关。引入均方差(MS)来反映变异程度：
组内均方差MSE=SSE/(n-g)
组间均方差MSB=SSB/(g-1)
构造F统计量： F=MSB/MSE~F(g-1,n-g)
H0：不存在处理效应
H1：存在处理效应
即：
H0:mu1=mu2=…=mug
H1:mu1、mu2、…、mug不全部相等
若H0成立则组内变异与组间变异都只反应随机误差大小，MSE与MSB相当，F接近1。
若H0不成立，即存在处理效应，则组内变异仍只反映随机误差的大小，但组间变异不仅反映随机误差大小还有处理效应，此时组内均方差MSB>组内均方差MSE。因此若F远大于1，即MSB远大于MSE即认为存在处理效应，H0不成立。
应用条件：（1）个观测值相互独立(这一条一般跟实验设计有关，数据处理时不做检验，只要实验设计得当，这一条一般认为满足）
（2）服从正态分布
（3）方差齐性
Example1:
为研究煤矿粉尘作业环境对尘肺的影响，将18只大鼠随机分到A、B、C 3个组，每组6只，分别在地面办公楼、煤炭仓库和矿井下染尘，12周后测量大鼠全肺湿重（g），数据如下表，问不同环境下大鼠全肺湿重有无差别？
Rcode:
####(1)one-way ANOVA
#创建数据
XA=c(4.2,3.3,3.7,4.3,4.1,3.3)
XB=c(4.5,4.4,3.5,4.2,4.5,4.2)
XC=c(5.6,3.6,4.5,5.1,4.9,4.7)
X=c(XA,XB,XC)
Treat=gl(3,6,label=c("A","B","C"))
#正态性检验
qqnorm(X)
qqline(X)
library(fBasics)
normalTest(X)
jarqueberaTest(X)
#方差齐性检验
bartlett.test(X~Treat)
#先画个箱线图看看各组均值是否大致在一个水平上
plot(X~Treat)
#方差分析
fit=aov(X~Treat)
summary(fit)
####模型诊断
layout(matrix(c(1,2,3,4),2,2)) # optional layout
plot(fit)
两种正态性检验的结果，p值均大于0.1，所以不能否定数据服从正态分布
Bartlett检验结果，p值大于0.1，不能否定各组数据方差相等。
箱线图显示各组均值存在较明显差异，有可能存在处理效应
ANOVA结果显示，P值小于0.05，，所以在alpha=0.05的显著性水平上，我们认为H0不成立，即有95%把握认为存在处理效应，不同环境下大鼠全肺湿重有差别。
模型诊断图：上左和上右表示残差的波动程度基本不随观测值而变化，下左看残差是否服从正态分布，下右各个因子水平上残差波动程度相当；上左、上右，下右是对方差齐性的诊断
，下左是对正态性的诊断。
（2）two-way ANOVA（双因素方差分析）
<1>不可重复（无交互效应）
平方和分解：
总离差平方和(SST)=随机误差平方和(SSE)+行因素误差平方和（SSR)+列因素误差平方和（SSC)
自由度关系：
kr-1= (k-1)(r-1) + k-1 + r-1
MSE=SSE/(k-1)(r-1) MSR=SSR/(k-1) MSC=SSC/(r-1)
行因素检验统计量： F=MSR/MSE~F(k-1，(k-1)(r-1))
列因素检验统计量： F=MSC/MSE~F(r-1，(k-1)(r-1))
Example2
研究A、B、C三种营养素对小白鼠体重增加的影响，已知窝别为影响因素。拟用6窝小白鼠，每窝3只，随机地安排喂养A、B、C三种营养素之一种，8周后观察小白鼠体重增加情况，数据如下表。问：（1）不同营养素之间小白鼠的体重增加是否不同？（2）不同窝别之间小白鼠的体重增加是否不同？
Rcode:
##############创建数据集
XA=c(64,53,71,41,50,42)
XB=c(65,54,68,46,58,40)
XC=c(73,59,79,38,65,46)
X=c(XA,XB,XC)
nut=gl(3,6,label=c("A","B","C"))
loc=gl(6,1,length=18)
#############正态性检验及方差齐性检验#################
normalTest(X)
##############画箱线图初步分析营养液和窝点两因素各自的效应
par(mfrow = c(1, 2))
plot(X~nut+loc)
#############ANOVA
fit2=aov(X~nut+loc)
summary(fit2)
##########ANOVA模型诊断图
par(mfrow=c(2,2))
plot(fit2)
P>0.1，正态性检验通过
行因素列因素方差齐性检验p值均大于0.1，方差齐性检验通过。
箱线图结果表明行因素列因素均对X产生影响（由图可知均值不在一个level上）
Anova结果显示行因素（loc)列因素（nut）F检验P值均小于0.05，故认为行因素和列因素对X均有显著影响。
==========================================================
关于我们：
创作团队：数据加工厂
淘宝店地址：http://shop114096406.taobao.com/?spm=a230r.7195193.1997079397.40.asBLlF
或者在淘宝首页搜店铺“数据加工厂”，我们的logo如下：
<2>可重复（交互效应）
若“重复”m次实验，行因素有r个level，列因素有k个level
SST = SSE + SSR + SSC + SSRC
总离差平方和 = 误差项平方和 + 行因素平方和 + 列因素平方和 + 交互作用平方和
自由度：
n-1 = kr(m-1) + k-1 + r-1 + (k-1)(r-1)
MSE=SSE/kr(m-1) MSR=SSR/(k-1) MSC=SSC/(r-1) MSRC=SSRC/(k-1)(r-1)
行因素检验统计量： F=MSR/MSE~F(k-1，kr(m-1))
列因素检验统计量： F=MSC/MSE~F(r-1，kr(m-1))
交互效应检验统计量：F=MSRC/MSE~F((k-1)(r-1),kr(m-1))
Example3
有一个关于检验毒品强弱的试验, 给48只老鼠注射I、 II、 III三种毒药(行因素), 同时有A、 B、 C、 D 4种治疗方案(列因素), 这样的试验在每一种因素组合下都重复四次测试老鼠的存活时间, 数据如下表. 试分析毒药和治疗方案以及它们的交互作用对老鼠存活时间有无显著影响。
Rcode:
#################创建数据集
rats<-data.frame(
Time=c(0.31, 0.45, 0.46, 0.43, 0.82, 1.10, 0.88, 0.72, 0.43, 0.45,
0.63, 0.76, 0.45, 0.71, 0.66, 0.62, 0.38, 0.29, 0.40, 0.23,
0.92, 0.61, 0.49, 1.24, 0.44, 0.35, 0.31, 0.40, 0.56, 1.02,
0.71, 0.38, 0.22, 0.21, 0.18, 0.23, 0.30, 0.37, 0.38, 0.29,
0.23, 0.25, 0.24, 0.22, 0.30, 0.36, 0.31, 0.33),
Toxi=gl(3, 16, 48, labels = c("I", "II", "III")),
Cure=gl(4, 4, 48, labels = c("A", "B", "C", "D"))
)
###############正态性检验
normalTest(rats$Time)
####################不满足正态进行box-cox变换
boxcox(rats$Time~rats$Toxi*rats$Cure)
boxcox(rats$Time~rats$Toxi*rats$Cure,lambda = seq(-1, 0, 0.05))
box_cox <- function(x, lambda=1) {
if (lambda==0) {
stopifnot(all(x)>0)
log(x)}
else {((x)^lambda - 1)/lambda}
}
Time.hat=box_cox(rats$Time,-0.8)
normalTest(Time.hat)
####检查每个格子里的数据方差是否大致相等（方差齐性检验）
leveneTest(Time.hat~rats$Toxi*rats$Cure)
#########将box-cox变换后的数据与行因素列因素重新整合到一个数据集中
Rats=data.frame(time=Time.hat,toxi=rats$Toxi,cure=rats$Cure)
############画箱线图初步分析毒药和治疗法案两因素各自的效应
par(mfrow=c(1, 2))
plot(time~toxi+cure, data=Rats)
####交互作用图，可初步判断，一般若相交则存在交互作用，若平行则不存在交互作用
par(mfrow=c(1, 2))
with(Rats,interaction.plot(toxi, cure, time, trace.label="cure"))
with(Rats,interaction.plot(cure, toxi, time, trace.label="toxi"))
######ANOVA
fit3<-aov(time~toxi*cure, data=Rats)
summary(fit3)
#####诊断图
par(mfrow=c(2,2))
plot(fit3)
P值小于0.1，正态性检验没通过，需进行box-cox变换
由极大似然函数图知在lambda=-0.8时似然函数值最大。取lambda=-08进行box-cox变换。
变换后数据p>0.1，通过正态性检验
变换后数据也通过了方差齐性检验（若通不过，一般进行平方根、对数、反正弦等变换）
箱线图表明行因素列因素均对结果有显著影响（均值不在一个level上）
交互作用图没有显示明显的平行于交叉，还需通过ANOVA进行更精确的分析。
ANOVA结果表明行因素（toxi)列因素（cure)均对结果有显著影响（P<0.1)
但交互作用（toxi:cure)不显著（p>0.1)
创作团队：数据加工厂
淘宝店地址：http://shop114096406.taobao.com/?spm=a230r.7195193.1997079397.40.asBLlF
或者在淘宝首页搜店铺“数据加工厂”，我们的logo如下：
<2>可重复（交互效应）
若“重复”m次实验，行因素有r个level，列因素有k个level
SST = SSE + SSR + SSC + SSRC
总离差平方和 = 误差项平方和 + 行因素平方和 + 列因素平方和 + 交互作用平方和
自由度：
n-1 = kr(m-1) + k-1 + r-1 + (k-1)(r-1)
MSE=SSE/kr(m-1) MSR=SSR/(k-1) MSC=SSC/(r-1) MSRC=SSRC/(k-1)(r-1)
行因素检验统计量： F=MSR/MSE~F(k-1，kr(m-1))
列因素检验统计量： F=MSC/MSE~F(r-1，kr(m-1))
交互效应检验统计量：F=MSRC/MSE~F((k-1)(r-1),kr(m-1))
Example3
有一个关于检验毒品强弱的试验, 给48只老鼠注射I、 II、 III三种毒药(行因素), 同时有A、 B、 C、 D 4种治疗方案(列因素), 这样的试验在每一种因素组合下都重复四次测试老鼠的存活时间, 数据如下表. 试分析毒药和治疗方案以及它们的交互作用对老鼠存活时间有无显著影响。
Rcode:
#################创建数据集
rats<-data.frame(
Time=c(0.31, 0.45, 0.46, 0.43, 0.82, 1.10, 0.88, 0.72, 0.43, 0.45,
0.63, 0.76, 0.45, 0.71, 0.66, 0.62, 0.38, 0.29, 0.40, 0.23,
0.92, 0.61, 0.49, 1.24, 0.44, 0.35, 0.31, 0.40, 0.56, 1.02,
0.71, 0.38, 0.22, 0.21, 0.18, 0.23, 0.30, 0.37, 0.38, 0.29,
0.23, 0.25, 0.24, 0.22, 0.30, 0.36, 0.31, 0.33),
Toxi=gl(3, 16, 48, labels = c("I", "II", "III")),
Cure=gl(4, 4, 48, labels = c("A", "B", "C", "D"))
)
###############正态性检验
normalTest(rats$Time)
####################不满足正态进行box-cox变换
boxcox(rats$Time~rats$Toxi*rats$Cure)
boxcox(rats$Time~rats$Toxi*rats$Cure,lambda = seq(-1, 0, 0.05))
box_cox <- function(x, lambda=1) {
if (lambda==0) {
stopifnot(all(x)>0)
log(x)}
else {((x)^lambda - 1)/lambda}
}
Time.hat=box_cox(rats$Time,-0.8)
normalTest(Time.hat)
####检查每个格子里的数据方差是否大致相等（方差齐性检验）
leveneTest(Time.hat~rats$Toxi*rats$Cure)
#########将box-cox变换后的数据与行因素列因素重新整合到一个数据集中
Rats=data.frame(time=Time.hat,toxi=rats$Toxi,cure=rats$Cure)
############画箱线图初步分析毒药和治疗法案两因素各自的效应
par(mfrow=c(1, 2))
plot(time~toxi+cure, data=Rats)
####交互作用图，可初步判断，一般若相交则存在交互作用，若平行则不存在交互作用
par(mfrow=c(1, 2))
with(Rats,interaction.plot(toxi, cure, time, trace.label="cure"))
with(Rats,interaction.plot(cure, toxi, time, trace.label="toxi"))
######ANOVA
fit3<-aov(time~toxi*cure, data=Rats)
summary(fit3)
#####诊断图
par(mfrow=c(2,2))
plot(fit3)
P值小于0.1，正态性检验没通过，需进行box-cox变换
由极大似然函数图知在lambda=-0.8时似然函数值最大。取lambda=-08进行box-cox变换。
变换后数据p>0.1，通过正态性检验
变换后数据也通过了方差齐性检验
箱线图表明行因素列因素均对结果有显著影响（均值不在一个level上）
交互作用图没有显示明显的平行于交叉，还需通过ANOVA进行更精确的分析。
ANOVA结果表明行因素（toxi)列因素（cure)均对结果有显著影响（P<0.1)
但交互作用（toxi:cure)不显著（p>0.1)
============================================================================================================================================
创作团队：数据加工厂
淘宝店地址：http://shop114096406.taobao.com/?spm=a230r.7195193.1997079397.40.asBLlF
（可复制打开）
或者在淘宝首页搜店铺“数据加工厂”，我们的logo如下：
扫码加好友，拉您进群
全部版块
我的主页
收藏
