Probability of getting certain score in Drop Dead

In the dice game of Drop Dead, https://en.wikipedia.org/wiki/Drop_Dead_(dice_game), what is the probability of getting a score of 0? 1? 10? 20?

The easiest score to work out is zero. To get zero, the player must roll at least one 2 or 5 on every roll. If we define a Markov process, with states 0, F, 5, 4, 3, 2, 1 (in that order), where 0 means a score of zero has been achieved, $F$ means a score greater than $0$ has been achieved, and $5$ through $1$ are the current number of dice being rolled, we have the following transition matrix:

$$P_0 = \pmatrix{1&0&0&0&0&0&0\\0&1&0&0&0&0&0\\\frac{1}{243}&\frac{32}{243}&0&\frac{80}{243}&\frac{80}{243}&\frac{40}{243}&\frac{10}{243} \\\frac{1}{81}&\frac{16}{81}&0&0&\frac{32}{81}&\frac{8}{27}&\frac{8}{81}\\\frac{1}{27}&\frac{8}{27}&0&0&0&\frac{4}{9}&\frac{2}{9}\\\frac{1}{9}&\frac{4}{9}&0&0&0&0&\frac{4}{9}\\\frac{1}{3}&\frac{2}{3}&0&0&0&0&0}$$

Why do we need the state $F$?

Solution 1:

Relabel the states according to the row and column indices of the transition matrix. Let $\ S_t\ $ be the state at time $\ t\ $, $\ T\ $ the first time that the chain enters one of the absorbing states $\ s=1\ $ or $\ s=2\ $, and $$ \pi_i=P\big(S_T=1\,\big|\,S_t=i\big)\ , $$ the probability, given that the chain is in state $\ i\ $ at time $\ t\ $, that it will end in the state where a final score of zero has been obtained. Then $\ \pi_3\ $ is the probability of obtaining a zero score from the starting state, and $\ \pi\ $ must satisfy \begin{align} \pi_1&=1\\ \pi_2&=0\\ \pi_i&=\sum_{j=1}^7p_{ij}\pi_j\ , \end{align} where $\ p_{ij}\ $ is the entry in row $\ i\ $ and column $\ j\ $ of the transition matrix $\ P_0\ $. Or, in matrix form, $$ \pmatrix{1\\0\\\pi_3\\\pi_4\\\pi_5\\\pi_6\\\pi_7} =\pmatrix{1&0&0&0&0&0&0\\0&1&0&0&0&0&0\\\frac{1}{243}&\frac{32}{243}&0&\frac{80}{243}&\frac{80}{243}&\frac{40}{243}&\frac{10}{243} \\\frac{1}{81}&\frac{16}{81}&0&0&\frac{32}{81}&\frac{8}{27}&\frac{8}{81}\\\frac{1}{27}&\frac{8}{27}&0&0&0&\frac{4}{9}&\frac{2}{9}\\\frac{1}{9}&\frac{4}{9}&0&0&0&0&\frac{4}{9}\\\frac{1}{3}&\frac{2}{3}&0&0&0&0&0}\pmatrix{1\\0\\\pi_3\\\pi_4\\\pi_5\\\pi_6\\\pi_7}\ .$$ The final row of this equation gives $\ \pi_7=\frac{1}{3}\ $, and then successively substituting back through the preceding rows gives \begin{align} \pi_6=&\frac{1}{9}+\frac{4\pi_7}{9}=\frac{7}{27}\\ \pi_5=&\frac{1}{27}+\frac{4\pi_6}{9}+\frac{2\pi_7}{9}=\frac{55}{243}\\ \pi_4=&\frac{1}{81}+\frac{32\pi_5}{81}+\frac{8\pi_6}{27}+\frac{8\pi_7}{81}=\frac{4163}{3^9}\\ \pi_3=&\frac{1+80\pi_4+80\pi_5+40\pi_6+10\pi_7}{243}=\frac{978853}{3^{14}}\approx0.205\ . \end{align} This last probability , $\ \pi_3=\frac{978853}{3^{14}}\ $, or about $\ \frac{1}{5}\ $ , is the probability of scoring zero in Drop Dead with $5$ dice.

Finding the probability of any other score could be handled in similar way, but the states of the appropriate Markov chain would need to incorporate the current score as one of its components, and it looks to me like the process will become quite unwieldy as the score increases. However, there's a recursive procedure by which the distributions of scores for $1$-, $2$-, $3$-, $4$- and $5$-dice versions of Drop Dead can be successively calculated without too much trouble. If you're playing $\ i$-dice Drop Dead, then on each throw you will either get a score $\ A_i\ $ of somewhere between $\ i\ $ and $\ 6i\ $, in which case you will continue playing with $\ i\ $ dice, or you will get a score $\ A_i\ $ of zero, in which case your remaining score $\ S\ $ will be obtained by playing the $\ J_i$-dice version of Drop Dead for some random $\ J_i=0,1,\dots,$$i-1\ $. Let $$ d_{is}=P_{DD_i}(S=s) $$ be the probability of achieving a score $\ s\ $ in Drop Dead with $\ i\ $ dice, \begin{align} a_{is}&=P\big(A_i=s\big)\ \ \text{ for }\ \ s>0\ \text{ and}\\ b_{ij}&=P\big(A_i=0, J_i=j\big)\ \ \text{ for }\ \ 0\le j\le i-1\ . \end{align} Calculation of $\ a_{is}\ $ and $\ b_{ij}\ $ are both straightforward, and \begin{align} d_{0s}&=\cases{1&if $\ s=0$\\ 0&otherwise}\ ,\\ d_{is}&=\sum_{k=i}^s a_{ik}d_{i\,s-k}+\sum_{j=0}^{i-1}b_{ij}d_{j\,s}\\ &=\sum_{k=0}^{s-i} a_{i\,s-k}d_{ik}+\sum_{j=0}^{i-1}b_{ij}d_{j\,s}\ \ \text{ for }\ \ i\ge1\ . \end{align} Given the known values of $\ d_{0s}\ $, this formula can be used to successively calculate those of $\ d_{1s}\ $, $\ d_{2s}\ $, $\ d_{3s}\ $, $\ d_{4s}\ $, and $\ d_{5s}\ $.

Calculation of $\ b_{ij}\ $ and $\ a_{is}\ $

The formula for $\ b_{ij}\ $ is \begin{align} b_{ij}&={i\choose i-j}\left(\frac{1}{3}\right)^{i-j}\left(\frac{2}{3}\right)^j\\ &={i\choose j}\frac{2^j}{3^i}\ , \end{align} which is just the probability that $\ i-j\ $ twos and fives appear when $\ i \ $ dice are thrown.

When $\ i\ $ dice are thrown, the conditional distribution of the face showing on each die, given that no twos or fives were thrown, is given by $$ \nu_1=\nu_3=\nu_4=\nu_6=\frac{1}{4}\\ \nu_i=0\ \ \text{ for }\ \ i\ne 1,3,4\ \text{ or }\ 6\ . $$ The conditional distribution of $\ A_i\ $ (i.e. the sum of the faces) given no twos or fives were thrown is the $\ i$-fold convolution $\ \nu^{*i}\ $ of this distribution. To get $\ a_{is}\ $ we have to multiply this convolution by the probability, $\ \left(\frac{2}{3}\right)^i\ $, that no twos or fives were thrown: $$ a_{is}=\left(\frac{2}{3}\right)^i\nu^{*i}_s\ . $$ Results

I've written a Magma script to perform the calculations described above, and run it in the online Magma calculator to tabulate the values of $\ d_{is}\ $ for $\ 1\le i\le5\ $ and $\ 0\le s\le20\ $. The results are given in the table below, and a copy of the script is given below the table. The entry in column $\ i\ $ and row $\ s\ $ of the table is the probability of achieving a score of $\ s\ $ in Drop Dead when you start with $\ i\ $ dice. Thus, the probabilities of scoring $\ 0,1,10,$ or $\ 20\ $ in $5$-dice Drop Dead are $$ 0.204654, 0.0212758,0.0321336\ \ \text{and}\ \ 0.0199990\ , $$ respectively.

Magma actually performs all the calculations in exact rational arithmetic. For $\ 5$-dice Drop Dead, however, the numerators and denominators of some of the exact probabilities are upwards of $20$ digits long. For the purposes of displaying the results, therefore, I've rounded them to $6$ significant digits. The exact values of the four probabilities quoted above are $$ \frac{978853}{3^{14}}, \frac{305285}{3^{15}}, \frac{2436213835}{2^63^{19}}, \frac{116576288911019013119}{2^{20}3^{33}}, $$ respectively. \begin{array}{c|c|c|c|c|c|} &1&2&3&4&5\\ \hline 0&0.333333&0.259259&0.226337&0.211502&0.204654\\ \hline 1&0.0555556&0.0246913&0.0233196&0.0220156&0.0212758\\ \hline 2&0.00925925&0.0113169&0.00708733&0.00706758&0.00690395\\ \hline 3&0.0570987&0.0260631&0.0253201&0.0233648&0.0226679\\ \hline 4&0.0743313&0.0477538&0.0378499&0.0366068&0.0354321\\ \hline 5&0.0231910&0.0268061&0.0202438&0.0182475&0.0180652 \\ \hline 6&0.0704803&0.0418532&0.0378484&0.0349727&0.0337666 \\ \hline 7&0.0429110&0.0513853&0.0361148&0.0344697&0.0335945 \\ \hline 8&0.0249487&0.0266114&0.0276699&0.0223774&0.0220318 \\ \hline 9&0.0292865&0.0356564&0.0286695&0.0276046&0.0258869 \\ \hline 10&0.0361682&0.0409402&0.0377166&0.0327515&0.0321336 \\ \hline 11&0.0212032&0.0243690&0.0285918&0.0250880&0.0232016 \\ \hline 12&0.0243198&0.0310250&0.0281388&0.0278927&0.0260151 \\ \hline 13&0.0221142&0.0255739&0.0306239&0.0266317&0.0257261 \\ \hline 14&0.0174057&0.0249209&0.0251908&0.0267092&0.0239210 \\ \hline 15&0.0153692&0.0210825&0.0230799&0.0229015&0.0224902 \\ \hline 16&0.0163286&0.0233715&0.0250409&0.0254577&0.0240678 \\ \hline 17&0.0128420&0.0207316&0.0210109&0.0232448&0.0223747 \\ \hline 18&0.0116561&0.0179636&0.0202223&0.0204600&0.0210895 \\ \hline 19&0.0109114&0.0184984&0.0196777&0.0213054&0.0208299 \\ \hline 20&0.00958127&0.0159901&0.0182795&0.0189537&0.0199990\\ \hline \end{array} $$

maxdice:=5;
maxscore:=20;

b:=ZeroMatrix(RationalField(),maxdice,maxdice);
for i in [1..maxdice] do
  for j in [0..i-1] do
     b[i,j+1]:=Binomial(i,j)*2^j/3^i;
  end for;
end for;

convolve:= function(arg1,arg2)
  result:=[];
  for i in [0..#arg1+#arg2-2] do
     r := 0;
     for j in [Max(0,i-#arg2+1)..Min(i,#arg1-1)] do
        r:=r+arg1[j+1]*arg2[i-j+1];
     end for;
     result := result cat [r] ;
   end for;
   return result;
 end function;

 dist :=[0,1/4,0,1/4,1/4,0,1/4];
 a:=ZeroMatrix(RationalField(),maxdice,Max(6*maxdice+1,maxscore+1));

ac:=dist;
for i in [2..maxdice] do
  for s in [0..#ac-1] do
    a[i-1, s+1] := ac[s+1]*(2/3)^(i-1);
  end for;
  ac:= convolve(dist,ac);
end for;
for s in [0..#ac-1] do
   a[maxdice, s+1] := ac[s+1]*(2/3)^maxdice;
end for;

d:= ZeroMatrix(RationalField(),maxdice+1, maxscore+1);
d[1,1]:=1;
for i in [1..maxdice] do
  for s in [0..i-1] do
     d[i+1,s+1]:=0;
     for j in [0..i-1] do
      d[i+1,s+1]:=d[i+1,s+1]+b[i,j+1]*d[j+1,s+1];
     end for;
   end for;
   for s in [i..maxscore] do
     d[i+1,s+1]:=0;
     for k in [i..s] do
        d[i+1,s+1]:=d[i+1,s+1]+a[i,k+1]*d[i+1,s-k+1];
     end for;
     for j in [0..i-1] do
        d[i+1,s+1]:=d[i+1,s+1]+b[i,j+1]*d[j+1,s+1];
     end for;
   end for;
end for;

dr:=ZeroMatrix(RealField(6),maxdice+1,maxscore+1);
for i in [0..maxdice] do
   for s in [0..maxscore] do
     dr[i+1,s+1]:=RealField(6)!d[i+1,s+1];
   end for;
end for;
print Transpose(dr);
print Transpose(d);