3.2两因素等重复试验下的方差分析
课 时 授 课 计 划
课次序号:15
一、课 题:§3.2 两因素等重复试验下的方差分析
二、课 型:新授课
三、目的要求:1.掌握两因素等重复实验下的方差分析理论与方法、模型的建立
与显著性检验;
2.掌握利用方差分析的SAS过程解决有关实际问题.
四、教学重点:方差分析方法的基本理论;利用方差分析的SAS过程解决有关
实际应用问题.
教学难点:方差分析方法的基本理论;利用方差分析的SAS过程解决有关
实际应用问题. 五、教学方法及手段:传统教学与上机实验相结合.
六、参考资料:
《应用多元统计分析》,高惠璇编,北京大学出版社,2005; 《使用统计方法与SAS系统》,高惠璇编,北京大学出版社,2001; 《多元统计分析》(二版),何晓群编,中国人民大学出版社,2008; 《应用回归分析》(二版),何晓群编,中国人民大学出版社,2007; 《统计建模与R软件》,薛毅编著,清华大学出版社,2007. 七、作业:3.6 八、授课记录:
九、授课效果分析:
复习 单因素方差分析 1.统计模型
因变量Y—因素A,水平A1,A2, ,Aa,Ai上观测值yi1,yi2, ,yin1
i
⎧
⎪yij=μ+δi+εij,j=1,2, ,ni,⎪⎪2
⎨εij~N(0,σ),且各εij相互独立⎪a
⎪∑nδ=0
ii
⎪⎩i=1
a
ni
2
a
ni
2
i=1,2, ,a
SST=∑∑(yij-)=∑∑(yij-i∙)+∑ni(i∙-)2=SSE+SSA
i=1j=1
i=1j=1
i=1
a
SSE
σ2
~χ(n-a),σ=
H0真
2
∧
2
SSE
n-a
SSA
σ2
~χ2(a-1)
2.显著检验
H0:μ1=μ2= =μa 或 H0:δ1=δ2= =δa=0
SSA/(a-1)MSAH0为真
F==~F(a-1,n-a)
SSE/(n-a)MSE
p=PH0(F≥f)=P(F(a-1,n-a)≥f)
3.置信区间
ˆi=i∙=∑yij/ni (i=1,2, ,a) μ
j=1
ni
μi置信度1-α的置信区间 i∙-t
⎝
⎛
1-
α(n-a)MSE/ni,2
i∙+t
⎫ ⎪,(n-a)MS/nαEi⎪1-
2⎭
μi-μj置信度1-α的置信区间为
⎛
i∙-j∙-tα(n-a)(1+1)MSE, 1-n1n2
2⎝
i∙-j∙+t
α(n-a)(2
1-
11
+)MSEn1n2
⎫⎪ ⎪⎭
m个μi-μj的置信度至少1-α的Bonferroni同时置信区间
⎛
i∙-j∙-tα(n-a)(1+1)MSE, 1-n1n2
2m⎝
i∙-j∙+t
(n-a)(
11
+)MSEn1n2
⎫⎪ ⎪⎭
1-
α2m
§3.2 两因素等重复实验下的方差分析
3.2.1 统计模型
设影响Y的因素有两个,分别记为A和B,其中A有a个不同水平A1,A2, ,Aa,B有b个水平B1,B2, ,Bb.在因素A和B的各水平下均做c(c>1)次实验,记yijk为水平组合(Ai,Bj)下第k次实验的Y的观测值,则两因素等重复试验下的方差分析数据可表示为表3.7的形式.
表3.7 两因素等重复方差分析数据
对于任一水平组合(Ai,Bj)(总体), Y观测值为yij1,yij2, ,yijc,
则各样本间是相互独立的.样本观察值yijk可看成是来自均值(i=1,2, ,a;j=1,2, ,b),为μij的总体,即
yijk~N(μij,σ2), k=1,2, ,c
令εijk=yijk-μij,εij为水平组合(Ai,Bj)下Y的随机误差, 则εijk~N(0,σ2), 这样yijk=μij+εijk,就是其均值μij与随机误差εijk迭加而产生的.
因此,两因素重复试验下方差分析的统计模型:
⎧⎪yijk=μij+εijk,
⎨2⎪⎩εijk~N(0,σ),
i=1,2, ,a,j=1,2, ,b,k=1,2, ,c
且诸εijk相互独立
为便于统计分析,我们需要对水平组合(Ai,Bj)上的样本均值作进一步分解,为此引入如下记号:
1ab
μ=∑∑μij
abi=1j=1
1b
μi∙=∑μij,αi=μi∙-μ,i=1,2, ,a
bj=1
μ∙j
1a
=∑μij,βj=μ∙j-μ,aj=1
j=1,2, ,b
其中μ为总平均, n=abc,μij是因素Ai水平与因素Bj水平在ij单元上所有观察值的平均,αi为因素A的水平Ai的效应,βj为因素B的水平Bj的效应.
γij=μij-μi∙-μ∙j+μ,i=1,2, ,a,
j=1,2, ,b
进一步有γij=μij-μ-(μi∙-μ)-(μ∙j-μ)=(μij-μ)-(αi+βj)其差γij称为Ai与Bj的交互效应.因此
其中μij-μ反映了水平组合(Ai,Bj)对Y的效应.一般情况,μij-μ≠αi+βj
μij=μ+αi+βj+γij,i=1,2, ,a,
容易验证:
j=1,2, ,b
ij
∑α
i=1
a
i
=0,
∑β
j=1
b
j
=0,
∑γ
i=1
a
ij
=0,
∑γ
j=1
b
=0
因此两因素等重复下的方差分析模型等价地改写为如下形式:
⎧
⎪y=μ+ε=μ+α+β+γ+ε,i=1,2, ,a,ijkijijkiiijijk⎪⎪2
且诸εijk相互独立⎨εijk~N(0,σ),
bab⎪a
⎪∑αi=0,∑βj=0,∑γij=0,∑γij=0⎪j=1i=1j=1⎩i=1
j=1,2, ,b
3.2.2交互效应及因素效应的显著性检验
一.偏差平方和分解
下面先对Y的观测的总平方和进行分解:
ij∙=∑yijk/c=μ+αi+βj+γij+ij∙,i=1,2, ,a,
k=1
c
j=1,2, ,b
i∙∙
1bc1b1c1b
=∑∑yijk=∑(∑yijk)=∑ij∙=μ+αi+i∙∙,i=1,2, ,a bcj=1k=1bj=1ck=1bj=11ac1a1c=yijk=∑(∑yijk)=μ+βj+∙j∙,j=1,2, ,b ∑∑aci=1k=1ai=1ck=1
∙j∙
观测数据的总(偏差)平方和为
SST=∑∑∑(yijk-)2
i=1j=1k=1c
abc
=∑∑∑[(i∙∙-)+(∙j∙-)+(ij∙-i∙∙-∙j∙+)+(yijk-ij∙)]2
i=1j=1k=1
a
ab
=bc∑(i∙∙-)+ac∑(∙j∙-)+c∑∑(ij∙-i∙∙-∙j∙+)+∑∑∑(yijk-ij∙)2
2
2
i=1
j=1
i=1j=1
i=1j=1k=1
bab
2
abc
=SSA+SSB+SSAB+SSE
其中
SSA=bc∑(i∙∙-)=bc∑(αi+i∙∙-)2——因素A的平方和
2
i=1
i=1
aa
由于E(i∙∙-)=αi,为αi的无偏估计,故SSA度量A的各水平效应的估计量的变化.
SSB=ac∑(∙j∙-)=ac∑(βj+∙j∙-)2——因素B的平方和
2
j=1
j=1
bb
由于E(∙j∙-)=βj,为βj的无偏估计,故SSB度量B的各水平效应的估计量的变化.
SSAB=c∑∑(ij∙-i∙∙-∙j∙+)=c∑∑(γij+ij∙-i∙∙-∙j∙+)2——交互效
2
i=1j=1
i=1j=1
abani
应的平方和
由于E(ij∙-i∙∙-∙j∙+)=γij,为γij的无偏估计,故SSAB度量A和B的交互效应的
估计量的变化.
SSE=∑∑∑(ijk-ij∙)=∑∑∑(εijk-ij∙)2——误差平方和
2
i=1j=1k=1
i=1j=1k=1
abcabc
度量了来自各总体的观测值与其样本均值的差异,反映了误差的变化. 由于εijk~N(0,σ2)且相互独立,可得
E(SSA)=(a-1)σ+bc∑αi2
2
i=1
a
E(SSB)=(b-1)σ+ac∑βj2
2
j=1
b
2 E(SSAB)=(a-1)(b-1)σ+c∑∑γij
2
i=1j=1
ab
E(SSE)=ab(c-1)σ2
a-1,b-1,(a-1)(b-1),ab(c-1)分别称为SSA,SSB,SSAB,SSE的自由度,abc-1
称为SST的自由度,为上述四个自由度的和.
aSSAbc2A令 MSA=——因素的均方, 则E(MSA)=σ+αi2 ∑a-1a-1i=1
bSSBac2
MSB=——因素B的均方, E(MSB)=σ+βj2 ∑b-1b-1j=1
MSAB
abSSABc22 ——交互效应均方,=E(MSAB)=σ+γij∑∑(a-1)(b-1)(a-1)(b-1)i=1j=1
MSE=
SSE2
——误差均方, E(MSE)=σ2, MSE为σ的无偏估计
ab(c-1)
二.假设检验
对两因素的情况,方差分析的主要目的除了考察因素A或B的各水平对因变量Y的影响有无显著差异外,还要考虑A和B之间是否存在交互作用,因为交互作用的存在会直接影响到对A和B影响检验结果的解释.涉及如下三个检验问题:
HA0:μ1∙=μ2∙= =μa∙ HA1:μ1∙,μ2∙, ,μa∙不全相等 HB0:μ∙1=μ∙2= =μ∙b HB1:μ∙1,μ∙2, ,μ∙b不全相等
HAB0:γij=0,i=1,2, ,a,
检验问题也可以改写成:
j=1,2, ,bHAB0:至少有一个γij≠0
HA0:α1=α2= =αa=0HA1:至少有某个αi≠0
HB0:β1=β2= =βb=0HB1:至少有某个βj≠0
HAB0:γij=0,i=1,2, ,a,
j=1,2, ,bHAB0:至少有一个γij≠0
利用上述结果,构造适当的统计量检验上述假设.
MSE为σ2的无偏估计,如果假设H0成立,MSA,MSB,MSAB取值接近σ2,如果假
设不成立,则MSA,MSB,MSAB有增大的趋势.因此,针对检验HA0,HB0,HAB0,分别构造统计量,分别有
MSAHA0真
FA=~F(a-1,ab(c-1))
MSEMSBHB0真
FB=~F(b-1,ab(c-1))
MSE
FAB
MSABHA0真=~F((a-1)(b-1),ab(c-1)) MSE
如果,各检验统计量FA,FB,FAB的值变大,则拒绝原假设.各检验的p值分别为 pA=PHA0(FA≥fA)=P(F(a-1,ab(c-1))≥fA) pB=PHA0(FB≥fB)=P(F(b-1,ab(c-1))≥fB)
pAB=PHAB0(FAB≥fAB)=P(F(a-1)(b-1),ab(c-1))≥fAB) 其中fA,fB,fAB为统计量观测值.
给定显著性水平α,如pA
注意:检验步骤:
1)先检验HAB0,如不拒绝HAB0,即交互作用不显著时,再考察A和B的效应的显著性.
因为当γij(i=1,2, ,a,任何水平Bj上
1
2
1
j=1,2, ,b)全为0时,即A与B无交互作用,则在B的
2
1
2
μij-μij= (αi-αi),
1
2
μij-μij在B的各水平上均相等且完全由A在水平Ai1和Ai2的效应之差(αi-αi)
确定,因此,检验假设的结论真实的反映了仅由A的各水平对Y的影响是否显著. 2)如拒绝HAB0,交互作用显著时,通过估计和比较因素A和B各水平组合(A,B)上的均值考察因素的联合影响
如果A和B存在交互作用,则γij不全为0,对于A的两个水平Ai1和Ai2在B的第j个水平上的两个组合(Ai1,Bj)和(Ai2,Bj)下均值差为
μij-μij=(μ+αi+βj+γij)-(μ+αi+βj+γij)
1
2
1
1
2
2
=(αi1-αi2)+(γi1j-γi2j)) 因此当γij(i=1,2, ,a,
j=1,2, ,b)不全为0时,μi1j-μi2j除与(αi1-αi2)有关
外,还可能与(γi1j-γi2j)有关,即B处在不同水平Bj,μi1j-μi2j有所不同.
如图3.1所示,假设A与B均有两个水平A1,A2和B1,B2
11
μ1∙-μ2∙=[(μ11-μ21)+(μ12-μ22)]=(α1-α2)+[(γ11-γ21)+(γ12-γ22)]
22
B1
B2
(a) (b) 图3.1 有交互效应时A的各水平均值在B的不同水平上的差异
(a) A1,A2在B1,B2下均值差μ11-μ21=0而μ12-μ22≠0认为μ1∙≠μ2∙(即α1与,差异主要表现在A1,A2在B2水平上的差异; α2不全为0)
(b) A1,A2在B1,B2下均值差μ11-μ21与μ12-μ22相反,综合有μ1∙=μ2∙,主要表现在A1,A2在B1和B2水平上的差异相互抵消,使得综合差异为零.
因此,在有交互效应时,尤其是交互效应显著,而因素A与B的效应不显著时,检验每个因素显著差异HA0,HB0实际意义不大,应慎重.此情况下,要进一步考察各因素对Y的影响的显著性,只能将一个因素的各个水平逐个给定,在给定的水平上考察另一因素的各水平均值之间的差异来了解该因素对Y的影响.
例3.5 某高校为了了解数学专业和计算机科学专业的低年级学生、高年级学生及研究
生在人文科学知识方面的差异,从不同专业和不同年级的学生中任选四名学生参加有关考试,其成绩如表3.9
所示,假设考试成绩服从两因素的方差分析模型,对其作方差分析.
B1
B2
解:因变量Y为成绩,有两个因素,记专业因素为A,有两个水平A1(数学),A2(计算机);学生级别为因素B,有三个水平B1(低年级),B2(高年级),B3(研究生).因此,
a=2,b=3,c=4.利用SAS系统proc anova过程作方差分析,程序及结果如下:
data examp3_5;
input majors $ classes $ grade @@; /*输入majors课程,年级classes,成绩grade*/ cards;
a1 b1 81 a1 b1 78 a1 b1 79 a1 b1 78 a1 b2 75 a1 b2 80 a1 b2 78 a1 b2 73 a1 b3 82 a1 b3 80 a1 b3 85 a1 b3 88 a2 b1 89 a2 b1 82 a2 b1 77 a2 b1 90 a2 b2 79 a2 b2 80 a2 b2 75 a2 b2 78 a2 b3 93 a2 b3 93 a2 b3 86 a2 b3 95 ; run;
proc anova data=examp3_5; /* 调用方差分析的anova过程 */
class majors classes; /* 因素变量名称为majors专业、classes级别 */ model grade=majors classes majors*classes; /* 因变量grade,因素变量专
业、级别、因素效应 */
run;
SAS 系统 11:29 Tuesday, October 21, 2008
The ANOVA Procedure
Class Level Information 因素 水平 因素变量 Class Levels Values 因素变量A majors a=2 a1 a2 因素变量B classes b=3 b1 b2 b3
Number of observations n=abc=24 c=4
注:SSR=SSA+SSB+SSAB,因素变量为自变量,个数为自由度p-1=ab-1=5.
H0真
先检验H0:μij=0,统计量F=MSR~F(ab-1,ab(c-1))=F(5,18),
MSE观测值f=MSR=9.34,检验p值p=PH0(F(5,18)≥f)=0.0002.
MSE
SAS 系统 11:29 Tuesday, October 21, 2008 2
The ANOVA Procedure
因变量:成绩 Dependent Variable: grade
表 3.10 考试成绩的方差分析结果
Sum of
Source DF Squares Mean Square F Value Pr > F 方差来源 自由度 平方和SS 均方 F值 p值 Model(模型) p-1=ab-1=5 SSR=637.0000000 MSR=127.400000 f=MSR=9.34 0.0002
MSEError(误差) n-p=ab(c-1)=18 SSE=245.5000000 MSE=13.6388889 Corrected Total n-1=abc-123 SST=882.5000000
R-Square Coeff Var Root MSE grade Mean R2=SSR=0.721813 4.490075 3.693087 y=82.25000
SST
Source DF Anova SS Mean Square F Value Pr > F
方差来源 自由度 平方和SS 均方 F值 p值
专业A majors a-1=1 SSA=150.000000 MSA=150.000000 f=MSA=11.00 pA=0.0038
A
MSE 级别B classes b-1=2 SSB=444.000000 MSB=222.000000 f=MSB=16.28 pB
B
MSE交互majors*classes (a-1)(b-1)=2 SSAB=43.000000 MSAB=21.500000 f=MSAB=1.58 pAB=0.2340
AB
MSE
检验HAB0:由结果知pAB=0.2340>0.05,接受HAB0,认为专业与级别的交互不显著,即数学专业中各级别学生人文知识水平的差异与计算机科学专业中相应级别知识水平差异认为相同.
检验HA0和HB0,检验p值分别为pA=0.0038和pB=0.0001,拒绝HA0和HB0.说明因素A和因素B对成绩的影响均显著,即两个专业学生的人文知识水平是有显著差异的,不同级别的学生的人文社科知识水平也有显著差异.
例3.6 某计算机修理公司为了解三名修理工对不同类型磁盘驱动系统修理的工作效率,每位修理工被随机地指定修理三种类型的磁盘驱动器系统各5个,其完成修理工作的时间如表3.11所示.假设修理时间服从两因素方差分析模型,试对此数据作方差分析.
解:因变量Y为修理时间为数量指标,因素A代表修理工,有三个水平A1,A2,A3代表三个修理工,因素B,有三个水平B1,B2,B3,代表三种类型的驱动系统,各修理工修理时间即为Y在各组合水平下的观测值,a=b=3,程,程序如下:
data examp3_6;
input repairer $ type $ time @@; /* 输入修理工、磁盘类型、修理时间 */ cards;
c=5,n=abc=45.利用proc anova过
a1 b1 62 a1 b1 48 a1 b1 63 a1 b1 57 a1 b1 69 a1 b2 57 a1 b2 45 a1 b2 39 a1 b2 54 a1 b2 44 a1 b3 59 a1 b3 53 a1 b3 67 a1 b3 66 a1 b3 47 a2 b1 51 a2 b1 57 a2 b1 45 a2 b1 50 a2 b1 39 a2 b2 61 a2 b2 58 a2 b2 70 a2 b2 66 a2 b2 51 a2 b3 55 a2 b3 58 a2 b3 50 a2 b3 69 a2 b3 49 a3 b1 59 a3 b1 65 a3 b1 55 a3 b1 52 a3 b1 70 a3 b2 58 a3 b2 63 a3 b2 70 a3 b2 53 a3 b2 60 a3 b3 47 a3 b3 56 a3 b3 51 a3 b3 44 a3 b3 50 ;
run;
proc anova data=examp3_6; class repairer type; /* 因素变量repairer修理工、type类型 */ model time=repairer type repairer*type; /* 模型因变量time,因素变量
repairer type 主效应及交互效应repairer*type */
means repairer type repairer*type; /* 列出主效应及交互效应的样本均值、方差*/ run;
SAS 系统 11:29 Tuesday, October 21, 2008 3
The ANOVA Procedure
Class Level Information
Class Levels Values
因素A repairer a=3 a1 a2 a3
因素B type b=3 b1 b2 b3
Number of observations n=45
SAS 系统 11:29 Tuesday, October 21, 2008 4
The ANOVA Procedure Dependent Variable: time 即y
表3.12 修理时间的方差分析表
Sum of
Source DF Squares Mean Square F Value Pr > F 模型 Model ab-1=8 SSR=1268.177778 MSR=158.522222 f=MSR=3.05
MSE
误差 Error ab(c-1)=36 SSE=1872.400000 MSE=52.011111
总和 Corrected Total abc-1=44 SST=3140.577778
R-Square Coeff Var Root MSE time Mean
0.403804 12.91936 7.211873 55.82222
Source DF Anova SS Mean Square F Value Pr > F A修理工repairer a-1=2 SSA=24.577778 MSA=12.288889 f=MSA= 0.24 pA=0.7908 AMSE
B系统类型type b-1=2 SSB=28.311111 MSB=14.155556 f=MSB=0.27 pB=0.7633 BMSEp=0.0101
交互repairer*type (a-1)(b-1)=4 SSAB=1215.288889 MSAB=303.822222 f=MSAB=5.84 pAB=0.0010 AMSE
由此可见,pAB=0.0010
单看因素A(修理工),pA=0.7908>0.05 ,接受HA0.说明不同的修理工对修理时间的影响不显著;因素B(磁盘系统)pB=0.7633>0.05 ,接受HB0说明不同类型的驱动器对修理时间的影响也不显著.
因此,关于因素A或因素B的效应的检验结果并无多大参考价值.为进一步了解交互效应的本质,对每一个组合水平(Ai,Bj)上的观测数据,求得样本均值作为每个组合水平上的总体均值的估计,结果如下:
SAS 系统 11:29 Tuesday, October 21, 2008 5
The ANOVA Procedure
Level of -------------time------------
因素A各水平上的样本均值
i∙∙135132=∑∑yijk si∙∙=(ij∙-i∙∙)2∑15j=1k=13-1j=1 repairer N Mean i∙∙ Std Dev si∙∙
a1 15 55.3333333 9.22470802
a2 15 55.2666667 8.75431866
a3 15 56.8666667 7.79987790
因素B各水平上的样本均值
∙j∙131352(ij∙-∙j∙)2=∑∑yijk s∙j∙=∑3-1i=115i=1k=1
Level of -------------time------------ type N Mean ∙j∙ Std Dev s∙j∙
b1 15 56.1333333 8.83876742
b2 15 56.6000000 9.17138407
b3 15 54.7333333 7.75026881
表3.13 因素A与B的组合水平上的样本均值 ij∙15152=∑yijk sij∙=(yijk-ij∙)2 ∑5k=15-1k=1
Level of Level of -------------time------------ repairer type N Mean ij∙ Std Dev sij∙
a1 b1 5 59.8000000 7.85493475
a1 b2 5 A1中12∙=47.8000000最短 7.46324326
a1 b3 5 58.4000000 8.53229160
a2 b1 5 A2中21∙=48.4000000最短 6.76756973
a2 b2 5 61.2000000 7.32802838
a2 b3 5 56.2000000 8.04363102
a3 b1 5 60.2000000 7.32802838
a3 b2 5 60.8000000 6.30079360
a3 b3 5 A3中33∙=49.6000000最短 4.50555213
由结果可见,不同的修理工修理不同的磁盘驱动系统所花费的时间确有较大差异,修理工A1修理B2类型的驱动器系统所花平均时间最短(12∙=47.8000),A2修理B1类型的驱动器系统所花平均时间最短(21∙=48.4000),A3修理B3类型的驱动器系统所花平均时间最短(33∙=49.6000).
12 3123图3.2 各水平组合(Ai,Bj)
上的样本均值
而由于i∙∙(i=1,2,3)之间的差异不大导致因素A的影响不显著的检验结果;而由于∙j∙(j=1,2,3)之间的差异不大导致因素B的影响不显著的检验结果;因此交互效应可能会掩盖各因素对因变量Y的某些本质影响.
3.2.3无交互效应的各因素均值的估计与比较
在给定的显著水平α下,当假设检验的结论是因素A和B之间的交互效应不显著,并且因素A和B至少有一个对Y有显著影响,可以进一步对影响显著的因素在其各水平下的均值作出估计,并给出其本身及任意两个之差的置信区间.
一.因素A均值的估计和比较
1.μi∙的无偏估计及区间估计
i=1,2, ,a)由
E所以μi∙的无偏估计为 μˆi∙=i∙∙=∑∑yijk/bc~N(μi∙,
j=1k=1bcσ2bc)
由于 i∙∙-μi∙
/bc
222=(i∙∙-μi∙)σ~N(0,1) ˆ=MSE为σ的无偏估计 ,且在方差分析模型下可证, 而σ
ab(c-1)
σ2MSE=SSE
σ2~χ2(ab(c-1)) 且i∙∙与MSE相互独立,从而可得
(i∙∙-μi∙)/σ
SSE=bc(i∙∙-μi∙)MSE~t(ab(c-1))
σ2/(ab(c-1)
⎛(-μ)⎫i∙∙i∙ 对给定的显著性水平α,由 P≤tα(ab(c-1))⎪=1-α ⎪1-MSE2⎝⎭
可得μi∙的置信度为1-α的置信区间为
(i∙∙-t1-α
2(ab(c-1)MSE/bc,i∙∙+t1-α
2(ab(c-1)MSE/bc),(i=1,2, ,a)
2.A的各水平差μi1∙-μi2∙的区间估计
同理可得μi1∙-μi2∙的置信度为1-α的置信区间为
(i1∙∙-i2∙∙-t1-α
2(ab(c-1)2MSE/bc,i∙∙-i2∙∙+t1-α
2(ab(c-1)2MSE/bc),
若有m个μi1∙-μi2∙作同时比较,则它们的置信度不小于1-α的Bonferroni同时置信区间为
(i1∙∙-i2∙∙-t1-α
2m(ab(c-1)2MSE/bc,i1∙∙-i2∙∙+t1-α2m(ab(c-1)2MSE/bc)
二.因素B均值的估计和比较
如因素B对Y的影响显著,μ∙j(j=1,2, ,b)的无偏估计为
ˆ∙j=∙j∙=∑∑yijk/ac μi=1k=1ac
μ∙j的置信度为1-α的置信区间为
(∙j∙-t1-α
2(ab(c-1)MSE/ac,∙j∙+t1-α
2(ab(c-1)MSE/ac)
μ∙j∙-μ∙j∙(j=1,2, ,b)的置信度为1-α的置信区间为 12
(∙j1∙-∙j2∙-t1-α
2(ab(c-1)2MSE/ac,∙j1∙-∙j2∙+t1-α
2(ab(c-1)2MSE/ac)
μ∙j1∙-μ∙j2∙的置信度不小于1-α的Bonferroni同时置信区间为
(∙j1∙-∙j2∙-t1-α
2m(ab(c-1)2MSE/ac,∙j1∙-∙j2∙+t1-α
2m(ab(c-1)2MSE/ac) 例3.7 (续例3.5)根据表3.9的数据,给出两专业之间学生成绩的均值之差和各级别学生之间成绩的均值之差的置信度不小于95%的Bonferroni同时置信区间. 解:由例3.5结果可知,专业A与学生级别B之间无显著的交互作用,因素A和B均对成绩影响显著,因此可进一步通过比较A的各水平均值差异了解各专业与各级别社科知识的差异;通过比较B的各水平均值差异了解各级别社科知识的差异. data examp3_7;
input majors $ classes $ grade @@; cards; a1 b1 81 a1 b1 78 a1 b1 79 a1 b1 78
a1 b2 75 a1 b2 80 a1 b2 78 a1 b2 73
a1 b3 82 a1 b3 80 a1 b3 85 a1 b3 88
a2 b1 89 a2 b1 82 a2 b1 77 a2 b1 90 a2 b2 79 a2 b2 80 a2 b2 75 a2 b2 78
a2 b3 93 a2 b3 93 a2 b3 86 a2 b3 95
;
run;
proc anova data=examp3_7; class majors classes; /* 因素变量专业、级别*/ model grade=majors classes majors*classes; /* 因变量成绩,因素变量专业、级别主效应及交互效应 */
means majors classes; /* 计算专业、级别对应的因变量样本均值和标准差 */ means majors classes/bon cldiff alpha=0.05; /* 因素变量专业、级别在不同水平上的均值进行Bonferroni同时两两比较的t检验,显著性水平0.05,输出不同水平上的两两均值差的置信度不小于1-0.05的置信区间.*/
run;
SAS 系统 11:29 Tuesday, October 21, 2008 12
The ANOVA Procedure
Class Level Information
Class Levels Values
因素A majors a= 2 a1 a2
因素B classes b=3 b1 b2 b3
Number of observations n=abc=24
以下为方差分析表,同3.5例表3.10
SAS 系统 11:29 Tuesday, October 21, 2008 13
The ANOVA Procedure
Dependent Variable: grade
Sum of
Source DF Squares Mean Square F Value Pr > F Model 5 637.0000000 127.4000000 9.34 0.0002 Error 18 245.5000000 13.6388889
Corrected Total 23 882.5000000
R-Square Coeff Var Root MSE grade Mean
0.721813 4.490075 3.693087 82.25000
Source DF Anova SS Mean Square F Value Pr > F majors 1 150.0000000 150.0000000 11.00 0.0038 classes 2 444.0000000 222.0000000 16.28
SAS 系统 11:29 Tuesday, October 21, 2008 14
The ANOVA Procedure
Level of ------------grade------------
majors N Mean Std Dev
专业因素水平 观测个数 样本均值i∙∙ 样本标准差si∙∙
i∙∙134132 =∑∑yijksi∙∙=(ij∙-i∙∙)2∑12j=1k=13-1j=1
a1 12 79.7500000 4.04800737
a2 12 84.7500000 7.08551660
Level of ------------grade------------
classes N Mean Std Dev
级别水平 观测个数 样本均值∙j∙ 样本标准差s∙j∙
∙j∙131342(ij∙-∙j∙)2=∑∑yijk s∙j∙=∑3-1i=112i=1k=1
b1 8 81.7500000 5.06387768
b2 8 77.2500000 2.60494036 b3 8 87.7500000 5.49675229
以下对专业因素A求均值差的Bonferroni同时置信区间
SAS 系统 11:29 Tuesday, October 21, 2008 15
The ANOVA Procedure
Bonferroni (Dunn) t Tests for grade
NOTE: This test controls the Type I experimentwise error rate, but it generally has a higher Type II error rate than Tukey's for all pairwise comparisons.
Alpha α=0.05
=13.63889 Error Degrees of Freedom ab(c-1)=18 Error Mean Square MSE
Critical Value of t t1-α
2m(ab(c-1)=t0.975(18)= 2.10092 m=1
Minimum Significant Difference 3.1676
Comparisons significant at the 0.05 level are indicated by ***.
Difference Simultaneous
majors Between 95% Confidence
Comparison Means Limits
μi∙-μj∙ 样本均值差i∙∙-j∙∙ μi∙-μj∙的置信区间
a2 - a1 5.000 1.832 8.168 ***
a1 - a2 -5.000 -8.168 -1.832 ***
对于专业因素A,由上述结果求得1∙∙=79.75,2∙∙=84.75,α=0.05,由于仅有一个差μ1∙-μ2∙, 故m=1,t1-α/2m(ab(c-1))=t0.975(18)= 2.10092,μ1∙-μ2∙的置信度为95%的置信区间为(-8.168 -1.832).因此认为μ1∙
类似可得级别因素B求均值差的Bonferroni同时置信区间
SAS 系统 11:29 Tuesday, October 21, 2008 16
The ANOVA Procedure
Bonferroni (Dunn) t Tests for grade
NOTE: This test controls the Type I experimentwise error rate, but it generally has a higher Type II error rate than Tukey's for all pairwise comparisons.
m=b(b-1)/2=3
Alpha α=0.05
Error Degrees of Freedom ab(c-1)=18
Error Mean Square MSE=13.63889 Critical Value of t t1-α
2m(ab(c-1)=t0.9917(18)= 2.63914
Minimum Significant Difference 4.8733
Comparisons significant at the 0.05 level are indicated by ***.
Difference
classes Between Simultaneous 95%
Comparison Means Confidence Limits
b3 - b1 6.000 1.127 10.873 *** 表示显著差异 b3 - b2 10.500 5.627 15.373 ***
b1 - b3 -6.000 -10.873 -1.127 ***
b1 - b2 4.500 -0.373 9.373
b2 - b3 -10.500 -15.373 -5.627 ***
b2 - b1 -4.500 -9.373 0.373
对于级别因素B,由上述结果求得∙j∙=81.75,∙2∙=77.25,∙3∙=87.75,共有m=3个均值差,故t1-α/2m(ab(c-1))=t0.9917(18)=2.639,均值差的置信度为95%的Bonferroni同时置信区间分别为
(∙j1∙-∙j2∙-t
μ∙1-μ∙2:(-0.373 9.373), μ∙1-μ∙3:(-10.873 -1.127), μ∙2-μ∙3:(-15.373 -5.627) 1-α2m(ab(c-1)2MSE/ac,∙j1∙-∙j2∙+t1-α
2m(ab(c-1)2MSE/ac)
由结果可知,μ∙1
社科知识无显著差异. μ∙3,μ∙2
3.2.4 有交互效应时因素各水平组合上的均值估计与比较
如果因素A与B之间有显著的交互效应,这时单独考虑A或B各水平上均值的差异 无多大实际意义.可通过直接比较各水平组合(Ai,Bj)上的均值μij来了解差异.
对给定的因素B的某个水平Bj(j=1,2, ,b),μij的无偏估计
μij=ij∙ 对任意 i=1,2, ,a
σ2
c), ∧则有 ij∙~N(μij,
可证 (ij∙-μij)
MSE~t(ab(c-1)
对给定的显著性水平α,由
⎛(ij∙-μij)⎫ P≤tα(ab(c-1)⎪=1-α ⎪1-MSE2⎝⎭
可得μij的置信度为1-α的置信区间为
(ij∙-t1-α
2(ab(c-1)MSE/c,ij∙+t1-α
2(ab(c-1)MSE/c),(i=1,2, ,a)
可证 c(i1j∙-i2j∙)-(μi1j-μi2j))
MSE~t(ab(c-1)
m个μi1j-μi2j同时比较,则它们置信度不小于1-α的Bonferroni同时置信区间为 (i1j∙-i2j∙-t1-α
2m(ab(c-1)2MSE/c,i1j∙-i2j∙+t1-α
2m(ab(c-1)2MSE/c)
类似可得,对给定因素A的某个水平Ai,关于B的各水平均值的估计和成对比较结果. 例3.8 (续3.6)根据表3.11数据,试对每一类型的磁盘驱动系统,给出三名修理工修理时间均值两两之差的置信度不小于90%的Bonferroni同时置信区间.
解:因变量Y-修理时间,因素A-修理工,有三个水平,因素B-磁盘驱动系统,有三个水平,各水平观测次数5,a=b=3,c=5.利用proc anova过程,程序见例3.6:
表3.12 修理时间的方差分析表(见例3.6)
Sum of
Source DF Squares Mean Square F Value Pr > F Model 8 1268.177778 158.522222 3.05 0.0101
误差 Error 36 1872.400000 MSE=52.011111
总和 Corrected Total 44 3140.577778
表3.13 因素A与B的组合水平上的样本均值 s2
ij∙15=(yijk-ij∙)2∑5-1k=1
Level of Level of -------------time------------ repairer type N Mean ij∙ Std Dev sij∙ a1 b1 5 59.8000000 7.85493475 a1 b2 5 A1中12∙=47.8000000最短 7.46324326 a1 b3 5 58.4000000 8.53229160 a2 b1 5 A2中21∙=48.4000000最短 6.76756973 a2 b2 5 61.2000000 7.32802838 a2 b3 5 56.2000000 8.04363102 a3 b1 5 60.2000000 7.32802838 a3 b2 5 60.8000000 6.30079360 a3 b3 5 A3中33∙=49.6000000最短 4.50555213
由方差分析表得,MSE=52.011111.
对B1类型的磁盘驱动系统,由表3.13可知,
11∙=59.8000,21∙=48.4000,31∙=60.2000,α=0.1
对于m=3,求得α=0.1=0.0167,而tα(ab(c-1)=t0.9833(36)=2.213,Bonferroni同1-2m62m
时置信区间
(i1j∙-i2j∙-t1-α2m(ab(c-1)2MSE/c,i1j∙-i2j∙+t1-α
2m(ab(c-1)2MSE/c)
在B1水平上,μi1-μj1的置信度不低于90%的Bonferroni同时置信区间分别为
(1.306 21.494); μ11-μ31:(-10.494 9.694); μ21-μ31:(-21.894 -1.706); μ11-μ21:
由此可知,对B1类型的磁盘驱动系统,可以至少90%的置信度断言:修理工A3与A1的平均修理时间无显著差异,而A2的平均修理时间显著地小于A3与A1的平均修理时间,即修理工A2最擅长修理B1类型的磁盘驱动系统.
对B2类型的磁盘驱动系统,由表3.13可知,
12∙=45.8000,22∙=61.2000,32∙=60.8000
可得B2水平上置信度不低于90%的Bonferroni同时置信区间分别为
(-23.494 -3.306); μ12-μ32:(-23.094 -2.096); μ22-μ32:(-9.694 10.494) μ12-μ22:
由此可知,对B2类型的磁盘驱动系统,可以至少90%的置信度断言:修理工A2与A3对修理B2类型的磁盘驱动系统的平均时间无显著差异,而修理工A1最擅长修理B2类型的磁盘驱动系统.
对B3类型的磁盘驱动系统,置信度不低于90%的Bonferroni同时置信区间分别为
(-7.894 12.294); μ12-μ32:(-1.294 18.894); μ22-μ32:(-3.494 16.694); μ13-μ23:
至少以90%置信度断言:三个修理工修理B3类型磁盘驱动系统的平均时间无显著差异.作业 3.6 3.7
内容总结
1.统计模型
Y—因素A和B,水平A1,A2, ,Aa,B1,B2, ,Bb,组合水平(Ai,Bj)观测值 yij1,yij2, ,yijc ⎧⎪y=μ+ε=μ+α+β+γ+ε,i=1,2, ,a,ijkijijkiiijijk⎪⎪εijk~N(0,σ2),且诸εijk相互独立⎨
abab⎪αi=0,∑βj=0,∑γij=0,∑γij=0⎪∑⎪i=1j=1i=1j=1⎩
SST=∑∑∑(yijk-)2
i=1j=1k=1abcj=1,2, ,b,k=1,2, ,c
=bc∑(i∙∙-)+ac∑(∙j∙-)+c∑∑(ij∙-i∙∙-∙j∙+)+∑∑∑(yijk-ij∙)222
i=1j=1i=1j=1i=1j=1k=1abab2abc=SSA+SSB+SSAB+SSE
SSE~χ(ab(c-1)),σ=2∧2
σ2SSE2 为σ无偏估计. ab(c-1)
假设成立时,SSASSBSSAB222分别服从χ(a-1),χ(b-1),χ[(a-1)(b-1)]分布.
222σσσ
2.显著检验
假设检验问题:HA0:α1=α2= =αa=0HA1:至少有某个αi≠0
HB0:β1=β2= =βb=0HB1:至少有某个βj≠0 HAB0:γij=0,i=1,2, ,a,j=1,2, ,bHAB0:至少γij≠0
MSBHB0真MSAHA0真统计量 FA=~F(b-1,ab(c-1)) ~F(a-1,ab(c-1)) FB=MSEMSE
FABMSABHA0真=~F((a-1)(b-1),ab(c-1)) MSE
检验p值:如p
pA=PHA0(FA≥fA)=P(F(a-1,ab(c-1))≥fA) pB=PHA0(FB≥fB)=P(F(b-1,ab(c-1))≥fB) pAB=PHAB0(FAB≥fAB)=P(F(a-1)(b-1),ab(c-1))≥fAB)
3.无交互效应的各因素均值的估计与比较
μi∙的无偏估计 μi∙=i∙∙=∑∑yijk/bc
j=1k=1∧bc
μi∙的置信度为1-α的置信区间
(i∙∙-t1-α2(ab(c-1)MSE/bc,i∙∙+t1-α
2(ab(c-1)MSE/bc)
μi∙-μi∙的置信度为1-α的置信区间 12
(i1∙∙-i2∙∙-t1-α
2(ab(c-1)MSE/bc,i∙∙-i2∙∙+t1-α
2(ab(c-1)2MSE/bc),
m个μi1∙-μi2∙同时比较,置信度不小于1-α的Bonferroni同时置信区间 (i1∙∙-i2∙∙-t1-α
2m(ab(c-1)2MSE/bc,
∧i1∙∙-i2∙∙+t1-α2m(ab(c-1)2MSE/bc) 4.有交互效应时因素各水平组合上的均值估计与比较
组合(Ai,Bj)上均值μij的无偏估计 μij=ij∙ μij的置信度为1-α的置信区间
(ij∙-t1-α
2(ab(c-1)MSE/c,ij∙+t1-α
2(ab(c-1)MSE/c),(i=1,2, ,a)
μij-μij置信度不小于1-α的Bonferroni同时置信区间 12
(i1j∙-i2j∙-t
1-α2m(ab(c-1)2MSE/c,i1j∙-i2j∙+t1-α2m(ab(c-1)2MSE/c)