Comment puis-je modéliser des flips jusqu'à N succès?

17

Vous et moi décidons de jouer à un jeu où nous lançons à tour de rôle une pièce. Le premier joueur à retourner 10 têtes au total remporte la partie. Naturellement, il y a un débat sur qui devrait aller en premier.

Les simulations de ce jeu montrent que le joueur qui retourne en premier gagne 6% de plus que le joueur qui retourne en deuxième (le premier joueur gagne environ 53% du temps). Je suis intéressé par la modélisation analytique.

Ce n'est pas une variable aléatoire binomiale, car il n'y a pas de nombre fixe d'essais (retournez jusqu'à ce que quelqu'un obtienne 10 têtes). Comment puis-je modéliser cela? Est-ce la distribution binomiale négative?


Afin de pouvoir recréer mes résultats, voici mon code 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)
Demetri Pananos
la source
3
Lorsque vous lancez une pièce jusqu'à ce que r échoue, puis regardez la distribution du nombre de succès qui se produisent avant de terminer une telle expérience, alors c'est par définition une distribution binomiale négative .
Tim
2
Je ne peux pas reproduire la valeur de 2%. Je trouve que le premier joueur gagne du temps. 53.290977425133892%
whuber
1
@whuber oui, je crois que vous avez raison. J'ai exécuté ma simulation moins de fois que je n'aurais dû. Mes résultats sont proportionnels aux vôtres.
Demetri Pananos
1
Si l'un gagne 53% du temps, l'autre devrait être de 47%, alors la description ne devrait-elle pas indiquer "le premier joueur gagne 6% de plus que le deuxième joueur" ou "3% plus de la moitié du temps"? Pas (comme il est dit actuellement) "3% de plus que le joueur qui retourne deuxième"
JesseM
3
Avez-vous reçu cette question du FiveThirtyEight Riddler Express ?
foutandabout

Réponses:

19

La distribution du nombre de queues avant d' atteindre la tête est binomiale négative avec les paramètres 10 et une / deux . Soit f la fonction de probabilité et G la fonction de survie: pour chaque n 0 , f ( n ) est la chance du joueur d'avoir n queues avant 10 têtes et G ( n10101/2fGn0f(n)n10 est la chance du joueur d'avoir n ou plusieurs queues avant 10 têtes.G(n)n10

Parce que les joueurs roulent indépendamment, la chance que le premier joueur gagne en roulant exactement queues est obtenue en multipliant cette chance par la chance que le deuxième joueur roule n ou plusieurs queues, égal à f ( n ) G ( n ) .nnf(n)G(n)

La somme de tous les possibles donne au premier joueur des chances de gagnern

n=0f(n)G(n)53.290977425133892%.

C'est environ plus que la moitié du temps.3%

En général, remplacer par tout entier positif m , la réponse peut être donnée en fonction d'une fonction hypergéométrique: elle est égale à10m

1/2+22m12F1(m,m,1,1/4).

Lorsque vous utilisez une pièce biaisée avec une chance de têtes, cela se généralise pourp

12+12(p2m)2F1(m,m,1,(1p)2).

Voici une Rsimulation d'un million de ces jeux. Il rapporte une estimation de . Un test d'hypothèse binomiale pour le comparer au résultat théorique a un score Z de - 0,8430.53250.843 , ce qui est une différence insignifiante.

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))
whuber
la source
1
Tout comme une note qui peut ne pas être évidente d'un coup d'œil, nos réponses s'accordent numériquement: (.53290977425133892 - .5) * 2 est essentiellement exactement la probabilité que j'ai donnée.
Dougal
1
@Dougal Merci de l'avoir signalé. J'ai regardé votre réponse, j'ai vu les , et sachant qu'elle n'était pas d'accord avec la forme de la réponse demandée dans la question, je n'ai pas reconnu que vous aviez calculé correctement. En général, c'est une bonne idée d'encadrer une réponse à n'importe quelle question dans le formulaire qui est demandé, si possible: cela permet de reconnaître facilement quand c'est correct et de comparer facilement les réponses. 6.6%
whuber
1
@whuber Je répondais à la phrase "Les simulations de ce jeu montrent que le joueur qui retourne en premier gagne 2% (EDIT: 3% de plus après avoir simulé plus de jeux) plus que le joueur qui retourne en second". J'interpréterais "gagne 2% de plus" comme ; la valeur correcte est en effet de 6,6%. Je ne suis pas sûr d'un moyen d'interpréter "gagne 2% de plus" signifie "gagne 52% du temps", même si apparemment c'est ce qui était prévu. Pr(A wins)Pr(B wins)=2%
Dougal
@Dougal Je suis d'accord que la description de l'OP est déroutante et même erronée. Cependant, le code et son résultat indiquaient clairement qu'il voulait dire "3% plus de la moitié du temps" plutôt que "3% de plus que l'autre joueur".
whuber
1
@whuber D'accord. Malheureusement, j'ai répondu à la question avant la publication du code et je n'ai pas exécuté de simulation moi-même. :)
Dougal
15

Nous pouvons modéliser le jeu comme ceci:

  • Le joueur A lance une pièce à plusieurs reprises, obtenant les résultats A1,A2, jusqu'à ce qu'ils obtiennent un total de 10 têtes. Soit l'indice de temps des 10èmes têtes la variable aléatoire X .
  • Le joueur B fait de même. Soit l'indice de temps des 10èmes têtes la variable aléatoire Y , qui est une copie iid de X .
  • Si XY , le joueur A gagne; sinon le joueur B gagne. Autrement dit,
    Pr(A wins)=Pr(XY)=Pr(X>Y)+Pr(X=Y)Pr(B wins)=Pr(Y>X)=Pr(X>Y).

L'écart dans les taux de gain est donc

Pr(X=Y)=kPr(X=k,Y=k)=kPr(X=k)2.

Comme vous le soupçonniez, X (et Y ) sont distribués essentiellement selon une distribution binomiale négative. Les notations varient, mais dans le paramétrage de Wikipedia , nous avons les têtes comme un "échec" et les queues comme un "succès"; nous avons besoin de r=10 "échecs" (têtes) avant l'arrêt de l'expérience, et probabilité de réussite p=12 . Alors le nombre de "succès", qui estX10, a

Pr(X10=k)=(k+9k)210k,
et la probabilité de collision est
Pr(X=Y)=k=0(k+9k)222k20,
ce que Mathematica nous dit utilement est7649952511622614676.6%.

Ainsi, le taux de victoire du joueur B est Pr(Y>X)46.7% , et le joueur A est 619380496116226146753.3%.

Dougal
la source
les têtes n'ont pas besoin d'être alignées, seulement 10 au total. Je suppose que c'est ce que vous corrigez.
Demetri Pananos
6
(+1) J'aime mieux cette approche que celle que j'ai publiée car elle est plus simple sur le plan informatique: elle ne nécessite que la fonction de probabilité, qui a une expression simple en termes de coefficients binomiaux.
whuber
1
J'ai soumis une modification remplaçant le dernier paragraphe remettant en question la différence de l'autre réponse par une explication de la façon dont leurs résultats sont réellement les mêmes.
Monty Harder
1

Let Eij be the event that the player on roll flips i heads before the other player flips j heads, and let X be the first two flips having sample space {hh,ht,th,tt} where h means heads and t tails, and let pijPr(Eij).

Then pij=Pr(Ei1j1|X=hh)Pr(X=hh)+Pr(Ei1j|X=ht)Pr(X=ht)+Pr(Eij1|X=th)Pr(X=th)+Pr(Eij|X=tt)Pr(X=tt)

Assuming a standard coin Pr(X=)=1/4 means that pij=1/4[pi1j1+pi1j+pij1+pij]

solving for pij, =1/3[pi1j1+pi1j+pij1]

But p0j=p00=1 and pi0=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(ij) 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 Ekl 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, pk,l=11/2[pl,k+1+pl,0] where pil=pii=1,pki=0

This is a linear system with i2 unknowns and one unique solution.

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

pk,l,n+1=1/(1+ϵ)[ϵpk,l,n+11/2(pl,k+1,n+pl,0,n)]

Choose ϵ and pk,l,0 wisely and run the iteration for a few steps and monitor the correction term.

John Rambo
la source