プログラミングスキルをつかってベイズ統計を理解するという面白そうな本があったので購入してみました。そのメモです。
この本ではthinkbayes.pyで定義されるクラスや関数を使っています。
以下のリンクからthinkbayes.pyというモジュールをダウンロードして使います。
http://thinkbayes.com/thinkbayes.py
・使い方
正六面体のサイコロの分布結果
from thinkbayes import Pmf
pmf = Pmf ()
for x in [1,2,3,4,5,6]:
pmf.Set(x, 1/6.0)
//Setメソッドが各値に確率1/6を設定している。
//各サイコロの目の確率の出力方法は以下
print pmf.Prob(1)
0.166666666667
print pmf.Prob(2)
0.166666666667
print pmf.Prob(3)
0.166666666667
・ベイズの定理
p(H|D) = p(H)p(D|H)/p(D)
p(H) : データを見る前の仮設の確率。事前確率（prior probability）
p(H|D) : データを見た後の仮設の確率。事後確率（posterior probability or posterior）
p(D|H) : 仮説の下でのデータの確率。尤度（likelihood）
p(D) : どのような仮説でもデータの得られる確率。正規化定数（normalizing constant）
・クッキー問題
「クッキーの入った2つのボウルがある。
ボウル1 – バニラクッキー 30枚 :チョコレートクッキー 10枚
ボウル2 – バニラクッキー 20枚 :チョコレートクッキー 20枚
どちらかのボウルをランダムに選びクッキーを1つランダムに選んだ場合、
それがボウル1だという確率を求めよ。」
ステップ1 事前確率の設定-仮説に確率を割り当てる
仮説B1 ボウル1からバニラクッキーを選ぶ確率 1/2
仮説B2 ボウル2からバニラクッキーを選ぶ確率 1/2
from thinkbayes import Pmf
pmf = Pmf()
pmf.Set('Bowl 1', 0.5)
pmf.Set('Bowl 2', 0.5)
ステップ2 分布の更新（各事前確率に対応する尤度を掛ける）
ボウル1からクッキーを取り出す確率 3/4
ボウル2からクッキーを取り出す確率 1/2
pmf.Mult('Bowl 1', 0.75)
pmf.Mult('Bowl 2', 0.5)
Multメソッドは指定された仮説の仮設を取り出し、そこに指定された尤度を掛ける。
ステップ3 最正規化しボウル1の事後確率を求める。
pmf.Normalize()
print pmf.Prob('Bowl 1')
答え 0.6
このモジュールを使ってほかの問題も解いてみます。
以下は、ネット上で見つけた問題です。
・20年前の早稲田大学の入試問題
「５回に１回の割合で帽子を忘れるくせのあるＫ君が、正月に Ａ、Ｂ、Ｃ ３軒を順に年始回りをして家に帰ったとき、
帽子を忘れてきたことに気がついた。２軒目の家 Ｂ に忘れてきた確率を求めよ。」
http://d.hatena.ne.jp/pashango_p/20090809/1249805193
ステップ1 事前確率の設定-
A Aに帽子を忘れてきた確率1/3
B Bに帽子を忘れてきた確率1/3
C Cに帽子を忘れてきた確率1/3
ステップ2
A Aに帽子を忘れて確率（1/5）
B Aで帽子を忘れない確率（4/5）× Bで帽子を忘れる確率（1/5）
C Aで帽子を忘れない確率（4/5）× Bで帽子を忘れない確率（4/5） × Cで帽子を忘れる確率（1/5）
from thinkbayes import Pmf
pmf = Pmf()
pmf.Set('A', 1/3.0)
pmf.Set('B', 1/3.0)
pmf.Set('C', 1/3.0)
pmf.Mult('A', 1/5.0)
pmf.Mult('B', 4/25.0)
pmf.Mult('C', 16/125.0)
pmf.Normalize()
print pmf.Prob('B')
答え 0.327868852459
・喫煙者の推定の問題
「男性１０人、女性７人が一室でパーティーを開いた。男子の喫煙者は５人、女性は３人である。
部屋に入ったら煙草の吸殻が１本、灰皿の上にあった。このとき、吸った人が女性である確率を求めなさい（煙草の吸い回しはしていないと仮定する）」
http://www.yasuienv.net/BayesExamples.htm
1.事前確率
Men :10/17
Women :7/17
2.分布の更新
Men :1/2
Women :3/7
pmf = Pmf()
pmf.Set('Men', 10/17.0)
pmf.Set('Women', 7/17.0)
pmf.Mult('Men', 1/2.0)
pmf.Mult('Women', 3/7.0)
pmf.Normalize()
print pmf.Prob('Women')
答え 0.375
・雨の日のお誘いの問題
「酒好きのＡさんはＢ氏をよくお酒に誘う。統計をとると、雨の降っていない日に誘うと、Ｂ氏は５回中４回誘いに応じ、
雨の降っている日に誘うと、５回中３回誘いに応じることが分かった。Ｂ氏がＡさんの誘いに応じたとき、雨が降っていない確率を求めよ。
雨が降った日と降らない日の割合は７：１とする」
http://www.yasuienv.net/BayesExamples.htm
1.事前確率
sunny 7/8
rainy 1/8
2.尤度
sunny 4/5
rainy 3/5
pmf = Pmf()
pmf.Set('sunny', 7/8.0)
pmf.Set('rainy', 1/8.0)
pmf.Mult('sunny', 4/5.0)
pmf.Mult('rainy', 3/5.0)
pmf.Normalize()
print pmf.Prob('sunny')
答え 0.903225806452
このモジュールを使えば簡単に確率の問題が解けることが分かりました。
以上