Open Side Menu Go to the Top
Register
First to... probability question First to... probability question

07-23-2020 , 09:57 AM
Two players take turns rolling an independent die.
Player A gets 1 point if he rolls a 1/2/3/4
Player B gets 1 point if he rolls a 1/2/3


Questions:
% chance Player A is ahead after 10 rolls/20 rolls/100 rolls etc
% chance Player A reaches 10 points/20 points/100 points etc first


Player A's success rate is 66%. Player B's success rate is 50%. Am I right in saying A has a 2/1 edge?

Player A should get 1.33x as many points as Player B.

1) both miss = .33*.5=.165
2) both score = .66*.5=.33
3) A score, B miss = .66*.5=.33
4) A misses, B scores = .33*.5=.165

edge comes from 3) and 4) as scores remain same in 1) & 2). 33/16.5 = 2/1


How does variance work in all this? Obviously Player A is in better shape if it's first to 1000 points or who is ahead after 1000 rolls vs 10 points/10 rolls.

Thanks for any answers.
First to... probability question Quote
07-23-2020 , 03:16 PM
Quote:
% chance Player A is ahead after 10 rolls/20 rolls/100 rolls etc
Let's start with 10 rolls.

Let p = (2/3)(1/2)
Let q = (1/3)(1/2)
Let r = 1-p-q

I'll set up a transition matrix such that element (j,k) is the probability of transitioning from a Player A lead of 7-j to a lead of 7-k.

[1 0 0 ...]
[p r q 0 0 ...]
[0 p r q 0 0 ...]
[0 0 p r q 0 0 ...]
.
.
.
[... 0 0 1]

That's 13x13: one row & column for each Player A lead from +6 to -6. There are really 21 possible leads, but a lead of 6 is insurmountable when there are only 10 rolls, so we simplify by lumping all leads >5 together as one "state". The probability of remaining in that state is 1.

When we raise that matrix to the power of 10, the sum of elements (7,1) through (7,6) is the probability that Player A is ahead after 10 rolls when starting tied.

It's a nicely organized shorthand for brute-force adding up each possibility's probability one by one. If you were to start writing out the matrix multiplications by hand like a sicko, you'd see what I mean.

Plugging it into Julia:
Code:
m = fill(0.0, 13,13)
m[1,1] = m[13,13] = 1.0
for j=2:12 m[j, (j-1):(j+1)]=[p r q] end
sum((m^10)[7, 1:6])
which spits out 0.7069290440989687

Generalizing to n rolls, the dimension will be 2*floor(n/2 + 1) + 1

Code:
function leadafter(p1::Float64, p2::Float64, rolls::Int64)
    p = p1*(1-p2)
    q = (1-p1)p2
    r = 1-p-q
    leads = floor(Int64, rolls/2 + 1)
    len = 2*leads+1
    m = fill(0.0, len,len)
    m[1,1] = m[len,len] = 1.0
    for j=2:(len-1) m[j, (j-1):(j+1)]=[p r q] end
    return sum((m^rolls)[leads+1, 1:leads])
end
(Filling in the values row-wise like that instead of column-wise costs a few microseconds, but I was lazy.)

Results:
leadafter(2/3, 1/2, 10) = 0.706929044098969
leadafter(2/3, 1/2, 20) = 0.8220946467880509
leadafter(2/3, 1/2,100) = 0.9903563123093899

In my next post I'll answer the 2nd version of the question with a matrix as well.
First to... probability question Quote
07-23-2020 , 04:30 PM
Quote:
Originally Posted by heehaww
Let's start with 10 rolls.

Let p = (2/3)(1/2)
Let q = (1/3)(1/2)
Let r = 1-p-q

I'll set up a transition matrix such that element (j,k) is the probability of transitioning from a Player A lead of 7-j to a lead of 7-k.

[1 0 0 ...]
[p r q 0 0 ...]
[0 p r q 0 0 ...]
[0 0 p r q 0 0 ...]
.
.
.
[... 0 0 1]

That's 13x13: one row & column for each Player A lead from +6 to -6. There are really 21 possible leads, but a lead of 6 is insurmountable when there are only 10 rolls, so we simplify by lumping all leads >5 together as one "state". The probability of remaining in that state is 1.

When we raise that matrix to the power of 10, the sum of elements (7,1) through (7,6) is the probability that Player A is ahead after 10 rolls when starting tied.

It's a nicely organized shorthand for brute-force adding up each possibility's probability one by one. If you were to start writing out the matrix multiplications by hand like a sicko, you'd see what I mean.

Plugging it into Julia:
Code:
m = fill(0.0, 13,13)
m[1,1] = m[13,13] = 1.0
for j=2:12 m[j, (j-1):(j+1)]=[p r q] end
sum((m^10)[7, 1:6])
which spits out 0.7069290440989687

Generalizing to n rolls, the dimension will be 2*floor(n/2 + 1) + 1

Code:
function leadafter(p1::Float64, p2::Float64, rolls::Int64)
    p = p1*(1-p2)
    q = (1-p1)p2
    r = 1-p-q
    leads = floor(Int64, rolls/2 + 1)
    len = 2*leads+1
    m = fill(0.0, len,len)
    m[1,1] = m[len,len] = 1.0
    for j=2:(len-1) m[j, (j-1):(j+1)]=[p r q] end
    return sum((m^rolls)[leads+1, 1:leads])
end
(Filling in the values row-wise like that instead of column-wise costs a few microseconds, but I was lazy.)

Results:
leadafter(2/3, 1/2, 10) = 0.706929044098969
leadafter(2/3, 1/2, 20) = 0.8220946467880509
leadafter(2/3, 1/2,100) = 0.9903563123093899

In my next post I'll answer the 2nd version of the question with a matrix as well.
very cool. Great explanation, thanks.
First to... probability question Quote
07-23-2020 , 05:03 PM
Actually there's a simpler way to solve that 1st problem. My bad, I should have considered my options before jumping to the matrix solution.

This is just a sum of a product of binomial distributions:

P(ahead after n games) = sum of [P(k ties) * P(win the majority of the non-tied games)]
from k=0 to n-1

P(k ties) = C(n,k) / 2^n
P(win majority of non-tied games) = 1 - binomcdf(n-k, floor((n-k)/2 + 1), 2/3)

We plug p=2/3 into the Binomial CDF because we're using the probability of Player A increasing their lead given that somebody increases their lead, which is 1/3 / (1/3 + 1/6).

The results of that agree with my matrix method.
First to... probability question Quote
07-23-2020 , 09:37 PM
Btw, above, the chance of Player B finishing ahead is not 1-P(A) because there is a chance of a tie. But the formula I showed can handle that. Or with the matrix, P(finish tied) is the element in the center of the exponentiated matrix.

Quote:
Originally Posted by Jam-Fly
% chance Player A reaches 10 points/20 points/100 points etc first
This, too, is a combination of binomial distributions.

In this variation of the problem, there is unlimited time, so the first thing to note is that the rolls where neither player scores are irrelevant. Therefore, we'll pretend they don't exist and adjust the probabilities accordingly.

Raw probabilities:
P(A score, B no) = 1/3
P(A no, B yes) = 1/6
P(both) = 1/3
P(neither) = 1/6

Adjusted:
P(A and !B) = (1/3) / (5/6) = 2/5
P(!A and B) = 1/5
P(both) = 1-2/5-1/5 = 2/5

Now that scoring ties are the only ties possible, there can be at most N ties (in which case I assume neither player would win).

P(Player A wins) = sum of [P(k ties) * P(A wins at least n-k out of 2(n-k)-1 of the non-tied rolls)]
from k=0 to n-1

For instance if there are 0 ties, there can be at most 19 rolls because after 18, the score would be 9-9 and the next roll wins. Even though the game can end before 19 rolls, the chance of reaching 10 pts first is the same as the chance of winning the majority of 19 rolls in any order. (That might not be obvious at first, but rather than ruin the fun by explaining it now, I'll give you a chance to ponder it until the "Aha!" moment comes.)

So P(A wins | 0 ties) = 1-binomcdf(19, 9, p)

The p there has to be adjusted a 2nd time: p = P(A and !B) / P(!(A and B)) = 2/3

Crunching the numbers:
P(A 1st to 10) = 0.8723365794771935
P(A 1st to 20) = 0.9490535206605235
P(A 1st to 100) = 0.9998936433783479

(I can share the code if anyone wants.)

Since I showed the matrix method for the other problem, and since it's qualitatively different for this problem, I'll show it next post. It will be a good way of checking the above results.

Last edited by heehaww; 07-23-2020 at 09:48 PM.
First to... probability question Quote
07-27-2020 , 10:13 PM
#2 the matrix way. There is infinite time but 3 absorbing states (as I interpret it), meaning 3 states where the probability of staying is 1. It is possible to find the limit of the transition matrix raised to the # of rolls as they approach infinity, known as the limiting matrix. Our answer will be encoded in the bottom-left corner of the limiting matrix. The top-left corner will be a 3x3 identity matrix and the other two corners will be all 0's because after infinite time there is no chance of being anywhere but one of the absorbing states.

To find this limit, we need to set up our transition matrix in a more particular way known as standard form: the absorbing states must all come before the non-absorbing states. In #1 it didn't matter and I had two absorbing states on opposite corners.

The matrix will be much larger (103x103 when playing to 10) because the probability of reaching an absorbing state is a function of not only the score difference like before, but also the current total score. Making matters worse, I don't think there's a way to organize it so that there's a single visual pattern like before.

However, we can still arrange it in a way that's easy to code. For instance after the absorbing states, let's group all the Player A scores together: 0-0, 0-1, 0-2...0-9, 1-0, 1-1, 1-2...1-9, and so on. This way it's easy to create a formula for the score represented by a given row:

A's score = floor[(row-4) / n]
B's score = row - 4 - A*n
where n = victory score

The probabilities we'll need to enter are 2/5, 1/5, 2/5, (2/5 + 2/5), and (2/5 + 1/5).

If the score is (a,b):
There is a 2/5 chance of reaching score (a+1, b), which is column (row+n).
There is a 1/5 chance of reaching score (a, b+1), which is column (row+1).
There is a 2/5 chance of reaching score (a+1, b+1), which is column (row+n+1).

However, we must check that we're not about to leave the matrix. If the score we're about to increment is n-1, then instead of the columns given by those formulas, the value belongs in one of the first 3 columns instead (transitioning to an absorbing state). Except in the bottom row, the value will be either 4/5 or 3/5 because when a player reaches a score of n-1, they can wiin with a 1-1 roll.

Once we have the matrix all set, there are two sub-matrices of interest:
R = all the elements from rows 4 to n and columns 1 to 3.
Q = all the elements to the right of R.

P(A win) = 1st element of (I-Q)-1R
P(B win) = (1,2)
P(tie) = (1,3)

Code:
import LinearAlgebra.I
function firsto(p::Float64, q::Float64, n::Int64)
    r = 1-p-q
    x = p+r
    y = q+r
    len = n^2 + 3
    m = fill(0.0, len, len)
    m[len,1:3] = [p q r]
    for k=4:(len-1)
        a = floor((k-4)/n)
        b = k-4-n*a
        if a < n-1
            m[k, k+n] = p
            if b < n-1
                m[k, k+n+1] = r
                m[k, k+1] = q
            else m[k,2]=y end
        else
            m[k,1] = x
            m[k, k+1] = q
        end
    end
    return (inv(I-m[4:len, 4:len])*m[4:len, 1:3])[1,1:3]
end
One line of code at the end for the actual math, and all the preceding code for setting the matrix. We don't need to set the absorbing states since we're only using the elements below them.

Results:
Code:
julia> firsto(.4, .2, 10)
3-element Array{Float64,1}:
 0.8088217206170584
 0.13836485459260212
 0.05281342479033959

julia> firsto(.4, .2, 20)
3-element Array{Float64,1}:
 0.9068162142710839
 0.07089393110886681
 0.022289854620049475
This disagrees with my binomial formula. 80.88% is correct; something I did in my last post was wrong. When I use that same binomial formula for playing to a score of 2, it gives an answer of 58.67% which is wrong.

We can calculate by hand the case of playing to a score of 2:

P(A win) = p²(1 + 2q) + pr = .4²(1 + 2*.2) + 2*.4² = 54.4%
which the matrix method agrees with.

So next I have to figure out what's wrong with my binomial formula.

Last edited by heehaww; 07-27-2020 at 10:21 PM.
First to... probability question Quote
07-29-2020 , 11:38 AM
Found my mistake:
Quote:
P(Player A wins) = sum of [P(k ties) * P(A wins at least n-k out of 2(n-k)-1 of the non-tied rolls)]
from k=0 to n-1
The problem is that the number of ties doesn't tell us the number of rolls. The Binomial PMF for # of ties assumes n rolls, but then the CDF for # of wins assumes more than n rolls. While that's a valid pretense for the chance of reaching n Bernoulli successes, it's problematic here because we need the adjusted probabilities to balance out. For instance P(non-tie)*P(A win roll | non-tie) = .6(2/3) = .4

If playing to a score of 2, the bad formula is equivalent to:

.6^2 * [2(1/3)(2/3)^2 + (2/3)^2] +
2*.6*.4 * 2/3

The 1/3 in the first line remains 1/3 because the .6 is only squared. That's bad because 1/3 is a probability conditioned on something else, so it needs to be multiplied by the probability of that condition.

Instead would could say: .6³•2(1/3)(2/3)² + .6²(2/3)² + 2•.6•.4•2/3

A general formula for this problem wouldn't be as simple as I led on, because it's not Negative Binomial but I suppose Negative Multinomial. Using a matrix wasn't as overkill as I thought for the "first to _" problem.

When playing to three: p³(1 + 3q + 6q²) + 3rp² + 3²rqp² + 3pr² = .61696

It's easy to think of each specific formula, but may be harder to put in general form, at least I'm too lazy to right now.

Quote:
P(A win) = 1st element of (I-Q)-1R
By the way, if you wanna know the average duration of the game, it's the sum of the 1st row of (I-Q)-1. However, the transition matrix needs to be adjusted: multiply all the elements by 5/6 (or start with the "raw" probabilities to begin with) and then put the value 1/6 (the chance of a 0-0 roll) everywhere along the diagonal from top-left to bottom-right.

Code:
import LinearAlgebra.I
function firsto(p::Float64, q::Float64, r::Float64, n::Int64)
    z = 1-p-q-r
    x = p+r
    y = q+r
    len = n^2 + 3
    m = fill(0.0, len, len)
    m[len,1:3] = [p q r]
    for k=4:(len-1)
        m[k,k] = z
        a = floor((k-4)/n)
        b = k-4-n*a
        if a < n-1
        m[k, k+n] = p
        if b < n-1
            m[k, k+n+1] = r
            m[k, k+1] = q
            else m[k,2]=y end
        else
            m[k,1] = x
            m[k, k+1] = q
        end
    end
    m[len,len] = z
    return (inv(I-m[4:len, 4:len])*m[4:len, 1:3])[1,1:3]
end
As hoped, this gives the same result as the previous matrix. Now for the avg duration we replace the final line with

Code:
return sum(inv(I-m[4:len, 4:len])[1,1:end])
When playing to 10, that says the avg duration is about 14.6 rolls.
First to... probability question Quote
07-29-2020 , 03:39 PM
Alright here's the formula. Let n = score played to, and j = # of ties.

Σj=0 to n-1{C(n,j)rj * p(n-j) * Σk=0 to n-j-1[C(n+k-1, k)*qk]}

Code:
function firsto(p::Float64, q::Float64, n::Int64)
    r = 1-p-q
    psum = 0.0
    for j=0:(n-1)
        prd = binomial(n,j)*r^j * p^(n-j)
        innersum = 1.0
        for k=1:(n-j-1) innersum += binomial(n+k-1, k)*q^k end
        psum += prd*innersum
    end
    return psum
end
Code:
julia> firsto(.4, .2, 2)
0.544

julia> firsto(.4, .2, 3)
0.6169600000000001

julia> firsto(.4, .2, 10)
0.8088217206170585

julia> firsto(.4, .2, 20)
0.9068162142710838

julia> firsto(.4, .2, 100)
0.9990420337374744
First to... probability question Quote

      
m