[R-sig-ME] NegBin mixed models and gamlss.mx

[David Atkins]

Hi all--

I got the second edition of Hilbe's book on negative binomial 
regression, which now includes R code in addition to Stata, and I was a 
bit surprised to see code to fit negative binomial mixed models using 
the gamlss.mx package.  To be clear, I just wasn't aware of it, not that 
I've got anything against that package.

Given the somewhat regular traffic here about NB mixed models (and lack 
of R software for fitting them), I am guessing this might be flying 
under the radar generally.  (Or, it could just be me!)

Below is a little code showing how to fit such models and comparison 
with glmer() using over-dispersed Poisson mixed model.  Note that you 
are restricted to "two level" models, in the multilevel lingo. 
(Possibly also uncorrelated random-effects, I haven't played with models 
too much yet...)

cheers, Dave

### RAPI
rapi.df <- 
read.csv("http://depts.washington.edu/cshrb/newweb/stats%20documents/RAPI.Final.csv", 
header = TRUE)

### Table 2 - Model 1
rapi.df$over <- 1:nrow(rapi.df)

library(lme4)
rapi.glmer <- glmer(rapi ~ gender + time + (1 | id) + (1 | over),
					data = rapi.df, verbose = TRUE,
					family = poisson)
summary(rapi.glmer)					

### Model with random-effect for time
rapi.glmer2 <- glmer(rapi ~ gender + time + (time|id) + (1|over),
					data = rapi.df, verbose = TRUE,
					family = poisson)
summary(rapi.glmer2)

### try gamlss
library(gamlss.mx)

rapi.nbmix <- gamlssNP(rapi ~ gender + time,
					data = rapi.df,
					random = ~ 1 | id,
					family = NBI,
					mixture = "gq", k = 20)
summary(rapi.nbmix)					
summary(rapi.glmer)

rapi.nbmix2 <- gamlssNP(rapi ~ gender + time,
					data = rapi.df,
					random = ~ time | id,
					family = NBI,
					mixture = "gq", k = 20)
summary(rapi.nbmix2)
summary(rapi.glmer2)

-- 
Dave Atkins, PhD
Research Associate Professor
Department of Psychiatry and Behavioral Science
University of Washington
datkins at u.washington.edu

Center for the Study of Health and Risk Behaviors (CSHRB)		
1100 NE 45th Street, Suite 300 	
Seattle, WA  98105 	
206-616-3879 	
http://depts.washington.edu/cshrb/
(Mon-Wed)	

Center for Healthcare Improvement, for Addictions, Mental Illness,
   Medically Vulnerable Populations (CHAMMP)
325 9th Avenue, 2HH-15
Box 359911
Seattle, WA 98104
http://www.chammp.org
(Thurs)