1

I'd like to create a reference sheet of common distributions for my statistical theory class, but I'm having some issues understanding R's implementation of the negative binomial distribution in the stats package. Here's the code I've been using to visualize the distribution.

ggplot(data.frame(x=c(0:30)), aes(x)) +
stat_function(geom="point", n=16, fun=dnbinom, 

              # Arguments to pass to dnbinom()
              args = c(p = 0.2, size = 2))

By my understanding of the NBD, p = 0.2 should be the probability of success (or failure, depending on your framing), and size = 2 should be the number of successes we are looking to observe. However, the output of the plot is this image.

I'm having trouble understanding how it's possible to have any probability mass at 0, before we've run any trials at all. I've read the documentation for dnbinom(), but I'm still not sure why this is the result, and the docs don't offer details on the exact equation that was implemented, or a sufficiently basic example that would show how the function works. Any advice/input would be much appreciated. Thanks!

conveniencesample
  • 33
  • 6
  • Hi Sycorax, my question concerns the built-in `stats` package of R, not the `MASS` package. This person's question is also highly specific to the user's dataset, while mine regards making a plot of the distribution itself. Hope that clears things up. – conveniencesample Sep 21 '18 at 22:23

1 Answers1

1

According to the documentation ?dnbinom, the distribution is that of the number of failures seen before the requested number of successes are seen. So it's possible that the number of failures is observed to be zero, which would occur exactly when each of the first size trials is a success. Such an event has probability $p^\text{size}$, which agrees with the formula below (with n=size):

 The negative binomial distribution with ‘size’ = n and ‘prob’ = p
 has density

                  Gamma(x+n)/(Gamma(n) x!) p^n (1-p)^x              

 for x = 0, 1, 2, ..., n > 0 and 0 < p <= 1.

 This represents the number of failures which occur in a sequence
 of Bernoulli trials before a target number of successes is
 reached.  The mean is mu = n(1-p)/p and variance n(1-p)/p^2.
grand_chat
  • 2,632
  • 1
  • 8
  • 11
  • Thank you for your quick reply, grand_chat, this cleared things up right away. To drive the point home, I calculated 0.2 * 0.2 = 0.04, the probability of observing no failures before observing 2 successes, and this aligns with the distribution's value at 0. – conveniencesample Sep 21 '18 at 19:53