群馬大学 | 医学部 | サイトトップ | 医学情報処理演習
2009年11月2日
（注）前回の課題の解答例は，http://phi.med.gunma-u.ac.jp/medstat/it2009-03r.htmlに示してあります。
下のコードは，10000人（男女5000人ずつ）の母集団データをファイルsample04.txtからデータフレームdatに読み込み，身長を意味する変数HTについて男女別に（変数SEXごとに）サイズ10とサイズ100の標本（サンプル）を抽出する数値実験である。この画面で範囲選択コピーしてから，R Consoleのプロンプトにペーストすれば実行できる。
なお，最初の2行は，ランダムな標本抽出を行うために，まずRNGkind("Mersenne-Twister")によって乱数生成アルゴリズムをメルセンヌツイスター法に指定し（注：乱数については，このページを参照されたい），set.seed(54321)によって乱数を初期化しているものである。このことによって，解答者全員がまったく同じサンプリングを行うことができる。
RNGkind("Mersenne-Twister") set.seed(54321) dat <- read.delim("http://phi.med.gunma-u.ac.jp/medstat/sample04.txt") malesHT <- dat$HT[dat$SEX=="M"] femalesHT <- dat$HT[dat$SEX=="F"] s10M <- sample(malesHT,10) s100M <- sample(malesHT,100) s10F <- sample(femalesHT,10) s100F <- sample(femalesHT,100)この数値実験で得られたデータに基づいて，以下の内容を実行せよ。1.と2.と3.については解答フォームの上の欄にRのコードを貼り付け，4.については解答フォームの下の欄に文章で考察せよ。学籍番号と氏名を忘れないこと。
- 男女それぞれについて，母集団での身長（malesHT及びfemalesHT）の平均値を求めよ。
- 男女それぞれについて，サイズ10の標本（s10Mとs10F）から推定される母集団の身長の平均値と95％信頼区間を求めよ。
- 男女それぞれについて，サイズ100の標本（s100Mとs100F）から推定される母集団の身長の平均値と95％信頼区間を求めよ。
- 母集団での身長の真の平均値と，サイズ10の標本から推定された平均値及び95％信頼区間とサイズ100の標本から推定された平均値及び95％信頼区間を比べて考察せよ。
ごく単純に考えるならば，次のようなコードと実行結果が得られる。
mean(malesHT) mean(femalesHT)
[1] 169.9103 [1] 157.9724
t.test(s10M) t.test(s10F)
95 percent confidence interval: 168.1843 172.9157 sample estimates: mean of x 170.55 95 percent confidence interval: 153.3266 159.6934 sample estimates: mean of x 156.51
t.test(s100M) t.test(s100F)
95 percent confidence interval: 169.0315 170.8685 sample estimates: mean of x 169.95 95 percent confidence interval: 156.4066 158.6314 sample estimates: mean of x 157.519
4. 上記の結果からまずわかることは，サンプルサイズ100の標本から推定される平均値の方が，サンプルサイズ10の標本から推定される平均値よりも真の母平均に近いことである。また，95％信頼区間は，男女とも，サイズ100の標本からの計算値の方がサイズ10の標本からの計算値より幅が狭いので，平均値の信頼性がより高いと考えられる。
結果を見やすく出力したければ，以下のコードのように，いったん計算結果を変数に付値しておき，それを行列の形に加工してから表示することもできる。
getmci <- function(Z) { as.numeric(c(Z$estimate,Z$conf.int)) }
mM <- mean(malesHT)
mF <- mean(femalesHT)
t10M <- getmci(t.test(s10M))
t100M <- getmci(t.test(s100M))
t10F <- getmci(t.test(s10F))
t100F <- getmci(t.test(s100F))
res <- rbind(cbind(mM,mF),cbind(t10M,t10F),cbind(t100M,t100F))
colnames(res) <- c("男性","女性")
rownames(res) <- c("母平均","サイズ10標本平均","95%信頼区間下限","95%信頼区間上限","サイズ100標本平均","95%信頼区間下限","95%信頼区間上限")
res計算結果は以下のように表示される。もちろん考察は変わらない。
男性 女性 母平均 169.9103 157.9724 サイズ10標本平均 170.5500 156.5100 95%信頼区間下限 168.1843 153.3266 95%信頼区間上限 172.9157 159.6934 サイズ100標本平均 169.9500 157.5190 95%信頼区間下限 169.0315 156.4066 95%信頼区間上限 170.8685 158.6314
サイズ10の標本から推定される95％信頼区間と同程度の幅は，サイズ100の標本では何％信頼区間になるかを計算した解答や，サイズ4000の標本では95％信頼区間の幅がより狭くなることを示した解答もあった。そのように自ら発展的な探求をすると，より早く深く統計解析が身につくと思う。