医学情報処理演習第４回「標本統計量と母数推定」
2005 年 10 月 31 日 中澤 港¶
前回の課題の回答例³

演習で説明した関数を，下枠内のように X に適用するだけで，その下の枠内の結果が得られる。なお，\はバックスラッシュと読むが，日本語フォントで言えば￥に相当する。従って，キーボードでは “￥ｎ”と入力すれば，\n が入力できる。\n は改行を意味する。

¶³

data(ChickWeight)
attach(ChickWeight)
X <- weight[Time==20]
cat("平均=",mean(X)," 不偏標準偏差=",sd(X)," 中央値=",median(X)," 四分位偏差=",(fivenum(X)[4]-fivenum(X)[2])/2,"\n")detach(ChickWeight)µ¶
平均= 209.7174´³
不偏標準偏差= 66.51171

中央値= 204µµ

四分位偏差= 49.5´´

標本統計量と母数
我々が統計解析で相手にするのは標本データの場合が多く，得られる統計量も標本統計量であるが，本当に知りたいのは母集団の統計量（母数）であるため，標本統計量から母数を推定しなければならない。その仕組みを考えるため，まず標本抽出をシミュレートしてみる。

標本抽出
100 万人の成人男性の身長からなる母集団が，表１のようになっていたとしよう∗1 。下のヒストグラムはその図示である。

∗1 こんなにきれいな正規分布に近い表が現実のデータとして存在するはずがなく，実は，これは R の正規乱数を使って得た人工データである。母平均 168，母標準偏差 5 の正規乱数を 100 万個発生させるためには，X <- rnorm(1000000,168,5) とすればよい。実行するたびに異なるデータになってしまうのを避けるためには，このコマンドを実行する前に乱数の初期値をRNGkind("Mersenne-Twister",normal.kind="Inversion"); set.seed(1) などとして固定すればよい。

1

表 1. 成人男子の身長からなる母集団と母平均 µ，母分散 σ 2 の計算
身長 (cm) (x)
度数 相対度数（p(x))
x · p(x) (x − µ)2 · p(x)

144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
合計

2
3
6
17
19
56
125
219
463
915
1609
2649
4550
7214
11005
16081
22098
29903
39048
48312
57703
66639
73332
78051
79829
77866
73767
66321
57993
48410
39081
29967
22055
15810
10875
7309
4596
2726
1519
939
462
224
128
50
3114
54
1000000

0.000002
0.000003
0.000006
0.000017
0.000019
0.000056
0.000125
0.000219
0.000463
0.000915
0.001609
0.002649
0.004550
0.007214
0.011005
0.016081
0.022098
0.029903
0.039048
0.048312
0.057703
0.066639
0.073332
0.078051
0.079829
0.077866
0.073767
0.066321
0.057993
0.048410
0.039081
0.029967
0.022055
0.015810
0.010875
0.007309
0.004596
0.002726
0.001519
0.000939
0.000462
0.000224
0.0001280.000050
0.000031
0.000014
0.000005
0.000004
1.00

表１のように計算すれば，X の平均 µ =P0.000288
0.000435
0.000876
0.002499
0.002812
0.008344
0.018750
0.033069
0.070376
0.139995
0.247786
0.410595
0.709800
1.132598
1.738790
2.556879
3.535680
4.814383
6.325776
7.874856
9.463292
10.995435
12.173112
13.034517
13.411272
13.159354
12.540390
11.340891
9.974796
8.374930
6.800094
5.244225
3.881680
2.798370
1.935750
1.308311
0.827280
0.4934060.276458
0.171837
0.085008
0.041440
0.023808
0.009350
0.005828
0.002646
0.000950
0.000764µ = 167.9998

1.151979e−03
1.586970e−03
2.903942e−03
7.496844e−03
7.599834e−03
2.021553e−02
4.049901e−02
6.328937e−021.185248e−01
2.058690e−01
3.153541e−01
4.476659e−01
6.551761e−01
8.728592e−01
1.100452e+00
1.302498e+00
1.414195e+00
1.465155e+00
1.405625e+00
1.207694e+00
9.231469e−01
5.996634e−01
2.932638e−01
7.801682e−02
3.828679e−09
7.790011e−02
2.951326e−01
5.969761e−01
9.279896e−01
1.210356e+00
1.407019e+00
1.468475e+00
1.411597e+00
1.280672e+00
1.087548e+00
8.844242e−01
6.618482e−01
4.607095e−01
2.977333e−01
2.112812e−01
1.182752e−01
6.473767e−02
4.147301e−02
1.805042e−021.240027e−02
6.174129e−03
2.420048e−03
2.116040e−03
σ 2 = 25.09521
σ = 5.009512

x · p(x) = 167.9998，X の分散 σ 2 = 25.09521 である（X は成人男性の身長からなる元の母集団を代表している）。ここから 10 個の標本を抽出することを考える。なるべく母集団全体の情報を偏りなく代表させるためには，身長が書かれた 100 万枚のカードを袋に入れてよくかきま2ぜ，10 枚取り出すことが考えられる。無作為に抽出される最初のカードが x である確率は，p(x) である。この 1 回目の抽出を X1 であらわす。2 枚目のカードについても同じである（厳密に言えば，1 枚とっては元に戻して再びよくかきまぜて抽出する「復元抽出」でないと，1 枚目に何が出たかに依存して相対度数が変化するため，条件付き確率はまったく同じにはならないが，100 万くらいの大きな母集団からサイズ 10 の標本を取り出すのであれば，「非復元抽出」（1 度サンプルを取り出したら元に戻さない）であっても，p(x) はほとん
ど変化しないし，1 枚目に何がでるかという場合を総当りして合計すれば，2 枚目に何が出るかという確率も厳密に p(x) となる）。このように，無作為に選び出されるカードは，すべて母集団における相対度数によって与えられる確率で，母集団の値のどれかをとりうる，確率変数である。
1サイズ 10 の標本の標本平均 X̄ は，X̄ = 10(X1 + X2 + ... + X10 ) なので，その期待値は，E(X̄) =1
1
1(E(X1 ) + E(X2 ) + ... + E(X10 )) =(µ + µ + ...µ) =· 10 · µ = µ101010となって母平均に一致する。標本平均の分散はvarX̄ = var{

1(X1 + X2 + ...X10 )}101 2) {var(X1 ) + var(X2 ) + ... + var(X10 )}1011= ( )2 (σ 2 + σ 2 + ... + σ 2 ) = ( )2 · 10 · σ 21010σ2=10σ√
したがって，標本平均の標準偏差は，σX̄ = 10 となる。これは，標本から推定される母集団の平均値の信頼=(性を示す値なので，平均の標準誤差と呼ばれる。

R で標本抽出をする関数は，sample() である。replace=オプションで復元，非復元を選択できる（デフォルトは非復元）。sample() 関数は乱数を使って標本抽出しているので，乱数を初期化して揃えない限り，実行するたびに結果は変わる。

¶it04-1.R³

HT <- 144:191
NUM <- c(2,3,6,17,19,56,125,219,463,915,1609,2649,4550,7214,11005,16081,22098,29903,39048,48312,57703,66639,73332,78051,79829,77866,73767,66321,57993,48410,39081,29967,22055,15810,10875,7309,4596,2726,1519,939,462,
224,128,50,31,14,5,4)
X <- rep(HT,NUM)
par(mfrow=c(2,2))
hist(s10 <- sample(X,10))
hist(s100 <- sample(X,100))hist(s1000 <- sample(X,1000))
hist(s10000 <- sample(X,10000))cat("mean(s10)=",mean(s10)," mean(s100)=",mean(s100)," mean(s1000)=",mean(s1000)," mean(s10000)=",mean(s10000),"\n")µ´サンプルサイズが大きくなるにつれて，標本分布は母集団の分布に近づく。しかし，標本平均の分布については，もっと多くのことがいえる。実行するたびに標本抽出されるサンプルは異なるのだが，その平均は一定の法則に従うのである。

3中心極限定理
つまり，ほとんどの任意の母集団から抽出された無作為標本の平均 X̄ の分布は，標本の大きさが増大する√につれて，平均 µ，標準偏差 σ/ n の正規分布に近づくことがわかっている。これを中心極限定理という。Rで試してみよう。

¶it04-2.R³

X <- c(runif(5000,0,100),rnorm(5000,100,10))par(mfrow=c(2,2))
hist(X,xlim=c(0,140),freq=F,main="元の分布")
Z5 <- rep(0,1000)for (i in 1:1000) { Z5[i] <- mean(sample(X,5)) }hist(Z5,xlim=c(0,140),freq=F)
Y2 <- dnorm(0:140,mean(X),sd(X)/sqrt(5))
lines(0:140,Y2,col="red")
Z30 <- rep(0,1000)for (i in 1:1000) { Z30[i] <- mean(sample(X,30)) }hist(Z30,xlim=c(0,140),freq=F)
Y3 <- dnorm(0:140,mean(X),sd(X)/sqrt(30))
lines(0:140,Y3,col="red")Z200 <- rep(0,1000)for (i in 1:1000) { Z200[i] <- mean(sample(X,200)) }hist(Z200,xlim=c(0,140),freq=F)
Y4 <- dnorm(0:140,mean(X),sd(X)/sqrt(200))lines(0:140,Y4,col="red")par(mfrow=c(1,1))µ´

サンプルサイズが大きくなるほど，正規分布に近づくと同時に，標本平均の標準偏差（つまり標準誤差）が小さくなっていき，ヒストグラムの幅が狭くなって，理論分布（赤い線）に近づくのが一目瞭然である。

信頼区間
標本平均の期待値は，母平均に一致するので，標本平均は，母平均のよい点推定量になっている。しかし，必ず一致するとは限らない。そこで，それがどれくらい確からしい推定かということを示すことを考える。通常，標本から計算される，ある区間の中に，正しく µ を含む割合が 95%であるような区間を推定したい（つまり，95%の信頼度をもって区間推定をしたい）と考える。
標本が大きい場合は，X̄ の正規分布において，ちょうど 95%の確率を囲むような最短の範囲を選択すればよい。そのために，正規分布の両側の裾 2.5%分を除く中央部分をとることができる。この部分は，母平均 µを中心として，大きい方と小さい方に，X̄ の標準偏差（つまり標準誤差）の 1.96 倍（1.96 は標準正規分布の97.5%点である）だけ動かした範囲を含む。つまり，
σσP r(µ − 1.96 √ < X̄ < µ + 1.96 √ ) = 0.95nnである。この式を変形すると，

σσP r(X̄ − 1.96 √ < µ < X̄ + 1.96 √ ) = 0.95nn4となる。σ は未知なので，通常，標本の不偏標準偏差 sd(X) を用いる。つまり，µ の 95%信頼区間は，√√
X̄ − 1.96sd(X)/ n から X̄ + 1.96sd(X)/ n までとなる。

標本が小さい場合は，未知の σ に対して sd(X) を代入することが無視できない誤差の原因となる。そこで，

sd(X) を使って定義される t という統計量，t=X̄ − µ√sd(X)/ nが自由度 n − 1 の t 分布に従うことから，t 分布の 2.5%点から 97.5%点までを 95%信頼区間とすることで，この誤差を回避することができる。自由度が無限大になれば t 分布は正規分布に一致するので，結局，常にこちらで計算すればよいことになる。
すなわち，標本サイズ n，標本平均 X̄ ，標本の不偏標準偏差 sd(X) のとき，母平均の 95%信頼区間は，√X̄ − t0.025 sd(X)/ nから√X̄ + t0.025 sd(X)/ nまでになる。R で自由度 n − 1 の t 分布の 97.5%点を与える関数は qt(0.975,n-1) なので，R のプログラム上で，例えば it04-1.R で定義した s100 から母平均の 95%信頼区間を推定するには，次のようにすればよい。

¶it04-3.R³barX <- mean(s100)
sdX <- sd(s100)
t975 <- qt(0.975,length(s100)-1)
rootn <- sqrt(length(s100))print(barX - t975*sdX/rootn)
print(barX + t975*sdX/rootn)µ´

この結果から，95%信頼区間を表示するときは，[166.96, 168.96] のように記載するのが普通である。ただし，もちろん，R にはもっと楽にこの計算をしてくれる関数が用意されていて，標本データが付値されている
数値型変数 X について，t.test(X) とすれば，母平均がゼロであるという帰無仮説の検定結果（検定については後日説明する）とともに，母平均の推定値と母平均の 95%信頼区間が表示される。t.test(s100) が表示する母平均の推定値と 95%信頼区間は，上のコードで得られる結果と一致するはずである。各自確認されたい。

自由度
自由度について直感的に理解するためには，次のように考えると良い。n 個の数値からなる標本には，はじめ n の自由度がある。しかし，自由度の１つは，X̄ を計算するときに使われてしまい，sd(X) を計算するための偏差 Xi − X̄ に対しては n − 1 の自由度しか残らない。逆に見れば，サイズ n の標本について，はじめのn − 1 個の偏差は自由だが，最後の偏差は，全偏差の合計がゼロにならねばならない（そのように平均という量はとられている）ために決定済みであり，結局自由に決められる個数は n − 1 となる。
一般には，データの数から，推定した母数の数を引いた値が，その統計量や分布の自由度になる。

5

課題次のいずれかを選んで回答せよ。

(1) 今日の講義の内容について問う問題を作成せよ。
(2) 2001 年に厚生科学研究で行われた「少子化の見通しに関する専門家調査」の結果の一部を見てみる。この調査は，「人口学，経済学，家族社会学，公衆衛生学を中心とした専門家を対象として少子化研究会のメンバーが対象候補者を抽出し，回答者の偏りや不足等について検討を加えた上で，748 名を対象として調査を実施した」もので，回収率は 44 ％であった。この調査では，2025 年の合計出生率がいくつになるかという推定値が尋ねられていて，生データを見ると，1.38 1.50 1.30 1.40 1.40 1.15 1.31 1.37 1.50 1.55 1.55 1.56
1.50 1.56 1.50 1.38 1.50 1.20 1.20 1.50 1.25 1.25 1.22 1.40
1.80 1.37 1.35 1.70 1.35 1.50 ...（後略）のようになっていて，回答数は 311，平均値は 1.385，不偏分散は 0.0252 であった。以上の情報から，専門家母集団における，2025 年の合計出生率の推定値の 95%信頼区間を求めよ。

(1)(2) いずれに回答する場合も，配布する紙に学籍番号，氏名と共に自筆して提出すること。提出をもって出席確認とする。

6