# Sample and rbinom generating disparate results?

I am sure I am just misunderstanding this but rbinom and sample seem to be generating disparate results and I am not sure why. I am trying to create a categorical factor with a balanced third in each levels (so .33 get 0, .33 get 1, .33 get 2).

If I attempt to sample evenly 1/3 from rbinom(), I get uneven distributions:

table(rbinom(100,2,c(.33,.33,.33)))
0  1  2
47 42 11
0  1  2
49 41 10
0  1  2
46 40 14
0  1  2
45 43 12
0  1  2
46 38 16


The third level is not getting many successes and I am not sure why.

If I do the same thing for the sample() function, I get more satisfactory results.

table(sample(c(1,2,3),100,replace = T,prob = c(.33,.33,.33)))
1  2  3
32 32 36
1  2  3
39 33 28
1  2  3
30 41 29
1  2  3
33 34 33
1  2  3
36 25 39
1  2  3
39 31 30


Here, each level (1, 2, and 3) is getting around 33% of successes which is what I wanted.

How do I get rbinom() to do the same thing, distributing successes evenly across 1, 2, and 3?

• You want the the multinomial, not the binomial distribution. Try rmultinom(5, 100, c(.33,.33,.33)). Commented Jun 21, 2022 at 6:06
• for rbinom, the probability of getting a 2 is from a size 2 with success of 0.33 is 0.1089. Hence 11% of the time you will get a 2. and that is exactly what you have. dbinom is not same as discrete uniform. You should use rmultinom as suggested or use sample as you are doing so far Commented Jun 21, 2022 at 6:52

As noticed in the comments, your code does not do what you want it to do. First, many of R's functions vectorized the passed parameter values, in case of your code rbinom(100,2,c(.33,.33,.33)) is the same as rbinom(100,2,.33), or rbinom(100,2,c(.33,.33)), or rbinom(100,2,rep(.33, 100)) because all the values provided to prob are the same, if they differed, it'd cycle through them until drawing the intended 100 samples.

The code itself produces draws from a binomial distribution with sample size=2, so you are independently tossing a biased coin with prob=0.33 two times. Your second snippet draws each of the values in c(1,2,3) with an equal probability of 0.33. So the first code would produce the result 2 with probability 0.33 * 0.33 = 0.1089 (two independent coin throws), while the second the result 3 with a probability of 0.33 (you roll a three-sided dice). In the first case, the probability distribution is

$$\Pr(X=k) = {2 \choose k} \, 0.33^k \,0.66^{2-k}$$

while in the second case, it's just $$\Pr(X=k) = 0.33$$. They are not the same.

This is what you observe if you run the code.

> n <- 1000000
> table(rbinom(n, 2, .33)) / n

0        1        2
0.448398 0.443002 0.108600
> table(rbinom(n, 2, c(.33,.33,.33))) / n

0        1        2
0.449160 0.441886 0.108954
> table(sample(c(1,2,3), n, replace = T, prob = c(.33,.33,.33))) / n

1        2        3
0.333505 0.332922 0.333573


You can achieve the same result as when using rbinom with sample if you use this unnecessarily slow code:

> table(replicate(n, sum(sample(c(0, 1), size=2, replace=T, prob=c(.66,.33))))) / n

0        1        2
0.444294 0.444688 0.111018


When calling

rbinom(100,2,c(a,b,c))


the outcome is equivalent to

c(rbinom(2,a),rbinom(2,b),rbinom(2,c),rbinom(2,a),rbinom(2,b),...)