[R-sig-ME] Does MCMCglmm allow offset?

[David Atkins]

Hi all--

I am wondering whether MCMCglmm will take an offset?  Specifically, we 
have a study with number of drinks as the outcome, but we also have an 
estimate of over how many hours those drinks were consumed.  Thus, might 
make sense to model it as a rate of drinks per hour.

As far as I can tell, glmer does allow an offset, but it seems like 
maybe MCMCglmm does not.  A search of MCMCglmm's course notes and JSS 
paper does not turn up any discussion.

Below is a real dataset, with a "toy" offset term that would seem to 
support this.

Any input appreciated.

cheers, Dave

[sessionInfo() at very end]

### MCMCglmm and offset
#
### read data
rapi.df <- 
read.csv("http://depts.washington.edu/cshrb/newweb/stats%20documents/RAPI.Final.csv")
head(rapi.df)

### "rapi" is the rutgers alcohol problems index and was measured
### over previous 30 days -- thus, constant exposure, but we can still
### use it as toy example
#
### over-dispersion term
rapi.df$over <- 1:nrow(rapi.df)

### make up "offset" term
rapi.df$off30 <- rpois(n = nrow(rapi.df), lambda = 30)

### GLM (no random-effects) just as baseline comparator
glm.nooff <- glm(rapi ~ gender*time,
						data = rapi.df, family = "poisson")
glm.off <- glm(rapi ~ gender*time + offset(log(off30)),
						data = rapi.df, family = "poisson")
summary(glm.nooff)						
summary(glm.off)						

### over-dispersed Poisson GLMM in glmer
library(lme4)

### without offset
glmer.nooff <- glmer(rapi ~ gender*time + (time | id) + (1 | over),
						data = rapi.df, verbose = TRUE, family = "poisson")

### with offset
glmer.off <- glmer(rapi ~ gender*time + offset(log(off30)) + (time | id) 
+ (1 | over),
						data = rapi.df, verbose = TRUE, family = "poisson")
print(glmer.nooff, corr = FALSE)	
print(glmer.off, corr = FALSE)	

### same models using MCMCglmm
library(MCMCglmm)

mcmc.nooff <- MCMCglmm(rapi ~ gender*time, data = rapi.df,
								family = "poisson", verbose = TRUE,
								random = ~ us(1 + time):id)
mcmc.off <- MCMCglmm(rapi ~ gender*time + offset(log(off30)), data = 
rapi.df,
								family = "poisson", verbose = TRUE,
								random = ~ us(1 + time):id)
summary(mcmc.nooff)					
summary(mcmc.off)	


 > sessionInfo()
R version 2.12.1 (2010-12-16)
Platform: i386-apple-darwin9.8.0/i386 (32-bit)

locale:
[1] en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods
[7] base

other attached packages:
[1] MCMCglmm_2.10      corpcor_1.5.7      ape_2.6-2
[4] coda_0.14-2        tensorA_0.36       lme4_0.999375-37
[7] Matrix_0.999375-47 lattice_0.19-16

loaded via a namespace (and not attached):
[1] gee_4.13-16   grid_2.12.1   nlme_3.1-97   stats4_2.12.1
[5] tools_2.12.1

-- 
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)