Follow

Keep Up to Date with the Most Important News

By pressing the Subscribe button, you confirm that you have read and are agreeing to our Privacy Policy and Terms of Use
Contact

How to generate a negative binomial distribution with different sample sizes for power analyses in R?

I try to do a power simulation with an outcome variable that is zero-inflated. So I use a negative binomial distribution. What I need is the following distribution:

library(tidyverse)
set.seed(123)
rt_random <- 
rnbinom(n = 2000, mu = 25, size = .9) 
qplot(rt_random)

desired distribution

Now, I tried to create a function where I can use different sample sizes for the power simulation. I use the pmapfunction for this. However this does not seem to work. The final distribution is not zero-inflated at all and the mean is not close to the defined mean:

MEDevel.com: Open-source for Healthcare and Education

Collecting and validating open-source software for healthcare, education, enterprise, development, medical imaging, medical records, and digital pathology.

Visit Medevel

generate_distribution <- function(n, mus, sizes){
  tibble(n = n, 
       t_mu = mus, 
       t_size = sizes)%>% 
  mutate(rt = pmap(list(n, t_mu, t_size),
                     ~rnbinom(..1, ..2, ..3))) %>% 
  unnest(rt) 
  
}
set.seed(123)
rt_df <- 
  generate_distribution(n = 2000,
                        mus = 25,
                        sizes = .9)
qplot(rt_df$rt)

wrong distribution

Is there an easy way to change my code so that I get the desired distribution?

>Solution :

The problem is that in your lambda function, you are not using the argument names for rnbinom. The default order for the arguments if you don’t name them is n, size, prob, mu, so you are passing 2000 to n, 25 to size and 0.9 to prob. Just name the arguments explicitly as you did in your first example, and your code will work.

generate_distribution <- function(n, mus, sizes){
  tibble(n = n, 
       t_mu = mus, 
       t_size = sizes)%>% 
  mutate(rt = pmap(list(n, t_mu, t_size),
                     ~rnbinom(n = ..1, mu = ..2, size = ..3))) %>% 
  unnest(rt) 
}

set.seed(123)
rt_df <- 
  generate_distribution(n = 2000,
                        mus = 25,
                        sizes = .9)
qplot(rt_df$rt) 

enter image description here

Add a comment

Leave a Reply

Keep Up to Date with the Most Important News

By pressing the Subscribe button, you confirm that you have read and are agreeing to our Privacy Policy and Terms of Use

Discover more from Dev solutions

Subscribe now to keep reading and get access to the full archive.

Continue reading