1

I'm trying to write a function in R which produces a set of n i.i.d Bernoulli random variables X1,X2,...,Xn ∼ B(1,π).

To generate M = 1000 random samples. When n = 10 and pi = 0.4.

I'm having trouble starting the code. If anyone could help out that would be greatly appreciated.

Here is what I have run to this point:

rbinom(n=10,prob=.4,size=1)

Now I want to assemble those variables into a matrix of 1000 of those samples run through the rbinom function.

pc8807
  • 47
  • 1
  • 6
  • Have you gotten anywhere at all? `?rbinom` and `?runif` could both be useful, as could `?replicate` (although the latter is not really necessary). Can you post an example of the structure of your desired output (e.g. list, row-wise matrix, column-wise matrix, ... ?), with small values of `M` and `n` for convenience? – Ben Bolker Sep 26 '16 at 21:58
  • I appologize, I'm in the process of learning R right now. I have run rbinom and recieved the output: [1] 26 47 38 48 37 42 42 44 39 49. I want that output stored in a column wise matrix, which is where I am running into the issues. Which should have 1000 total samples created through the rbinom command – pc8807 Sep 26 '16 at 22:05
  • no need to apologize. Can you *edit your question* to specify *how* you ran `rbinom` (i.e., what arguments did you use? Fill in the blanks: `rbinom(n=___,prob=____,size=___)` You can use `matrix()` to convert results to a matrix, but your results don't look like Bernoulli draws to me (hint, `size=1`)! – Ben Bolker Sep 26 '16 at 22:06
  • Thank you for the advice. I had the size all wrong. I'm learning how to post appropriate questions to the site, and since I'm new with R I don't know how to frame the question appropriately yet. – pc8807 Sep 26 '16 at 22:16

2 Answers2

2

Taking what you've done so far:

set.seed(101) ## for reproducibility
rr <- replicate(1000,rbinom(n=10,prob=.04,size=1))

creates 1000 samples.

dim(rr)
## [1] 10 1000
mean(rr)
## [1] 0.0411

i.e. this is 10 rows x 1000 columns; you can use t() to transpose it if you like. If you specify simplify=FALSE in the replicate() call you'd get a list of vectors of length 10 instead.

(Your question is inconsistent as to whether you want probability of 0.04 or 0.4, but you should be able to sort that out for yourself.)

For this particular question you could also use

matrix(rbinom(10000,prob=0.04,size=1),nrow=10)

or

matrix(as.numeric(runif(10000)<0.04),nrow=10)
Ben Bolker
  • 211,554
  • 25
  • 370
  • 453
  • @ Ben, Can you help to this https://stackoverflow.com/questions/59595292/simulation-bernoulli-data-in-r –  Jan 04 '20 at 22:35
0

A fast solution is to use sample instead of runif or rbinom

n <- 10L
m <- 1000L
p <- .04

set.seed(1)
microbenchmark::microbenchmark(
  replicate               = replicate(m, rbinom(n = n, prob = p, size = 1)), 
  `matrix rbinom`         = matrix(rbinom(m * n, prob = p, size = 1), nrow = n),
  `matrix runif`          = matrix(as.numeric(runif(m * n) < p), nrow = n), 
  `matrix sample`         = matrix(
    sample(1:0, size = m * n, replace = TRUE, prob = c(p, 1 - p)), nrow = n),
  `matrix sample logical` = matrix(
    sample(c(T, F), size = m * n, replace = TRUE, prob = c(p, 1 - p)), nrow = n))
#R Unit: microseconds
#R                   expr      min        lq       mean    median        uq       max neval
#R              replicate 2453.795 2499.4770 3397.60451 2539.1285 2583.8825 70546.023   100
#R          matrix rbinom  347.826  349.6820  359.90301  357.5655  359.4210   424.812   100
#R           matrix runif  288.464  290.7830  301.35238  297.0440  299.1310   406.262   100
#R          matrix sample  105.276  111.0735  116.48064  113.1600  116.4065   192.928   100
#R  matrix sample logical   83.943   86.0295   89.31318   87.1890   88.5800   135.885   100
Benjamin Christoffersen
  • 4,703
  • 15
  • 37