Proof-Check: Beta-distribution sample generation through variable transformation

133 Views Asked by At

Purpose of this thread: I want your feedback on my proof for the following problem and correct it where necessary.

Problem: For $\alpha > 0$ and $\beta > 0$, consider the following accept–reject algorithm:

  1. Generate $U_1$ and $U_2$ iid uniform$(0,1)$ random variables. Set $V_1 = U_1^{1/\alpha}$ and $V_2 =U_2^{1/\beta}$.

  2. Set $W = V_1 + V_2$. If $W \leq 1$, set $X = V_1/W$; else go to step 1.

  3. Deliver X.

Show that X has a beta distribution with parameters $\alpha$ and $\beta$.

Related threads: This question as a complete proof here.

My attempt: Using the so called probability integral transform, we get that the cdf of $V_1$ and $V_2$ are $F_{V_1}(v_1) = v_1^\alpha$ and $F_{V_2}(v_2) = v_2^\beta$. We want $$ \begin{array}{ccl} F_X(x) &= P(X \leq x) &= P\left(V_1/(V_1+V2 \leq x) \bigr\rvert V_1 + V_2 \leq 1 \right) \\ &&= \cfrac{P \left( \frac{V_1}{V_1+V_2} \leq x, V_1+V_2 \leq 1 \right)}{P\left( V_1 + V_2 \leq 1\right)} \\ \end{array} $$

One approach is to simply manage with the joint pdf of $V_1$ and $V_2$ as was done here. My approach was to define certain relevant transformations which could directly lead us to the exact expression for the cdf of the beta distribution.

Consider the transformations $Y = V_1/(V_1+V_2)$ and the given $W = V_1+V_2$. Then we have $V_1 = YW$ and $V_2 = W-YW$, the Jacobian would be $W$ and so the joint pdf ends up being \begin{array}{ll} f_{Y,W}(y,w) &= \alpha(yw)^{\alpha-1}\beta(w-yw)^{\beta-1} \times w \\ &= \alpha \beta w^{\alpha + \beta -1} y^{\alpha -1 }(1-y)^{\beta - 1} \end{array}

Check if joint pdf integrates to $1$:

As $v_1$ and $v_2$ change from $0$ to $1$ ($0$ and $1$ are not a part of the support as question specifically mentions "uniform$(0,1)$"), the transformed variable $y$ changes from $0$ to $1$ while $w$ changes from $0$ to $2$. However, as $v_1 < 1$, $yw < 1$. Also, $v_2 < 1$ implies $w-yw < 1$ which leads to additional constraints on the integration limits as $(w-1)/w < y < 1/w$ (this is valid as $0 < w < 2$).

These constraints when expressed for $w$ in terms of $y$ end up being

$$ \begin{array}{ll} 0 < w < \frac{1}{1-y} & \text{, if } 0 < y < 1/2 \text{,} \\ 0 < w < \frac{1}{y} & \text{, if } 1/2 < y < 1 \text{.} \end{array} $$

The above defines the region of integration for the pdf which needs to integrate to unity if derived correctly. With that, we need to show that $$ \int_0^{1/2} \int_{0}^{1/(1-y)} f_{Y,W}(y,w) dw dy + \int_{1/2}^{1} \int_{0}^{1/y} f_{Y,W}(y,w) dw dy = 1. $$

Now

$$ \begin{array}{ll} \int_0^{1/2} \int_{0}^{1/(1-y)} \alpha \beta w^{\alpha + \beta -1} y^{\alpha -1 }(1-y)^{\beta - 1} dw dy + \int_{1/2}^{1} \int_{0}^{1/y} \alpha \beta w^{\alpha + \beta -1} y^{\alpha -1 }(1-y)^{\beta - 1} \\ = \frac{\alpha \beta}{\alpha + \beta} \left\{ \int_0^{1/2} \frac{y^{\alpha -1 }}{(1-y)^{\alpha + 1}} dy + \int_{1/2}^{1} \frac{(1-y)^{\beta - 1}}{y^{\beta +1 }} dy \right\}\\ = \frac{\alpha \beta}{\alpha + \beta} \left\{ \frac{1}{\alpha} + \frac{1}{\beta} \right\} \\ =1. \end{array} $$

This establishes that (at least :)) the derived joint pdf is valid!

Derivation of $F_X(x):$

Now that we have established that the joint pdf is valid, we evaluate the numerator and the denominator of the expression for $F_X(x)$. As the conditional probability involves the region where $w \leq 1$, we just have to integrate the joint pdf in the region $0 < y < 1$ and $0 < w \leq 1$.

$$ \begin{array}{ll} P \left( \frac{V_1}{V_1+V_2} \leq x, V_1+V_2 \leq 1 \right) &= P \left( Y \leq x, W \leq 1 \right)\\ &= \int_0^{x} \int_{0}^{1} \alpha \beta w^{\alpha + \beta -1} y^{\alpha -1 }(1-y)^{\beta - 1} dw dy \\ &= \frac{\alpha \beta}{\alpha + \beta} \frac{\Gamma(\alpha) \Gamma(\beta)}{\Gamma(\alpha + \beta)} I(x,\alpha,\beta) \end{array} $$

where $I(x,\alpha,\beta)$ is the incomplete beta function which is the standard beta distribution cdf according to Wikipedia. Next, the denominator in the expression for $F_X(x)$ is

$$ \begin{array}{ll} P \left( V_1+V_2 \leq 1 \right) &= P \left( W \leq 1 \right)\\ &= \int_0^{1} \int_{0}^{1} \alpha \beta w^{\alpha + \beta -1} y^{\alpha -1 }(1-y)^{\beta - 1} dw dy \\ &= \frac{\alpha \beta}{\alpha + \beta} \frac{\Gamma(\alpha) \Gamma(\beta)}{\Gamma(\alpha + \beta)} \end{array} $$

and thus we can finally conclude that $F_X(x) = I(x,\alpha,\beta)$.

Please let me know if there are any errors in the proof.