1
$\begingroup$

Although it is unknown dice, the symmetry of the evidence tells us, that we can treat the dice as fair, so the chance should be exactly 50%.

But if we simulate it by hand, the result is less then 50%:

  • Mathematica code

In[1]:= nsym = 200000 (*Some very large number of simulations*)
In[46]:= sample = RandomReal[DirichletDistribution[Table[100, {6}]], nsym]; 
In[51]:= counts = If[#1 > 1/6, 1, 0] & /@ #1 & /@ Transpose[sample]; 
In[52]:= N[Plus @@ Transpose[counts]/nsym] 
Out[52]= {0.49387, 0.48934, 0.487965, 0.488245, 0.49123}
  • R code

require(gtools)
nsym=200000  #Some very large number of simulations
sample<-rdirichlet(nsym,rep(100,6))
apply(sample,2,function(x,testval){sum(x>testval)},testval=1/6)/nsym
# [1] 0.489715 0.490490 0.489800 0.490315 0.491900 0.488165

What am I doing wrong??


update:

I've found out empirically, that the mean of the difference between the above probabilities and the "expected" 0.5 equals approximately $0.096/\sqrt{100}$, where "100" is the number of prior trials.

$\endgroup$
3
  • $\begingroup$ What justifies drawing from the Dirichlet distribution to simulate outcomes of the rolls of a die? $\endgroup$
    – whuber
    May 6, 2014 at 16:59
  • $\begingroup$ @whuber I do Bayesian reasoning here: I don't simulate outcomes of the rolls of die - for that I'd use categorial distribution. I simulate the probabilities with which the die rolls, given the already observed historical outcomes, i.e. the posterior distribution of the parameters $\pi_{1\ldots6}$ that define the die. $\endgroup$ May 14, 2014 at 15:29
  • 1
    $\begingroup$ Thank you. I think it would help to edit your question to make that approach clearer, rather than relying on the code to convey what you are doing. Consider improving the title, too: because it refers to "frequency" it suggests something other than the (posterior) probabilities which you seem to be estimating. $\endgroup$
    – whuber
    May 14, 2014 at 16:52

2 Answers 2

1
$\begingroup$

Seems like 2% of the time you probably get exactly 100 rolls which is neither greater or less than 1/6. It is exactly 1/6.

$\endgroup$
0
$\begingroup$

Dirichlet distribution is not multivariate normal, so it makes no sense to expect that mean ($1/6$) equals median (i.e. one observes values below mean with frequency 0.5).

The bigger is evidence (i.e. the number "100" is greater), the more the dirichlet distribution converges to normal distribution.

Why the expected 0.5 frequency is bigger then the observed?

It so happens, that (similar to the beta distribution) dirichlet distribution is positively skewed when the shape parameter is < 0.5 (in our case it equals 1/6). Hence the mean (0.5) is greater then the median.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.