CODEWITHSUNDEEP
ralgorithmrandomstatisticssimulation
user9802913
user9802913

Reputation: 245

Simulation of multiple binomial random numbers in R

I have the following algorithm

Step 1. Generate X1=x1~Bin(6,1/3)

Step 2. Generate X2|X1=x1~Bin(6-x1,(1/3)/(1-1/3))

Step 3. Generate X3|X1=x1,X2=x2~Bin(6-x1-x2,(1/3)/(1-1/3-1/3))

Step 4. Repeat step 1-3 N times.

Here is my approach to implement this algorithm in R:

mult_binom<-function(n) #n=6
{
  n=1000
  random_vectors<-Matrix(0,n,3)

  for(i in 1:n){

    X1<-rbinom(n,3,1/3) 

    X2<-rbinom(n-X1,3,(1/3)/(1-(1/3)))

    X3<-rbinom(n-X1-X2,3,(1/3)/(1-(1/3)-(1-3)))


    arr<-c(X1,X2,X3)


  }
  for(j in 1:n){

    random_vectors[j]<-arr[j]
  }
  return(random_vectors)
}

Calling the function as mult_bin(6) yields a similar matrix as the below one

1000 x 3 sparse Matrix of class "dgCMatrix"

  [1,] 1 . .
  [2,] 1 . .
  [3,] 1 . .
  [4,] 2 . .
  [5,] 1 . .
  [6,] 1 . .
  [7,] 1 . .
  [8,] . 3 .

and continues up to [1000,]

I was not expecting this result.

Why are dots there?

What did I do wrong?

Upvotes: 1

Views: 480

Answers (1)

Cettt
Cettt

Reputation: 11981

there are several errors in your implementation. The most important one is that the first argument of rbinom is not the parameter n from the binomial distribution but the number of random numbers you want to generate.

Here is my solution. My function only returns on experiment. I then use replicate to return the result of multiple (in my case 5) experiments:

myfun <- function(){

  x1 <- rbinom(1, 6, 1/3) 
  x2 <- rbinom(1, 6 - x1, (1/3)/(1-(1/3)))
  x3 <- rbinom(1, 6 - x1 - x2, (1/3)/(1-(1/3)-(1/3)))

  return(c(X1 = x1, X2 = x2, X3 = x3))  
}

    set.seed(1)
    replicate(5, myfun())
   [,1] [,2] [,3] [,4] [,5]
X1    1    4    4    0    3
X2    2    0    1    2    1
X3    3    2    1    4    2

In this output every column is the result of one experiment. You can see that the numbers always add up to 6. Note also, that I set a random seed with set.seed. This ensures that your results are reproducible.

In your output the dots appear because you use the Matrix package to create a Matrix object instead of working with "normal" matrices. Usually you create a matrix with matrix not Matrix.

Upvotes: 3

Related Questions