Como posso modelar flips até N ter sucesso?

Você e eu decidimos jogar um jogo em que revezamos o lançamento de uma moeda. O primeiro jogador a virar 10 cabeças no total ganha o jogo. Naturalmente, há uma discussão sobre quem deve ir primeiro.

Simulações deste jogo mostram que o jogador que vira primeiro ganha 6% a mais do que o jogador que vira o segundo (o primeiro jogador vence aproximadamente 53% do tempo). Estou interessado em modelar isso analiticamente.

Esta não é uma variável aleatória binomial, pois não há um número fixo de tentativas (vire até que alguém ganhe 10 cabeças). Como posso modelar isso? É a distribuição binomial negativa?

Para poder recriar meus resultados, aqui está o meu código python:

import numpy as np
from numba import jit


@jit
def sim(N):

    P1_wins = 0
    P2_wins = 0

    for i in range(N):

        P1_heads = 0
        P2_heads = 0
        while True:

            P1_heads += np.random.randint(0,2)

            if P1_heads == 10:
                P1_wins+=1
                break

            P2_heads+= np.random.randint(0,2)
            if P2_heads==10:
                P2_wins+=1
                break
    return P1_wins/N, P2_wins/N


a,b = sim(1000000)

A distribuição do número de caudas antes de atingir cabeças é negativo binomial com parâmetros 10 e 1 / 2 . Let $10$$10$$1/2$$f$ a função de probabilidade e a função de sobrevivência: para cada n ≥ 0 , f ( n ) é a chance de n caudas antes de 10 cabeças e G ( n ) é a chance de n ou mais caudas antes de 10 cabeças.$G$$n\ge 0$$f(n)$$n$$10$$G(n)$$n$$10$

Como os jogadores rolam independentemente, a chance do primeiro jogador vencer rolando exatamente caudas é obtida multiplicando essa chance pela chance de o segundo jogador rolar n ou mais caudas, igual a f ( n ) G ( n ) .$n$$n$$f(n)G(n)$

A soma de todos os possíveis fornece as chances de vitória do primeiro jogador como$n$

Isso é cerca de mais da metade do tempo.$3\%$

Em geral, substituindo por qualquer número inteiro positivo m , a resposta pode ser dada em termos de uma função hipergeométrica: é igual a$10$$m$

Ao usar uma moeda tendenciosa com uma chance de cabeças, isso generaliza para$p$

Aqui está uma Rsimulação de um milhão desses jogos. Ele relata uma estimativa de . Um teste de hipótese binomial para compará-lo com o resultado teórico tem um escore Z de - 0,843 , que é uma diferença insignificante.$0.5325$$-0.843$

n.sim <- 1e6
set.seed(17)
xy <- matrix(rnbinom(2*n.sim, 10, 1/2), nrow=2)
p <- mean(xy[1,] <= xy[2,])
cat("Estimate:", signif(p, 4), 
    "Z-score:", signif((p - 0.532909774) / sqrt(p*(1-p)) * sqrt(n.sim), 3))

Podemos modelar o jogo assim:

• O jogador A joga uma moeda repetidamente, obtendo os resultados $A_1, A_2, \dots$ até que eles tenham um total de 10 cabeças. Deixe o índice de tempo do 10.º cabeças ser a variável aleatória $X$ .
• O jogador B faz o mesmo. Deixe o índice de tempo dos 10o cabeças ser a variável aleatória $Y$ , que é uma cópia de iid $X$ .
• Se $X \le Y$ , o jogador A vence; caso contrário, o jogador B vence. Ou seja, $$\begin{align} \Pr(A\text{ wins})&= \Pr(X \ge Y) = \Pr(X > Y) + \Pr(X = Y)\\ \Pr(B\text{ wins})&= \Pr(Y > X) = \Pr(X > Y). \end{align}$$

A diferença nas taxas de ganho é, assim,

Como você suspeitava, $X$ (e $Y$ ) são distribuídos essencialmente de acordo com uma distribuição binomial negativa. As notações para isso variam, mas na parametrização da Wikipedia , temos cabeças como "falha" e coroas como "sucesso"; precisamos de $r = 10$ "falhas" (cabeças) antes que o experimento seja interrompido e a probabilidade de sucesso $p = \tfrac12$ . Então o número de "sucessos", que é$X - 10$, tem

Assim, a taxa de vitória do jogador B é $\Pr(Y > X) \approx 46.7\%$ e a do jogador A é $\frac{619\,380\,496}{1\,162\,261\,467} \approx 53.3\%$.

$E_{ij}$$X$$\{ hh,ht,th,tt\}$$p_{ij} \equiv Pr(E_{ij})$

$p_{ij}=Pr(E_{i-1j-1}|X=hh)*Pr(X=hh)+Pr(E_{i-1j}|X=ht)*Pr(X=ht)+Pr(E_{ij-1}|X=th)*Pr(X=th)+Pr(E_{ij}|X=tt)*Pr(X=tt)$

Assuming a standard coin $Pr(X=*)=1/4$ means that $p_{ij}=1/4*[p_{i-1j-1}+p_{i-1j}+p_{ij-1}+p_{ij}]$

solving for $p_{ij}$, $= 1/3*[p_{i-1j-1}+p_{i-1j}+p_{ij-1}]$

But $p_{0j}=p_{00}=1$ and $p_{i0}=0$, implying that the recursion fully terminates. However, a direct naive recursive implementation will yield poor performance because the branches intersect.

An efficient implementation will have complexity $O(i*j)$ and memory complexity $O(min(i,j))$. Here's a simple fold implemented in Haskell:

Prelude> let p i j = last. head. drop j $ iterate ((1:).(f 1)) start where
  start = 1 : replicate i 0;
  f c v = case v of (a:[]) -> [];
                    (a:b:rest) -> sum : f sum (b:rest) where
                     sum = (a+b+c)/3 
Prelude> p 0 0
1.0
Prelude> p 1 0
0.0
Prelude> p 10 10
0.5329097742513388
Prelude> 

UPDATE: Someone in the comments above asked whether one was suppose to roll 10 heads in a row or not. So let $E_{kl}$ be the event that the player on roll flips i heads in a row before the other player flips i heads in a row, given that they already flipped k and l consecutive heads respectively.

Proceeding as before above, but this time conditioning on the first flip only, $p_{k,l} = 1-1/2*[p_{l,k+1}+p_{l,0}]$ where $p_{il}=p_{ii}=1, p_{ki}=0$

This is a linear system with $i^2$ unknowns and one unique solution.

To convert it into an iterative scheme, simply add an iterate number $n$ and a sensitivity factor $\epsilon$:

$p_{k,l,n+1} = 1/(1+\epsilon)*[\epsilon*p_{k,l,n} +1-1/2*(p_{l,k+1,n}+p_{l,0,n})]$

Choose $\epsilon$ and $p_{k,l,0}$ wisely and run the iteration for a few steps and monitor the correction term.