Let $0 \leq \gamma < 1$, $X_i \sim \text{Beta}(\alpha, \beta)$, and

$$Y \sim \sum_{i = 0}^\infty \gamma^i X_i$$

What is the distribution of $Y$? Does it have a closed form? Can it be sampled efficiently? We know that \begin{align} \operatorname{supp}(Y) &= \left[0, \frac{1}{1-\gamma}\right) \\ \operatorname{E}(Y) &= \frac{\operatorname{E}(X)}{1-\gamma} \\ \operatorname{Var}(Y) &= \frac{\operatorname{Var}(X)}{1-\gamma^2} \\ \end{align}

Also, note that \begin{align} Y = X + \gamma Y' \end{align} where $Y$ and $Y'$ have the same distribution. Therefore, denoting the characteristic function of a random variable $A$ as $\phi_A$, \begin{align} \phi_Y(t) &= \phi_{X + \gamma Y'}(t) \\ &= (\phi_X \cdot \phi_{\gamma Y'})(t) \\ &= \phi_X(t) \cdot \phi_{\gamma Y'}(t) \\ &= \phi_X(t) \cdot \phi_{Y'}(\gamma t) \\ &= \phi_X(t) \cdot \phi_Y(\gamma t) \\ \frac{\phi_Y(t)}{\phi_Y(\gamma t)} &= \phi_X(t) \\ &= {}_1F_1(\alpha; \alpha + \beta; i t) \end{align}

where the last RHS is the characteristic function of the beta distribution. Is there a way to solve for $\phi_Y(t)$ or at least sample efficiently from $Y$'s distribution?

  • 2
    $\begingroup$ This is another case of the situation considered by Vincent Granville at stats.stackexchange.com/questions/432136. The connection is this: let $\phi:\mathbb{R}_{+}\to\mathbb{R}_{+}$ be $\phi(x)=\sqrt{x+1}-1$ for him and $\phi(x)=\gamma x$ for you, so that in both cases $0$ is the sole attractor of $\phi.$ Take a sequence of iid positive random variables $X_1,X_2,X_3,\ldots$ and consider the distribution of $\phi(X_1+\phi(X_2+\phi(X_3+\cdots))).$ These two situations are thereby quite alike, indicating an analysis of one will suggest an effective analysis of the other. $\endgroup$
    – whuber
    Dec 4 '19 at 4:55
  • $\begingroup$ I noticed, thanks. I had also worked out the relationship between the cfs, but stopped there because it doesn't appear to get us any closer to a solution. You could pursue it a little further to deduce the cumulants of $Y$--but then, when you explore the possibility of writing a closed formula even for the cumulants of a Beta distribution, you realize that's not going to be useful either. $\endgroup$
    – whuber
    Dec 4 '19 at 15:24

Note I've used slightly different notation to you.

Firstly, let's consider the underlying random variables we are dealing with. We have that all


where all random variables are independent of each other. If we now consider some scaled form of these random variables:

$$Y=a X$$

where we have that $a<1$ because $\gamma<1$. Now, for a continuous random variable

$$\begin{align} f_{Y}(y)&=\frac{f_{X}(y/a)}{a}\\ &=\frac{(y/a)^{\alpha-1}(1-y/a)^{\beta-1}}{a\cdot\text{B}(\alpha,\beta)} \end{align}$$

Notice that the support has changed from $$0<x<1$$ to $$0<y<a$$

If we define


as a variant on the usual Beta function (this function equals the Beta function for $a=1$), then we can show that


This means that

$$\begin{align} f_{Y}(y)&=\frac{(y/a)^{\alpha-1}(1-y/a)^{\beta-1}}{a\cdot\text{B}(\alpha,\beta)}\\ &=\frac{y^{\alpha-1}(a-y)^{\beta-1}}{\text{B}'(\alpha,\beta,a)}\cdot\frac{\text{B}'(\alpha,\beta,a)}{a^{\alpha+\beta-1}\cdot\text{B}(\alpha,\beta)}\\ &=\frac{y^{\alpha-1}(a-y)^{\beta-1}}{\text{B}'(\alpha,\beta,a)} \end{align}$$

For the purposes of this answer, let's denote this random variable as

$$Y\sim\text{Beta}'(\alpha,\beta,a)$$ This is essentially a scaled Beta distribution, where the support has been limited to $[0,a]$. The important point here is that random values generated from this random variable will be bound between $[0,a]$.

Now, in your summation, each Beta distribution $X_{i}$ is scaled by $\gamma^{i}$ for increasing $i$. This means that each summand will have an increasingly reduced support $[0,\gamma^{i}]$ and each subsequent summand's probability mass will be increasingly narrow.

We know this series will converge for $0\leq\gamma<1$ because, for each $i$, the maximum value that each subsequent random value can take is $\gamma^{i}$. So, at most, the series is


which is a simple geometric series converging to


Interestingly, this tells us the upper bound of the support of $Y$ for a given $\gamma$. If we define


then the support for $Z$ would then be $[0,1/(1-\gamma)]$.

Now, you can derive the moment generating function of $Y$ as

$$M_{Y_{i}}(t)=1+\sum_{k=1}^{\infty}\Bigg(\prod_{r=0}^{k-1}\bigg(\frac{\alpha+r}{\alpha+\beta+r}\bigg)\frac{(\gamma^{i}\cdot t)^{k}}{k!}\Bigg)$$

Then the moment generating function of $Z$ is

$$\begin{align} M_{Z}(t)&=\prod_{i=0}^{\infty}M_{Y_{i}}(t)\\ &=\Bigg[1+\sum_{k=1}^{\infty}\Bigg(\prod_{r=0}^{k-1}\bigg(\frac{\alpha+r}{\alpha+\beta+r}\bigg)\frac{t^{k}}{k!}\Bigg)\Bigg]\cdot \\ &\Bigg[1+\sum_{k=1}^{\infty}\Bigg(\prod_{r=0}^{k-1}\bigg(\frac{\alpha+r}{\alpha+\beta+r}\bigg)\frac{(\gamma\cdot t)^{k}}{k!}\Bigg)\Bigg]\cdots\\ \end{align}$$

I'm not sure how useful this is. But it's useful to note that for increasing $i$ the subsequent terms will decay in size so that $M_{Y_{i}}(t)\rightarrow 1$ for large $i$.

If we look at some simulations of this distribution, there are some interesting results. First we'll define some functions for $Y$.

    dsbeta=function(x, alpha, beta, gamma) {
                alpha+beta-1)*beta(alpha, beta)))
    psbeta=function(q, alpha, beta, gamma) {
      return(pbeta(q/gamma, alpha, beta))
    qsbeta=function(p, alpha, beta, gamma) {
      return(gamma*qbeta(p, alpha, beta))
    rsbeta=function(n, alpha, beta, gamma) {
      return(gamma*rbeta(n, alpha, beta))

    #Perform some simulations:
    gamma=0.2 #Scaling of each X
    N=10000   #Number of samples of Z
    M=100     #Number of Y in each Z
    #Generate random values of Z:
    r=as.vector(matrix(rbeta(N*M, a, b), N, M) %*%  
    #Final simulations:
    hist(r, freq = FALSE, breaks = 50)
    #Perform fit (here we provide the upper bound gamma as the 
    #             theoretical upper limit):
    f=fitdistrplus::fitdist(r, "sbeta", start = list(alpha = 2, 
         beta = 3), fix.arg = list(gamma = 1/(1-gamma)))
    #Visualize the fit:
    x=seq(min(r), max(r), length.out = 1000)
    lines(x, dsbeta(x, f$estimate[1], f$estimate[2], 
             1/(1-gamma)), lwd = 2, col = "blue")

    #GOF tests:
    ks.test(r, "psbeta", f$estimate[1], f$estimate[2],  
    goftest::cvm.test(r, "psbeta", f$estimate[1], f$estimate[2], 
goftest::ad.test(r, "psbeta", f$estimate[1], f$estimate[2], 
    #qq plot:
    plot(sort(r), qsbeta(ppoints(length(r)), f$estimate[1], 
     f$estimate[2], f$fix.arg$gamma), ylab = "Theoretical 
         quantiles", xlab = "Empirical quantiles", col = "blue", 
                     type = "l", lwd = 2)
    abline(0, 1, lty = 2, lwd = 2, col = "red")

It seems from the simulations that $Z$ can be well approximated by the scaled beta


but as per @whuber's comments, for small $\alpha,\beta$ and $\gamma\rightarrow 0$ this breaks down. I think it has something to do with the bi-modality for $\alpha,\beta<1$ and the fewer 'non-degenerate' random variables due to $\gamma$ small.

  • $\begingroup$ I think you will get more insight into the underlying issues by considering cases where at least one of $\alpha,\beta$ is close to zero (preferably both). This will more clearly show the limitations of approximating $Z$ as a scaled Beta distribution. Focus on relatively small $\gamma,$ too: for large $\gamma \approx 1$ there's a CLT effect. $\endgroup$
    – whuber
    Dec 5 '19 at 14:15
  • $\begingroup$ @whuber I'm only seeing a poor scaled Beta fit for $\alpha\approx\beta\approx 1$ and $\gamma\rightarrow 0$ where the simulated data looks quite uniform. $\endgroup$
    – epp
    Dec 5 '19 at 22:28
  • $\begingroup$ @whuber I wonder whether it can be proven that it approaches a scaled beta as $\gamma \rightarrow 1$ (or as $\alpha, \beta \rightarrow \infty$). $\endgroup$
    – user76284
    Dec 8 '19 at 20:15
  • $\begingroup$ Under those conditions, properly standardized, the sum ought to approach a Normal distribution. $\endgroup$
    – whuber
    Dec 8 '19 at 21:48

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.