Dice Stats

Jim Mahoney, September 2003

Setup

Load various packages.
See Add-ons / Standard Packages / Statistics  and / Graphics in the Help Browser under the Help menu.

<< Statistics`DescriptiveStatistics`;

    << Statistics`ConfidenceIntervals`;

    << Statistics`NormalDistribution`;

    << Statistics`DataManipulation`;

    << Statistics`HypothesisTests`;

    << Statistics`ANOVA`;

    << Graphics`MultipleListPlot` ; << 
        Graphics`Graphics` ;

One Die

Here are the possible outcomes from the roll of a single die.  Each outcome is equally likely.  This is our population.

oneDice  =  {1, 2, 3, 4, 5, 6};

What is the mean of this population?

μ1  =  Mean[oneDice] // N

Tr[oneDice]/Length[oneDice]

σ1  =  StandardDeviationMLE[oneDice] // N

s1 = StandardDeviation[oneDice] // N

Histogram[oneDice]

Is that last graph lying with statistics or what?  Here's a more complicated way of asking for that, specifying that the counting is to be done in bins that are 1 wide from 0.5 to 6.5.

Histogram[BinCounts[oneDice, {0.5, 6.5, 1}], 
      FrequencyData -> True]

TwoDice

The sum of two dice can be anwhere from 2 to 12.  Since there are 6 ways the first dice can come up, and 6 ways for the second, there are 36 equally likely outcomes.  They are

twoDice  =  Flatten[Table[i + j, {i, 6}, {j, 6}]]

Histogram[twoDice]

Histogram[BinCounts[twoDice, {1.5, 12.5, 1}], 
      FrequencyData -> True]

μ2  =  Mean[twoDice] // N

σ2 = StandardDeviationMLE[twoDice] // N

s2 = StandardDeviation[twoDice] // N

Three Dice

The number of possibilities get bigger, but the idea is the same.  Now there are (6^3) = 216equally likely outcomes.  The smallest is 1+1+1=3 and the largest is 6+6+6=18.  We can let the computer generate the list.

threeDice  =  
      Flatten[Table[i + j + k, {i, 6}, {j, 6}, {k, 6}]]

Histogram[threeDice]

This histogram is starting to look like a bell curve - which is exactly what its heading towards.

μ3 = Mean[threeDice] // N

σ3 = StandardDeviationMLE[threeDice] // N

s3 = StandardDeviation[threeDice] // N

Eight Dice

Way too big to do by hand, but we can still (barely) have the computer do the whole population for us.

dice = 0;
Do[
        dice = Table[i + dice, {i, 1, 6}], {8}];
dice = Flatten[dice];

diceFreq  =  BinCounts[dice, {6, 48, 1}]

Now there's a true Normal distribution.

Histogram[diceFreq, FrequencyData -> True]

The total population size should be

(6^8)

which is (just checking) the same as the number of entries in our computer generated population.

Length[dice]

And the population parameters are

μ8 = Mean[dice] // N

σ8 = StandardDeviationMLE[dice] // N

s8 = StandardDeviation[dice] // N

Result Summary

\({{"", "one", "two", "three", "...", "eight"}, 
{"μ", μ1, μ2, μ3, "...", μ8}, 
{"σ", σ1, σ2, σ3, "...", σ8}, 
{"s", s1, s2, s3, "...", 
          s8}} // MatrixForm\)

There are several trends to notice here:

First, the difference between σ and s becomes less noticeable as your sample size gets larger.

Second, the mean of a sum of randrom numbers gets bigger much faster than the standard deviation does.  That means that typically these histograms get narrower as the number of items in the population gets bigger.  

Third, for say 20 dice even a good computer will crap out: (6^20)is just too big.  So either we find a way to get what we want analytically (which is pretty straightforward for a situation as simple as dice goes), or we sample the population, and estimate things from that sample.

Estimation of Population Parameters (μ,σ) from Sample Statistics (m,s)

For example, suppose we didn't know the result for 8 dice, and wanted to find out what it was.  Well then, roll eight dice a few times, call that your sample, and then use those numbers to estimate (i.e. guess) what the real results are.

rollDice[n_] := 
        Sum[Random[Integer, {1, 6}], {n}];
diceSample = Table[rollDice[8], {10}]

What that just did was find the sum of 6 dice (random integers from 1 to 6) ten times.

Mean[diceSample] // N

StandardDeviation[diceSample] // N

The true values which we found before are 28.0 and 4.83, so you see that we came pretty close even with only 10 randomly chosen members of the entire population of 1.7 million .

σ vs sComparison

σ[x_] := StandardDeviationMLE[x];
s[x_] := StandardDeviation[x];

t  =  {90. ,  80. ,  85. , 85. };

Now find σ and s for the list t.  Make sure you understand the difference: the first is the standard deviation of that population of values.  The second is the best guess for what the standard deviation is for some large population of which t is a sample.  The second one is bigger because we had to guess the mean of the population from the mean of the sample, which makes the result more uncertain.

{σ[t], s[t]}