Definite integral over a simplex

Question

Let $T^d:=\{(x_1,..,x_d):x_i \geq 0, \sum_{i=1}^{d}x_i \leq 1\}$ be the standard simplex in $\mathbb{R}^d$. Compute the integral $$\int_{T^d} x_1^{\nu_1-1}x_2^{\nu_2-1}...x_d^{\nu_d-1}(1-x_1-...-x_d)^{\nu_0-1}$$ where $\nu_i>0$.

Remark: I know the answer is $$\frac{\prod_{i=0}^{d}\Gamma(\nu_i)}{\Gamma(\sum_{i=0}^{d}\nu_i)}.$$ I evaluated for the case $d=2$ by using the transformation $(p-1)\iiint\limits_{T^{3}} x^{m-1}y^{n-1}z^{p-2} \mathrm{d}z\mathrm{d}y\mathrm{d}x= \iint\limits_{T^{2}} x^{m-1}y^{n-1}(1-x-y)^{p-1}\mathrm{d}y\mathrm{d}x$

and the substitutions $ \left\{\begin{matrix}x=u^2& &\\y=v^2& &\\z=w^2& &\end{matrix}\right.$and $ \left\{\begin{matrix}u=r\sin\varphi\cos\theta& &\\v=r\sin\varphi\sin\theta& &\\w=r\cos\varphi& &\end{matrix}\right.,$ but this method is complex for computing the general case.

Answer 1

For any $a \in \mathbb{R}-\{0\}$ and $m,n \in (0,+\infty)$ one has \begin{align*} \frac{1}{a^{m+n-1}}\int_{0}^{a}y^{m-1}(a-y)^{n-1}\mathrm{d}y&=\frac{1}{a^{m+n-1}}\int_{0}^{a}a^{m+n-2}\left(\frac{y}{a}\right)^{m-1}\left(1-\frac{y}{a}\right)^{n-1}\mathrm{d}y\\&=\frac{1}{a}\int_{0}^{a}\left(\frac{y}{a}\right)^{m-1}\left(1-\frac{y}{a}\right)^{n-1}\mathrm{d}y\\&=\int_{0}^{1}x^{m-1}(1-x)^{n-1}\mathrm{d}x \end{align*} Thus $\int_{0}^{a}y^{m-1}(a-y)^{n-1}\mathrm{d}y=a^{m+n-1}B(m,n)$. With above observation we can integrate $x_1^{\nu_1-1}\cdots x_d^{\nu_d-1}(1-x_1-...-x_d)^{\nu_{0}-1}$ over $T^{d}$ by integrating out variables one at each step \begin{align*} &\mathrel{\phantom{=}} \int_{T^{d}}x_1^{\nu_1-1}...x_d^{\nu_d-1}(1-x_1-...-x_d)^{\nu_{0}-1}\mathrm{d}\boldsymbol x\\ &=\int_{0}^{1}\int_{0}^{1-x_1}...\int_{0}^{1-x_1-...-x_{d-1}}x_1^{\nu_{1}-1}x_2^{\nu_{2}-1}...x_d^{\nu_{d}-1}(1-x_1-...-x_d)^{\nu_{0}-1}\mathrm{d}x_d\cdots\mathrm{d}x_2\mathrm{d}x_1\\ &=\int_{0}^{1}x_1^{\nu_{1}-1}\int_{0}^{1-x_1}x_2^{\nu_{2}-1}...\int_{0}^{1-x_1-...-x_{d-1}}x_d^{\nu_{d}-1}(1-x_1-...-x_d)^{\nu_{0}-1}\mathrm{d}x_d\cdots\mathrm{d}x_2\mathrm{d}x_1\\&=B(\nu_0,\nu_d)\int_{0}^{1}x_1^{\nu_{1}-1}\int_{0}^{1-x_1}x_2^{\nu_{2}-1}...\int_{0}^{1-x_1-...-x_{d-2}}x_{d-1}^{\nu_{d-1}-1}(1-x_1-...-x_{d-1})^{\nu_{0}+\nu_{d}-1}\mathrm{d}x_{d-1}\cdots\mathrm{d}x_2\mathrm{d}x_1\\&=B(\nu_{d-1},\nu_{0}+\nu_{d})B(\nu_0,\nu_d)\int_{0}^{1}x_1^{\nu_{1}-1}...\int_{0}^{1-x_1-...-x_{d-2}}x_{d-2}^{\nu_{d-2}-1}(1-x_1-...-x_{d-2})^{\nu_{0}+\nu_{d-1}+\nu_{d}-1}\,\mathrm{d}x_{d-2}\cdots\mathrm{d}x_1\\&=\cdots\\&=B(\nu_1,\nu_0+\nu_d+\nu_{d-1}+\cdots+\nu_{2})B(\nu_2,\nu_0+\nu_d+\nu_{d-1}+\cdots+\nu_{3})\cdots B(\nu_{d-1},\nu_0+\nu_d)B(\nu_{0},\nu_d)\\&=\frac{\Gamma(\nu_{0})\Gamma(\nu_{1})\cdots\Gamma(\nu_{d})}{\Gamma(\nu_{0}+\nu_{1}+\cdots+\nu_{d})}\\&=\frac{\Gamma(\boldsymbol {\nu})}{\Gamma(|\boldsymbol {\nu}|)} \end{align*}

Answer 2

If one wishes to use a change of variables, an appropriate transformation would be

$(x_1,x_2,...,x_d)\to(y_1,y_2,...,y_d)$ such that

$\qquad y_1=x_1+x_2+\cdots+x_{d-1}+x_d$

$\qquad y_1y_2=x_1+x_2+\cdots+x_{d-2}+x_{d-1}$

$\qquad y_1y_2y_3=x_1+x_2+\cdots+x_{d-3}+x_{d-2}$

$\qquad \vdots$

$\qquad \prod_{i=1}^{d-1} y_i=x_1+x_2$

$\qquad \prod_{i=1}^{d} y_i=x_1$

$\implies x_1=y_1y_2y_3\,...\,y_{d-1}y_d$

$\qquad x_2=y_1y_2y_3\,...\,y_{d-1}(1-y_d)$

$\qquad x_3=y_1y_2y_3\,...\,y_{d-2}(1-y_{d-1})$

$\qquad \vdots$

$\qquad x_{d-1}=y_1y_2(1-y_3)$

$\qquad x_d=y_1(1-y_2)$ $\qquad\qquad$, so that $0\le y_i\le 1$ for all $i=1,2,...,d$.

The Jacobian of the transformation can be calculated as $\det(J)=-\prod_{i=1}^{d-1} y_i^{d-i}$

The integral thus reduces to

$$\displaystyle I=\int_0^1 y_1^{({\sum_{i=1}^{d}\nu_i-1})}(1-y_1)^{\nu_0-1}\mathrm{d}y_1\int_0^1 y_2^{({\sum_{i=1}^{d-1}\nu_i-1})}(1-y_2)^{\nu_d-1}\mathrm{d}y_2\cdots\int_0^1 y_{d-1}^{(\nu_1+\nu_2-1)}(1-y_{d-1})^{\nu_3-1}\mathrm{d}y_{d-1}\int_0^1 y_{d}^{(\nu_1-1)}(1-y_{d})^{\nu_2-1}\mathrm{d}y_{d}$$

$$=B\left(\sum_{i=1}^d\nu_i,\nu_0\right)B\left(\sum_{i=1}^{d-1}\nu_i,\nu_d\right)\cdots B\left(\nu_1+\nu_2,\nu_3\right)B\left(\nu_1,\nu_2\right)$$

$\qquad\qquad\qquad\quad=\displaystyle\frac{\prod_{i=0}^d\Gamma(\nu_i)}{\Gamma\left(\sum_{i=0}^d\nu_i\right)}$

Answer 3

Another, elegant way to perform your integration is explained by Jaynes in §18.10 of "Probability Theory: The Logic of Science". It uses a Laplace transform:

Slightly modifying your notation, your integral can be rewritten as

$$I(r) := \int_{[0,\infty)^n} x_1^{a_1} \dotsm x_n^{a_n} \,\mathrm{\delta}(r-x_1-\dotsb-x_n)\, \mathrm{d}x_1\dotsm\mathrm{d}x_n\qquad\text{with $r=1$}.$$

Note that the delta allows us to extend the integration to $\infty$, and renders the integration bounds independent from one another.

Now take the Laplace transform of $I(r)$ and eliminate the delta from it:

$$\int_0^{\infty} \exp(-zr) I(r)\,\mathrm{d}r = \int_{[0,\infty)^n} x_1^{a_1} \dotsm x_n^{a_n} \,\exp[-z(x_1+\dotsb+x_n)]\, \mathrm{d}x_1\dotsm\mathrm{d}x_n \\= \prod_{i=1}^n\int_0^\infty x_i^{a_i}\,\exp(-z x_i)\,\mathrm{d}x_i =\prod_{i=1}^n\frac{a_i!}{z^{a_i+1}}.$$

Now take the inverse Laplace transform of the last expression to obtain back $I(r)$:

$$I(r) = \frac{1}{2\pi\mathrm{i}}\int_{-\mathrm{i}\infty}^{+\mathrm{i}\infty} \exp(rz)\,\prod_{i=1}^n\frac{a_i!}{z^{a_i+1}}\,\mathrm{d}z \\=\frac{\prod_{i=1}^n a_i!}{(n+\sum_ia_i-1)!} r^{n+\sum_ia_i-1}$$

and evaluate it at $r=1$. The path of integration of the last integral "passes to the right of the origin, and is closed by an infinite semicircle over the left half-plane" (Jaynes).

Hope this helps.

Answer 4

Put $(n_0,n_1,\ldots, n_d)=:n$, $\sum_{k=0}^d n_k=:|n|$ and for $\lambda\geq0$ define $$Q(n,\lambda):=\int\nolimits_{\lambda T^d}\prod_{1\leq k\leq d} x_k^{n_k-1} \ (\lambda -x_1-\ldots -x_d)^{n_0-1}\ {\rm d}(x)\ .$$ The substitution $x:=\lambda y\ (y\in T^d)$ gives $$Q(n,\lambda)=\lambda^{d+|n|-(d+1)}\int\nolimits_{T^d}\prod_{1\leq k\leq d} y_k^{n_k-1} \ (1 -y_1-\ldots -y_d)^{n_0-1}\ {\rm d}(y)=\lambda^{|n|-1} Q(n,1)\ .$$ Letting the "outer integration" be with respect to the last variable we now obtain $$\eqalign{Q(n,1)&=\int_0^1 x^{n_d-1}\int\nolimits_{(1-x)T^{d-1}}\prod_{1\leq k\leq d-1} x_k^{n_k-1} (1-x_1-\ldots-x_{d-1}-x)^{n_0-1}\ {\rm d}(x')\ dx\cr &= \int_0^1 x^{n_d-1} Q(n',1-x)\ dx=Q(n',1)\int_0^1 x^{n_d-1}(1-x)^{|n'|-1}\ dx\ . \cr}$$ Here the last integral evaluates to $B(|n'|,n_d)$, so that we end up with the recursion $$Q(n,1)=Q(n',1) B(|n'|,n_d)\ .$$ From here it should not be too difficult to arrive at the desired result.

Answer 5

Another way is the following: we can rewrite the integral as $$ I_n(a)=\int_{\mathbb R_+^n} x_1^{\nu_1-1}\cdots x_n^{\nu_n-1} \delta(x_1+\cdots+x_n-a) dx_1 \cdots dx_n $$ with $n=d+1$ and $\nu_0\equiv \nu_n$. By rescaling each $x_i$ by $a$, we see that $$ I_n(a)=a^{\nu_1+\cdots + \nu_n-1} I_n(1). $$ Then, $$\begin{aligned} I_n(1) &=\int_0^1 dx\, x^{\nu_n-1} I_{n-1}(1-x) \\ &= I_{n-1}(1)\int_0^1 x^{\nu_n-1} (1-x)^{\nu_1+\cdots+\nu_{n-1}-1}dx\\ &=I_{n-1}(1) B(\nu_n,\nu_1+\cdots+\nu_{n-1}). \end{aligned}$$ Proceeding by induction, $$ I_n(1)=I_1(1) B(\nu_2, \nu_1) B(\nu_3,\nu_1+\nu_2)\cdots B(\nu_{n-1},\nu_1+\cdots + \nu_{n-2}) B(\nu_n, \nu_1+\cdots+\nu_{n-1}) $$ and after a few cancellations between numerator and denominator, $$ I_n(1)=\frac{\Gamma(\nu_1)\cdots \Gamma(\nu_n)}{\Gamma(\nu_1+\cdots+\nu_n)}\,. $$

Answer 6

Here is the solution form my professor at TU Darmstadt: We start by realizing \begin{align*} \int_{[0,\infty]^{n}}\prod_{i=1}^{n}\exp(-x_{i})x_{i}^{\alpha_{i}-1}\mathrm{d}\mathbf{x} & =\prod_{i=1}^{n}\int_{0}^{\infty}\exp(-x_{i})x_{i}^{\alpha_{i}-1}\mathrm{d}x_{i}=\prod_{i=1}^{n}\Gamma(\alpha_{i}). \end{align*} We simplify \begin{align*} & \exp\underbrace{\left(-\sum_{i=1}^{n}x_{i}\right)}_{=-t}\prod_{i=1}^{n}x_{i}^{\alpha_{i}-1}\\ & =\exp\left(-t\right)\prod_{i=1}^{n-1}\left(tp_{i}\right)^{\alpha_{i}-1}\left(t\left(1-\sum_{j=1}^{n-1}p_{j}\right)\right)^{\alpha_{n}-1}\\ & =\left(\exp\left(-t\right)t^{\left(\sum_{i=1}^{n}\alpha_{i}\right)-n}\right)\left(\prod_{i=1}^{n-1}p_{i}^{\alpha_{i}-1}\left(1-\sum_{j=1}^{n-1}p_{j}\right)^{\alpha_{n}-1}\right) \end{align*} using the transformations \begin{align*} t & =\sum_{j=1}^{n}x_{j},\\ x_{i} & =tp_{i},\\ x_{n} & =t\left(1-\sum_{i=1}^{n-1}p_{i}\right), \end{align*} which yields a Jacobian \begin{align*} \mathbf{J} & =\begin{bmatrix}t\mathbf{I} & \mathbf{p}\\ -t\mathbf{1}_{n-1}^{\intercal} & \left(1-\mathbf{1}_{n-1}^{\intercal}\mathbf{p}\right) \end{bmatrix}. \end{align*} The blockdiagonal structure of the upper $(n-1)\times(n-1)$ matrix makes the determinant

\begin{align*} \det\mathbf{J}=\left(\left(1-\mathbf{1}_{n-1}^{\intercal}\mathbf{p}\right)-\left(-t\mathbf{1}_{n-1}^{\intercal}\right)\left(t\mathbf{I}\right)^{-1}\mathbf{p}\right)\det\left(t\mathbf{I}\right)=t^{n-1}. \end{align*}

Now, we use the substitution formula \begin{align*} & \int_{[0,\infty]^{n}}\prod_{i=1}^{n}\exp(-x_{i})x_{i}^{\alpha_{i}-1}\mathrm{d}\mathbf{x}\\ & =\int_{0}^{\infty}\left(\exp\left(-t\right)t^{\left(\sum_{i=1}^{n}\alpha_{i}\right)-n}\right)\int_{[0,1]^{n-1}}\prod_{i=1}^{n-1}p_{i}^{\alpha_{i}-1}\left(1-\sum_{j=1}^{n-1}p_{j}\right)^{\alpha_{n}-1}t^{n-1}\mathrm{d}\mathbf{p}~\mathrm{d}t\\ & =\int_{0}^{\infty}\left(\exp\left(-t\right)t^{\left(\sum_{i=1}^{n}\alpha_{i}\right)-1}\right)\mathrm{d}t\int_{[0,1]^{n-1}}\prod_{i=1}^{n-1}p_{i}^{\alpha_{i}-1}\left(1-\sum_{j=1}^{n-1}p_{j}\right)^{\alpha_{n}-1}\mathrm{d}\mathbf{p}\\ & =\Gamma\left(\sum_{i=1}^{n}\alpha_{i}\right)\int_{[0,1]^{n-1}}\prod_{i=1}^{n-1}p_{i}^{\alpha_{i}-1}\left(1-\sum_{j=1}^{n-1}p_{j}\right)^{\alpha_{n}-1}\mathrm{d}\mathbf{p} \end{align*} and, thus, we have

\begin{align*} \int_{[0,1]^{n-1}}\prod_{i=1}^{n-1}p_{i}^{\alpha_{i}-1}\left(1-\sum_{j=1}^{n-1}p_{j}\right)^{\alpha_{n}-1}\mathrm{d}\mathbf{p}=\frac{\prod_{i=1}^{n}\Gamma(\alpha_{i})}{\Gamma\left(\sum_{i=1}^{n}\alpha_{i}\right)}. \end{align*}