globals [
numdice
observedfrequency
grandtotal
pmf
standarddev
plottheoretical?
maxtheoreticalprobability
degreesoffreedom
chisquared
totalroll
rollsplotted
]
breed [dice]
diceown [
spots
]
to startup
setup
end
to setup
clearall
setupglobals
setupdice
resetticks
setupplot
end
to setupglobals
set numdice numberofdice
set observedfrequency (nvalues (5 * numdice + 1) [0])
set grandtotal 0
set pmf theoreticalpmf numdice
set standarddev (sqrt (numdice * 35 / 12))
set maxtheoreticalprobability (max pmf)
set plottheoretical? showtheoretical?
set degreesoffreedom (5 * numdice)
end
to setupdice
setdefaultshape dice "die1"
createdice numdice [
set xcor who
set color white
rollone
]
end
to setupplot
setplotxrange (numdice  0.5) (0.5 + numberofdice * 6)
setplotyrange 0 ((ceiling (100 * maxtheoreticalprobability)) / 100)
if plottheoretical? [
plottheoreticaldistribution
plotmoments
]
end
to go
roll
tick
end
to roll
ask dice [
rollone
]
updateobservations (sum [spots] of dice)
end
to rollone
set spots (1 + random 6)
set shape (word "die" spots)
end
to updateobservations [value]
let head (sublist observedfrequency 0 (value  numdice + 1))
let tail
(sublist observedfrequency (value  numdice + 1) (5 * numdice + 1))
let newfrequency (1 + last head)
set observedfrequency (sentence (butlast head) newfrequency tail)
set totalroll value
set grandtotal (grandtotal + value)
end
to updateplotrange
let peakcount (max observedfrequency)
let ylimit (max (list maxtheoreticalprobability (peakcount / ticks)))
setplotyrange 0 ((ceiling (50 * ylimit)) / 50)
end
to updateplotaverage
let average (grandtotal / ticks)
setcurrentplotpen "average"
plotpenreset
plotxy average 0.1
plotxy average 1
end
to updatehistogram
setcurrentplotpen "observed"
plotpenreset
(foreach (nvalues (5 * numdice + 1) [numdice + ?]) observedfrequency [
plotxy (?1  0.5) (?2 / ticks)
])
set rollsplotted ticks
end
to updateplot
if (ticks > 0) [
setcurrentplot "Relative Frequencies"
updateplotrange
updateplotaverage
updatehistogram
updatestatistics
]
end
to updatestatistics
set chisquared
(sum (map [(?1 * ?1) / (?2 * ticks)] observedfrequency pmf)  ticks)
end
to plottheoreticaldistribution
setcurrentplot "Relative Frequencies"
setcurrentplotpen "pmf"
plotpenreset
(foreach (nvalues (5 * numdice + 1) [numdice + ?]) pmf [
plotxy (?1  0.5) ?2
plotxy (?1 + 0.5) ?2
])
end
to cleartheoreticaldistribution
setcurrentplotpen "pmf"
plotpenreset
end
to plotmoments
setcurrentplotpen "moments"
plotpenreset
plotxy (3.5 * numdice  standarddev) 0.1
plotxy (3.5 * numdice  standarddev) 1
plotxy (3.5 * numdice) 1
plotxy (3.5 * numdice) 0.1
plotxy (3.5 * numdice + standarddev) 0.1
plotxy (3.5 * numdice + standarddev) 1
end
to clearmoments
setcurrentplotpen "moments"
plotpenreset
end
toreport theoreticalpmf [n]
let result []
let possiblevalues (5 * n + 1)
let setcardinality (6 ^ n)
foreach (nvalues (ceiling (possiblevalues / 2)) [?]) [
let eventcardinality (combinations n (n + ?))
set result (lput (eventcardinality / setcardinality) result)
]
set result (sentence result reverse
ifelsevalue ((possiblevalues mod 2) = 1) [butlast result] [result])
report result
end
toreport combinations [n x]
report enumerate n x []
end
toreport enumerate [n x dicelist]
let result 0
ifelse ((length dicelist) = n) [
set result (factorial n)
foreach (removeduplicates dicelist) [
let currentvalue ?
let occurrences (length (filter [? = currentvalue] dicelist))
if (occurrences > 1) [
set result (result / factorial occurrences)
]
]
]
[
let diceremaining (n  (length dicelist))
let runningsum (sum dicelist)
let maxpossible (min (list
(ifelsevalue (empty? dicelist) [6] [last dicelist])
(x  runningsum  diceremaining + 1)))
let minpossible
(max (list 1 (ceiling ((x  runningsum) / diceremaining))))
if (maxpossible >= minpossible) [
foreach (nvalues (maxpossible  minpossible + 1) [maxpossible  ?]) [
set result (result + enumerate n x lput ? dicelist)
]
]
]
report result
end
toreport factorial [n]
let product 1
foreach (nvalues n [1 + ?]) [
set product (product * ?)
]
report product
end
toreport plotdaemon?
let result? false
if (showtheoretical? xor plottheoretical?) [
set plottheoretical? showtheoretical?
ifelse (plottheoretical?) [
plottheoreticaldistribution
plotmoments
set result? true
]
[
cleartheoreticaldistribution
clearmoments
]
]
if (ticks > rollsplotted) [
updateplot
set result? true
]
report result?
end
;; Copyright (c) 2013, Nicholas Bennett. All rights reserved.
;; For more informationon on permitted uses, see the Information tab.
@#$#@#$#@
GRAPHICSWINDOW
10
10
740
101
1
0
60.0
1
10
1
1
1
0
1
1
1
0
11
0
0
1
1
1
Rolls
30.0
SLIDER
10
110
160
143
numberofdice
numberofdice
1
12
1
1
1
NIL
HORIZONTAL
MONITOR
190
290
265
335
NIL
plotdaemon?
17
1
11
BUTTON
35
270
135
303
Roll Forever!
go
T
1
T
OBSERVER
NIL
R
NIL
NIL
0
BUTTON
35
150
135
184
Setup
setup
NIL
1
T
OBSERVER
NIL
S
NIL
NIL
1
SWITCH
10
310
160
343
showtheoretical?
showtheoretical?
1
1
1000
BUTTON
35
190
135
223
Roll Once
go
NIL
1
T
OBSERVER
NIL
O
NIL
NIL
0
BUTTON
35
230
135
263
Roll 100
repeat 100 [\n go\n]
NIL
1
T
OBSERVER
NIL
1
NIL
NIL
0
PLOT
170
110
740
345
Relative Frequencies
NIL
NIL
0.0
10.0
0.0
10.0
false
false
"" ""
PENS
"pmf" 1.0 0 2382653 true "" ""
"moments" 1.0 0 6759204 true "" ""
"average" 1.0 0 4539718 true "" ""
"observed" 1.0 1 16777216 true "" ""
MONITOR
680
110
740
155
Average
ifelsevalue (ticks > 0) \n [grandtotal / ticks]\n [\"N/A\"]
3
1
11
@#$#@#$#@
## WHAT IS IT?
This model repeatedly rolls 1 to 12 dice, computes the sum, and displays a histogram of the empirical distribution of the resulting values. Optionally, it also displays the theoretical distribution of the values.
## HOW IT WORKS
This is a very simple Monte Carlo simulation (an experimental method where random or pseudorandom numbers are generated and used to simulate a process with a high degree of randomness). Here, each trial consists of generating one or more uniform discrete pseudorandom variables, from the set {1, 2, 3, 4, 5, 6} – i.e. modeling the roll of dice – and summing them. This process is illustrated by creating one agent per die, and dynamically changing the die shape, based on the random number generated. After each trial, the result is used to update a tally – which is in turn used to construct a histogram of relative frequencies.
## HOW TO USE IT
Use the controls to configure and run the simulation as follows.
* The **numberofdice** slider controls the number of sixsided dice used in each roll.
* The **Setup** button creates agents for the dice display, resets the cell counts used for the histogram, and clears the plot.
* **Roll Once** executes a single roll of the dice.
* The **Roll 100** button rolls the dice 100 times.
* The **Roll Forever!** button rolls repeatedly until the user presses the button again. As it does so, the **Relative Frequencies** plot is automatically updated periodically.
* **showtheoretical?** toggles the display of the probability mass function (PMF) of the theoretical distribution, and the μ and μ ± σ lines in the **Relative Frequencies** plot (see below).
The **Average** monitor displays the current average value, computed across all the rolls (since the last time **Setup** was clicked).
The **Relative Frequencies** plot includes the following data.
* The relative frequency observed for each value is plotted in black as a histogram.
* The sample mean (average) of the observed values is shown with a gray vertical line.
Optionally (using the **showtheoretical?** switch), additional information is displayed in **Relative Frequencies**:
* The probability mass function (the theoretical relative frequencies for the dice totals) is shown in magenta, as the upper outline of a histogram.
* μ (the mean, or theoretical average) and μ ± σ (the mean plus or minus one standard deviation) are displayed with cyan vertical lines.
(In the histogram of relative frequencies, each bar begins at the relevant value, less 0.5; for example, with 1 die, the first bar start at an X value of 0.5 and ends at an X value of 1.5. This doesn't follow the standard convention for frequency histograms for quantitative sample values – in which each bar begins with the smallest value in the interval – but it can be less confusing to the casual user.)
## THINGS TO NOTICE
With a single die, the histogram of dice rolls becomes roughly flat, after a large number of rolls. However, with two or more dice, that's no longer the case. Why not?
Examining the possible values produced by the sum of two dice shows us that even though each of the numbers 16 is equally likely on a single die, the sums 212 obtained from a pair of dice aren't equally likely.

Die #2 

(1 + 2) 
1 
2 
3 
4 
5 
6 
Die #1
 1 
2 
3 
4 
5 
6 
7 
2 
3 
4 
5 
6 
7 
8 
3 
4 
5 
6 
7 
8 
9 
4 
5 
6 
7 
8 
9 
10 
5 
6 
7 
8 
9 
10 
11 
6 
7 
8 
9 
10 
11 
12 
As the table above shows, there are 36 different combinations of rolls obtained with a pair of dice, and some combinations give the same sum as other combinations. The proportion of occurrence of a given sum in the table corresponds to its longterm frequency of occurrence in the actual dice rolls. For example, 6 of the 36 combinations result in a sum of 7, while only 1 of the 36 gives a sum of 12. Thus, over a large number of rolls, we should expect a sum of 7 to appear approximately 6/36^{ths}, or 1/6^{th} of the time; on the other hand, 12 should only appear about 1/36^{th} of the time.
As the number of dice increases further, the histogram of relative frequencies becomes more tightly grouped around the average value, which also tends to be very close to the midpoint of the range of possible values. This is seen especially clearly when **showtheoretical?** is selected: as the number of dice increase, the band of possible values between the standard deviation lines is a smaller fraction of the total range of possible values. Notice also that the line showing the average of the observed values moves closer and closer to the line for the mean, as the number of rolls increases.
Both of the above are examples of the phenomenon first described mathematically by Bernoulli, then called the "law of large numbers" by Poisson. The law of large numbers is a theorem showing that the sample mean converges to the true mean as the sample size goes to infinity, regardless of the underlying distribution.
You might also notice that the histograms of observed and theoretical relative frequencies become more and more bellshaped, as the number of dice increases. This phenomenon is explained by the "central limit theorem," originally stated by de Moivre, refined significantly by Laplace, and rigorously generalized by Lyapunov. (The "central limit" name was first applied by Polya.) This theorem shows that the distribution of the sum (or average) of independent, identically distributed random values, with finite mean and variance, converges to the normal distribution.
The law of large numbers and the central limit theorem are the first two fundamental theorems of probability and statistics.
## THINGS TO TRY
To get a taste of how the computational capabilities of modern computers have revolutionized Monte Carlo methods, turn the **view updates** checkbox off, and click the **Run Forever!** button, to run experiments as fast as possible. Depending on the processing speed of your computer, the model should be able to roll the dice hundreds of thousands – even millions – of times within a minute or two. Imagine trying to perform the same number of trials by hand – or even with the primitive computing facilities available to the scientists working at Los Alamos National Laboratory in the 1940s. (Monte Carlo methods were first used in the Manhattan Project, and were subsequently studied and extended at LANL and elsewhere.)
Of course, realworld Monte Carlo simulations are rarely as simple as this one, but the increase in computing power has had a dramatic impact across the spectrum of model types and magnitudes.
## EXTENDING THE MODEL
The model could easily be extended to include a Pearson's chisquared goodness of fit test, for use by students studying statistical hypothesis testing. In that case, the null hypothesis would be that the generated outcomes are distributed according to the relevant theoretical distributions. (In fact, the model already declares and updates two global variables, `chisquared` and `degreesoffreedom`, that could be used for this purpose.)
## NETLOGO FEATURES
Although the histogram primitive can plot a histogram from a list of data with minimal code, performance degrades as the list of data gets longer and longer. For that reason, this model keeps a list of cell counts, instead of the raw data, and it uses plotxy to construct the histogram (as well as the other plotted displays) explicitly. This results in constant plotting performance, even after millions of trials.
For performance reasons, plotting isn't driven by the NetLogo v5+ tickbased `updateplots` mechanism, and thus isn't automatically updated on every tick. Instead, it's driven by a hidden monitor (located under the **Relative Frequencies** plot) that displays the value of a reporter procedure. That reporter procedure checks to see if the plot needs to be updated (i.e. if the number of ticks is higher than the number of rolls plotted), and does so if necessary. This means that at high speed, several ticks (sometimes thousands) can go by between plot updates, but the updates still occur often enough to appear smooth.
Monitoring the **showtheoretical?** switch (which is done even when the dice aren't being rolled) is done with the same hidden monitor and reporter procedure. That procedure watches for changes in the switch value and updates the plot accordingly.
For some models, this monitordriven behavior is a useful technique for implementing dæmons that run asynchronously to the main simulation loop (e.g. a forever button) – and that can even run when no button is pressed.
## RELATED MODELS AND ACTIVITIES
This model was inspired in part by the Project GUTS "Dice & Data" activity, in which students form pairs and roll one die and then two dice, tallying the outcomes in each case. These tallies are aggregated and used to construct histograms, illustrating the effect on central tendency produced by summing random variables.
The "Dice Stalagmite" model in the standard NetLogo models library is somewhat similar to this one, in that histograms are drawn for one and twodice rolls. However, the emphasis in that model is on using animation to construct the histograms, rather than comparing the observed and theoretical distributions over a large number of trials; it is also limited to rolls of two dice. (No code from the "Dice Stalagmite" model is used in this model.)
## REFERENCES
Mlodinow, L. _The Drunkard's Walk: How Randomness Rules Our Lives_. Random House Inc., New York, NY, 2009.
Boslaugh, S. and Watters, P. A. _Statistics in a Nutshell_. O'Reilly Media, Inc., Sebasopol, CA, 2008.
"Central limit theorem". http://en.wikipedia.org/wiki/Central_limit_theorem. Wikipedia, Dec. 2012. (Accessed: 26 Jan. 2013.)
"Law of large numbers". http://en.wikipedia.org/wiki/Law_of_large_numbers. Wikipedia, Jan. 2013. (Accessed: 26 Jan. 2013.)
"Monte carlo method". http://en.wikipedia.org/wiki/Monte_Carlo_method. Wikipedia, Jan. 2013. (Accessed: 26 Jan. 2013.)
Lee, I. "Dice & Data to Random Walks". http://www.projectguts.org/files/D1Dice&Data_0.doc. Project GUTS, 2010. (Accessed: 28 Oct. 2010.)
Abrahamson, D. and Wilensky, U. NetLogo Dice Stalagmite model. http://ccl.northwestern.edu/netlogo/models/DiceStalagmite. Center for Connected Learning and ComputerBased Modeling, Northwestern University, Evanston, IL, 2005.
## COPYRIGHT
Copyright © 2013, Nicholas Bennett.
All rights reserved.
Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
THIS SOFTWARE IS PROVIDED BY NICHOLAS BENNETT "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL NICHOLAS BENNETT BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
## ACKNOWLEDGEMENTS
Development of this model was funded in part by [Project GUTS](http://www.projectguts.org) and the [Santa Fe Alliance for Science](http://www.sfafs.org/). As noted above, the concepts of this model were inspired in part by the Project GUTS "Dice & Data" activity.
@#$#@#$#@
default
true
0
Polygon 7500403 true true 150 5 40 250 150 205 260 250
die1
false
0
Rectangle 7500403 true true 15 60 285 240
Rectangle 7500403 true true 60 15 240 285
Circle 7500403 true true 15 15 88
Circle 7500403 true true 195 15 90
Circle 7500403 true true 15 195 90
Circle 7500403 true true 195 195 90
Circle 16777216 true false 120 120 60
die2
false
0
Rectangle 7500403 true true 15 60 285 240
Rectangle 7500403 true true 60 15 240 285
Circle 7500403 true true 15 15 88
Circle 7500403 true true 195 15 90
Circle 7500403 true true 15 195 90
Circle 7500403 true true 195 195 90
Circle 16777216 true false 195 45 60
Circle 16777216 true false 45 195 60
die3
false
0
Rectangle 7500403 true true 15 60 285 240
Rectangle 7500403 true true 60 15 240 285
Circle 7500403 true true 15 15 88
Circle 7500403 true true 195 15 90
Circle 7500403 true true 15 195 90
Circle 7500403 true true 195 195 90
Circle 16777216 true false 195 195 60
Circle 16777216 true false 45 45 60
Circle 16777216 true false 120 120 60
die4
false
0
Rectangle 7500403 true true 15 60 285 240
Rectangle 7500403 true true 60 15 240 285
Circle 7500403 true true 15 15 88
Circle 7500403 true true 195 15 90
Circle 7500403 true true 15 195 90
Circle 7500403 true true 195 195 90
Circle 16777216 true false 195 195 60
Circle 16777216 true false 45 45 60
Circle 16777216 true false 45 195 60
Circle 16777216 true false 195 45 60
die5
false
0
Rectangle 7500403 true true 15 60 285 240
Rectangle 7500403 true true 60 15 240 285
Circle 7500403 true true 15 15 88
Circle 7500403 true true 195 15 90
Circle 7500403 true true 15 195 90
Circle 7500403 true true 195 195 90
Circle 16777216 true false 195 195 60
Circle 16777216 true false 45 45 60
Circle 16777216 true false 45 195 60
Circle 16777216 true false 195 45 60
Circle 16777216 true false 120 120 60
die6
false
0
Rectangle 7500403 true true 15 60 285 240
Rectangle 7500403 true true 60 15 240 285
Circle 7500403 true true 15 15 88
Circle 7500403 true true 195 15 90
Circle 7500403 true true 15 195 90
Circle 7500403 true true 195 195 90
Circle 16777216 true false 195 195 60
Circle 16777216 true false 45 45 60
Circle 16777216 true false 45 195 60
Circle 16777216 true false 195 45 60
Circle 16777216 true false 45 120 60
Circle 16777216 true false 195 120 60
dot
false
0
Circle 7500403 true true 90 90 120
turtle
true
0
Polygon 10899396 true false 215 204 240 233 246 254 228 266 215 252 193 210
Polygon 10899396 true false 195 90 225 75 245 75 260 89 269 108 261 124 240 105 225 105 210 105
Polygon 10899396 true false 105 90 75 75 55 75 40 89 31 108 39 124 60 105 75 105 90 105
Polygon 10899396 true false 132 85 134 64 107 51 108 17 150 2 192 18 192 52 169 65 172 87
Polygon 10899396 true false 85 204 60 233 54 254 72 266 85 252 107 210
Polygon 7500403 true true 119 75 179 75 209 101 224 135 220 225 175 261 128 261 81 224 74 135 88 99
@#$#@#$#@
NetLogo 5.0.3
@#$#@#$#@
@#$#@#$#@
@#$#@#$#@
@#$#@#$#@
@#$#@#$#@
default
0.0
0.2 0 1.0 0.0
0.0 1 1.0 0.0
0.2 0 1.0 0.0
link direction
true
0
Line 7500403 true 150 150 90 180
Line 7500403 true 150 150 210 180
@#$#@#$#@
1
@#$#@#$#@