"Dad? How many times do I have to roll a die until all six sides appear?"
I stopped what I was doing to consider my son's question. Although I could figure out the answer mathematically, sometimes experiments are more powerful than math equations for showing how probability works.
"Why don't we find out?" I said. "You roll and I'll keep track of the numbers that appear."
Each trial took at least six rolls, but, eventually, we completed 20 trials. The following SAS DATA step contains the number of rolls required until all six faces appeared on the die. The UNIVARIATE procedure displays a histogram of the data:
data RollsUntilAllFaces; input rolls @@; datalines; 6 7 7 7 8 9 9 9 10 10 11 12 12 12 15 16 17 19 20 26 ; proc univariate data=RollsUntilAllFaces; histogram rolls / endpoints=6 to 27 by 3; run; |
The Expected Value
It's not hard to write down the expected number of rolls for a single die. You need one roll to see the first face. After that, the probability of rolling a different number is 5/6. Therefore, on average, you expect the second face after 6/5 rolls. After that value appears, the probability of rolling a new face is 4/6, and therefore you expect the third face after 6/4 rolls. Continuing this process leads to the conclusion that the expected number of rolls before all six faces appear is
6/6 + 6/5 + 6/4 + 6/3 + 6/2 + 6/1 = 14.7 rolls.When I told my son that "you should expect to roll the die 14.7 times before all six faces appear," he was unimpressed.
"That's too high," he retorted. "Most of the time we only need 10–12 rolls."
The Median Value
He and I had previously discussed the mean, median, and mode, so I knew what he meant when he told me that I was "using the wrong average." And I agree with him. The median value in our experiment is 10.5, but the mean value (which is strongly influenced by the trials with values 19, 20, and 26) is 12.1. The long tail in the histogram makes it apparent that the median and mean values will be substantially different.
So what is the expected median number of rolls before all six faces appear? I don't know. The only way that I know to estimate the median is through simulation. If you know a reference or a better way, leave a comment.
The General Problem
This problem is called the "coupon collector's problem." If you roll two dice, it is harder (but still possible) to find the expected number of rolls because the outcomes 2–12 do not have equal probability.
The general formulation is "If you are trying to collect N unique items, and the probability of getting the ith item is pi, what is the expected number of items that you need to obtain before you complete the collection?"
Now I understand why I never completed that 1977 Star Wars trading card set! I kept buying and buying until my jaws hurt from chewing all the gum, but with 66 cards in the set, my quest was essentially futile.
Alas, the coveted "light saber" card shall never be mine: the expected number of cards required to complete the set is 315, assuming all cards have an equal probability.
10 Comments
If you actually simulate it, it seems like the median number of rolls is actually 13 and the lower and upper quartiles are 10 and 18. This is with almost 3.4 million samples (from 50 million rolls)! Anyway, it's closer to the mean (14.703 in the simulation) than I would have thought Also, it's pretty close to the lognormal distribution.
Wow, that was fast! Here's my simulation, along with a discussion of simulating "time-to-event outcomes": http://blogs.sas.com/content/iml/2011/07/22/simulating-the-coupon-collectors-problem/
Interesting problem! I'm sure that the median number of rolls to observe all 6 values has a solution through combinatorics. Consider a simpler problem with a 4-sided die (equilateral triangles on the sides). What is the median number of trials necessary to observe all 4 values?
If we roll the die 4 times, then we have 4*3*2*1 different ways to observe all 4 values out of 4^4 possible rolls. The probability of observing all 4 in 4 rolls is 24/256 which is less than 10%. Thus, we need to keep going.
If we roll the die again, then all of the permutations which resulted in observation of every value in 4 rolls also results in observation of every value in 5 rolls. There are 4*3*2*1*4 such permutations. In addition, we might have observed 3 of the 4 outcomes in the first 4 trials, and observe the final outcome in the 5th trial. This can be observed with repetitions in rolls (1,2), (1,3), (1,4), (2,3), (2,4), (3,4). That is, there are 6 ways that we could observe the 4 outcomes in 5 trials, with the final outcome occurring in the last trial. Note that there are still 4*3*2*1 orders in which we could observe the 4 different outcomes with 6 different combinations of repeated values. So, there are 4*3*2*1*6 permutations which result in observation of all 4 values at the 5-th roll. This is in addition to the 4*3*2*1*4 permutations which produced all 4 values by the 4th roll, with any additional value on the 5th roll. So, with 5 rolls, we have
4*3*2*1*4 + 4*3*2*1*6 = 4*3*2*1*(10) = 240
ways that we could observe all values out of 4^5=1024 possible rolls. The probability of observing all 4 patterns through 5 rolls is, then, just a little under 25%. 5 trials gets us almost to the first quartile.
If we go to 6 rolls, then every success through 4 rolls will be a success after rolls 5 and 6 regardless of what was rolled in the last two trials. This is 4*3*2*1*4*4 ways to achieve success through the first 4 rolls with 2 additional rolls. There are also 4*3*2*1*6*4 ways to achieve success through 5 rolls, with anything allowed on the 6th roll. Alternatively, we could observe one value repeated 3 times, with the other three values observed in the other three trials. Enumerating these possible combinations, we have (1,2,3), (1,2,4), (1,2,5), (1,3,4), (1,3,5), (1,4,5), (2,3,4), (2,3,5), (2,4,5), (3,4,5). Thus, there are 10 combinations of 3 repeat positions out of 5. We could also have pairs of repeats like the first and second rolls result in the same value and the third and fourth rolls result in another value. Again, enumerating all of the possibilities, we have: (1,2;3,4), (1,2;3,5), (1,2;4,5), (1,3;2,4), (1,3;2,5), (1,3;4,5), (1,4;2,3), (1,4;2,5), (1,4;3,5), (1,5;2,3), (1,5;2,4), (1,5;3,5), (2,3;4,5), (2,4;3,5), (2,5;3,4). There are 15 combinations of two repeated values (5 ways to choose the first pair times 3 ways to choose the second pair) for every one of the 4*3*2*1 orders in which we could observe all possible values. Adding these up, through 6 rolls, we could have
4*3*2*1*4*4 + 4*3*2*1*6*4 + 4*3*2*1*10 + 4*3*2*1*15 = 1560
Through 6 rolls, the probability of observing all 4 values is 1560/4096, about 38%.
Through 7 rolls, we have to consider: 1) observing all outcomes in 4 trials (with 3 extra trials which do not alter the outcome, but produce 4^3 observations for every permutation of the first 4), 2) observing all outcomes in 5 trials (with 2 extra trials), 3) observing all outcomes in 6 trials (with 1 extra trial), 4) observing all outcomes in 7 trials with various repetition patterns. The repetition patterns we could observe are: i) one value repeated 4 times in the first 6 trials, ii) one value repeated 3 times, a second value repeated 2 times, and iii) three values observed in the first 6 trials, all of which are observed twice.
Enumerating all of the different permutations, we have
4*3*2*1*4*4*4 + /** observed in the first 4 trials **/
4*3*2*1*6*4*4 + /** observed in the first 5 trials **/
4*3*2*1*25*4 + /** observed in the first 6 trials **/
4*3*2*1*(15 + /** observed at 7th trial; 1 value observed 4 times **/
60 + /** observed at 7th trial; 1 val observed 3 times **/
15) /** observed at 7th trial; each val observed 2 times **/
This produces 8400 outcomes which result in all 4 values of our 4-sided die being observed out of 4^7=16,384 possible outcomes. At this point, a little more than 51% of the outcomes result in observing all values. Therefore, for our 4-sided die, the median number of rolls necessary to observe all 4 outcomes would be 7 (assuming a fair die with equal probability for all possible outcomes).
Combinatorics will produce an exact answer to the question "What is the median number of trials necessary to observe all outcomes?" However, it is simpler to write a program which returns the result. For our 4-sided die, we can code the following:
data sim;
do i=1 to 100000;
outcome = cat(rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25),
rantbl(12928377,0.25,0.25,0.25,0.25));
success=0;
do trial=4 to 12;
if indexc(substr(outcome,1,trial),'1') &
indexc(substr(outcome,1,trial),'2') &
indexc(substr(outcome,1,trial),'3') &
indexc(substr(outcome,1,trial),'4') then leave;
end;
output;
end;
run;
proc format;
value trials
13 = ">=13";
run;
proc freq data=sim;
tables trial;
format trial trials.;
run;
This produces the following cumulative probabilities:
trial CumProb Combinatorics
4 0.0936 24/256=0.0938
5 0.2333 240/1024=0.2344
6 0.3820 1560/4096=0.3809
7 0.5131 8400/16384=0.5127
... ...
The simulation compares quite favorably with the combinatorics. And simulation is so much easier!
I remember doing this exact experiment with my daughter. You are so right that kids interest in probability or any other kind of math can be sparked far less with an answer than with that simple phrase "let's find out ....
Pingback: Simulating the Coupon Collector's Problem - The DO Loop
I used the Markov Chain approach to solve this problem. There are six states corresponding to the number of sides of the die already seen and the transition matrix is bidiagonal, with only the current state of the next highest one possible. It is then possible to create the set of difference equations linking the probabilty vector after N rolls of the die to the vector after N-1 rolls. These equations solve to give values of p1(N),p2(N),...,p6(n) which are the probabilities of having seen 1,2,..6 faces after n rolls, respectively. It is more convenient to write the probabilities are (n+1) throws so:
p1(n+1) = (1/6)^n
p2(n+1) = 5 * [ (2/6)^n - (1/6)^n ]
p3(n+1) = 10 * [ (3/6)^n - 2(2/6)^n + (1/6)^n ]
p4(n+1) = 10 * [ (4/6)^n - 3(3/6)^n + 3(2/6)^n - (1/6)^n ]
p5(n+1) = 5 * [ (5/6)^n - 4(4/6)^n + 6(3/6)^n - 4(2/6)^n +(1/6)^n ]
p6(n+1) = 1 - 5(5/6)^n + 10(4/6)^n - 10(3/6)^n + 5(2/6)^n - (1/6)^n ]
The pattern here, for anyone familiar with Pascal's Trinagle and combinatorics is clear, with the sign of the combinations alternating between + and -, and the proportion being summed reducing by 1/6.
The probability that it will need N throws is just the (1/6) times the probability that 5 sides had been seen after N-1 throws, i.e.
P(N rolls before all 6 sides are seen) =
(5/6) * [ (5/6)^(N-2) - 4*(4/6)^(N-2) + 6*(3/6)^(N-2) - 4*(2/6)^(N-2) + (1/6)^(N-2) ]
for N=2,3,4,...
and 0, otherwise.
By the way, it can be quicly seen that for N=2,3,4,5, this formula gives 0, as it should.
The mean of a random varable with thie PDF can be worked out to be 14.7 using basic results of sums of infinite series.
The Median is solved by seeing at which value of N the cumulative probabilities reach 0.5.
The CDF, i.e. P(number of rolls <=N) is, by summing the series:
1 - 5*(5/6)^(N-1) + 10*(4/6)^(N-1) - 10*(3/6)^(N-1) + 5*(2/6)^(N-1) - (1/6)^(N-1)
This cannot be solved explicitly, but by iteration gives a solution of the median as 12.8110473
It's not exactly SAS, but I hope you enjoy anyway!
p1(n+1)
I do! Thanks for typing up the math.
Thanks for sharing this explanation, it was a huge help! In a long blog post, I used it as a starting point for an explanation involving several additional steps (e.g., diagonalizing the matrix to work out the equations). If anyone's interested: http://www.untrammeledmind.com/2020/04/dice-probability-median-rolls-to-see-all-6-sides-coupon-collectors-problem/
Martin Robert's directed me to this StackExchange discussion, which states that the median is approximately (asymptotically) n! / n^m * S2(m,n), where S2(m,n) is the Stirling number of the second kind.
That is the exact probability actually, so setting one of the values and finding the other such that the expression evaluates to 0.5 will get the exact median.