Rolling a dice 10 times and getting exactly r same roll. I am interested in the probability of rolling a six-sided dice 10 times and getting exactly r same roll. I wrote a simple brute-force Monte Carlo simulation program to get the probability, and I also calculated the probability using binomial distribution mass function.
My Monte Carlo simulation flows like this:
-> roll a dice 10 times (an example roll: 2,1,3,4,3,6,5,3,3,2)
-> count the number of same rolls (for the example above, 2:2, 3:4)
-> extract the highest number of same rolls (for the example above, extract 4)
-> record in a list accumulating the counts (for the example above, the key "4" will increment by 1 since it's 1 case where the roll result in 4 of the same) 
-> repeat 100,000 times
I used binomial mass function like this:
for example, getting exactly 8 of the same 6*10C8*(1/6)^8*(5/6)^2 = 0.0111%
The Monte Carlo result and the binomial function result converge when the r is larger than 4, but diverge when r <= 4 . See below list of results:
M.C. Result: {2: 6.7621, 3: 52.9124, 4: 31.0536, 5: 7.8071, 6: 1.3046, 7: 0.1483, 8: 0.0111, 9: 0.00046, 10: 1e-05} %
Binomial Result: {2: 174.426, 3: 93.0272, 4: 32.5595, 5: 7.8143, 6: 1.3024, 7: 0.1488, 8: 0.0111, 9: 0.00050, 10: 9.922e-06} %
I'm just wondering why it works up until 3 or 4...
 A: Consider this roll: 1, 1, 1, 2, 3, 6, 6, 6, 6, 6. It is a valid case for $r=3$ in your function $6\cdot \binom{10}{r}\cdot \left( \frac{1}{6}\right) ^r \cdot \left( \frac{5}{6}\right) ^{10-r}$. However, we can see the problem - this roll should be considered only for $r=5$ because of five same rolls. The Monte Carlo counts this for $r=5$, but it is also one of rolls involved in calculation for $r=3$.
This is obviously caused by the factor $\left( \frac{5}{6}\right) ^{10-r}$. The initial 6 in your formula chooses, which number will be that, which is $r$-times repeated, but this doesn't prevent possibility, that some (maybe more than $r$) of the other will be same to each other. The problem is mostly visible on small $r$'s, because there is the highest probability (actually, on $r>5$ it is rather impossible), that this situation happens.
You will need to find a correct distribution. However, it will be much more complicated.
A: The issue is that the probability is actually more complicated than the binomial formula you have.  Let's look at the case of $r = 2$: you're calculating the probability that one number comes up exactly twice and then multiplying by six.  What this quantity is, is the expected amount of numbers that come up exactly twice.  This is (significantly) different than the probability that $2$ is the largest frequency.  In particular, you need to account for two facts:


*

*That $2$ will only be recorded if nothing comes up three times

*That even if two (or more) numbers come up twice then you only count $2$ once.
