唯物是真 @Scaled_Wurm

プログラミング(主にPython2.7)とか機械学習とか

ソシャゲで確率を推定するにはどれくらい試行が必要か

ソシャゲをやっていてドロップ確率などを調べるには何回ぐらい試せばよさそうか気になったので調べてみた

一定確率\(p\)で成功、\(1-p\)で失敗するという単純なモデルを仮定する
つまりある確率\(p\)でアイテムがドロップするとしたときに、試行回数に対する成功回数の割合で成功確率\(p\)を推定する
何回ぐらいやれば確率\(p\)の推定結果がどれぐらいばらつくかを調べる
運営の闇の魔術によって久しぶりにプレイするとレアが出たり連続では当たりが出なかったりという操作があるかもしれないがそういうのは考えない

二項分布の信頼区間を求めたりしてみてもよいがWikipediaをみるとわりと複雑そうなので、実際にシミュレーションして値のばらつきを調べてみる(正規分布近似してもよいかもだけど

バイオリンプロット

バイオリンプロットという分布の可視化の方法を見かけたのでせっかくなので使ってみる(わかりやすいかはさておき
こんな感じの図でバイオリン(？)みたいな形になっている左右の2つの弧は、データの推定された確率分布の形を表している

f:id:sucrose:20161203222027p:plain

弧の中の四角形は箱ひげ図を表している(箱ひげ図以外にもいろいろな可視化が使われる

箱ひげ図

分布の可視化だと箱ひげ図はよく出てくる方法で、以下のような形をしている
f:id:sucrose:20161203223431p:plain

箱ひげ図はデータの分布を可視化するためのもので、中央の箱の上下が四分位点(データを順番に並べたときに上と下から25%までのデータが含まれる位置)で箱の真ん中の線が中央値(データを順番に並べたときに片側から半分までのデータが含まれる位置)
つまり箱の中には分布の50%が含まれる

箱から伸びている線(ひげ)は、四分位点の間の範囲かける1.5の領域を箱の端からとったときにその中に含まれるデータのうち最も遠いものまでの線になっている(ひげの書き方はいろいろ種類があるらしい
ひげの外側の点は上記の範囲に含まれなかったデータで、外れ値的なものを表すとされる
ひげの部分には箱に含まれなかった50%の大部分が入っていると考えられる(もちろん分布の形による)

上の箱ひげ図はバイオリンプロットの方と同じように2つ山がある分布を可視化しているが、箱ひげ図だとその情報は失われてしまう

結果

というわけでシミュレーション結果の分布をバイオリンプロット(中は箱ひげ図)で作図してみた

20, 100, 1000, 10000回のそれぞれの回数試行して確率を推定するのを1セットとし、十分に大きなセット数(今回は100000)やったときの推定された確率の分布を調べた
具体的な実装としては、0から1の乱数を生成して与えられた確率値以下のものの割合を調べた
ちなみに一番小さな試行回数が10回ではなく20回なのは、試行回数が少なすぎるとグラフの形が残念になることが多かったので少し増やしています

\(p=0.02\)

というわけで成功確率が0.02のときのプロット結果を見てみる(見づらかったので横向きにしました)
f:id:sucrose:20161203225454p:plain
真ん中の太い黒い部分が箱ひげ図の箱の部分(上下の四分位点の間)に当たる
つまり試行回数が20回の時には0.05から0の間に50%以上のデータが含まれることになる
この箱の部分に注目すると、試行回数100回でもそれなりに幅があり、1000回になると明らかに狭くなってくる感じがする

確率を変えてもあまり傾向に変化はないので以下の図の説明は省略する

\(p=0.05\)

f:id:sucrose:20161203225940p:plain

\(p=0.1\)

f:id:sucrose:20161203230007p:plain

\(p=0.3\)

f:id:sucrose:20161203230023p:plain

まとめ

2桁や3桁ぐらいのオーダーの1人の人間がちょっとやる程度の試行回数だと結構大きな差が出そう
というか、バイオリンプロットを見てもよくわからないですねorz

ソースコード

バイオリンプロットにはSeabornを使いました

import random
import numpy as np
import numpy.random as npr
import seaborn as sns

N = [20, 100, 1000, 10000]
prob = [0.02, 0.05, 0.1, 0.3]
M = 100000

for p in prob:
    temp = []
    for n in N:
        result = []
        for i in xrange(M):
            result.append(np.count_nonzero(npr.random(n) <= p) / float(n))
        temp += result
    
    y = []
    for n in N:
        y += [n] * M
    
    sns.violinplot(x=temp, y=y, bw=1, scale='width', cut=0, orient='h')
    sns.plt.title('p = ' + str(p))
    sns.plt.ylabel('Number of Trials')
    sns.plt.xlabel('Distribution of Estimated Probability')
    sns.plt.show()
'); jQuery.noConflict(true);