The Dirichlet function is the conjugate prior of the multinomial. So the posterior is also Dirichlet given some observations.
If e.g. I observe the counts $X=(10,3,4)$ from 17 trials (10 for class 1, 3 for class 2 and 4 for class 3).
Assuming a non-informative prior, with concentration parameter $\alpha$ equal to $(1.,1.,1.)$, the posterior Dirichlet should have the concentration parameter $\alpha_{\text{post}}=(11.,4.,5.)$ (according to my knowledge).

If I draw lots of samples from a multinomial with such an posterior $\alpha_{\text{post}}$ the resulting maximum likelihood estimator of the probabilities of the samples differs from the mle of the observations mle. (see code below).

I don't understand why: Is the posterior Dirichlet alpha given by: $$ (\alpha_0 + X_0 - 1, \alpha_1 + X_1 - 1, , \alpha_2 + X_2 - 1 ) $$ In my simulation that seems to be the correct answer.

import numpy as np
import scipy.stats
x = np.array([ 10., 3., 4.])

# posterior concentration parameters for Dirichlet
# with uninformative prior, i.e. concentration parameters [1,1,1]
a = x + 1
a = x

# simulate 
nb_samples = 200000
d = scipy.stats.dirichlet.rvs(a, nb_samples)
c_ = np.ndarray(nb_samples)
# without a loop possible??
for i in range(nb_samples):
    c_[i] = np.random.choice(3, 1, p=d[i])

nbs = float(nb_samples) 
counts = np.bincount(c_.astype(int))
print "concentration parameter of the Dirichlet:", a
print "simulation sample probabilities (mle)", counts /nbs
print "observation mle probabilities:", x/float(x.sum())

>concentration parameter of the Dirichlet: [ 11.   4.   5.]
>simulation sample probabilities (mle) [ 0.550255  0.198655  0.25109 ]
>observation mle probabilities: [ 0.58823529  0.17647059  0.23529412]

>concentration parameter of the Dirichlet: [ 10.   3.   4.]
>simulation sample probabililities (mle) [ 0.58836  0.17733  0.23431]
>observation mle probabilities: [ 0.58823529  0.17647059  0.23529412]

1 Answer 1


I do not think this has anything to do with a wrong definition of the Dirichlet prior or posterior: simply, when $$(x_1,\ldots,x_k)\sim\mathcal{D}(\alpha_1,...,\alpha_k)$$ the mean is given by $$\mathbb{E}[(x_1,\ldots,x_k)]=\dfrac{(\alpha_1,...,\alpha_k)}{\sum_{i=1}^k\alpha_i}$$ which explains for the discrepancy with the MLE.

[1] 0.55 0.20 0.25
> c(10.,3.,4.)/sum(c(10.,3.,4.))
[1] 0.5882353 0.1764706 0.2352941

If instead you use the mode of the Dirichlet distribution, $$(x_1^\text{mode},\ldots,x_k^\text{mode})=\dfrac{(\alpha_1-1,...,\alpha_k-1)}{\sum_{i=1}^k\alpha_i-k}$$ which recovers the MLE. And this makes complete sense because the MAP is the MLE under a flat prior.

  • 1
    $\begingroup$ So I computed the mean from the samples (see code above), which differs from the mode in case of the dirichlet distribution. I didn't know that the mode is different from the mean for the dirichlet, so I came to the wrong conclusion. $\endgroup$ Nov 23, 2015 at 14:41

Your Answer

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

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