ナード戦隊データマン

データサイエンスを用いて悪と戦うぞ

シミュレーションによって確率予測する

単純な条件で、厳密に定義された問題では、機械学習を用いるよりシミュレーションをしたほうが理論値の予測に役立てられる場合があります。ここでは、ギャンブラーの破産問題という確率論の問題を通じて、Rによるシミュレーションを試します。

ギャンブラーの破産問題

10円持っているAliceと、5円持っているBobがいます。50%の確率で表裏が出るコインを投げ、表が出ればAliceはBobから1円もらい、裏が出ればBobはAliceから1円もらいます。どちらか一方の所持金がなくなるまでゲームを繰り返すとき、両者の勝率はどのぐらいでしょうか?また、初期資金を一般に変数A, Bとしたときにはどうなるでしょうか。

R言語のコード

シミュレーション用の3つの関数を作成します。

coin_gamble <- function(M, N){
  tmp_M = M; tmp_N = N;
  while(T){
    coin = rbinom(n=1, size=1, prob=0.5)
    if(coin==1){tmp_M = tmp_M-1; tmp_N = tmp_N+1}
    else{tmp_N = tmp_N-1; tmp_M = tmp_M+1}
    if(tmp_M == 0){return(0)}
    if(tmp_N == 0){return(1)}
  }
}

coin_gamble_repeat <- function(M, N, size){
  out = c()
  for(i in 1:size){
    out[i] = coin_gamble(M, N)
  }
  return(out)
}

coin_gamble_stat <- function(M, N, repeat_size, mean_size){
   out = c()
   for(i in 1:mean_size){
     out[i] = mean(coin_gamble_repeat(M,N,repeat_size))
   }
   return(out)
}

これらの関数は、ギャンブラーの破産問題を実行する関数です。coin_gambleは、Aliceが勝てば1, Bobが勝てば0をリターンします。coin_gamble_repeatによって、指定回数ゲームを続けます。coin_gamble_statは、指定回数のゲームを、指定回数続けます。言い換えれば、二重ループです。内側のループはゲーム回数であるのに対し、外側のループは平均を求める回数です。

シミュレーションによって求めた勝率をプロットする

シミュレーションによって勝率を求めてからプロットしてみましょう。

results = coin_gamble_stat(10, 5, 100, 1000)
hist(results, xlim=c(0,1), breaks=100)

Screenshot from 2017-10-25 13-31-11.png

0.7あたりを中心にした正規分布っぽくなっていることがわかります。これは、中心極限定理によるものです。

平均と標準偏差を求める

> mean(results)
[1] 0.66595
> sd(results)
[1] 0.04701391

このように、2/3が平均になっています。つまり、Aliceが10円、Bobが5円持っているときのAliceの勝率は2/3ということになります。

他の値で試して一般化の予想をする

AliceとBobの資金がなんであっても成立する公式を導きたいので、他の値で試してみます。

for(i in 1:10){
  out[i] = mean(coin_gamble_stat(i, 5, 10, 100))
}

Bob = rep(5,10)
Alice = seq(10)

results.df = data.frame(y=out, A=Alice, B=Bob)
results.df$A2AB <- results.df$A/(results.df$A+results.df$B)
>results.df
       y  A B      A2AB
1  0.177  1 5 0.1666667
2  0.304  2 5 0.2857143
3  0.378  3 5 0.3750000
4  0.444  4 5 0.4444444
5  0.503  5 5 0.5000000
6  0.540  6 5 0.5454545
7  0.579  7 5 0.5833333
8  0.637  8 5 0.6153846
9  0.627  9 5 0.6428571
10 0.683 10 5 0.6666667

> mean((results.df$y - results.df$A2AB)^2)
[1] 0.0001493324

このように、Aliceの勝率はA/(A+B)によって求めると良い精度が出ることがわかります。また、Bobの勝率はB/(A+B)だと考えることができるでしょう。その最小二乗誤差はわずか0.0001です。

参考文献

なし