Simplex Distributions
The simplex probabilities have support on the unit \(K\)-simplex for a specified \(K\). A \(K\)-dimensional vector \(\theta\) is a unit \(K\)-simplex if \(\theta_k \geq 0\) for \(k \in \{1,\ldots,K\}\) and \(\sum_{k = 1}^K \theta_k = 1\).
Dirichlet distribution
Probability density function
If \(K \in \mathbb{N}\) and \(\alpha \in (\mathbb{R}^+)^{K}\), then for \(\theta \in \text{$K$-simplex}\),
\[\begin{equation*} \text{Dirichlet}(\theta|\alpha) = \frac{\Gamma \! \left( \sum_{k=1}^K \alpha_k \right)}{\prod_{k=1}^K \Gamma(\alpha_k)} \prod_{k=1}^K \theta_k^{\alpha_k -1} \end{equation*}\]
Warning: If any of the components of \(\theta\) satisfies \(\theta_i = 0\) or \(\theta_i = 1\), then the probability is \(0\) and the log probability is \(-\infty\). Similarly, the distribution requires strictly positive parameters, with \(\alpha_i > 0\) for each \(i\).
Meaning of Dirichlet parameters
A symmetric Dirichlet prior is \([\alpha, \ldots, \alpha]^{\top}\). To code this in Stan,
data {
int<lower=1> K;
real<lower=0> alpha;
}
generated quantities {
vector[K] theta = dirichlet_rng(rep_vector(alpha, K));
}Taking \(K = 10\), here are the first five draws for \(\alpha = 1\). For \(\alpha = 1\), the distribution is uniform over simplexes.
1) 0.17 0.05 0.07 0.17 0.03 0.13 0.03 0.03 0.27 0.05
2) 0.08 0.02 0.12 0.07 0.52 0.01 0.07 0.04 0.01 0.06
3) 0.02 0.03 0.22 0.29 0.17 0.10 0.09 0.00 0.05 0.03
4) 0.04 0.03 0.21 0.13 0.04 0.01 0.10 0.04 0.22 0.18
5) 0.11 0.22 0.02 0.01 0.06 0.18 0.33 0.04 0.01 0.01That does not mean it’s uniform over the marginal probabilities of each element. As the size of the simplex grows, the marginal draws become more and more concentrated below (not around) \(1/K\). When one component of the simplex is large, the others must all be relatively small to compensate. For example, in a uniform distribution on \(10\)-simplexes, the probability that a component is greater than the mean of \(1/10\) is only 39%. Most of the posterior marginal probability mass for each component is in the interval \((0, 0.1)\).
When the \(\alpha\) value is small, the draws gravitate to the corners of the simplex. Here are the first five draws for \(\alpha = 0.001\).
1) 3e-203 0e+00 2e-298 9e-106 1e+000 0e+00 0e+000 1e-047 0e+00 4e-279
2) 1e+000 0e+00 5e-279 2e-014 1e-275 0e+00 3e-285 9e-147 0e+00 0e+000
3) 1e-308 0e+00 1e-213 0e+000 0e+000 8e-75 0e+000 1e+000 4e-58 7e-112
4) 6e-166 5e-65 3e-068 3e-147 0e+000 1e+00 3e-249 0e+000 0e+00 0e+000
5) 2e-091 0e+00 0e+000 0e+000 1e-060 0e+00 4e-312 1e+000 0e+00 0e+000Each row denotes a draw. Each draw has a single value that rounds to one and other values that are very close to zero or rounded down to zero.
As \(\alpha\) increases, the draws become increasingly uniform. For \(\alpha = 1000\),
1) 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10
2) 0.10 0.10 0.09 0.10 0.10 0.10 0.11 0.10 0.10 0.10
3) 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10
4) 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10
5) 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10 0.10Sampling statement
theta ~ dirichlet(alpha)
Increment target log probability density with dirichlet_lupdf(theta | alpha).
Stan functions
The Dirichlet probability functions are overloaded to allow the simplex \(\theta\) and prior counts (plus one) \(\alpha\) to be vectors or row vectors (or to mix the two types). The density functions are also vectorized, so they allow arrays of row vectors or vectors as arguments; see section vectorized function signatures for a description of vectorization.
real dirichlet_lpdf(vectors theta | vectors alpha)
The log of the Dirichlet density for simplex(es) theta given prior counts (plus one) alpha
real dirichlet_lupdf(vectors theta | vectors alpha)
The log of the Dirichlet density for simplex(es) theta given prior counts (plus one) alpha dropping constant additive terms
vector dirichlet_rng(vector alpha)
Generate a Dirichlet variate with prior counts (plus one) alpha; may only be used in transformed data and generated quantities blocks