(引き続きハンバーガー統計学にようこそ!の第2章をRでやってみたものです。)
2.1 平均的ポテトを推定する
前回も説明したとおり、標準のvar()
、sd()
でそれぞれ不偏分散、標本標準偏差が計算される。
wakupote <- c(47, 51, 49, 50, 49, 46, 51, 48, 52, 49) mean(wakupote); var(wakupote); sd(wakupote)
2.3 信頼区間/区間推定
テキスト通りにやるならこのように。
# t値を求める qt(0.025, 9) # 自由度9の場合の95%信頼区間の下側t値 qはquantileから qt(0.975, 9) # 同上側t値 qt(0.005, 9) # 自由度9の場合の99%信頼区間の下側t値 qt(0.995, 9) # 同上側t値 # 信頼区間を求める se <- function(x){ sd(x) / sqrt(length(x)) } # 標本標準誤差を求める関数 # 95%信頼区間 mean(wakupote) + qt(0.025, 9) * se(wakupote) # 下側 mean(wakupote) + qt(0.975, 9) * se(wakupote) # 上側 # 99%信頼区間 mean(wakupote) + qt(0.005, 9) * se(wakupote) # 下側 mean(wakupote) + qt(0.995, 9) * se(wakupote) # 上側
もちろん、Rを使っているのであればt.test()
に突っ込むだけという姿勢が正しい。必要であれば信頼区間だけを取り出すこともできる。
> t.test(wakupote, conf.level = 0.95)$conf.int [1] 47.85957 50.54043 attr(,"conf.level") [1] 0.95 > t.test(wakupote, conf.level = 0.99)$conf.int [1] 47.27432 51.12568 attr(,"conf.level") [1] 0.99
t値が計算できるのだから実用的な意味は特に無いがt分布表も作ってみよう。
t_table <- data.frame( df = (df <- c(1:30, 40, 60, 120, Inf)), p95 = qt(0.975, df), p99 = qt(0.995, df) ) knitr::kable(t_table)
df | p95 | p99 |
---|---|---|
1 | 12.706205 | 63.656741 |
2 | 4.302653 | 9.924843 |
3 | 3.182446 | 5.840909 |
4 | 2.776445 | 4.604095 |
5 | 2.570582 | 4.032143 |
6 | 2.446912 | 3.707428 |
7 | 2.364624 | 3.499483 |
8 | 2.306004 | 3.355387 |
9 | 2.262157 | 3.249836 |
10 | 2.228139 | 3.169273 |
11 | 2.200985 | 3.105806 |
12 | 2.178813 | 3.054540 |
13 | 2.160369 | 3.012276 |
14 | 2.144787 | 2.976843 |
15 | 2.131449 | 2.946713 |
16 | 2.119905 | 2.920782 |
17 | 2.109816 | 2.898230 |
18 | 2.100922 | 2.878440 |
19 | 2.093024 | 2.860935 |
20 | 2.085963 | 2.845340 |
21 | 2.079614 | 2.831360 |
22 | 2.073873 | 2.818756 |
23 | 2.068658 | 2.807336 |
24 | 2.063899 | 2.796940 |
25 | 2.059539 | 2.787436 |
26 | 2.055529 | 2.778714 |
27 | 2.051830 | 2.770683 |
28 | 2.048407 | 2.763263 |
29 | 2.045230 | 2.756386 |
30 | 2.042273 | 2.749996 |
40 | 2.021075 | 2.704459 |
60 | 2.000298 | 2.660283 |
120 | 1.979930 | 2.617421 |
Inf | 1.959964 | 2.575829 |
2.4 平均的チキンを推定する
🐔
chicken <- c(568, 530, 581, 554, 536, 518, 564, 552) mean(chicken) var(chicken) sd(chicken) se(chicken) t.test(chicken, conf.level = 0.95)$conf.int t.test(chicken, conf.level = 0.99)$conf.int
2.9 通過テスト
- (1)
- 母集団
- 標本
- 無作為抽出
- 不明である
dについて実際どのような解答が想定されているか分からないが、問題文の情報量だけで度数分布の形状を推定することは不可能である。例えば、テストの総問題数が1問であった場合を考えてみれば想像できるだろう。実際にはどのような分布も想定すべきでないという姿勢で分析にかかるべきではないのか。
- (2)
もっとうまいやり方があるような気がして仕方がない。
# 95%信頼区間 65 + qt(0.025, 500) * sqrt(60)/sqrt(500) 65 + qt(0.975, 500) * sqrt(60)/sqrt(500) # 99%信頼区間 65 + qt(0.005, 500) * sqrt(60)/sqrt(500) 65 + qt(0.995, 500) * sqrt(60)/sqrt(500)
- (3)
信頼区間の計算式においては確率pのみが変数であるので、逆に計算された信頼区間を変数として確率pが計算されるような関数を考えることができる。このとき、pの範囲が2.5%〜97.5%となるような区間が95%信頼区間、同様に0.5%〜99.5%となるような区間が99%信頼区間である。
「信頼区間に母平均が含まれる確率が95%または99%である」という説明は信頼区間の説明としては正しくない(例えばRで楽しむ統計 (Wonderful R 1)に詳しい解説がある。上述の方法とは異なる信頼区間の構成方法も紹介されている。)。