少し前に参加した続・わかりやすいパターン認識の読書会でモンティ・ホール問題について知ったので”司会が正解を知っているか否かによる結果の違い”について、ちょうど触ってみたかったGoで確かめてみました。(Go触りたかっただけ)

モンティ・ホール問題って？

(1) 3つのドア (A, B, C) に（景品、ヤギ、ヤギ）がランダムに入っている。
(2) プレイヤーはドアを1つ選ぶ。
(3) モンティは残りのドアのうち1つを必ず開ける。
(4) モンティの開けるドアは、必ずヤギの入っているドアである。
(5) モンティはプレーヤーにドアを選びなおしてよいと必ず言う。

モンティ・ホール問題 - Wikipedia

(5)で選び直した方が良いのかどうかという問題です。

ここでは、プレイヤーが最初に選んだドアをA、司会のモンティがあけるドア(はずれ)をB、残っている開けられていないドアをCとします。

司会が正解を知っている場合、Aが正解の確率は1/3、Cが正解の確率は2/3となり、選び直した方が良い事になります。確率はベイズの定理で求まります。

f:id:nihma:20141229004918j:plain

ここで司会が正解を知らなかった場合、確率が変わります。 Aが正解の確率は1/2、Cが正解の確率は1/2でどちらを選んでも同じ確率です。こちらもベイズの定理で求まります。

f:id:nihma:20141229005058j:plain

P(Bが選ばれる|Bが当たり)とP(Bが選ばれる|Cが当たり)が違います。

確かめてみました

とりあえずGoで10000000回の試行による確率を求める計算をそれぞれ５回ずつ行ってみました。大数の法則を信じて10000000回もやれば大丈夫と信じました。

コードは次の通りです。
問題の(4)の条件から、Bが必ずはずれであるのでAorCが当たりとなり、司会が正解を知らない場合にBを選んで当たってしまったケースは試行に含まれないです。
(tryA関数は無い方が速いですがgoroutineしてみたかったのでしました。)

package main

import (
  "fmt"
  "math/rand"
  "time"
  "sync"
  "runtime"
)

func main() {
  // CPU全コアを使いたい
  runtime.GOMAXPROCS(runtime.NumCPU())

  // ランダムの初期化
  rand.Seed(time.Now().UnixNano())

  // それぞれ5回ずつprintしてみる
  label_know := map[string] bool {
    "司会が正解を知っている場合": true,
    "司会が正解を知らない場合":   false,
  }
  n := 5
  printP := func(p_a float64, p_c float64) {
    fmt.Printf("P(A) = %.3f, P(C) = %.3f\n", p_a, p_c)
  }
  for label, know_answer := range label_know {
    fmt.Println(label)
    tryNtimes(know_answer, n, printP)
  }
}

/* n回やってみる
 * 引数1：know_answer 司会が当たりを知っているか？
 * 引数2：n やってみる回数
 * 引数3：f 結果を渡す関数
 */
func tryNtimes(know_answer bool, n int,
               f func(float64, float64)) {
  var wg sync.WaitGroup
  for i := 0; i < n; i++ {
    wg.Add(1)
    go func() {
      var p_a float64 = probabilityA(know_answer)
      var p_c float64 = 1 - p_a
      f(p_a, p_c)
      wg.Done()
    }()
  }
  wg.Wait()
}

/* ドアAが当たりの確率
 * 引数：know_answer 司会が当たりを知っているか？
 */
func probabilityA(know_answer bool) float64 {
  try_num := 10000000  // 試行回数(このくらいやれば収束する？)
  a_num := 0 // ドアAが当たりの回数

  ch := tryA(know_answer)
  for i := 0; i < try_num; i++ {
// こっちの方が速いですが->    if isA(know_answer) {
    if <- ch {
      a_num++
    }
  }
  return float64(a_num) / float64(try_num)
}

/* ドアAが正解であるかの試行を生成
 * 引数：know_answer 司会が当たりを知っているか？
 */
func tryA(know_answer bool) chan bool {
  ch := make(chan bool)
  go func() {
    for {
      ch <- isA(know_answer)
    }
  }()
  return ch
}

/* ドアAが正解であるか？
 * 引数：know_answer 司会が当たりを知っているか？
 */
func isA(know_answer bool) bool {
  /* [答えの候補]
   * 0: ドアA <- 最初に選択
   * 1: ドアB <- 司会が選択(はずれ)
   * 2: ドアC <- こちらに変更できる
   */
  // 最初に答えを決める
  answer := rand.Intn(3)  // 0〜2

  if ! know_answer { // 答えを知らない場合
                     // ドアBが当たりになる試行は無効、捨てる
    for answer == 1 {
      answer = rand.Intn(3)
    }
  }
  return answer == 0 // ドアAが当たりか？ドアCがはずれか？
}

結果は理論通りになったようです。

$ go build MontyHall.go
$ ./MontyHall
司会が正解を知っている場合
P(A) = 0.333, P(C) = 0.667
P(A) = 0.334, P(C) = 0.666
P(A) = 0.334, P(C) = 0.666
P(A) = 0.333, P(C) = 0.667
P(A) = 0.333, P(C) = 0.667
司会が正解を知らない場合
P(A) = 0.500, P(C) = 0.500
P(A) = 0.500, P(C) = 0.500
P(A) = 0.500, P(C) = 0.500
P(A) = 0.500, P(C) = 0.500
P(A) = 0.500, P(C) = 0.500