[R] How to generate a matrix of Beta or Binomial distribution

[R. Michael Weylandt]

On Mon, Aug 27, 2012 at 9:49 AM, JeffND <Zuofeng.Shang.5 at nd.edu> wrote:
> Hi folks,
>
> I have a question about how to efficiently produce random numbers from Beta
> and Binomial distributions.
>
> For Beta distribution, suppose we have two shape vectors shape1 and shape2.
> I hope to generate a 10000 x 2 matrix X whose i th rwo is a sample from
> reta(2,shape1[i]mshape2[i]). Of course this can be done via loops:
>
> for(i in 1:10000)
> {
> X[i,]=rbeta(2,shape1[i],shape2[i])
> }
>
> However, the above code is time consuming. It would be better to directly
> generate X without any loops. I tried the following code
>
> X<- rbeta(2,shape1,shape2)

No, not quite: this only makes 2 random variates, while before you
made 2*10,000 = 20k. You can make them all in one fell swoop like
this:

X <- rbeta(20000, shape1, shape2)

Now you might ask, "shape1" and "shape2" are only 10k long: how does
this work? R will recycle shape1 and shape2 as needed, so you get a
call as if you actually used

rbeta(20000, c(shape1, shape1), c(shape2, shape2))

or more generally

rbeta(20000, rep(shape1, length.out = 20000), rep(shape2, length.out = 20000))

If I understand your code correctly, you want each row to use the same
shape1/shape2 parameters, so we will need to reshape the result of
rbeta() into a matrix: it turns out that this is actually pretty easy
to do for your case:

matrix(rbeta(20000, shape1, shape2), ncol = 2)

R will get the correct number of rows by default (it does the 20000
elements / 2 columns = 10000 rows division internally) and because R
uses column-major ordering, you get the parameter arrangement by happy
coincidence. (Basically, R fills up the first column first, then the
second; this aligns nicely with the recycling behavior in this
problem, so we get back to shape1[1] as soon as we start filling the
second column)

This is a little bit of happy hackery, and you might instead prefer
something more explicit like this:

matrix(rbeta(20000, rep(shape1, each = 2), rep(shape2, each = 2)),
ncol = 2, byrow = TRUE)

which will repeat shape1 and shape2 elementwise two times (so c(1,2,3)
becomes c(1,1,2,2,3,3)) and then fills rowwise. It's a little clearer,
at the expense of verbosity.

>
> but it looks not right. So I was wondering if there is an R code doing this.
>
> For Binomial distribution, I have a similar question.  I hope to generate a
> 10000 x n matrix Y whose i th row is a sample from rbinom(n,2,p[i]), where p
> is a vector of binomial probabilities. We can also do it by loops
>
> for(i in 1:10000)
> {
> Y[i,]=rbinom(n,2,p[i])
> }

This should follow pretty easily from the above.

Hope this helps,
Michael

>
> But it would be nice to do it without using any loops. Is this possible in
> R?
>
> Many thanks for great help and suggestions!
> Jeff
>
>
>
>
>
>
>
>
>
> --
> View this message in context: http://r.789695.n4.nabble.com/How-to-generate-a-matrix-of-Beta-or-Binomial-distribution-tp4641422.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.