このページは、令和２年４月２５日に一部更新しました。
| 1.6.4.1 SAS による単純な乱塊２要因デザイン ANOVA プログラ ムとその概要 |
| 1.6.4.2 SAS による単純な乱塊２要因デザイン ANOVA プログラ ムの出力結果の概要 |
さらに、この節にはつぎの SAS プログラムのダウンロードコーナーを用意 してある：
| 単純な乱塊２要因デザインのプログラム例 |
この節では、1.6.1 節で述べた乱塊２要因デザインすなわち RBF-pq デザインの 原型とでも言うべき反復測度でない乱塊２要因デザインデータのための SAS プロ グラムを示す。ただし、1.5.4 節の RB-p デザインの場合と同様、反復測度でない 乱塊２要因デザインデータの良い例が見あたらないので、ここでも 1.5.4 節と 同様、RBF-pq デザインの原型データではないデータを、そうみなしてそのための SAS プログラムを示す。
データは、Huynh and Mandeville (1979) の２要因反復測度デザインデータ である。このデータは、1.6 節の最初で紹介した 大局的球形仮定が満たされている データであり、したがって 加算的モデルによる分析 が必要なデータなのである。 しかし、１要因の場合に既に示した glm プロシジャの repeated 文 を使う SAS による反復測度 ANOVA 分析では、実は非加算的モデル による分析しかできない。つまり、repeated 文を使うと、加算的モデルが妥当な データに対しても、Table 1.15 におけるプールされた F-比による検定ではなく、 一律に非加算的モデルにおける (1.181) 式による主効果や交互作用 の検定を行うのでまずい。
この意味では、Huynh らのデータは、もともとつぎのようにして RBF-pq デザ インの原型としての反復測度でない乱塊２要因デザインデータとして、なおかつ 加算的モデルとして分析する必要があるデータなのである。
プログラムは、以下に示すように：
の６項から成る。単に、glm プロシジャによる乱塊２要因デザインデータの加算的 モデルによる分散分析のみを行うには、これらの内の (1) 及び (2-1) だけを実行 すればよいが、平均値のグラフ出力やモデルの誤差項の正規性検定まで行うには一連 のプログラムを実行するとよい。
*---------------------------------------------------------------------* | March 4, 1999 | | sas program--rbf3-3hm.sas | | An example of sas programs for analyzing RBF-pq design data by | | an additive model. The data is that of Table 2. by Huynh and Man- | | deville (1979). | | | *---------------------------------------------------------------------*; |
この項は、以下のプログラムに関する冒頭でのコメントを記述している。
options ps=60 ls=80;
/* (1) data input with a transformation */
data HuMa1979;
input no 2. (b1_c1 b1_c2 b1_c3 b2_c1 b2_c2 b2_c3
b3_c1 b3_c2 b3_c3) (3.);
label no='sample number'
b1_c1='combination of b1 and c1'
b1_c2='combination of b1 and c2'
b1_c3='combination of b1 and c3'
b2_c1='combination of b2 and c1'
b2_c2='combination of b2 and c2'
b2_c3='combination of b2 and c3'
b3_c1='combination of b3 and c1'
b3_c2='combination of b3 and c2'
b3_c3='combination of b3 and c3';
array rt(3,3) b1_c1 b1_c2 b1_c3 b2_c1 b2_c2 b2_c3
b3_c1 b3_c2 b3_c3;
block+1;
do factB=1 to 3;
do factC=1 to 3;
rtime=rt(factB,factC);
output;
end;
end;
cards;
1 53 20 12 14 42 30 10 5 63
2 23 55 77 10 2 30 56 30 50
3 20 30 50 43 12 30 53 21 20
4 3 20 77 45 53 32 65 30 20
5 23 22 21 12 32 30 3 54 33
6 33 89 53 65 45 42 2 10 23
7 30 33 55 42 87 30 2 10 30
8 36 56 32 3 65 86 54 23 30
9 23 3 78 63 68 68 54 12 39
10 98 65 63 32 45 75 86 63 21
11 53 65 86 96 63 32 12 45 25
12 33 22 42 21 3 35 63 32 54
13 10 30 53 65 43 20 32 45 65
14 32 35 63 65 66 33 53 63 32
15 12 30 56 22 30 56 42 30 12
16 56 33 89 65 78 99 63 30 24
17 53 30 36 65 22 33 22 54 33
18 32 30 36 65 63 32 30 36 65
19 30 36 65 66 33 22 12 30 32
20 98 65 63 65 45 12 2 36 30
21 00 22 11 42 53 63 32 53 32
22 30 35 63 66 33 22 56 52 21
;
|
この項では、cards 文の後に用意してある Huynh and Mandeville (1979) の２要
因反復測定デザインデータを、まず input 文で指示された形で読み込ませ、それら
を単純な乱塊２要因デザインデータとみなし、さらに加算的モデルで分析するため
に array 文や do 文を用いて組み直すためのプログラムを示している。
もちろん、ここで紹介している分散分析は、原則はあくまでも（独立測定要
因の場合の）単純な乱塊２要因デザインデータの分散分析であるので、この例のよ
うな反復測定要因データは特別なケースではあるが、この例のように反復測定
要因データであっても大局的球形仮説が満たされているようなデータの場合には、
むしろここでのやり方でデータを加算的モデルで分析するために変形するのがよい。
/* (2-1) repeated measures RBF-3.3 design data by an additive model */ title 'additive model analysis for a repeated measure RBF-IJ design'; proc glm data=HuMa1979; class block factB factC; model rtime=block factB factC factB*factC; output out=temp residual=resid; run; |
この項では、通常の乱塊２要因デザインデータの glm プロシジャによる加算的
モデルによる分散分析を行わせ、さらにモデルの残差の推定値を resid と名づけた
変数名で temp なる一時ファイルに保存するためのプログラムを示す。
この例のように、反復測定乱塊２要因デザインデータの場合でも、大局的球形
仮説が満たされている場合には、これを用いるのが適切であろう。
/* (2-2) compute means of repeated measures */
title 'plot the means for each cell';
proc summary data=HuMa1979 print n mean std;
var b1_c1 b1_c2 b1_c3 b2_c1 b2_c2 b2_c3 b3_c1 b3_c2 b3_c3;
output out=new mean=b1c1 b2c1 b3c1 b1c2 b2c2 b3c2 b1c3 b2c3 b3c3;
run;
data xyplot;
set new;
array av(3,3) b1c1 b2c1 b3c1 b1c2 b2c2 b3c2 b1c3 b2c3 b3c3;
do c=1 to 3;
do b=1 to 3;
aver=av(c,b);
output;
end;
end;
drop b1c1 b2c1 b3c1 b1c2 b2c2 b3c2 b1c3 b2c3 b3c3;
run;
options ps=40 ls=80;
proc plot data=xyplot;
plot aver*b=c;
run;
|
/* (3-1) posterior analysis (1): histogram plot of error estimates */ title 'histogram of the error estimates '; proc chart data=temp; vbar resid; run; |
/* (3-2) posterior analysis (2): test for normality of error estimates */ options ps=60 ls=80; title 'test for global normality of the dependent variable'; proc univariate data=temp normal; var resid; run; |
| rbf3-3hm.sas |
うえのプログラムを実行すると、以下のような出力結果が得られる：
additive model analysis for a repeated measure RBF-IJ design 8
General Linear Models Procedure
Dependent Variable: RTIME
Sum of Mean
Source DF Squares Square F Value Pr > F
Model 29 23807.772727 820.957680 1.82 0.0104
Error 168 75818.090909 451.298160
Corrected Total 197 99625.863636
R-Square C.V. Root MSE RTIME Mean
0.238972 52.33630 21.243779 40.590909
Source DF Type I SS Mean Square F Value Pr > F
BLOCK 21 16005.863636 762.183983 1.69 0.0367
FACTB 2 3042.939394 1521.469697 3.37 0.0367
FACTC 2 656.121212 328.060606 0.73 0.4849
FACTB*FACTC 4 4102.848485 1025.712121 2.27 0.0635
Source DF Type III SS Mean Square F Value Pr > F
BLOCK 21 16005.863636 762.183983 1.69 0.0367
FACTB 2 3042.939394 1521.469697 3.37 0.0367
FACTC 2 656.121212 328.060606 0.73 0.4849
FACTB*FACTC 4 4102.848485 1025.712121 2.27 0.0635
|
この項では、プログラムで指定した加算的モデルによる乱塊２要因分散分
析の出力結果を示す。1.6.1 節の表 1.15 から明らかなように、このデザインに
おける誤差項は加算的モデルの場合、どの効果の検定の場合にも通常の単一誤差項
SSE を選ぶ必要があるが、上記のプログラムを用いると確かにその方式
で誤差項が選ばれていることが上の結果の簡単なチェックによりわかる。もち
ろん、出力結果のうち Corrected Total で Sum of Squares の項の 99625.863636
がそれである。
このデータは、釣り合い型と見れるので、SAS の平方和はすべて等しく、ここでも
Type I 平方和も Type III 平方和の場合も検定結果は変わらない。これらを見ると、
このデータでは、ブロック因子の効果及び最初の要因 (FACTB) の主効果のみが５パ
ーセント水準で統計的に有意であるが２つ目の要因（FACTC) の主効果は有意でない
ことや、（全体的）交互作用は傾向程度であることがわかる。
もっとも、この節の最初に指摘したように、このデータはもともと独立測定要因
デザインデータではなく反復測定デザインデータであり、かつ大局的球形仮説が満
たされているケースであるが、SAS で repeated 文を用いて同一データを分析
すると常に非加算的モデルによる分析を行っているので、このような場合、
repeated 文を用いた分析結果は、必ずしも上の加算的モデルによる結果とは異なっ
たものになる可能性があろう。
実際、このデータで repeated 文による分析を行うと、ε修正を
しない場合で見ても（この場合、加算的モデルなので修正しない方がよい）、検定
結果はうえの結果と一部異なり、最初の要因の主効果は５パーセント水準では有意
でなくなる。この違いは、このデータでは加算的モデルの方がより適切なのに、
SAS では反復測定要因データの場合、repeated 文を用いると常に非加算的モデルを
選んでしまう、すなわちモデルの誤差項に SSEを用いないことに起因
する。
もっとも、通常の独立測定要因の場合の乱塊２要因デザインの場合には、反復測定
要因の場合のような（大局的球形仮説の検定による）加算性の有無の検討はできない
ので、何らかの他の方法なり過去の知見なりから加算性の適否の検討が必要になろ
う。
plot the means for each cell 9
Variable Label N Mean Std Dev
-------------------------------------------------------------------
B1_C1 combination of b1 and c1 198 35.5000000 24.8413749
B1_C2 combination of b1 and c2 198 37.5454545 19.6240155
B1_C3 combination of b1 and c3 198 53.6818182 21.7205316
B2_C1 combination of b2 and c1 198 46.9090909 24.2403745
B2_C2 combination of b2 and c2 198 44.6818182 22.3385300
B2_C3 combination of b2 and c3 198 41.4545455 22.4127288
B3_C1 combination of b3 and c1 198 36.5454545 24.8376547
B3_C2 combination of b3 and c2 198 34.7272727 16.9063219
B3_C3 combination of b3 and c3 198 34.2727273 15.1445275
-------------------------------------------------------------------
plot the means for each cell 10
プロット : AVER*B. 使用するプロット文字 : C の値．
AVER |
|
55 +
| 1
|
|
|
50 +
|
|
|2
|
45 + 2
|
|
|
| 2
40 +
|
| 1
|3
|1
35 + 3
| 3
|
|
|
30 +
|
-+-----------------------------------+-----------------------------------+-
1 2 3
B
histogram of the error estimates 11
|
histogram of the error estimates 11
Frequency
| *****
| *****
| ***** *****
40 + ***** *****
| ***** *****
| ***** ***** ***** *****
| ***** ***** ***** *****
| ***** ***** ***** *****
30 + ***** ***** ***** *****
| ***** ***** ***** *****
| ***** ***** ***** *****
| ***** ***** ***** *****
| ***** ***** ***** *****
20 + ***** ***** ***** *****
| ***** ***** ***** *****
| ***** ***** ***** *****
| ***** ***** ***** ***** ***** *****
| ***** ***** ***** ***** ***** *****
10 + ***** ***** ***** ***** ***** *****
| ***** ***** ***** ***** ***** *****
| ***** ***** ***** ***** ***** ***** ***** *****
| ***** ***** ***** ***** ***** ***** ***** *****
| ***** ***** ***** ***** ***** ***** ***** ***** *****
----------------------------------------------------------------------------
-42 -30 -18 -6 6 18 30 42 54
RESID 中間点
|
この項では、プログラム (1-1) の項の glm プロシジャであらかじめ計算し一時 ファイルに出力しておいたモデルの誤差項の推定値を chart プロシジャによりヒスト グラムとして書かせた結果を示す。
test for global normality of the dependent variable 12
Univariate Procedure
Variable=RESID
Moments
N 198 Sum Wgts 198
Mean 0 Sum 0
Std Dev 19.61794 Variance 384.8634
Skewness 0.175913 Kurtosis -0.10674
USS 75818.09 CSS 75818.09
CV . Std Mean 1.394186
T:Mean=0 0 Pr>|T| 1.0000
Num ^= 0 198 Num > 0 97
M(Sign) -2 Pr>=|M| 0.8312
Sgn Rank -167 Pr>=|S| 0.8367
W:Normal 0.981173 Pr < W 0.4032
Quantiles(Def=5)
100% Max 56.86869 99% 51.82323
75% Q3 14.24242 95% 33.15657
.........................
|
この項では、モデルの誤差項の正規性を univariate プロシジャで検定した結果 を示す。もっとも、このプロシジャでは上に見るように、各種の統計量も同時に 出力されるが、中程の W:NOrmal の項の右端の p-値を見ればよい。これが、0.05 以下ならば分布の正規性が棄却されるが、この場合誤差項は正規分布に従うと見な してよい。
