4
$\begingroup$

I need to generate random numbers from Beta distribution using random variables from Uniform distribution.
If I have two random variables $Y_1=U_1^{1/\alpha}$ and $Y_2=U_1^{1/\beta}$,
and If $Y_1+Y_2<= 1$ then $X=Y_1/(Y_1+Y_2)$ is from Beta distribution with parameters $\alpha\ and\ \beta$.

This is what I've done so far:
I've generated 10,000 random variables from a Uniform distribution.
I then used those rv to generate $Y_1 and\ Y_2 $ using $Y_1=U_1^{1/\alpha}$ and $Y_2=U_1^{1/\beta}$.
Then I generated $X$ using the formula $X=Y_1/(Y_1+Y_2)$.

Here is the R code I used:

library(dplyr)
# my plot
alpha <- 2
beta <- 4

set.seed(10)
u <- runif(10000,0,1)

y1 <- u^(1/alpha)
y2 <- u^(1/beta)


x <- data.frame(y1,y2 ) %>% 
  filter(y1+y2<=1) %>% # check that y1+y1 <=1
  mutate(x = y1/(y1+y2)) %>% # generate x random variable
  select(x) 

# plot x
hist(x$x)

# plot from beta distribution
hist(rbeta(10000,shape = 2,shape2 =4)) 

But when I plot this, the histogram seems to be going in the opposite direction to when I just generate rv from Beta distribution. Here is my plot: enter image description here

And here is the plot if I generate rv from actual beta distribution:
enter image description here

Why does my curve look different to the one from Beta distribution?
This question seems very similar, but the solution uses qbeta() function. I think I need to use my $X=Y_1/(Y_1+Y_2)$

$\endgroup$
  • 3
    $\begingroup$ Hint: You need the $Y_i$ to be independent. $\endgroup$ – whuber May 26 '18 at 10:53
  • 4
    $\begingroup$ I believe it is a bit rude to downvote this question. The question is well explained and developed. The answer may be simple, but, it is a minor mistake anyone of us can do. $\endgroup$ – Jon Nagra May 27 '18 at 17:33
  • $\begingroup$ Please add the [self-study] tag $\endgroup$ – kjetil b halvorsen May 27 '18 at 19:13
  • 1
    $\begingroup$ @whuber I think I got it. I generated only 1 sample from the Uniform distribution, where as I need to get generate 2 If I change my code to this, it looks like it's working: y1 <- runif(10000,0,1)^(1/alpha) y2 <- runif(10000,0,1)^(1/beta) $\endgroup$ – jmich738 May 28 '18 at 9:40

Your Answer

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

Browse other questions tagged or ask your own question.