### Example 19-3 The 2D6 Simulation

### This Monte Carlo simulation simulates rolling 2D6 10000 times

# Generate a vector of random values of 1-6 for each of 2 dice
# Parameter replace=TRUE allows the generation of more than 6 elements

dice1 <- sample(1:6, 10000, replace = TRUE)
dice2 <- sample(1:6, 10000, replace = TRUE)

# Bind into a dataframe with sum of the two die rolls

results <- data.frame(sum=dice1+dice2,D1=dice1,D2=dice2)

# Generate the expected probabilities for each 2D6 roll

expected_probs <- data.frame(
  val = 2:12,
  prob = c(1, 2, 3, 4, 5, 6, 5, 4, 3, 2, 1) / 36
)

# Generate the expected frequencies for each bin

expected_probs$count <- expected_probs$prob * 10000

### Create a basic histogram of the sum using 11 columns

hist(results$sum, breaks=11)