Как получить подробные результаты (вероятность дерева) выполнений Бернулли эксперимент большого числа раз

Question

, предполагающих следующий эксперимент: Выполнить же испытание Бернулли (с вероятностью успеха P) N количества раз

I нужна следующая информация: все возможные последовательности успеха/неудачи с его вероятностью.

Пример: Бернулли эксперимент с вероятностью успеха P = 40% выполненный в 3 раза бы дают следующие результаты (S имеет успех, Р сбой):

FFF 0,216

SFF 0,144

FSF 0,144

SSF 0,096

ФФС 0,144

SFS 0,096

ФСС 0,096

SSS 0,064

Я пытался BruteForce его, чтобы получить результаты, но он быстро дроссели только с 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); 
     } 
    } 
} 

мне нужно, чтобы быть в состоянии поддерживать до N = 3000 своевременно (получение результата менее чем за 3 секунды для любого P)

Есть математический способ сделать это оптимально?

Answer 1

Вот другой подход, который является точными и достаточно быстро быть только квадратичными. Ожидаемое значение самой длинной выигрышной полосы равно

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

Мы рассуждаем о вероятности следующим образом. Либо запись начинается с выигрышной полосы k (вероятность pwin**k), либо открывается j выигрышей за j in 0..k-1, за которыми следует потеря (вероятность pwin**j * (1 - pwin)), на которой условие равно вероятности длины - подряд в n - (j + 1) пытается. Мы используем memoization для оценки повторяемости, которую эта логика подразумевает в 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)) 

Answer 2

Поскольку вы заинтересованы в длине самой длинной победной полосы, в любой момент испытаний есть только два важных факта об истории: (1) как долго самая длинная победная полоса до сих пор (м) и (2), как долго текущая выигрышная полоса (k). Начальное состояние (0, 0). Переходы на выигрыш (m, k) -> (max (m, k + 1), k + 1) и (m, k) -> (m, 0) на убыль. Как только мы узнаем о вероятности всех конечных состояний, мы можем усреднить.

EDIT: в пересмотренной версии этого кода есть оптимизация, чтобы сократить вычисления, необходимые для очень маловероятных событий. В частности, мы игнорируем все состояния с полосой длины, большей или равной некоторым cutoff. Учитывая абсолютный бюджет ошибки abserr, мы определяем, что мы можем исключить массу вероятности до abserr/n из расчета конечных ожидаемых значений, так как самая длинная полоса может быть равна n. Минимальное количество для начало строки в столбцах: n и где бы то ни было, есть вероятность, что длина строки cutoff равна pwin**cutoff. Использование сырого союза, связанный, нам нужно

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

где последнее неравенство переворачивает направление, потому что log(pwin) < 0.

Этот пересмотренный код работает менее чем за три секунды, несмотря на плохой оценщик качества (т. Е. Интерпретатор Python 3 на оборудовании эпохи 2014 года).

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))