# Marginal distribution of the difference of two elements of a Dirichlet distributed vector

Let $(X_1, X_2, \ldots X_n) \sim \text{Dirichlet}(\alpha_1, \alpha_2, \ldots \alpha_n)$, what is the marginal distribution of $X_i - X_j$?

• It has a piecewise formula in terms of multidimensional hypergeometric functions. It looks messy, especially for $n\gt 2$. – whuber Jan 20 '17 at 19:53
• You can collapse the distribution to $n = 3$, for what good that does, by combining all the elements of the vector other than $i,j$ into one. – jbowman Jan 20 '17 at 21:55

### Analysis

The thread at Construction of Dirichlet distribution with Gamma distribution shows that the Dirichlet distribution with parameters $$(\alpha_1, \alpha_2, \ldots, \alpha_{n+1})$$ arises as the distribution of the ratios

$$X_i=\frac{Y_i}{Y_1+Y_2+\cdots + Y_{n+1}},$$

$$i=1, 2, \ldots, n$$ where the $$Y_j$$ are independently distributed with Gamma$$(\alpha_j)$$ distributions. This permits two simplifications, because (a) we may choose the order of the $$\alpha_i,$$ reducing the question to the difference $$X_2 - X_1$$ and (b) since the sum $$Y_3 + \cdots +Y_{n+1}$$ has a Gamma$$(\alpha_3+\cdots+\alpha_{n+1})$$ distribution, we only have to consider the case $$n=2.$$

### Solution

As usual for computing marginal distributions, change the variables in the Dirichlet distribution from $$(X_1,X_2)$$ to $$(X,X+Y)$$ where $$Y$$ represents the difference $$X_2-X_1.$$ This transformation has unit Jacobian, so the integrand remains otherwise unchanged. The (unnormalized) integrand is proportional to

$$X^{\alpha_1-1} (X + Y)^{\alpha_2-1} (1 - X - (X+Y))^{\alpha_3-1}$$

and we have to integrate out $$X$$ to find the marginal distribution of $$Y$$.

Since none of the three terms can be negative, the integral breaks into two parts: one from $$x=0$$ to $$x=(1-y)/2$$ when $$y \ge 0$$ and the other from $$x=-y$$ to $$x=(1-y)/2$$ when $$y \lt 0.$$ These are well-known integrals--Mathematica or Maple or even tables of integrals should provide answers for you. The results are

$$f_Y(y) = \frac{\Gamma (\alpha_1+\alpha_2+\alpha_3)2^{-\alpha_1} }{\Gamma (\alpha_2)} y^{\alpha_2-1} (1-y)^{\alpha_1+\alpha_3-1} \, _2\tilde{F}_1\left(\alpha_1,1-\alpha_2;\alpha_1+\alpha_3;\frac{y-1}{2 y}\right)$$

for $$0 \le y \lt 1$$ and

\eqalign{f_Y(y) &= \csc \left(\pi \left(\alpha _1+\alpha _2\right)\right) (1-y)^{\alpha _3-2} \\ &\left(\frac{(y-1) \sin \left(\pi \alpha _1\right) \Gamma \left(\alpha\right) (-y)^{\alpha _1+\alpha _2} }{y \Gamma \left(\alpha _3\right)} \,_2\tilde{F}_1\left(\alpha _1,1-\alpha _3;\alpha _1+\alpha _2;\frac{2 y}{y-1}\right)\\ -\frac{\pi 2^{-\alpha_1-\alpha_2+1} \left(\alpha-1\right) (1-y)^{\alpha_1+\alpha_2}}{\Gamma \left(\alpha _1\right) \Gamma \left(\alpha _2\right)} \,_2\tilde{F}_1\left(1-\alpha _2,2-\alpha;-\alpha _1-\alpha _2+2;\frac{2 y}{y-1}\right)\right) }

for $$-1 \lt y \lt 0,$$ where $$\alpha=\alpha_1+\alpha_2+\alpha_3$$ and $$\, _2\tilde{F}_1$$ is the regularized Hypergeometric function.

For integral values of $$\alpha_1+\alpha_2$$ and negative values of $$y$$ you also have to take a limit (because the cosecant blows up). For small integral values of the parameters the function is algebraic (because the Hypergeometric functions that are involved reduce to polynomials). You can compute these by elementary means if you wish.

### Verification

Here are histograms of independent simulations of $$X_2-X_1$$ using 50,000 iterations each for three combinations of $$(\alpha_1,\alpha_2, \alpha_3+\cdots+\alpha_{n+1}).$$ On each is superimposed the graph of $$f_Y.$$ All show close agreement.

• Wow thank you for the detailed reply almost 2 years later! I ended up just sampling from the posterior in my application which was sufficient. The solution is pretty messy but it's cool to know it exists. – fouronnes Dec 12 '18 at 10:03