From: Bayesian Models for Astrophysical Data, Cambridge Univ. Press

(c) 2017,  Joseph M. Hilbe, Rafael S. de Souza and Emille E. O. Ishida  

 

you are kindly asked to include the complete citation if you used this material in a publication

Code 6.9 Synthetic negative binomial data and model in R

====================================================

library(MASS)

set.seed(141)
nobs <- 2500

x1 <- rbinom(nobs,size=1, prob=0.6)
x2 <- runif(nobs)
xb <- 1 + 2.0*x1 - 1.5*x2


a <- 3.3

theta <- 0.303                                                                 # 1/a


exb <- exp(xb)

nby <- rnegbin(n=nobs, mu=exb, theta=theta)
negbml <- data.frame(nby, x1, x2)

nb2 <- glm.nb(nby ~ x1 + x2, data=negbml)

summary(nb2)

====================================================

Output on screen:

Call:

glm.nb(formula = nby ~ x1 + x2, data = negbml, init.theta = 0.295684269,

             link = log)

 

Deviance Residuals:

        Min                     1Q        Median                   3Q Max

-1.58697            -1.12176      -0.78956           0.06661 3.05031

 

Coefficients:

                      Estimate Std.        Error       z value        Pr(>|z|)

(Intercept)             0.98613     0.08845          11.15       <2e-16 ***

x1                          2.03995     0.08081          25.24       <2e-16 ***

x2                        -1.61330     0.13544          -11.91       <2e-16 ***

---

Signif.    codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

 

(Dispersion parameter for Negative Binomial(0.2957) family taken to be 1)

Null deviance:            3030.9     on      2499     degrees     of     freedom

Residual deviance:     2355.5     on      2497     degrees     of     freedom

AIC: 11553

 

Number of Fisher Scoring iterations: 1

 

                   Theta:       0.2957

               Std. Err.:      0.0108

 

2 x log-likelihood:     -11544.6130