# How can I model flips until N successes?

You and I decide to play a game where we take turns flipping a coin. The first player to flip 10 heads in total wins the game. Naturally, there is an argument about who should go first.

Simulations of this game show that the player to flips first wins 6% more than the player who flips second (the first player wins approx 53% of the time). I'm interested in modelling this analytically.

This isn't a binomial random variable, as there are no fixed number of trials (flip until someone gets 10 heads). How can I model this? Is it the negative binomial distribution?

So as to be able to recreate my results, here is my python code:

import numpy as np
from numba import jit

@jit
def sim(N):

P1_wins = 0
P2_wins = 0

for i in range(N):

while True:

P1_wins+=1
break

P2_wins+=1
break
return P1_wins/N, P2_wins/N

a,b = sim(1000000)

• When you toss a coin until $r$ failures and then look at the distribution of the number of successes that happen before finishing such experiment, then this is by definition negative binomial distribution.
– Tim
Jan 20, 2017 at 17:26
• I cannot reproduce the 2% value. I find that the first player wins $53.290977425133892\ldots\%$ of the time.
– whuber
Jan 20, 2017 at 17:34
• @whuber yes, I believe you are right. I ran my simulation fewer times than I should. My results are commensurate with yours. Jan 20, 2017 at 17:35
• If one wins 53% of the time, the other should be 47%, so shouldn't the description read "the first player wins 6% more than the second player," or "3% more than half the time" ? Not (as it currently says) "3% more than the player who flips second" Jan 20, 2017 at 21:49
• Did you get this question from the FiveThirtyEight Riddler Express? Jan 20, 2017 at 23:11

The distribution of the number of tails before achieving $10$ heads is Negative Binomial with parameters $10$ and $1/2$. Let $f$ be the probability function and $G$ the survival function: for each $n\ge 0$, $f(n)$ is the player's chance of $n$ tails before $10$ heads and $G(n)$ is the player's chance of $n$ or more tails before $10$ heads.

Because the players roll independently, the chance the first player wins with rolling exactly $n$ tails is obtained by multiplying that chance by the chance the second player rolls $n$ or more tails, equal to $f(n)G(n)$.

Summing over all possible $n$ gives the first player's winning chances as

$$\sum_{n=0}^\infty f(n)G(n) \approx 53.290977425133892\ldots\%.$$

That is about $3\%$ more than half the time.

In general, replacing $10$ by any positive integer $m$, the answer can be given in terms of a Hypergeometric function: it is equal to

$$1/2 + 2^{-2m-1} {_2F_1}(m,m,1,1/4).$$

When using a biased coin with a chance $p$ of heads, this generalizes to

$$\frac{1}{2} + \frac{1}{2}(p^{2m}) {_2F_1}(m, m, 1, (1 - p)^2).$$

Here is an R simulation of a million such games. It reports an estimate of $0.5325$. A binomial hypothesis test to compare it to the theoretical result has a Z-score of $-0.843$, which is an insignificant difference.

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

• Just as a note that may not be obvious at a glance, our answers agree numerically: (.53290977425133892 - .5) * 2 is essentially exactly the probability I gave. Jan 20, 2017 at 17:48
• @Dougal Thank you for pointing that out. I looked at your answer, saw the $6.6\%$, and knowing that it did not agree with the form of the answer requested in the question, I did not recognize that you had computed correctly. In general it's a good idea to frame an answer to any question in the form that is requested, if possible: that makes it easy to recognize when it's correct and easy to compare answers.
– whuber
Jan 20, 2017 at 17:56
• @whuber I was responding to the phrase "Simulations of this game show that the player to flips first wins 2% (EDIT: 3% more after simulating more games) more than the player who flips second". I'd interpret "wins 2% more" as $\Pr(A\text{ wins}) - \Pr(B\text{ wins}) = 2\%$; the correct value is indeed 6.6%. I'm not sure of a way to interpret "wins 2% more" means "wins 52% of the time", though apparently that is what was intended. Jan 20, 2017 at 17:59
• @Dougal I agree that the OP's description is confusing and even wrong. However, the code and his result made it clear he meant "3% more than half the time" rather than "3% more than the other player."
– whuber
Jan 20, 2017 at 18:06
• @whuber Agreed. Unfortunately, I answered the question before the code was posted, and didn't run a simulation myself. :) Jan 20, 2017 at 18:06

We can model the game like this:

• Player A flips a coin repeatedly, getting results $$A_1, A_2, \dots$$ until they get a total of 10 heads. Let the time index of the 10th heads be the random variable $$X$$.
• Player B does the same. Let the time index of the 10th heads be the random variable $$Y$$, which is an iid copy of $$X$$.
• If $$X \le Y$$, Player A wins; otherwise Player B wins. That is, \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}

The gap in the win rates is thus $$\Pr(X = Y) = \sum_k \Pr(X = k, Y = k) = \sum_k \Pr(X = k)^2 .$$

As you suspected, $$X$$ (and $$Y$$) are distributed essentially according to a negative binomial distribution. Notations for this vary, but in Wikipedia's parameterization, we have heads as a "failure" and tails as a "success"; we need $$r = 10$$ "failures" (heads) before the experiment is stopped, and success probability $$p = \tfrac12$$. Then the number of "successes," which is $$X - 10$$, has $$\Pr(X - 10 = k) = \binom{k + 9}{k} 2^{-10 - k},$$ and the collision probability is $$\Pr(X = Y) = \sum_{k=0}^\infty \binom{k + 9}{k}^2 2^{-2k - 20} ,$$ which Mathematica helpfully tells us is $$\frac{76\,499\,525}{1\,162\,261\,467} \approx 6.6\%$$.

Thus Player B's win rate is $$\Pr(Y > X) \approx 46.7\%$$, and Player A's is $$\frac{619\,380\,496}{1\,162\,261\,467} \approx 53.3\%$$.

• the heads need not be in a row, just 10 total. I assume that is what you are fixing. Jan 20, 2017 at 17:08
• (+1) I like this approach better than the one I posted because it is computationally simpler: it requires only the probability function, which has a simple expression in terms of binomial coefficients.
– whuber
Jan 20, 2017 at 18:04
• I've submitted an edit replacing the last paragraph questioning the difference from the other answer with an explanation of how their results are actually the same. Jan 20, 2017 at 21:49

Let $$E_{ij}$$ 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 $$p_{ij} \equiv Pr(E_{ij})$$.

Then $$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.