SAS讲义 第二十五课方差分析

第二十五课 方差分析

当影响观察结果的影响因素（原因变量或分组变量）的水平数大于2或原因变量的个数大于1个，一元时常用F检验（也称一元方差分析），多元时用多元方差分析（最常用Wilks’∧检验）。

一、 方差分析概述

方差分析（analysis of variance）又称变异数分析，可简记为ANOVA，主要用于检验计量资料中的两个或两个以上均值间差别显著性的方法。当欲比较几组均值时，理论上抽得的几个样本，都假定来自正态总体，且有一个相同的方差，仅仅均值可以不相同。还需假定每一个观察值都由若干部分累加而成，也即总的效果可分成若干部分，而每一部分都有一个特定的含义，称之谓效应的可加性。所谓的方差是离均差平方和除以自由度，在方差分析中常简称为均方MS（mean square）。 1. 方差分析的基本思想

根据效应的可加性，将总的离均差平方和分解成若干部分，每一部分都与某一种效应相对应，总自由度也被分成相应的各个部分，各部分的离均差平方除以相应部分的自由度得出各部分的均方，然后列出方差分析表算出F值，作出统计推断。

方差分析的关键是总离均差平方和的分解，分解越细致，各部分的含义就越明确，对各种效应的作用就越了解，统计推断就越准确。方差分析表的一般形式见表25.1所示：

表25.1 方差分析表形式

变异来源 source 效应S1 效应S2 ?? 效应Sm 误差Se 总变异ST 离差平方和 SS SS1 SS2 ?? SSm SSe + SSm+ SSe 自由度 df df1 df2 ?? dfm dfe + dfm + dfe 均方 MS MS1= SS1/df1 MS2= SS2/df2 ?? MSm= SSm/dfm MSe= SSe/dfe MST= SST/dfT F统计量 F F1(df1, dfe)= MS1/ MSe F2(df2, dfe)= MS2/ MSe ?? Fm(dfm, dfe)= MSm/ MSe FT(dfT, dfe)= MST/ MSe P概率值 P P1 P2 Pm PT SST= SS1+ SS2+?dfT=df1+ df2+?

表中变异来源一栏，可分为总变异（total），误差（residual）,各个效应（effect）相对应的项。效应项与试验设计或统计分析的目的有关，一般有：主效应（包括各种因素），交互影响项（因素间的多级交互影响），协变量（来自回归的变异项），等等。

当分析和确定了各个效应项S后，根据原始观察资料可计算出各个离均差平方和SS，再根据相应的自由度df，由公式MS=SS/df，求出均方MS，最后由相应的均方，求出各个变异项的F值，F值实际上是两个均方之比值，通常情况下，分母的均方是误差项的均方。根据F值的分子、分母均方的自由度f1和f2，在确定显著性水平为?情况下，由F(f1,f2)临界值表查得单侧F?界限值。当F?F?时，则P??，不拒绝原假设H0，说明不拒绝这个效应项的效应为0的原假设，也即这个效应项是可能对总变异没有实质影响的；如果F?F?,则

上海财经大学经济信息管理系IS/SHUFE

Page 1 of 30

P??，拒绝原假设H0，说明拒绝这个效应项的效应为0的原假设，也即这个效应项是很

可能对总变异有实质影响的。

2. 方差分析的试验设计

为了确定方差分析表中各个有关效应项，需要在试验设计阶段就作出安排，再根据设计要求进行试验，得出原始观察值，按原来设计方案算出方差分析表中的各项。在试验设计阶段常需要作主要四个方面的考虑： 1) 研究的主要变量

方差分析的主要变量，也称响应变量或因变量（dependent variable），它是我们试验所要观察的主要指标。一次试验时可以有多个观察指标，方差分析时也可以同时对多个因变量进行分析。

2) 因素和水平

试验的因素（factor）可以是品种、人员、方法、时间、地区等等，因素所处的状态叫水平（level）。在每一个因素下面可以分成若干水平。例如，某工厂的原料来自四个不同地区，那么用不同地区的原料生产的产品质量是否一致呢？所要比较的地区就是因素，四个地区便是地区这一因素的四个水平。当某个主要因素的各个水平间的主要因变量的均值呈现统计显著性时，必要时可作两两水平间的比较，称为均值间的两两比较。 3) 因素间的交互影响

多因素的试验设计，有时需要分析因素间的交互影响（interaction），2个因素间的交互影响称为一级交互影响，例如因素A与因素B的一级交互影响可记为A×B，3个因素间的交互影响称为二级交互影响，例如因素A与因素B与因素C的二级交互影响可记为A×B×C。当交互影响项呈现统计不显著时，表明各个因素独立，当呈现统计显著时，就需要列出这个交互影响项的效应，以助于作出正确的统计推断。

二、 单因素方差分析

单因素方差分析（one factor ANOVA 或one-way ANOVA）或称为完全随机设计的方差分析（completely random design ANOVA）。试验设计时按受试对象的抽取或分组的随机程度不同可细分为以下两类：

? 完全随机设计——从符合条件的总体中完全随机地抽取所需数目的受试对象，再

将全部受试对象完全随机地分配到k组中去。此时，受试对象与试验因素间无直接联系。

? 组内完全随机设计——按试验因素的k个水平将全部受试对象划分成k个子总体，

再分别从k个子总体中完全随机地抽取所需数目的受试对象。此时，试验因素的各水平决定了受试对象各自应该归属的组别。

设因素A有k个水平A1,A2,?,Ak，在每一个水平下考察的指标可以看成一个总体，现有k个水平，故有k个总体，并假定：

① 每一总体均服从正态分布； ② 每一总体的方差相同；

③ 从每一总体中抽取的样本相互独立。

我们要比较各个总体的均值是否一致，就是要检验各总体的均值是否相同，设第i个总体的均值为?i，那么就是要检验如下原假设：

上海财经大学经济信息管理系IS/SHUFE

Page 2 of 30

H0:?1??2????k

其备选假设为：

H1:?1,?2,?,?k不全相同。

设从第i个总体获得容量为ni的样本观察值为yi1,yi2,?,yini，i?1,2,?,k，各样本间还是相互独立的。样本观察值yij可看成是来自均值为?i的总体，这样yij就是其均值?i与随机误差?ij迭加而产生的。上面我们已经假定在Ai水平下的yij服从N(?i,?2)分布，则有

?ij~N(0,?2)。因此，我们有单因素方差分析的统计模型：

j?1,2,?,ni??yij??i??ij,i?1,2,?,k, ?2??各?ij相互独立,且都服从N(0,?)(25.1)

为了能更仔细地描述数据，常在方差分析模型中引人一般平均与效应的概念。称各个?i的加权平均

1k???ni?i

ni?1为总平均，其中n?(25.2)

?ni?1ki。称

ai??i??,i?1,2,?,k

(25.3)

为因素A在第i水平的主效应，也简称为Ai的效应，同时也表明第i个总体的均值是一般平均与其效应的迭加。容易看出效应间有如下关系式：

?naii?1ki?0

(25.4)

此时，单因素方差分析的统计模型可改写成包含效应的形式：

?yij??i?ai??ij,i?1,2,?,k,?k???niai?0?i?1?各?相互独立,且都服从N(0,?2)ij?所要检验的原假设也可改写成：

j?1,2,?,ni

(25.5)

H0:a1?a2???ak?0

现在，我们知道造成各yij间差异的原因可能有两个：一个可能是假设H0不真，即各水

上海财经大学经济信息管理系IS/SHUFE

Page 3 of 30

平下总体均值?i（或水平效应ai）不同，因此从各总体中获得的样本观察值也就有差异了；另一可能是H0为真，差异是由于随机误差引起的。为了进一步定量分析这些差异，我们需要把这些差异表达出来。由（25.1）可推导出：

yi???i??i?

其中yi??(25.6)

?yj?1niij/ni，?i????ij/ni。即组内样本观察值的平均值等于组内总体均值加上

j?1ni组内随机误差的平均值。还可由（25.5）推导出：

y???? 其中y?(25.7)

??yi?1j?1kniij/n，?????ij。即所有样本观察值的平均值等于总平均（各组均值的

i?1j?1kni加权平均）加上所有随机误差的平均值。这样，每一个观察值yij与总平均y的偏差可以分解成两部分：

yij?y?(yij?yi?)?(yi??y)

其中yij?yi?称为组内偏差，由（25.1）和（25.6）代入得到：

(25.8)

yij?yi??(?i??ij)?(?i??i?)??ij??i?

(3.2.9)

说明组内偏差仅仅反映了随机误差。而yi??y称为组间偏差，由（25.6）、（25.7）和（25.3）代入得到：

yi??y?(?i??i?)?(???)?ai??i???

说明第i组间偏差除了反映随机误差外还反映了第i个水平的效应ai。

各yij间总的差异大小可用总偏差平方和ST表示：

(25.9)

ST???(yij?y)2

i?1j?1kni(25.10)

由（25.9）随机误差引起的数据间的差异可以用组内偏差平方和表示，也称误差偏差平方和Se：

Se???(yij?yi?)2

i?1j?1kni(25.11)

由于组间偏差除了随机误差外，还反映了效应的差异，故由于效应不同引起的数据差异可以用组间偏差平方和表示，也称因素A的偏差平方和SA：

上海财经大学经济信息管理系IS/SHUFE

Page 4 of 30

SA??ni(yi??y)2

i?1k(25.12)

将表示总差异的平方和进行分解：

ST???(yij?y)???(yij?yi??yi??y)22i?1j?1knii?1j?1kknikni???(yij?yi?)???(yi??y)?2??(yij?yi?)(yi??y)22i?1j?1knii?1j?1ki?1j?1nikni (25.13)

???(yij?yi?)??ni(yi??y)22i?1j?1i?1?Se?SA其中

?(yj?1niij?yi?)?0。证明了：总的差异=组内差异+组间差异。由于

11?2?(yj?1niij?yi?)?2?2?(?j?1niij??i?)2~?2(ni?1)

(25.14)

又由?2分布的可加性可知

?1????2i?1??2Sekk?2(yij?yi?)??~?(?(ni?1))??2(n?k) ?j?1i?1?ni2(25.15)

还可证明，在H0为真时，即各组效应ai都为0

SA?因此可采用统计量

2~?2(k?1)

(25.16)

F?来假设检验。

SA/(k?1)~F(k?1,n?k)

Se/(n?k)(25.17)

三、 多重比较

当k组均值比较，如果经过F检验拒绝原假设，表明因素A是显著的，即k个水平对应的指标均值不全相等，但不一定两两之间都有差异。在一些实际问题中，当方差分析的结论是因素A显著时，还需要我们进一步去确认哪些水平间是确有差异的，哪些水平间无显著差异。同时比较任意两个水平均值间有无显著性差异的问题称为多重比较，即要以显著性水平

?，同时检验以下Ck2个假设：

上海财经大学经济信息管理系IS/SHUFE

Page 5 of 30

ijH0:?i??ji?j,i,j?1,2,?,k

(25.18)

均值间的多重比较的方法从形式上可分为几类：临界值相对固定的两两比较、临界值不

固定的多级检验、全部处理组均值与一个对照组均值比较。每一种类型中，根据所控制误差的类型和大小不同，又有许多不同的具体方法。如T（成组比较t检验法）、Bon（Bonforroni t检验法）、Dunnett（与对照组均数比较）、SNK（Student-Newman-Keuls或称q检验法）、Tukey（学生化极差HSD或称最大显著差）、Duncan（新多极差检验法）、LSD（最小显著差）、SIDAK（Sidak不等式进行校正t检验法）、SCHEFFE（Scheffe的多重对比检验）、Waller-Duncan（k比率t检验）、GT2或SMM（学生化最大模数和Sidak不等式进行校正t检验法）、REGWF（多重F检验）、REGWQ（多重极差检验）。

在多重比较时，选用什么样的检验方法，首先要注意每种方法适用的试验设计条件，其次要关心所要控制的误差类型和大小。例如，某因素有10个水平，若采用通常的t检验进行

2多重比较，共需要比较的次数为C10?45次，即使每次比较时都把第一类错误?控制在0.05

水平上，但经过45次多重比较后，犯第一类错误的概率上升到：1?(1?0.05)45?0.90。从中我们可以看到选用t检验法进行多重比较，仅仅控制了每次比较的显著水平，但却大大增加了整体的显著水平。

下面是所要控制的几种误差类型和选用的检验方法： ? 第一类误差率——即犯第一类错误的概率?。 ? 比较误差率——即每一次单独比较时，所犯第一类错误的概率。可使用T法、LSD

法、DUNCAN法。

? 试验误差率——即完成全部比较后，整体所犯第一类错误的概率。

? 完全无效假设下的试验误差率——即在H0假设完全无效下的试验误差率。可使用

SNK法。

? 部分无效假设下的试验误差率——即在H0假设部分无效下的试验误差率。 ? 最大试验误差率——即在在H0假设完全或部分无效下，完成全部比较后所犯第一类错误的最大概率。可使用BON法、SIDAK法、SCHEFFE法、TUKEY法、GT2/SMM

法、GABRIEL法、REGWQ法、REGWF法、DUNNETT法。

1) T检验和Bonforroni检验

当因素有k个水平时，对任意两个水平均值间的差异的显著性检验，可用 t统计量

tij?yi??yj??Se??1?1??n?k??ninj?~t(n?k)

(25.19)

2两两比较的次数共有m?Ck=k(k?1)/2，因此，共有m个置信水平，每次比较的显著水平：

T检验的方法取?。完成所有比较后的整体显著水平等于

1?(1??)m

(25.21)

当比较次数m越大，试验误差就越大。而Bonforroni检验的方法取?/m。完成所有比较后

上海财经大学经济信息管理系IS/SHUFE

Page 6 of 30

的整体显著水平等于

1?(1??/m)m??

(25.22)

即最大试验误差率小于?。 2) LSD检验

既可以通过两两比较的显著水平的特定限制来控制最终的试验误差率，也可以通过两两比较的绝对差异界限来判别显著性。最容易想到的这个界限就是在两两比较中采用的t检验法而得到Fisher最小显著差（LSD）为

Se?11??LSDij?t?(n?k)??

?n?k?ninj?2?当yi??yj??LSDij时，则P??。

(25.23)

3) SNK检验和Duncan检验

SNK法和Duncan法都属于多级检验法中的一种，使用多级检验可以获得同时检验的更高效率。多级检验分为步长增加法和步长减少法，SAS系统采用步长减少法。当因素有k个水平时，即有k个均值需要比较，检验步骤为：

① 将均值由大到小排队，即y1??y2??,?,?yk?。

② 比较y1?与yk?是否有显著差异。此时跨度a?k。若两者之间无显著差异，说明

其他均值之差比它小的任何两个水平均值之间的差别也无显著性，所以停止一切比较；反之，则继续进行下一步。

③ 比较y1?与yk?1?，比较y2?与yk?是否有显著差异。此时这2个比较的跨度

a?k?1。若两者之间的比较无显著差异，则停止一切比较。如果每一步都有不

满足停止比较的对比组存在，最后应到达跨度为2的所有需要比较的相邻两水平均值间都作完比较时为止。

多级检验在作每一级比较时，通过控制比较误差率?a的显著水平来实现其最终要控制的试验误差率。要注意的是?a在每一级比较时可能是不同的，它是跨度a和整体试验误差率?的函数，即?a?f(a,?)。另外，要注意的是?a其实就是每一级比较时特定统计量分布的显著水平。常用的两种方法是SNK检验和Duncan检验。它们的检验统计量为q（也称学生化极差统计量），如下

qij?yi??yj??Se??1?1??2(n?k)??ninj?~q(a,n?k)

(25.24)

其中a是yi?和yj?之间的跨度值，q分布的自由度是a和n?k，显著水平为?a。SNK检验和Duncan检验的区别主要在于?a取值

上海财经大学经济信息管理系IS/SHUFE

Page 7 of 30

? SNK检验：?a??。注意，当比较次数很大时，最大试验误差率将趋向于1。 ? Duncan检验：?a?1?(1??)a?1。

四、 随机单位组设计的方差分析

随机单位组设计（randomized block design）又称随机区组设计或随机配伍组设计，它是两样本配对试验的扩大。欲比较因素A中的k个水平的各个均值，试验设计时，先将受试对象按性质相同或相近者组成单位组，每个单位组有k个受试对象，分别随机分配到因素A的k个水平上。这时每个水平的受试对象，不仅数量相同，而且性质也相同或相近，就能缩小误差，提高试验效率。这样的设计可将单位组看作一个因素，就成为两个因素的设计（因素与单位组），由于两个因素的各水平仅仅交叉1次，所以重复数为1，在这样的意义下，随机单位组设计可看作为两因素重复数为1的设计，一般这种设计不考虑交互影响。

设有因素A具有k个水平，受试对象按性质相同或相近者分成b个单位组，每个单位组有k个受试对象，分别随机分配到因素A的k个水平上。那么，随机单位组设计的方差分析表见表25.2所示：

表25.2 方差分析表形式

变异来源 source 因素A 单位组 误差Se 总变异ST 离差平方和 SS SSA SS单 SSe 自由度 df 均方 MS MSA= SSA/( k－1) MS单= SS单/( b－1) F统计量 F FA= MSA/ MSe F单= MS单/ MSe P概率值 P PA P单 PT k－1 b－1 bk－k－b+1 MSe= SSe/( bk－k－b+1) MST= SST/( bk－1) FT= MST/ MSe SST= SSA+ SS单+SSe bk－1

五、 析因设计的方差分析

析因设计（factorial design）是一种多因的设计。各因素在试验中所处的地位基本平等，而且因素之间存在一级（即2个因素之间）、二级（即3个因素之间）乃至更复杂的交互作用。例如，两个因素时，第1个因素有3个水平，第2个因素有2个水平，全部水平组合共有3×2=6种组合，每种组合都作试验时就是析因试验设计，也可称为3×2析因试验设计。同样3×4×2析因试验设计，则代表3个因素，分别有3，4，2个水平，全部试验后的水平组合为3×4×2=24种。在每一种组合下，适当重复几次，称为重复数。重复数可以不相等，一般地说，重复数相等时，效率最高。

析因设计能够检验每个因素的各水平间主要变量的平均值的统计差异，也能检验因素间的交互影响。当存在交互影响时，表示一个因素各水平间的差异会随着另一个因素的水平改变而不同；当不存在交互影响时，则各个因素独立，即一个因素的水平改变时不影响另一个因素的各个水平之效应。析因设计的方差分析因为能研究交互影响，所以能提供较多信息。但是，当有较高级（二级以上）的交互影响时，由于涉及多个因素，各有多个水平，情况将错综复杂，可能会引起解释上的困难。

析因设计的方差分析同样是从数据差异的总平方和开始分解。例如，对于A×B双因素方差分析，这个总差异能分解成：A因素的各个水平之间的差异，B因素的各个水平之间的差

上海财经大学经济信息管理系IS/SHUFE

Page 8 of 30

异，A与B的各种不同组合之间的差异，以及观察数据必然会产生的随机误差这四部分。方差分析的主要目的就是要将这四部分从总平方和中分离出来，再以各个平方和与误差平方和作比较。假设A因素有r个水平，B因素有c个水平，每一种水平下的重复数为m，那么总的观察数据有n=r×c×m个，方差分析表见表25.3所示：

表25.3 双因素（r×c）重复数m的方差分析表形式

变异来源 source 因素A 因素B A×B 误差Se 总变异ST 离差平方和 SS SSA SSB SSAB SSe SST= SSA+ SSB+ SSAB +SSe 自由度 df 均方 MS MSA= SSA/( r－1) MSB= SSB/( c－1) F统计量 F FA= MSA/ MSe FB= MSB/ MSe FT= MST/ MSe P概率值 P PA PB PAB PT r－1 c－1 (r－1)(c－1) MSAB= SSAB/(( r－1)( c－1)) FAB= MSAB/ MSe r×c×(m－1) MSe= SSe/( rc(m－1)) r×c×m－1 MST= SST/( rcm－1)

六、 拉丁方设计的方差分析

若试验中涉及到3个因素，当它们之间不存在交互作用或交互作用可以忽略不计，且各因素均取相同水平时，适合于选择拉丁方设计。用K个拉丁字母排成K行K列的方阵，使每行每列中每个字母仅出现1次，这样的方阵称为拉丁方（latin square）。然后将3个因素分别放置到拉丁方的行、列及字母上面。例如，三个4×4的拉丁方为：

A B C D A B C D A B C D B A D C B A D C D C B A D C B A C D A B B A D C C D A B D C B A C D A B

四个5×5的拉丁方为：

A B C D E A B C D E A B C D E A B C D E B C D E A C D E A B D E A B C E A B C D C D E A B E A B C D E A B C D D E A B C D E A B C B C D E A B C D E A C D E A B E A B C D D E A B C C D E A B B C D E A

使用时可选择其中一个。拉丁方试验设计的关键是这3个因素之间不存在交互作用或者交互作用可以忽略不计，一般情况是仅涉及到1个试验因素，因此就不存在交互作用。试验因素有K个水平（如A、B、C、D、E），还有2个是非处理因素，或者说是2个区组因素，让这2个区组因素也正好取K个水平，同时把这2个区组因素放在K×K拉丁方阵的横向和纵向上，构成了K×K个区组水平组合，每种组合下伴有试验因素K个水平中的1个水平。

七、 proc anova和proc glm过程

SAS系统的STAT软件提供了anova过程和glm过程等几个过程进行方差分析。anova过程主要处理均衡数据，所谓均衡数据是指自变量（或称分类变量）的每种组合中的观察数是相等的，如果不相等则称为非均衡数据。。虽然glm过程能够处理均衡和不均衡的两种数据，但是anova过程考虑到均衡设计的特殊构造，对于均衡数据使用anova比使用glm计算快且占用存储少，还可以处理拉丁方设计、若干不完全的均衡区组设计等等。因此，无论何

上海财经大学经济信息管理系IS/SHUFE

Page 9 of 30

时作方差分析，一旦可能都应该用anova过程来完成。如果试验设计不均衡，也不是上述的几种特殊情况之一，那么应该使用glm过程。 1. anova过程的语句格式

anova过程的主要控制语句如下：

proc anova 输入数据集名 <选项列表> ;

class 变量列表 ;

model 因变量列表=自变量列表 ; means 效应列表 ; test E=效应列表; run ;

其中class语句、model语句是必需的，而且class语句必须出现model语句之前。test语句必须放在model语句之后。

1) proc means 语句中的<选项列表>。

? manova——按多元方式删除那些含有丢失值的观察，也即在因变量中有丢失值就从这次分析中删除这个观察。

? outstat=输出数据集名——生成一个输出数据集，它包含模型中每个效应的平方和、F统计量和概率水平。 2) class语句。

在anova过程中要使用的分类变量、区组变量必须首先在class语句的变量列表中说明。Class语句是必需的，且必须放在model语句前面。Class变量可以是数值型，也可以是字符型。

3) model 语句。

该语句用来规定因变量和自变量效应。如果没有规定自变量的效应，则只拟合截距，假设检验为因变量的均值是否为0。Model语句的主要形式有四种：

① 主效应模型 Model y=a b c;

② 含有交叉因素的模型

Model y=a b c a*b a*c b*c a*b*c; ③ 嵌套模型

Model y=a b c(a b);

④ 包含嵌套、交叉和主效应的模型 Model y=a b(a) c(a) b*c(a); Model语句的选项列表有：

int——打印与截距有关的假设检验结果。anova过程总是把截距作为模型的一个效应进行处理，缺省时，不打印结果。

? nouni——不打印单变量分析结果。 4) means 语句。

该语句是用来计算在means语句后列出的每个效应所对应的因变量均值。Anova过程可以对出现在model语句等号右边的任一效应计算因变量的均值。不过这些均值没有针对模型中的效应进行修正。如果需要修正的均值，应该调用glm过程，使用其中的lsmenas语句。在anova过程里可以使用任意多个means语句，它们放在model语句后面。

Means语句的选项列表主要有两个内容，一是选择多重比较的检验方法，二是规定这些

上海财经大学经济信息管理系IS/SHUFE

Page 10 of 30

检验的细节，注意这些细节选项只能用于主效应。

① 多重比较的检验方法

? bon——对所有主效应均值之差进行Bonferroni的t检验。 ? duncan——对所有主效应均值进行Duncan的多重极差检验。

? dunnett<(‘格式化对照值’)>——进行Dunnett的双尾t检验。用以检验对所有主效应均值的某个水平作为对照，处理有无显著差异。为了规定这个对照效应的水平，在括号内用单引号把这个水平的格式化值括起来。缺省时，效应的第一个水平作为对照。

? dunnettl<(‘格式化对照值’)>——进行Dunnett的单尾t检验。它检验是否任一个处理显著地小于这个对照。

? dunnettu<(‘格式化对照值’)>——进行Dunnett的单尾t检验。它检验是否任一个处理显著地大于这个对照。

? gabriel——对所有主效应均值进行Gabriel的多重对比检验。

? regwf——对所有主效应均值进行Ryan-Einot-Gabriel-Welsch的多重F检验。 ? regwq——对所有主效应均值进行Ryan-Einot-Gabriel-Welsch的多重极差检验。 ? scheffe——对所有主效应均值进行Scheffe的多重对比检验。 ? sidak——对所有主效应均值水平依据Sidak不等式进行调整后，对其均值之差两两进行t检验。

? Smm|gt2——当样本量不等时，基于学生化最大模和Sidak不相关t不等式，等到Hochberg的GT2方法，对主效应均值进行两两对比检验。

? snk——对所有主效应均值进行Student-Newman-Keuls的多重极差检验。

? t|lsd——对所有主效应均值进行两两t检验，它相当于在单元观察数相等时Fisher的最小显著差（Fisher’s least-significant-difference）检验。

? tukey——对所有主效应均值进行Tukey的学生化极差检验。

? waller——对所有主效应均值进行Waller-Duncan的k比率（k-ratio）检验。 ② 多重比较的检验细节

? alpha=p——给出均值间对比检验的显著性水平。缺省值是0.05。 ? cldiff——要求把两两均值之差的结果用置信区间的形式输出。 ? clm——对变量的每个水平的均值按置信区间形式输出。

? e=效应——指定在多重对比检验中所使用的误差均方。如果缺省，使用残差均方（MS）。指定的效应必须是在model语句中出现过的效应。

? kratio=值——给出Waller-Duncan检验的类型1/类型2的误差限制比例。Kratio的合理值为50、100、500，大约相当于两水平时alpha值为0.1、0.05、0.01。缺省值为100。

? lines——按下降次序列出所有检验方法产生的均值，并用一条线段在均值旁指出非显著的子集。

? hovtest——要求输出组间方差齐性的Levene检验。 5) test语句

? 在分析中，如果这个语句缺省，仍然使用残差均方（MS）作为误差项对所有平方和（SS）计算F值。但用户可以使用本语句要求使用其他效应作为误差项，得到另外的F检验。可以使用多个test语句，把它们放在model语句后面。Test语句的选项为：

? h=效应——规定模型里哪些效应用来作为假设的效应。

? e=效应——规定一个而且只能是一个效应用来作为误差项，这个说明项是必须的。 2. glm过程的语句格式

proc glm是分析符合一般线性模型（General Linear Models）的数据，因此取名GLM。它能被用在许多不同的分析中，如简单回归、多元回归、方差分析、协方差分析、加权回归、

上海财经大学经济信息管理系IS/SHUFE

Page 11 of 30

多项式回归、偏相关分析、多元方差分析等。

在glm过程中的大多数方差分析的语句和选项与anova过程中基本相同。用anova过程编写的程序几乎不用修改就可在glm过程中运行。glm过程仅仅是附加了三条语句：contrast、estimate和lsmeans。contrast和estimate语句允许你测试和估计均值的某种功能。lsmeans语句允许你计算调整后的均值。

glm过程的主要控制语句如下：

proc glm 输入数据集名 <选项列表> ;

class 变量列表 ;

model 因变量列表=自变量列表 ; contrast ‘标签’ 效应 值表 ; estimate ‘标签’ 效应 值表 ; lsmeans 效应列表 ; means 效应列表 ;

output <统计量关键字=变量名列表>; test E=效应列表; run ;

其中class语句、model语句是必需的，而且class语句必须出现model语句之前。其他语句必须放在model语句之后。下面主要介绍与anova过程相比不同的语句和新增加的语句。 1) model语句。

在glm过程的model语句中可以使用几种不同效应，下面是使用这些效应的几个例子，a、b和c代表分类变量；y1、y2、x1和x2代表连续变量。

Model y=x1; (简单回归) Model y=x1 x2; (多重回归) Model y=x1 x1*x1; (多项式回归) Model y1 y2=x1 x2; (多元回归)

Model y=a; (单因素方差分析) Model y=a b c; (主效应模型) Model y=a b a*b; (因素模型) Model y=a b(a) c(b a); (嵌套模型)

Model y1 y2=a b; (多元方差分析模型) Model y=a x1 (协方差分析模型)

Model语句的主要选项有（与anova过程中的model语句选项相同不再列出）： ? solution——打印正规方程的解，即参数估计值。

e1/e2/e3/e4——打印模型中每一效应的类型1/类型2/类型3/类型4的可估函数，并计算相应的平方和。

ss1/ss2/ss3/ss4——对每个效应，打印与类型1/类型2/类型3/类型4的可估函数相关的平方和。

alpha=0.01/0.05/0.1——指定置信区间的?水平。缺省值为0.05。

cli/clm——打印每一观察的预测值/预测均值的置信限，两者不能同时使用。 p——打印自变量没有缺失值的每一观察值、预测值和残差值。同时还打印Durbin-Waston统计量。

xpx——打印叉积矩阵X?X。

i——打印矩阵X?X的逆矩阵或广义逆矩阵。

上海财经大学经济信息管理系IS/SHUFE

Page 12 of 30

2) contrast语句。

提供一种获得一般假设检验的技巧。其中，效应可以是截距，用字符intercept表示。通过规定L向量或M矩阵来构造一元假设检验L??0或多元假设检验L?M?0。例如，当发现某两个因素的交互作用项有显著性时，我们可用本语句来实现一个因素被控制在某水平上，对另一个因素的各水平间进行两两比较的目的。

设M因素有三个水平a、b、c，V因素有两个水平1、2，且M?V有显著性。如果我们要比较

?a?(?b??c)

的差异，那么有几种不同的比较方法：

① 在因素V的每一个水平上，分别比较因素M的三个水平a、b、c均值的之间的线

性关系假设是否显著。也即

12H0:?a1?0.5?b1?0.5?c1?0和H0:?a2?0.5?b2?0.5?c2?0。

② 在因素V平均的所有水平上，比较因素M的三个水平a、b、c均值的之间的线性

关系假设是否显著。也即

H0:0.5(?a1?0.5?b1?0.5?c1)?0.5(?a2?0.5?b2?0.5?c2)?0。

③ 在因素V平均的子集上，比较因素M的三个水平a、b、c均值的之间的线性关系

假设是否显著。也即

H0:(?a1?0.5?b1?0.5?c1)?(?a2?0.5?b2?0.5?c2)?0

glm模型为双因素试验设计的方差分析指定了下面的效应公式

?ij????i??j?(??)ij

(25.25)

其中?ij是因素Mi水平与因素Vj水平在ij单元上所有观察值的平均。?为总平均。?i是因素M在i水平上的主效应，?j是因素V在j水平上的主效应，(??)ij为因素M和因素V在ij水平上的交互效应。因此，对任一观察值有

yijk??ij??ijk????i??i?(??)ij??ijk (25.26)

因此，根据单元均值给出的线性组合可以转换成效应模型的合并参数形式，即L??0，如

?a1?0.5?b1?0.5?c1????a??1?(??)a1?0.5??0.5?b?0.5?1?0.5(??)b1 ?0.5??0.5?c?0.5?1?0.5(??)c1??a?0.5?b?0.5?c?(??)a1?0.5(??)b1?0.5(??)c1同理

?a2?0.5?b2?0.5?c2??a?0.5?b?0.5?c?(??)a2?0.5(??)b2?0.5(??)c2

上海财经大学经济信息管理系IS/SHUFE

Page 13 of 30

相应的glm过程的语句为

proc glm ;

class M V ; model Y=M V M*V;

contrast ‘a vs b,c in v1’M 1 -0.5 -0.5 M*V 1 0 -0.5 0 -0.5 0; contrast ‘a vs b,c in v1’M 1 -0.5 -0.5 M*V 0 1 0 -0.5 0 -0.5; run ;

Contrast语句中的可选项： e——打印整个L向量。

e=效应——规定模型中的某个效应作为误差项。过程将把这一效应作为单变量F检验的分母。如果缺省，过程把均方误差（MSE）作为误差项。

etype=n——指明e=效应的类型（1、2、3、4）。如果指明e=而没有指明etype=，则使用最高类型。

3) Estimate语句

可用来估计参数的线性函数，通过用参数的估计b乘以向量L来得到Lb。其中

b?(X?X)?X?Y。Estimate语句的使用格式同contrast语句。

estimate语句中的可选项： e——打印整个L向量。

divisor=数字——为简便地输入效应的系数而规定的一个值，用该值除以所有系数使得分数系数可以作为整数输入。例如

estimate ‘1/3(a+b)－2/3c’ M 1 1 -2 /divisor=3; 可替代

estimate ‘1/3(a+b)－2/3c’ M 0.33333 0.33333 -0.66667; 4) Lsmeans语句

计算列在语句中的每一效应的最小二乘均值（LSM）。最小二乘均值估计是针对非均衡数据设计的，而类和子类的算术平均值是针对均衡数据设计的。

lsmeans语句中的可选项：

cov——在选项out=指明的输出数据集中输出协方差。 e——打印用以计算最小二乘均值的可估函数。 e=效应——规定模型中的某个效应作为误差项。 etype=n——指明e=效应的类型（1、2、3、4）。

out=输出数据集名——产生一个包含LSM值、标准差及协方差的输出数据集。

pdiff——打印假设检验H0:LSM(i)?LSM(j)的所有可能的概率值。 stderr——打印LSM的标准差和H0:LSM?0的概率值。

tdiff——打印假设检验H0:LSM(i)?LSM(j)的t值和相应的概率值。

pdiff=all/control/conroll/controlu——打印最小二乘均值之差的概率值。

adjust=bon/dunnett/scheffe/sidak/smm/gt2/tukey/t——要求多重比较对最小二乘均值之差的概率值和置信限进行调整。缺省值为t。

slice=效应——通过规定的这个效应来分开交叉的LSM效应。例如，假定交叉项A*B是显著的，如果想对B的每个效应检验A的效应，使用下面语句：

上海财经大学经济信息管理系IS/SHUFE

Page 14 of 30

lsmeans A*B /slice=B;

八、 实例分析

1． 单因素试验设计的均值比较

例25.1 考虑在5种不同品牌的人工合成胶合板材料上进行磨损时间测试，每种品牌的材料做四次试验，且都是采用的同一种磨损措施，所有的试验都是在完全随机的顺序下在相同的机器上完成的。

程序如下：

data study.veneer;

input brand $ wear @@; cards;

ACME 2.3 ACME 2.1 ACME 2.4 ACME 2.5 CHAMP 2.2 CHAMP 2.3 CHAMP 2.4 CHAMP 2.6 AJAX 2.2 AJAX 2.0 AJAX 1.9 AJAX 2.1 TUFFY 2.4 TUFFY 2.7 TUFFY 2.6 TUFFY 2.7 XTRA 2.3 XTRA 2.5 XTRA 2.3 XTRA 2.4 ;

proc anova data=study.veneer;

class brand; model wear=brand; means brand;

means brand /hovtest;

run;

程序说明：因为数据仅仅是按照brand值分类，所以在class语句中这是仅有的一个变量。变量wear是被分析的因变量，故wear出现在model语句等号的左边。在方差分析表中，除了总方差和误差外，方差的来源仅仅是由于各种不同brand值的变异造成的，因此brand出现在model语句等号的右边。Means语句计算主效应brand不同水平所对应的因变量均值，选项hovtest计算不同品牌组方差齐性的假设检验。

输出的结果见表25.4所示：

上海财经大学经济信息管理系IS/SHUFE

Page 15 of 30

The SAS System Analysis of Variance Procedure Class Level Information Class Levels Values BRAND 5 ACME AJAX CHAMP TUFFY XTRA Number of observations in data set = 20 Analysis of Variance Procedure Dependent Variable: WEAR Source DF Sum of Squares Mean Square F Value Pr > F Model 4 0.61700000 0.15425000 7.40 0.0017 Error 15 0.31250000 0.02083333 Corrected Total 19 0.92950000 R-Square C.V. Root MSE WEAR Mean 0.663798 6.155120 0.14433757 2.34500000 Source DF Anova SS Mean Square F Value Pr > F BRAND 4 0.61700000 0.15425000 7.40 0.0017 Levene's Test for Equality of WEAR Variance ANOVA of Squared Deviations from Group Means Sum of Mean Source DF Squares Square F Value Pr > F BRAND 4 0.000659 0.000165 0.5310 0.7149 Error 15 0.00466 0.00031 Analysis of Variance Procedure 表25.4 单因素设计的方差分析结果

结果分析：anova过程总是输出两个基本的方差分析表。一个是总体模型的方差分析表，一个是包含模型中各个变量的方差分析。首先输出class语句中规定的每个变量（brand）、分类变量的取值数（5）、具体取值（ACME AJAX CHAMP TUFFY XTRA）以及数据集中的观察个数（20）。

接着anova过程对model语句中每个因变量输出方差分析表。包括：因变量的总平方和（0.9295）、属于模型部分的平方和（0.6170）、属于误差部分的平方和（0.3125）、自由度DF（4、5、19）、模型的均方MS（0.15425=0.617/4）、误差的均方MSE（0.02083333=0.3125/15）、模型的F值（7.40=0.15425000/0.02083333）、分布大于7.40的概率（0.0017）、R2（0.663798=0.617/0.9295）、变异系数CV（6.155120=100×0.0208333、因变量的/2.345）标准差（0.14433757=0.0208333）、因变量均值（2.345）。

对模型中的每个效应，anova过程还输出方差分析表。brand自由度DF（4）、平方和（0.617）、

上海财经大学经济信息管理系IS/SHUFE

Page 16 of 30

均方MS（0.15425=0.617/4）、F值（7.40=0.15425000/0.02083333）、分布大于7.40的概率（0.0017）。

总体F检验是显著的（0.0017<0.05），表明模型是有意义的。品牌brand的F检验也是显著的（0.0017<0.05），表明不同品牌的均值不全相等。这里两个F检验是完全相同的，这仅仅是因为在模型中只有一项brand。注意，我们可以用glm过程替代这个anova过程，能得到相同的方差分析结果。最大区别是glm过程将计算每个效应的类型1和类型3平方和，而anova只计算类型1的平方和。对于单因素和多因素平衡数据来说，anova过程的SS1、glm过程的SS1和SS3都相同。

Levene的 方差齐性检验结果表明：不能拒绝（0.7149>0.05）不同品牌组里观察值的方差是相等的原假设。

最后输出的是每种品牌的观察数、均值和标准差。例如，ACME品牌的观察数为4，均值为2.32500000，标准差为0.17078251。 2． 均值的多重比较和置信区间

例25.2 继续上例的分析。由于品牌brand的F检验是显著的（0.0017<0.05），表明5种不同品牌的均值不全相等，但可能存在某2个或某3个或某4个品牌的均值相同。因此，常需要进一步的均值多重比较和置信区间分析。

程序如下：

proc anova data=study.veneer;

class brand; model wear=brand; means brand /ducan;

means brand /lsd clm cldiff;

run;

程序说明：第一个means语句选用了ducan选项，要求计算输出组间均值比较的新多重极差检验，结果见表25.5（a）。第二个means语句选用了lsd clm选项，对所有组均值进行两两t检验，输出各组均值的置信区间，结果见表25.5（b）。第二个means语句还选用了lsd cldiff选项，将对各组间均值之差采用最小显著差检验，输出各组间均值之差的置信区间，结果见表25.5（c）。

表25.5(a) Duncan的新多重极差检验

上海财经大学经济信息管理系IS/SHUFE

Page 17 of 30

The SAS System Analysis of Variance Procedure Duncan's Multiple Range Test for variable: WEAR NOTE: This test controls the type I comparisonwise error rate, not the experimentwise error rate Alpha= 0.05 df= 15 MSE= 0.020833 Number of Means 2 3 4 5 Critical Range .2175 .2280 .2346 .2390 Means with the same letter are not significantly different. Duncan Grouping Mean N BRAND A 2.6000 4 TUFFY A B A 2.3750 4 XTRA B A 表25.5（a）中结果分析：注意到各组均值按大到小排列（2.60，2.375，2.375，2.3250，2.0500），在标题“Duncan Grouping”下是一系列字母A、B、C，如果均值间差异不显著标上相同的字母，否则标上不同的字母。对于Duncan多重极差检验来说，5个均值之间的比较，只要看最大的均值与最小的均值之差的是否大于临界值0.239，因为2.600－2.050=0.55>0.239，则为显著，所以品牌TUFFY的均值不同与AJAX，应该标识不同的字母。因为存在5个均值之间最大差的显著性，接下来就需要比较4个均值之间差的显著性，临界值为0.2346。2.600－2.325=0.275>0.2346，显著，2.375－2.050=0.325>0.2346，显著，只要存在一个显著性，就需要继续比较3个均值之间差的显著性。虽然，均值2.600、2.375和2.375之间的差小于0.2280，均值2.375、2.375和2.325之间的差也小于0.2280，但由于存在2.375－2.050=0.325>0.2280，显著，继续比较2个均值之间差的显著性。2.600－2.375=0.225>0.2175，显著，2.325－2.050=0.275>0.2175，显著，其他相邻两均值比较不显著。

表25.5(b) 各组均值的t检验置信区间

T Confidence Intervals for variable: WEAR Alpha= 0.05 Confidence= 0.95 df= 15 MSE= 0.020833 Critical Value of T= 2.13 Half Width of Confidence Interval= 0.153824 Lower Upper BRAND N Confidence Mean Confidence Limit Limit TUFFY 4 2.44618 2.60000 2.75382 XTRA 4 2.22118 2.37500 2.52882 表25.5（b）中结果分析：均值t分布的95%置信区间的一半宽度为0.153824，因此TUFFY品牌均值置信区间的下限为2.600－0.153824=2.44618，上限为2.600＋0.153824=2.75382。其他品牌均值的置信区间计算，同样是均值加减0.153824而得到的。

表25.5(c) lsd最小显著差检验

上海财经大学经济信息管理系IS/SHUFE

Page 18 of 30

The SAS System T tests (LSD) for variable: WEAR NOTE: This test controls the type I comparisonwise error rate not the experimentwise error rate. Alpha= 0.05 Confidence= 0.95 df= 15 MSE= 0.020833 Critical Value of T= 2.13145 Least Significant Difference= 0.2175 Comparisons significant at the 0.05 level are indicated by '***'. Lower Difference Upper BRAND Confidence Between Confidence Comparison Limit Means Limit TUFFY - XTRA 0.0075 0.2250 0.4425 *** TUFFY - CHAMP 0.0075 0.2250 0.4425 *** TUFFY - ACME 0.0575 0.2750 0.4925 *** TUFFY - AJAX 0.3325 0.5500 0.7675 *** XTRA - TUFFY -0.4425 -0.2250 -0.0075 *** XTRA - CHAMP -0.2175 0.0000 0.2175 XTRA - ACME -0.1675 0.0500 0.2675 XTRA - AJAX 0.1075 0.3250 0.5425 *** CHAMP - TUFFY -0.4425 -0.2250 -0.0075 *** CHAMP - XTRA -0.2175 0.0000 0.2175 CHAMP - ACME -0.1675 0.0500 0.2675 CHAMP - AJAX 0.1075 0.3250 0.5425 *** ACME - TUFFY -0.4925 -0.2750 -0.0575 *** 表25.5（c）中结果分析：注意在显著水平为0.05上，两两比较的最小显著差为0.2175，如果显著则被标上“***”。例如，TUFFY均值减XTRA均值=2.600－2.375=0.225>0.2175，显著。综合分析的结果表明，AJAX品牌均值显著与其他品牌均值不同，且为最小的均值；TUFFY品牌均值也显著与其他品牌均值不同，且为最大的均值；XTRA、CHAMP、ACME三个品牌均值之间无显著差异。 3． 有计划的均值比较和参数估计

例25.3 继续上例的分析。有时在实际情况中，多重比较要按某种分类标准来进行，例如，假设我们知道5种品牌的制造商情况，品牌ACMX、AXAX和CHAMP来自美国U.S.制造商，而品牌TUFFY和XTRA来自非美国non-U.S.制造商。如果我们有兴趣比较美国品牌的均值与非美国品牌的均值是否有差异。

程序如下：

proc glm data=study.veneer;

class brand; model wear=brand;

contrast 'US vs NON-U.S.' brand 2 2 2 -3 -3;

上海财经大学经济信息管理系IS/SHUFE

Page 19 of 30

estimate 'US vs NON-U.S.' brand 2 2 2 -3 -3;

run;

程序说明：使用contrast语句来产生有计划的均值比较分析和使用estimate语句进行参数估计。注意在anova过程中没有这两条语句，必须使用glm过程。使用contrast语句前，应该首先表达出所关心的均值线性组合的原假设，如

H0:1/3(?ACME??AJAX??CHAMP)?1/2(?TUFFY??XTRA) 等价于H0:2(?ACME??AJAX??CHAMP)?3(?TUFFY??XTRA)?0contrast语句的三个基本参数，一是标签（'US vs NON-U.S.'），二是效应名（brand），三是效应的数字系数表（2 2 2 -3 -3）。应特别注意的是，数字系数的次序是匹配分类变量按字母数字次序的水平值。事实上，均值线性组合的系数同样是model语句中效应参数组合的系数，这是因为，?i????i，将它们分别代入均值线性组合后，可得到

2(?ACME??AJAX??CHAMP)?3(?TUFFY??XTRA)?2(???ACME????AJAX??CHAMP??XTRA)?3(???TUFFY????XTRA) ?2(?ACME??AJAX??XTRA)?3(?TUFFY??XTRA)所以，estimate语句的使用格式与contrast语句非常类同。

输出的主要结果见表25.6所示：

表25.6 有计划的均值比较和参数估计

The SAS System Contrast DF Contrast SS Mean Square F Value Pr > F US vs NON-U.S. 1 0.27075000 0.27075000 13.00 0.0026 T for H0: Pr > |T| Std Error of Parameter Estimate Parameter=0 Estimate 表3.16中结果分析：显示了美国品牌均值与非美国品牌均值比较的平方和为0.27075，F值为13=0.27075/0.020833，这个F（1，15）分布F值大于13的概率为0.0026小于0.05，因此原假设是显著的，拒绝接受，即美国品牌均值与非美国品牌均值是不同的。效应组合的参数估计为-1.425=3×（2.325+2.050+2.375）－2×（2.600+2.375），对于原假设参数是否为0的t检验，t统计量为-3.60，概率为0.0026小于0.05，拒绝接受。注意到t检验的p值为0.0026，与对比分析的F检验的p值相同，这是因为两种检验是相同的，F值等于t的平方。 4． 随机单位组试验设计的方差分析

例25.4 某食品公司对一种食品设计了四种包装。为了考察哪种包装最受欢迎，选了十个有近似相同销售量的商店作试验，其中两种包装各指定两个商店，另两种包装各指定三个商店销售。在试验期中各商店的货架排放位置、空间都尽量一致，营业员的促销方法也基本相同。观察在一定时期的销售量，数据见表25.7所示。试比较四种包装的销售量是否一致。

表25.7 四种包装在10个商店中的销售量

包装类型 （treat） A1 A2 商店（block） 1 12 14 2 18 12 3 13 商店数 n 2 3 Page 20 of 30

上海财经大学经济信息管理系IS/SHUFE

A3 A4 程序如下：

data pack;

19 24 17 30 21 3 2 input treat $ n ; do block=1 to n; input y @@; output; end; cards; A1 2 12 18 A2 3 14 12 13 A3 3 19 17 21 A4 2 24 30 ;

proc glm data=pack;

class block treat; model y=block treat; means block treat/snk;

means block treat/dunnett('1'); means block treat; run;

程序说明：由于包装类型A1和A4在商店3里没有进行试验，所以产生的是一种不平衡的试验数据，请注意随机单位区组中有不平衡数据集时，程序应该如何设计。本例有两个分组变量，一是处理组treat，包含四个水平A1、A2、A3、A4，二是单位组block，包含三个水平1、2、3，所以在class语句及model语句中都要写上这两个变量名。Means语句中的分组变量名可以是两个，也可以是一个。如果对单位组间的两两比较或各组均值不感兴趣，可去掉block，只用treat。

第一条means语句的snk选项，要求对均值的多级比较采用多极差检验法，也称Student-Newman-Keuls法或q法。第二条means语句的dunnett(‘1’)选项，要求所有分组均值分别与对照组均值进行比较，采用dunnett的双尾t检验，这个t检验是一种多对一t检验，也可用dunnetl（单尾t检验，分组的均值是否显著地小于这个对照组的均值）或dunnetu（单尾t检验，分组的均值是否显著地大于这个对照组的均值）。对照组在括号内规定为’1’，即分组变量的第一个水平分组，第1家商店和A1包装。第三个means语句，仅仅要求输出各个分组的均值和标准差。输出的主要结果见表25.8（a）和（b）所示。

表25.8（a） 随机单位区组试验设计的snk检验

上海财经大学经济信息管理系IS/SHUFE

Page 21 of 30

The SAS System General Linear Models Procedure Class Level Information Class Levels Values BLOCK 3 1 2 3 TREAT 4 A1 A2 A3 A4 Number of observations in data set = 10 Dependent Variable: Y Source DF Sum of Squares Mean Square F Value Pr > F Model 5 269.00000000 53.80000000 6.15 0.0515 Error 4 35.00000000 8.75000000 Corrected Total 9 304.00000000 R-Square C.V. Root MSE Y Mean 0.884868 16.43355 2.95803989 18.00000000 Source DF Type I SS Mean Square F Value Pr > F BLOCK 2 10.50000000 5.25000000 0.60 0.5917 TREAT 3 258.50000000 86.16666667 9.85 0.0256 Source DF Type III SS Mean Square F Value Pr > F BLOCK 2 11.00000000 5.50000000 0.63 0.5789 TREAT 3 258.50000000 86.16666667 9.85 0.0256 Student-Newman-Keuls test for variable: Y NOTE: This test controls the type I experimentwise error rate under the complete null hypothesis but not under partial null hypotheses. Alpha= 0.05 df= 4 MSE= 8.75 WARNING: Cell sizes are not equal. Harmonic Mean of cell sizes= 3 Number of Means 2 3 Critical Range 6.7057809 8.6078784 Means with the same letter are not significantly different. SNK Grouping Mean N BLOCK A 19.250 4 2 A A 17.250 4 1 A A 17.000 2 3 Student-Newman-Keuls test for variable: Y NOTE: This test controls the type I experimentwise error rate under the complete null hypothesis but not under partial null hypotheses. 表25.8（a）中结果分析：总的模型的方差分析结果p=0.0515，基本上是有显著意义的。模型中的变异（269）基本上反映了总的变异（304），判定系数R=0.884868=269/304。对于

上海财经大学经济信息管理系IS/SHUFE

Page 22 of 30

2单因素不平衡数据的方差分析，类型1和类型3的平方和就不相同了，分组变量的变异计算应该采用类型3的平方和。分组变量block的方差分析结果p=0.5789，不具有显著意义，说明食品在3家不同商店进行销售时，销售量的均值没有显著差异。分组变量treat的方差分析结果p=0.0256，具有显著意义，说明4种不同包装食品的销售量的均值具有显著差异，但没有指出具体哪几种包装之间有显著差异。

对block组进行snk多极差检验，3个组比较时，大均值与小均值之差的临界值为8.6078784，而2个组比较时，临界值为6.7057809，结果显示在“SNK Grouping”标题下，3个商店标有相同的字母“A”，说明了3个商店的销售量均值没有显著差异。对treat组进行snk多极差检验，在“SNK Grouping”标题下，出现了标有不相同的字母“A”和“B”，snk的检验方法简单地说，将不同包装类型的食品销售量均值从大排到小，然后从两个最大的跨度组开始比较，例如A4比较A2为27－13=14>10.99260，有显著意义。Snk的检验方法最后给出4种不同包装食品的销售量的均值有显著差异，并指出了A4包装的销售量均值最高，其他三种包装具有相同的效果。单元观察数的调和均数2.4=4/(1/2+1/3+1/3+1/2)。

表25.8（b） 随机单位区组试验设计的dunnett检验

上海财经大学经济信息管理系IS/SHUFE

Page 23 of 30

The SAS System General Linear Models Procedure Dunnett's T tests for variable: Y NOTE: This tests controls the type I experimentwise error for comparisons of all treatments against a control. Alpha= 0.05 Confidence= 0.95 df= 4 MSE= 8.75 Critical Value of Dunnett's T= 3.336 Comparisons significant at the 0.05 level are indicated by '***'. Simultaneous Simultaneous Lower Difference Upper BLOCK Confidence Between Confidence Comparison Limit Means Limit 2 - 1 -4.977 2.000 8.977 3 - 1 -8.795 -0.250 8.295 Dunnett's T tests for variable: Y NOTE: This tests controls the type I experimentwise error for comparisons of all treatments against a control. Alpha= 0.05 Confidence= 0.95 df= 4 MSE= 8.75 Critical Value of Dunnett's T= 3.579 Comparisons significant at the 0.05 level are indicated by '***'. Simultaneous Simultaneous Lower Difference Upper TREAT Confidence Between Confidence Comparison Limit Means Limit A4 - A1 1.413 12.000 22.587 *** A3 - A1 -5.665 4.000 13.665 A2 - A1 -11.665 -2.000 7.665 General Linear Models Procedure Level of --------------Y-------------- BLOCK N Mean SD 表25.8（b）中结果分析：用 Dunnett双侧检验求得的临界t统计量为3.579，A4组与A1组均值之差为27－15=12>2×3.579，有显著意义，在旁边标上“***”，置信区间的下限计算为12－3.579×8.75(1/2?1/2)=1.413，上限计算为12+3.579×8.75(1/2?1/2)=22.587。而A3组与A1组之差为19－15=4<2×3.579，无显著意义，置信区间的下限计算为4－3.579×8.75(1/3?1/2)=-5.665，上限计算为4+3.579×8.75(1/3?1/2)=13.665。最后给出每个分组的均值与标准差，其中A4包装的均值为27.0000000，标准差为4.24264069。 5． 双因素实验设计的方差分析

例25.5 研究饮食和健美操对减肥的作用。我们认为饮食对减肥肯定有一定作用，适当的健美操对减肥也有效果。如果饮食加上健美操作为减肥的手段，就存在哪一种饮食配上哪一样健美操最为有效的问题，因为饮食与饮食这两种减肥手段之间存在着交互作为，会强化减肥的作用。现有三套饮食方案称为a、b、c，五种不同的健美操标记为1、2、3、4、5。构成成了3×5=15种水平组合，选择了情况基本相同的90个肥胖人进行试验，将他们随机地指派

上海财经大学经济信息管理系IS/SHUFE

Page 24 of 30

到这15个组中且每组6人。经过一段时间后，体重的下降结果见表25.9所示。

表25.9 3×5双因素设计的试验结果

饮食方案 food a b c

程序如下：

健美操train 1 22.1 24.1 19.1 22.1 25.1 18.1 13.5 14.5 11.5 6.0 27.0 18.0 19.0 22.0 20.0 14.5 19.0 16.0 2 27.1 15.1 20.6 28.6 15.1 24.6 16.9 17.4 10.4 19.4 11.9 15.4 20.0 22.0 25.5 16.5 18.0 17.5 3 22.3 25.8 22.8 28.3 21.3 18.3 15.7 10.2 16.7 19.7 18.2 12.2 16.4 14.4 21.4 19.9 10.4 21.4 4 19.8 28.3 26.8 27.3 26.8 26.8 15.1 6.5 17.1 7.6 13.6 21.1 24.5 16.0 11.0 7.5 14.5 15.5 5 20.0 17.0 24.0 22.5 28.0 22.5 21.8 22.8 18.8 21.3 16.3 14.3 11.8 14.3 21.3 6.3 7.8 13.8 data fatness;

do i=1 to 3; Input food $ ; do train=1 to 5; do j=1 to 6;

input y @@; output;

end; end; end; cards; a

22.1 24.1 19.1 22.1 25.1 18.1 27.1 15.1 20.6 28.6 15.1 24.6 22.3 25.8 22.8 28.3 21.3 18.3 19.8 28.3 26.8 27.3 26.8 26.8 20.0 17.0 24.0 22.5 28.0 22.5 b

13.5 14.5 11.5 6.0 27.0 18.0 16.9 17.4 10.4 19.4 11.9 15.4

上海财经大学经济信息管理系IS/SHUFE

Page 25 of 30

15.7 10.2 16.7 19.7 18.2 12.2 15.1 6.5 17.1 7.6 13.6 21.1 21.8 22.8 18.8 21.3 16.3 14.3 c

19.0 22.0 20.0 14.5 19.0 16.0 20.0 22.0 25.5 16.5 18.0 17.5 16.4 14.4 21.4 19.9 10.4 21.4 24.5 16.0 11.0 7.5 14.5 15.5 11.8 14.3 21.3 6.3 7.8 13.8 ;

proc glm data=fatness;

class food train;

model y=food train food*train; lsmeans food train food*train;

lsmeans food*train/slice=food slice=train;

Contrast ‘t1 vs t4 in f1’train 1 0 0 -1 0 food*train 1 0 0 -1 0 ;

Contrast ‘t2 vs t4 in f1’train 0 1 0 -1 0 food*train 0 1 0 -1 0 ;

Contrast ‘t3 vs t4 in f1’train 0 0 1 -1 0 food*train 0 0 1 -1 0 ;

Contrast ‘t4 vs t5 in f1’train 0 0 0 1 -1 food*train 0 0 0 1 -1 ;

Contrast ‘t2 vs t5 in f3’train 0 1 0 0 -1 food*train 0 0 0 0 0 0 0 0 0 0 0 1 0 0 -1 ; run;

程序说明：本例中有两个因素food和train，因此在class语句中要有这两个分组变量名。由于除了要考察这两个因素的主效应外，还要考察这两个因素的交互效应，表示为food*train，所以需要在model语句的后面加上这个交互效应。用lsmeans语句替代means语句的主要原因是，对于非均衡的试验数据需要计算最小二乘均值，它是一种调整后的均值。第二条lsmeans语句的作用，考虑到交叉项food*train是显著情况时，通过Slice选项规定的food效应和train效应来分开交叉的food*train效应。Contrast语句是作更进一步的对比，前四条Contrast语句是把因素food固定在第一个水平a上，然后对food因素有显著交互作用的train因素的某两个水平之间进行比较，最后一条contrast语句是把因素food固定在第三个水平c上，对train因素的第二个水平均值和第五个水平均值进行比较。要注意food*train交叉效应的参数化形式的规则为，交叉组合下标里最右边的变量水平比最左边的变量水平变化快，即f1*t1、f1*t2、f1*t3、f1*t4、 f1*t5、f2*t1、f2*t2、f2*t3、f2*t4、 f2*t5、f3*t1、f3*t2、f3*t3、f3*t4、 f3*t5。程序输出的主要结果见表25.10（a）和（b）所示。

表25.10（a） 双因素试验设计的方差分析表

上海财经大学经济信息管理系IS/SHUFE

Page 26 of 30

The SAS System General Linear Models Procedure Class Level Information Class Levels Values FOOD 3 a b c TRAIN 5 1 2 3 4 5 Number of observations in data set = 90 Dependent Variable: Y Source DF Sum of Squares Mean Square F Value Pr > F Model 14 1339.02488889 95.64463492 4.87 0.0001 Error 75 1473.76666667 19.65022222 Corrected Total 89 2812.79155556 R-Square C.V. Root MSE Y Mean 0.476048 24.04225 4.43285712 18.43777778 Source DF Type I SS Mean Square F Value Pr > F FOOD 2 953.15622222 476.57811111 24.25 0.0001 TRAIN 4 11.38044444 2.84511111 0.14 0.9648 表25.10（a）中结果分析：总的模型方差分析结果表明，F=4.87，p=0.0001，模型效应是显著的。模型中有三个效应：两个主效应food和train及一个交互效应food*train，其中主效应food和交互效应food*train是显著的，而主效应train，F=0.14，p=0.9648，是不显著的。所以我们可以得出的基本结论为：饮食控制和健美操对减肥是有作用的，3种不同的饮食控制方案对减肥效果是有区别的，而5种不同的健美操对减肥效果是没有区别的，同时饮食方案和健美操的不同组合对减肥效果也是有区别的。

表25.10（b） 最小二乘均值和对比分析

上海财经大学经济信息管理系IS/SHUFE

Page 27 of 30

Least Squares Means FOOD Y LSMEAN a 23.0100000 b 15.6966667 c 16.6066667 TRAIN Y LSMEAN 1 18.4222222 2 19.0000000 3 18.6333333 4 18.1000000 5 18.0333333 FOOD TRAIN Y LSMEAN a 1 21.7666667 a 2 21.8500000 a 3 23.1333333 a 4 25.9666667 a 5 22.3333333 b 1 15.0833333 b 2 15.2333333 b 3 15.4500000 b 4 13.5000000 b 5 19.2166667 c 1 18.4166667 c 2 19.9166667 c 3 17.3166667 c 4 14.8333333 c 5 12.5500000 FOOD*TRAIN Effect Sliced by FOOD for Y Sum of Mean FOOD DF Squares Square F Value Pr > F a 4 72.638667 18.159667 0.9241 0.4546 b 4 107.204667 26.801167 1.3639 0.2546 c 4 206.025333 51.506333 2.6212 0.0414 FOOD*TRAIN Effect Sliced by TRAIN for Y Sum of Mean TRAIN DF Squares Square F Value Pr > F

表25.10（b）中结果分析：由于主效应food是显著的，说明三种饮食方案对减肥的效果是不同的，再通过查看三种饮食方案减肥体重的最小二乘均值均值，可以得出a方案最好，c方案最差，且a方案和c方案的差异应该是显著的，至于a与b的比较及b与c比较，可以

上海财经大学经济信息管理系IS/SHUFE

Page 28 of 30

采用多重比较的方法进一步分析。为了知道交互效应food*train显著的具体原因，我们需要比较两因素在各种组合时的均值差异，以便寻找最好的组合方案。先作切片（slice）分析，分别固定food因素在三个水平a、b、c上，再对交互效应food*train的五种不同组合均值进行分析，其中在切片a、b上无显著性（0.4546>0.05，0.2546>0.05），而在c上有显著性（0.0414<0.05），即只有当选择c饮食方案时，选择不同的健美操才会存在减肥效果区别。对train的五种水平作切片分析，结果都是显著的（0.0383、0.0341、0.0100、0.0001、0.0010），说明无论采用哪种健美操，选择不同的饮食方案对减肥效果都存在区别。进一步的分析，我们把因素food固定在减肥效果最好的第一个水平a上，然后把train因素的每个水平与第四个水平进行比较，结果显示都是无显著性（0.1050、0.1119、0.2718、0.1599），与前面的切片分析是一致的。最后把因素food固定在第三个水平c上，对train因素的最大水平均值和最小水平均值进行比较，结果显示是显著的，同样证实了前面的切片分析。综上所述，最佳效果的减肥措施是选择a饮食方案搭配5种健美操中的任何一种都可以。 6． 拉丁方试验设计的方差分析

例25.6 研究5种不同的防护服（A、B、C、D、E）对脉搏数的影响。采用5×5拉丁方试验设计，选用5个受试者，在5个不同日期进行试验，在行、列与字母上分别安排3个因素（日期、受试者、防护服），得到试验结果数据见表25.11所示。

表25.11 5×5拉丁方试验设计的数据

日期 date 1 2 3 4 5 程序如下：

受试者person 1 A 129.8 B 144.4 C 143.0 D 133.4 E 142.8 data pulse;

do date=1 to 5; do person=1 to 5;

input cloth $ y @@; output;

end; end; cards;

A 129.8 B 116.2 C 114.8 D 104.0 E 100.6 B 144.4 C 119.2 D 113.2 E 132.8 A 115.2 C 143.0 D 118.0 E 115.8 A 123.0 B 103.8 D 133.4 E 110.8 A 114.0 B 98.0 C 110.6 E 142.8 A 110.6 B 105.8 C 120.0 D 109.8 ;

proc anova data=pulse;

class date person cloth; model y=date person cloth; run;

2 B 116.2 C 119.2 D 118.0 E 110.8 A 110.6 3 C 114.8 D 113.2 E 115.8 A 114.0 B 105.8 4 D 104.0 E 132.8 A 123.0 B 98.0 C 120.0 5 E 100.6 A 115.2 B 103.8 C 110.6 D 109.8 程序说明：对于拉丁方试验设计应该要用anova过程。实际上拉丁方试验设计是一种特

上海财经大学经济信息管理系IS/SHUFE

Page 29 of 30

殊类型的3个因素试验设计，其水平数必须相同，因此在class语句中有3个分类变量名（date、person、cloth）。在3个水平交叉的单元上只有一次试验，且不存在3个分类变量的交互效应，所以在model语句等号的右边也只有这3个分类变量名。所作的三个原假设为：①各种防护服的平均脉搏数相同；②各个受试者的平均脉搏数相同；③不同日期的平均脉搏数相同。如果欲进一步比较某个因素的任两个水平的平均脉数是否相同，可增加means或contrast语句。程序输出的主要结果见表25.12所示。

表25.12 5×5拉丁方试验设计的方差分析

The SAS System Analysis of Variance Procedure Class Level Information Class Levels Values DATE 5 1 2 3 4 5 PERSON 5 1 2 3 4 5 CLOTH 5 A B C D E Number of observations in data set = 25 Dependent Variable: Y Source DF Sum of Squares Mean Square F Value Pr > F Model 12 3579.77280000 298.31440000 6.80 0.0011 Error 12 526.14080000 43.84506667 Corrected Total 24 4105.91360000 R-Square C.V. Root MSE Y Mean 表25.12中结果分析：5种防护服对脉搏数作用的差异，统计上显示无显著意义（0.3445>0.05）。而受试者之间的差异是有统计意义的（0.0001<0.05）,说明试验的差异主要来自受试者个体的情况，如身体状况、心理状况等，而与5种不同的防护服无关，基本上也与试验在哪一天发生无关。

上海财经大学经济信息管理系IS/SHUFE

Page 30 of 30

本文来源：https://www.bwwdw.com/article/drgo.html