Comments
Description
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