means irrig/dunnett alpha=0.1; /* 选项dunnett表示与对照组进行比较,缺省的对照组为以字母顺序排列的第一个水平,此处即basin,也可指定,如下例 */ means irrig/dunnett ('flood') alpha=0.1; run;
lsd: Least significant difference. comparisonwise error rate; duncan:duncan's new multiple range. error rate comparable to k-1 orthogonal comparisons tested simultaneously;
waller:waller-duncan. error rate dependent on value of analysis-of-variance F-test;此检验是基于Bayesian principle,利用
kratio?type I error来确定类似于犯第一型错误的概率. type II errortukey: tukey's honest significant difference. experimentwise error rate. lsd 犯第一错误的概率最大,tukey犯第一类错误的概率最小,duncan waller则介于二者之间.
3.6 拉丁方设计
两个响应变量的拉丁方设计(wtloss,shrink为两个响应变量) (相当于有两个区组因素)
data garments;
input run pos mat $ wtloss shrink; cards; 2 4 A 251 50 2 2 B 241 48 2 1 D 227 45 2 3 C 229 45 3 4 D 234 46 3 2 C 273 54 3 1 A 274 55 3 3 B 226 43 1 4 C 235 45 1 2 D 236 46 1 1 B 218 43 1 3 A 268 51 4 4 B 195 39 4 2 A 270 52 4 1 C 230 48 4 3 D 225 44 ;
proc print data=garments; run;
proc glm data=garments; class run pos mat;
model wtloss shrink = run pos mat; run;
注:/*在(2,4,A)这个搭配中,实验测量两个指标wtloss,shrink的值,从而看出这个拉丁方设计*/
3.7 两因素方差分析
A two-way factorial experiment(两因素方差分析)
3.7.1 图形识别交互效应
data grasses;
input method $ 1 variety 2 y1-y6 trt $ 1-2;/*列方式输入的好处*/ y1=y1/10; y2=y2/10; y3=y3/10; y4=y4/10; y5=y5/10; y6=y6/10; cards;
a1 221 241 191 221 251 181 a2 271 151 206 286 151 246 a3 223 258 228 283 213 183 a4 198 283 268 273 268 268 a5 200 170 240 225 280 225 b1 135 145 115 60 270 180 b2 169 174 104 194 119 154 b3 157 102 167 197 182 122 b4 151 65 171 76 136 211 b5 218 228 188 213 163 143 c1 190 220 200 145 190 160 c2 200 220 255 165 180 175 c3 164 144 214 199 104 214 c4 245 160 110 75 145 155 c5 118 143 213 63 78 138 ;
*proc print; *run;
data fctorial; set grasses; drop y1-y6; yield=y1; output; yield=y2; output; yield=y3; output; yield=y4; output;
yield=y5; output; yield=y6; output; run;
/*以上数据步将行变为列的方式*/ proc sort;
by method variety;/*按两因素的不同搭配排序*/ proc means data=fctorial noprint; by method variety;
output out=factmean mean=yldmean;/*按不同搭配输出均值*/ proc print data=factmean; run;
axis1 value=(font=swiss2 h=2) label=(f=swiss h=2 'yield mean'); axis2 value=(font=swiss h=2 )label=(f=swiss h=2 'Variety'); legend1 value=(font=swiss h=2 ) label=(f=swiss h=2 'Method'); symbol1 color=black interpol=join line=1 value='A' font=swiss; symbol2 color=black interpol=join line=2 value='B' font=swiss; symbol3 color=black interpol=join line=20 value='C' font=swiss; /*这些均为全程语句,不在数据步中,也不在过程步中*/ proc gplot data=factmean;
plot yldmean*variety=method/caxis=black ctext=black vaxis=axis1 haxis=axis2 legend=legend1; run;
3.7.2 带交互作用的书写,每个搭配均值的计算,均值的多重比较
proc glm data=fctorial; class method variety;
model yield=method|variety;/*交互效应的写法一*/
means method variety method*variety;/*不同因素的各水平,及交互效应的均值*/ lsmeans method*variety/slice=(variety method);/*计算cell的均值,选项slice规定分组变量*/
lsmeans method*variety/cl pdiff adjust=tukey;/*均值的多重比较,cl要求计算置信限,pdiff要求计算均值差的概率,adjust规定计算概率时用的方法,缺少为lsd*/ run;
3.7.3 用混合效应mixed过程更自然
proc mixed data=fctorial; class variety method;
model yield=method variety method*variety; lsmeans method*variety/diff; ods output diffs=cld; run;
data smpleff; set cld;
if variety=_variety;
proc print data=smpleff;
var variety _variety method _method estimate stderr df tvalue probt; run; run;
3.7.4 自定义的线性检验和估计
write estimate and contrast statement as following: 1) Write the linear combination you want to test or estimate in terms of means; 2) Convert means into model parameters; 3) Gather like terms; 例子见书p83-84,89;
proc glm data=fctorial; class method variety;
model yield=method variety method*variety / ss1; estimate 'A vs B,C in V1' method 1 -.5 -.5
method*variety 1 0 0 0 0 -.5 0 0 0 0 -.5 0 0 0 0; estimate 'A vs B,C in V2' method 1 -.5 -.5
method*variety 0 1 0 0 0 0 -.5 0 0 0 0 -.5 0 0 0; estimate 'A vs B,C in V3' method 1 -.5 -.5
method*variety 0 0 1 0 0 0 0 -.5 0 0 0 0 -.5 0 0; estimate 'A vs B,C in V4' method 1 -.5 -.5
method*variety 0 0 0 1 0 0 0 0 -.5 0 0 0 0 -.5 0; estimate 'A vs B,C in V5' method 1 -.5 -.5
method*variety 0 0 0 0 1 0 0 0 0 -.5 0 0 0 0 -.5; estimate 'A vs B,C Overall' method 1 -.5 -.5
method*variety .2 .2 .2 .2 .2 -.1 -.1 -.1 -.1 -.1 -.1 -.1 -.1 -.1 -.1; /*估计0.2[?A1?0.5(?B1??C1)]?...?0.2[?A5?0.5(?B5??C5)]*/ contrast 'A vs BC * Varieties'
method * variety 1 0 0 0 -1 -.5 0 0 0 .5 -.5 0 0 0 .5, method * variety 0 1 0 0 -1 0 -.5 0 0 .5 0 -.5 0 0 .5, method * variety 0 0 1 0 -1 0 0 -.5 0 .5 0 0 -.5 0 .5, method * variety 0 0 0 1 -1 0 0 0 -.5 .5 0 0 0 -.5 .5; /*检验[?A1?0.5(?B1??C1)]?...?[?A5?0.5(?B5??C5)]；即
[?A1?0.5(?B1??C1)]?[?A5?0.5(?B5??C5)]和[?A2?0.5(?B2??C2)]?[?A5?0.5(?B5??C5)]和 [?A3?0.5(?B3??C3)]?[?A5?0.5(?B5??C5)]和
[?A4?0.5(?B4??C4)]?[?A5?0.5(?B5??C5)]四个等式同时成立*/
contrast 'A vs B,C * V1,V2'
method*variety 1 -1 0 0 0 -.5 .5 0 0 0 -.5 .5 0 0 0; /*两组内的比较*/
contrast 'A vs B,C * V3,V4'
method*variety 0 0 1 -1 0 0 0 -.5 .5 0 0 0 -.5 .5 0; estimate 'A vs B,C in V1,V2' method 1 -.5 -.5
method*variety .5 .5 0 0 0 -.25 -.25 0 0 0 -.25 -.25 0 0 0; /*两组合并*/ run;
options ls=78;
3.8 用回归方式进行方差分析
可将其看成是一般线性回归模型,只不过要作适当的量化,即将分类变量进行量化,设置所谓的哑变量,可用第一种方法,这可通过调用SAS/GLMMOD过程来完成,其调用格式与GLM过程完全一样,下面给出一个例子:
下面是单因素方差分析用常规回归模型做的例子（见王松桂—线性统计模型，p144）
data a;
input type x1 x2 x3 y;