...

問題 - So-net

by user

on
Category: Documents
12

views

Report

Comments

Transcript

問題 - So-net
【問題】t 検定・分散分析とデザイン行列
H22 年度 BioS 継続勉強会:第1回
土居正明
【はじめに】
本勉強会の目的と注意
本勉強会では【演習問題】として「自分でデータセットを作り、解析せよ」という問題を毎回出題する。本勉強会の目的
の1つは、
「勉強会で扱わないテーマ・プロシジャ・解析に対しても、自分で Help を調べたりテストデータやシミュレー
ションを用いたプログラミングを行い、確認できるようになる」ことである。言いかえれば、「解答が丁寧に書かれていな
い状況で、自力で正解を探り当てる」能力をつけることである。
そのための訓練なので、【演習問題】は必ず全問解答してから次の問題へ進むこと。解答がついていないことに意味があ
るので、「解答がないからやらない」という風には絶対にしないこと。
一方【補足問題】は補足なので、時間に余裕がある場合だけでよい。
記号
本稿では、自由度 2 の t 分布の下側 97.5 %点を t(0.975, 2) と書く。また、自由度 4 の χ2 分布の下側 2.5 %点は
χ2 (0.975, 4) と表す。ただし、SAS の Help では t0.975,2 や χ20.975,4 のように記載されているので、Help からの引用はその
まま書くこととする*1 。
基礎知識:パーセント点の求め方
SAS で、標準正規分布の 下側97.5 %、つまり z(0.975) を求めたい場合、データセットで
data d1;
z=quantile(’norm’, 0.975);
put z;
run;
と入力すればよい。
【準備の問題】
標準正規分布でない分布・他の分布では、上のプログラムをどのように変更すればよいか調べよ。また、quantile 関数の
親戚として pdf 関数と cdf 関数がある。これらは何を計算する関数か調べ、t 分布について何でもよいのでこれらを使って
求めよ。
*1
確率や自由度を全て添え字にするのは、字が小さくなるため個人的な好みに合わないので、引用以外の部分では用いない。
1
問題1:t 検定と一元配置分散分析
以下のデータに対して、以下の各問いに答えよ。なお、データは全て独立に正規分布に従い、群ごとに平均が一定である
こと、分散は両群で等しいことは前提とする。また、データの数値は小さい方が「薬効がある」とする。
プラセボ群
1, 11, 15
実薬群
2, 3, 10
(0) 投与群ごとにデータをプロットせよ。
(i) 投与群間の平均値を比較する 両側 検定を行いたい。適切に文字をおいて、帰無仮説 H0 ・対立仮説 H1 を述べよ。
(ii) (i) の検定を、有意水準 0.05 で SAS の proc ttest を用いて行え。
(iii) proc ttest が平均の信頼区間をどのような数式で求めているか、Help で調べて述べよ。
(iv) (ii) で得られた出力のうち、
「標準偏差の信頼限界」と、t 検定の「Satterthwaite」の欄と等分散性の検定の欄 以外 の全
てを手計算*2 で求めよ。
【演習問題(必須)
】
以上を参考にして自分で同様のデータセットを1つ作り、手計算で上の (ii) と同様の計算をせよ。その後、proc ttest で
計算を行い、結果が一致するかどうか確認せよ。
なお、作成するデータは両群の例数が異なるようにせよ。
【補足問題(任意)
】
(v) Help で標準偏差の信頼区間の上限・下限をどのように計算しているか調べよ。
(vi) 今のデータを用いて、(v) で調べた信頼区間を手計算で求めよ。
(vii) Help で Satterthwaite と等分散の検定をどのように計算しているか調べよ。
(viii) 今のデータに対して、(vii) の各値を手計算で求めよ。
問題2:一元配置分散分析(2 群)
以下のデータは「問題1」と全く同じデータである。以下の各問いに答えよ。
プラセボ群
1, 11, 15
実薬群
2, 3, 10
2−1:SAS での実行と t 検定との比較
(i) 投与群のみを要因とした 1 元配置分散分析で、群間の平均値が等しいかどうかをみる 両側 検定を行いたい。適切に文字
をおいて、帰無仮説 H0 ・対立仮説 H1 を述べよ。
(ii) (i) の検定を、有意水準 0.05 で SAS の proc glm を用いて行え。なお、モデルステートメントで/ ss3 を入れて、
Type I の出力をせず、Type III だけの出力をせよ。
*2
SAS のデータステップにある関数は使ってよいが、その場合は式を全てノートに書いた上で関数を使うこと(たとえばプラセボ群の平均の場合、
「(1 + 11 + 15)/3 =」と書いた上で SAS を使うようにし、
「平均=」とは書かないこと)
。プロシジャの使用は不可とする。以下「手計算」とある
場合は同様。
2
(iii) (ii) の結果を、「問題1」の検定の結果と比較せよ。特に、(a) p 値、 (b)「問題1」の t 統計量と本問の F 統計量の
関係、(c)「問題1」のプールした分散の値については必ず述べよ。また、(d) 表中の F 値から p 値を算出せよ。
2−2:デザイン行列と最小二乗推定値
(iv) (ii) のプログラムの model ステートメントに / solution を付け加え、各パラメータの推定値を求めよ。
また、「問題1」で得られた各群の平均値を、この推定値から求めよ。
(v) (iv) の結果を参考にしつつ、今回の解析の線形モデルを述べよ。また、proc glm の用いるデザイン行列 X を
(制約条件なしの形で)述べよ。
(vi) (v) に、(iv) の出力と合致するような制約条件を考え、制約条件を考慮したデザイン行列 X1 を求めよ*3 。
(vii) (vi) のデザイン行列 X1 に対して、
(X01 X1 )−1 X01 y
を求めよ。結果の数値が、(ii) の出力のどこに出てくるかを調べよ。
(viii) 制約条件を入れる前のデザイン行列 X に (iii) とは別の制約条件を入れ、その際できるデザイン行列 X2 を求めよ。
(ix) (viii) の X2 に対して
(X02 X2 )−1 X02 y
を求め、結果(推定値)が (vii) と異なることを確認せよ。この結果をもとに、(iv) の推定値の出力中の’B’ と表の下
に出る NOTE の意味を考えよ。
(x) (ix) のデザイン行列 X2 をもとにした分散分析を、proc glm で行え。なお、データセットを変更してよいものとする。
この出力推定値(solution の出力)と、(iv) の出力の推定値を比べよ。
2−3:予測値と残差
両群プールした全体の平均の推定値を µ
b とする。また、群ごとの平均の推定値を、プラセボ群 µ
b1 、実薬群 µ
b2 とする。こ
こで、




ȳ = 



µ
b
µ
b
µ
b
µ
b
µ
b
µ
b








 , ŷ = 






µ
b1
µ
b1
µ
b1
µ
b2
µ
b2
µ
b2








とおく。このとき、以下の問いに答えよ。
(xi) (iv), (x) のプログラムの model ステートメントに / p を付け加え、予測値と残差を出力せよ。この場合、
予測値・残差とはそれぞれ何か?また、(iv) と (x) の予測値・残差を比較せよ。
(xii) model ステートメントの下に、アウトプットステートメント”output out=out1 p=pred r=resid;”と書くことで、
予測値と残差をデータセットに出力せよ*4 。さらに、残差を群ごとにプロットせよ。
(xiii) ベクトル x の大きさを ||x|| と書くこととする。今、
||y − ȳ||2 = ||y − ŷ||2 + ||ŷ − ȳ||2
が成り立つことを確認せよ。また、||y − ȳ||2 , ||y − ŷ||2 , ||ŷ − ȳ||2 がそれぞれ分散分析表のどこに出てくるかを
確認せよ。
*3
厳密には、proc glm は制約条件を考えているわけではなく「一般化逆行列」というものを利用している。しかし、本質的には制約条件を考えるの
と全く同じことである。
*4 pred は予測値”predicted value”の最初の部分。resid は残差”residual”の最初の部分。
3
【演習問題(必須)
】
以上を参考にして自分で同様のデータセットを1つ作り、上の問題と同様に手計算と SAS の結果を照合せよ。
4
Fly UP