ベイズの考え方をRで簡単に実験してみよう。

ベイズの勉強をしたいと思ったが、何をみてもいまいち実感が湧かなかったため、Rで実験したくなった。(初めの段階の理解は以下の記事参照)

頻度論とベイズは同じ式の形をしてても（同じ確率密度関数の形）解釈が違う。例えば以下のHPで説明されている

しかし私自身が実際やってみないと何を言っているのかわからなかったためRで手を動かして実験してみる¹。

まずは頻度論の考え方

目の前にコインがあるとする。一般的にはこのコインが表になる確率は1/2だろうと考える(仮定する)。5回くらいコインを投げて全部表が出た場合に、そんな確率は0.03125(=1/32：3%ちょっと)しかないんだからこのコインはいかさまだ！！と考えるのが頻度論の考えである。

頻度論ではこのコインはまずは「理想的なコインである」と仮定して（考えて）おり、理想的なコインとしてふさわしくない行動を取るとそれはいかさまコインだと考える。

分布を考えて(plotを書いて)、それより低い所の総和(面積)が0.05%以下となっている所なのか否かを検定する事になる。

目の前に得体の知れないコインがあるとする。このコインがどのくらいの確率で表が出るのか(表が出る真の確率であるpの値がいくらと考えるのが妥当か)を実際やってみて(データから)推定したい。

まずは事前情報がないので、事前分布(今回は共役事前分布²としてベータ分布を使う)を一様分布：uniとする

pには「このコインが本当はどのくらいの確率で表が出やすいんだろうか」という確率が入っており、縦軸の高さはわかりやすく言うと「そのp値(表の出やすさ)が妥当な値(真の値)と決めてしまっても良いよーという度合い」が入っている³：自分は理解の上でここが一番難しかった

p <- seq(from=0,to=1,length=100)
uni <- dbeta(p,1,1)
plot(p,uni)

初めの状態は、何の情報もないので、均一(どのp値も同じくらい真の値でありうる)

次に、このコインを5回回すと、表が3回・裏が2回でたとする。

aft <- dbeta(p,1+3,1+2)
plot(p,aft)

表が出やすいかもしれないという程度のグラフになる。(細かく言うと、今ある情報だけから考えると、pは0.5よりちょっと大きめなコインなのかもしれないな、でも他の値である確率も十分高いな)

例えばここで、本当はp=0.2じゃないか？(イカサマコインだから表は5回に1回しか出ない,というのがp=0.2の意味)とか言われると、まだその可能性も十分あるよね、もう少しやってみないとわからないよね、と返事する事になる。

さっきは表が3回・裏が2回 (3,2)としたが、自分でこの値を(4,0)とか(4,1)とか(5,1)とか自分で実験していくと、だいぶイメージがついてくると思う⁴。

そしてそんな数回ではなくコインをどんどん回してみて、最終的には例えば表が1000回、裏が1000回出たとするとこんなグラフに成長する。

matomo <- dbeta(p,1+1000,1+1000)
plot(p,matomo,type="l") #グラフを直線表記にした

だいぶ尖っている。

これをみてさっきのグラフと比較すると、p=0.5なのだろう(表が出る確率が0.5なのだろう)という考えに近づいていく事が可視化される。

つまりこれはp=0.5のまともなコインなのだろう、という結論に近づく。

ここでp=0.2じゃないか？とか言われると、さすがにそれはないんじゃないかな、といえる。

どっちが良いかという結論はないだろうけど、

気がするので、今後置き換わってくるかもしれない。ビジネスに使うなら人間の考え方に近いベイズなのではないだろうか。