This is a follow-up to the question Square root of a Beta(1,1) random variable, which received two great answers.

If $XY \sim \text{Beta}(\alpha, \beta)$, and $X$ and $Y$ are two independent identically distributed random variables, do we know their distribution? If there is no simple general solution, can we figure out the distribution of $X$ and $Y$ for the case $XY \sim \text{Beta}(1, 1)$?

From my search so far I guess this is a deconvolution problem for which no readily worked-out solution, such as the Kumaraswamy distribution in the previous question, exists. However, I would appreciate to be proven wrong in this regard.


1 Answer 1


In the case of $Z=XY\sim \text{Beta}(\alpha,1)$, the moment generating function (mgf) of $-\ln(XY)=-\ln X-\ln Y$ is \begin{align} M_{-\ln(XY)}(t) &= E(e^{-t\ln Z}) \\ &=E(Z^{-t}) \\ &=\int_0^1 z^{-t}\alpha z^{\alpha-1}dz \\ &= \frac{1}{1-t/\alpha} \end{align} which is the mgf of an exponentially distributed random variable with rate parameter $\alpha$.

Thus, if $X$ and $Y$ are iid, this means that $-\ln X$ and $-\ln Y$ must be Gamma distributed with shape parameter $1/2$ and rate parameter $\alpha$ (see wikipedia). Backtransforming, the pdf of $X$ and $Y$ is $$ f(x)=\sqrt{-\frac\alpha{\pi\ln x}}x^{\alpha-1} $$ for $0\le x\le 1$.

A similar calculation in the general Beta$(\alpha,\beta)$-case leads to $$ M_{-\ln(XY)}(t)=\frac{\Gamma(\alpha-t)\Gamma(\alpha+\beta)}{\Gamma(\alpha+\beta-t)\Gamma(\alpha)} $$ and the mgf of $-\ln X$ and $-\ln Y$ would need to be the square root of this. But that is not the mgf of any known distribution and such a distribution (that when convolved with itself produces the target density) may not even exist.

The pdf of $-\ln X$ can be computed via numerical integration of the inversion formula for the characteristic function $\varphi(t)=(M_{-\ln(XY}(it))^{1/2}$. This seems to produce valid pdfs for $\beta\le 2$ but for $\beta > 2$, the resulting function no longer appear to be a valid pdf:

enter image description here

R code:

f <- function(x,phi,...) {
  integrand <- function(t) {
f <- Vectorize(f,vectorize.args = "x")
phi <- function(t, alpha, beta) {
  sqrt(gammaz(alpha - 1i*t)*gammaz(alpha + beta)/(gammaz(alpha + beta - 1i*t)*gamma(alpha)))
x <- seq(-1,3,by=.01)
fx <- f(x, phi, alpha=3,beta=2.5)
  • $\begingroup$ This looks extremely useful. Can you maybe add a hint, how you derived the relationships in the first paragraph? Or how they are called? $\endgroup$
    – LuckyPal
    May 4, 2021 at 9:08
  • $\begingroup$ @LuckyPal I added some more explanation and made the result slightly more general. $\endgroup$ May 4, 2021 at 9:37
  • $\begingroup$ "but for β>2, the resulting function no longer appear to be a valid pdf" Is it not valid because of the negative values? But, those negative values only appear outside of the support $0<x<1$. $\endgroup$ Sep 20, 2021 at 3:35

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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