大数の法則を確かめよう〜「数学ガール 確率の冒険」研究問題1-1
数学ガール
結城浩先生の「数学ガールの秘密ノート 確率の冒険」を読んでいます。
もう、さすがとしかいいようがありません。第1章は、わたしが書こうとしている確率変数の本の第1章の内容が、これでもかというくらいに、丁寧に、そして読みやすく、書かれている。ああ、もう、私の仕事はないな。そう思えてしまいます。
といいつつ、私は私で、確率変数と確率分布について、そもそもそれはどのようにイメージすれば良いのかということ、それから、確率変数を「足す」という説明が統計学の教科書にはたびたび登場しますが、それはいったい何をしているのか、どうイメージすれば良いのか、ということ。などについて、私の方法で書いてみたいと思っています。
といいつつ。私は統計の専門家でもなんでもなく、ただの横好きなので、できるだけ間違いのないように、でも、これまでの統計の本とは違う味わいをもたせた本にしたいと思っているのです。
ところで。結城先生の本には練習問題があり、巻末には研究問題もあります。これがまた。私が別の本で勉強したことや、私の本に書こうと思っていたことが、もうほとんどそのまま形になっているという、優れた問題です。問題を引用しつつ、私の回答を書いて、私の本のための準備にしたいと思っております。
研究問題1-1
研究問題1-1はこんな問題です。
これは要するに、大数の法則を目で見て確かめよう、という趣旨です。実際にコイン投げをするとわかりますが、10回くらい投げたところで、表が5回、裏が5回などというふうには、めったになりません。(もちろん、そうなることだってあるんですけど。)
けれど、投げる回数を増やしていくと、だいたい表と裏が半々くらい出てるなあ、という結果に近づいていくよね、というのが大数の法則です。
(こんなざっくりしすぎた説明は許さん!という方は、どうぞお引き取りください。続きは読まないでください。)
本当にサイコロ投げをしてもいいのですが、確率0.5で表が出るフェアなコインを仮定して、つまり、コンピュータプログラムとして用意して、何回もコインを投げてもらうことにします。では、R Studioの登場です。
先に申し上げておきますが、以下、ゆるっと、ふわっとした説明です。
「変数」とか、「代入」とか「forループ」とか、ちゃんと正しい用語で説明せんかい!という方も、どうぞお引き取りください。
# 数学ガール研究問題
# 研究問題1-1 確率と相対度数
# 結果を逐一画面に表示する
OUTPUT <- FALSE
# フェアなコイン:表が出る確率を0.5とする
prob <- 0.5
# コイン投げ10回から10000回まで試す
trial <- seq(10, 10000, 10)
# 相対度数を記録する配列
result <- numeric(length(trial))
# コイン投げ実験開始
for (t in 1:length(trial)) {
# 記録のための配列
n <- numeric(trial[t])
# コイン投げ開始
for (i in 1:trial[t]) {
n[i] <- rbinom(1, 1, prob)
}
result[t] <- mean(n)
if (OUTPUT) print(paste("trial=",trial[t],", success=", sum(n), ", prob = ", mean(n), sep=""))
}
# 記録
result
# グラフ
plot(trial, result, type="l", ylim=c(0,1), main=paste("Probability =",prob))
abline(h=prob, col="red", lty=2)
実際にコイン投げをしているのは「#コイン投げ開始」という部分です。再掲します。
# コイン投げ開始
for (i in 1:trial[t]) {
n[i] <- rbinom(1, 1, prob)
}
真ん中あたりに見える「rbinom」が乱数を発生させる関数です。括弧に入っている(1, 1, prob)は、先頭から順に、「乱数は1個でいいよ」「コイン投げ1回だけやってね」「表が出る確率はprobとしてね」という意味です。「prob」って何なん? というのは、元のコードの上の方に書いてあります。
# フェアなコイン:表が出る確率を0.5とする
prob <- 0.5
ここで、「prob」という名前が出てきたら、それは「0.5」ということなんで、よろしく!と書いているのです。
で、元に戻って、「for」というのが書いてあって、その後になんやら書いてありますが、要するに「rbinom」でコイン投げをして、表なのか裏なのかを記録するという作業を、必要な回数繰り返してちょうだい!というのが「for」の役割です。
またしても書いておきますが、ここはもっと簡単に書けます。「for」なぞ使わなくても、1行で済みます。でも、やっていることをできるだけそのまま形にしようと思って、「わざわざ」冗長に書いています。ご理解ください。(「わざわざ」と鉤括弧で書くあたり、「本当はちゃんと1行で書けるんだぞ」という気持ちの表れです。・・・何を説明しているんだか。はあ。)
プログラムの実行結果
実行結果は次のようになります。(疑似乱数を使っているので、実行結果は実行するごとに少しずつ異なりますが、以下に説明する特徴は共通しています。)
問題の指定通り、横軸はコイン投げをした回数、縦軸は表が出た相対度数(表が出た回数÷コイン投げをした回数)です。相対度数は「割合」と言い換えても同じです。
最初の方こそ、0.5からずいぶん遠い値が出ていますが、1000回近く繰り返すと、ほぼ0.5に近い値しか出ていません。どうして0.5に近い値になってくるかというと、コイン投げをして表が出る確率を、最初に、0.5と仮定したからです。
インチキなコインならどうなん? というのも試せます。さきほどの「prob<-0.5」のところを、「prob<-0.8」とかにすれば、「やたらに表ばかり出るコイン」(そんなコインは現実にはないと思います)で試すことができます。最初に0.8で設定すれば、ほとんど0.8に近い値しか出なくなります。大数の法則すごい。直感的には「当たり前じゃん?」と思われるのですが、これをちゃんと研究したヤコブ・ベルヌーイ先生は、やっぱりすごい。
このグラフだと細かい数値は確認できないので、Rが実行できる環境をお持ちの方は、上のコードをコピーして実行してみてください。グラフを描く前の、相対度数の値を、いったんぜんぶ表示する設定になっています。
え? Rなんて使えない?
それは困った。
Javascriptでシミュレーションしよう
はい、「Rなんて使えない」という方、「いや、使えないことはないんだけどさ…」という方のために、ブラウザ上で同じようなプログラムの実行ができるサイトを用意しています。ブラウザがあればOKです。「ブラウザ」というのは、インターネットするときの、ほら、Google Chromeとか、Safariとか、Firefoxとか。そういうのです。Edgeというのもありますね。
というわけで、次の記事で、そのことについて書こうと思います。
カバー写真
カバー写真はUnsplashのAndre Taissinが撮影した写真でした。「coin」で検索すると、ビットコインぽい画像ばかり出てきて面白くないので、「coin toy」で検索したら、なんと可愛い! ブタの貯金箱くんがコインを食べようとしているではないですか! 気に入った。