2017-02-09 13 views
3

多数の実行の詳細な結果(確率ツリー)を取得する方法:回ベルヌーイ実験に以下の実験を想定回

IのN数を(成功Pの確率で)同じベルヌーイ試行を実行します次の情報が必要です:すべての可能な成功/失敗のシーケンスが起こる確率。

例:成功P = 40%の確率で Aベルヌーイ実験を3回実行し、以下の結果をもたらすであろう(Sは成功であり、Fは失敗である):

FFF 0.216

SFF 0.144

FSF 0.144

SSF 0.096

FFS 0.144

SFS 0.096

FSS 0.096

SSS 0.064

Iの結果を得ることをブルートフォースしようとしたが、それは唯一のN = 25と急速にチョーク、私はOutOfMemoryExceptionを取得します...

using System; 
using System.Linq; 
using System.Collections.Generic; 
using System.Text.RegularExpressions; 

namespace ConsoleApplication 
{ 
    class Program 
    { 
     static Dictionary<string, double> finalResultProbabilities = new Dictionary<string, double>(); 

     static void Main(string[] args) 
     { 
      // OutOfMemoryException if I set it to 25 :(
      //var nbGames = 25; 
      var nbGames = 3; 
      var probabilityToWin = 0.4d; 

      CalculateAverageWinningStreak(string.Empty, 1d, nbGames, probabilityToWin); 

      // Do something with the finalResultProbabilities data... 
     } 

     static void CalculateAverageWinningStreak(string currentResult, double currentProbability, int nbGamesRemaining, double probabilityToWin) 
     { 
      if (nbGamesRemaining == 0) 
      { 
       finalResultProbabilities.Add(currentResult, currentProbability); 
       return; 
      } 

      CalculateAverageWinningStreak(currentResult + "S", currentProbability * probabilityToWin, nbGamesRemaining - 1, probabilityToWin); 
      CalculateAverageWinningStreak(currentResult + "F", currentProbability * (1 - probabilityToWin), nbGamesRemaining - 1, probabilityToWin); 
     } 
    } 
} 

私は

(すべてのPのために3秒未満で結果を得る)タイムリーに= 3000 Nまでサポートできるようにする必要がありますが最適にこれを行うための数学的な方法はありますか?

+0

なぜあなたは*すべて*の結果を保存したいですか? 'P(F ... S ... F ... S)== P(S)**(Sの数)* P(F)**(Fの数)' –

+0

@DmitryBychenko一番長い勝利ストリークの平均(例えば0.216 * 0 + 0.144 * 1 + 0.144 * 1 + 0.096 * 2 + 0.144 * 1 + 0.096 * 1 + 0.096 * 2 + 0.064 * 3 = 1.104) – ibiza

答えて

1

ここ正確かつ唯一の二次であること十分な速さの異なるアプローチが、です。最長最長ストリークの期待値は、

n 
sum Pr(there exists a win streak of length at least k). 
k=1 

となります。確率は次のようになります。レコードが長さk勝利(確率pwin**k)で開くか、またはで始まり、何らかのj in 0..k-1で勝利した後、損失(確率pwin**j * (1 - pwin))が発生し、確率は長さの確率に等しい - k勝ちn - (j + 1)のストリークが試みます。私たちはこの論理が暗示している再発を評価するためにメモを使用します。pwinstreak; fastpwinstreakのより速いバージョンは、繰り返しの合計を避けるために代数を使用します。ことができます

def avglongwinstreak(n, pwin): 
    return sum(fastpwinstreak(n, pwin, k) for k in range(1, n + 1)) 


def pwinstreak(n, pwin, k): 
    memo = [0] * (n + 1) 
    for m in range(k, n + 1): 
     memo[m] = pwin**k + sum(pwin**j * (1 - pwin) * memo[m - (j + 1)] 
           for j in range(k)) 
    return memo[n] 


def fastpwinstreak(n, pwin, k): 
    pwink = pwin**k 
    memo = [0] * (n + 1) 
    windowsum = 0 
    for m in range(k, n + 1): 
     memo[m] = pwink + windowsum 
     windowsum = pwin * windowsum + (1 - pwin) * (memo[m] - pwink * 
                memo[m - k]) 
    return memo[n] 


print(avglongwinstreak(3000, 0.4)) 

バージョンエラー:

def avglongwinstreak(n, pwin, abserr=0): 
    avg = 0 
    for k in range(1, n + 1): 
     p = fastpwinstreak(n, pwin, k) 
     avg += p 
     if (n - k) * p < abserr: 
      break 
    return avg 


def fastpwinstreak(n, pwin, k): 
    pwink = pwin**k 
    memo = [0] * (n + 1) 
    windowsum = 0 
    for m in range(k, n + 1): 
     memo[m] = pwink + windowsum 
     windowsum = pwin * windowsum + (1 - pwin) * (memo[m] - pwink * 
                memo[m - k]) 
    return memo[n] 


print(avglongwinstreak(3000, 0.4, 1e-6)) 
+0

もう一度ありがとう!私はそれを見て試してみる。 Q:最初のバージョンでは、 'pwinstreak'関数の目的は何ですか、私は間違っていると思って呼び出されませんか? – ibiza

+1

@ibizaあなたは正しいです。それは決して呼び出されません。それはもっと明らかに正しいが、遅い。高速版は同等であることがわかります。 –

+0

これはすごくうまくいきます!私はあなたの深い数学的なスキルを羨望:)あなたの知識を共有してくれてありがとう – ibiza

2

最長のストリークの長さに興味があるので、トライアルのどの時点でも、歴史に関する2つの関連する事実があります:(1)最長のストリークがどれだけ長く(m) (2)現在のストリークが(k)である期間。初期状態は(0、0)です。勝者は(m、k)→(max(m、k + 1)、k + 1)、損失は(m、k)→(m、0)となります。すべての最終状態の確率を知ると、平均をとることができます。

EDIT:このコードの改訂版では、非常に起こりそうなイベントに必要な計算を削減する最適化が行われています。具体的には、ある程度の長さのストライクがすべてcutoff以上のすべての状態を無視します。絶対誤差予算abserrが与えられた場合、可能な最長ストリークはnなので、最終予想値計算からabserr/nまでの確率質量を除外することができると判断します。ストリークを開始する場所がn未満であり、各場所で、長さがcutoffのスジリックの確率はpwin**cutoffです。バインドされた粗組合を使用して、我々はlog(pwin) < 0ので、最後の不等式は方向を反転し

n * pwin**cutoff <= abserr/n 
pwin**cutoff <= abserr/n**2 
log(pwin) * cutoff <= log(abserr/n**2) 
cutoff >= log(abserr/n**2, pwin), 

を必要としています。

この改訂版のコードは、品質の低い評価者(2014年のハードウェアのPython 3インタプリタ)にもかかわらず、3秒未満で実行されます。

import math 


def avglongwinstreak(n, pwin, abserr=0): 
    if n > 0 and pwin < 1 and abserr > 0: 
     cutoff = math.log(abserr/n**2, pwin) 
    else: 
     cutoff = n + 1 
    dist = {(0, 0): 1} 
    for i in range(n): 
     nextdist = {} 
     for (m, k), pmk in dist.items(): 
      winkey = (max(m, k + 1), k + 1) 
      if winkey[0] < cutoff: 
       nextdist[winkey] = nextdist.get(winkey, 0) + pmk * pwin 
      losskey = (m, 0) 
      nextdist[losskey] = nextdist.get(losskey, 0) + pmk * (1 - pwin) 
     dist = nextdist 
    return sum(m * pmk for ((m, k), pmk) in dist.items()) 


print(avglongwinstreak(3000, 0.4, 1e-6)) 
+0

うわー、ありがとうこの素晴らしい答えのために! – ibiza

+0

これは私のものより優れていますが、まだN = 350周りには遅すぎています...これは難しい問題です。 – ibiza

+0

おそらく、完全な木のいくつかの枝をスキップしてその結果を推測するために使用できるoptmization技法があります。 – ibiza