A quetion

June 21, 2010
By

(This article was first published on Shige's Research Blog, and kindly contributed to R-bloggers)

I posted a question to the R-Mixed-Model mailing list but have not had any responses yet:
---------------------------------------------------------------------------------------
Dear All,

I wonder if anybody has tried to make glmmBUGS work with JAGS. My
attempt was not successful.

Here is the simple example I copied from the glmmBUGS Vignettes:

------------------------
library(MASS)
data(bacteria)

bacterianew <- bacteria
bacterianew$yInt = as.integer(bacterianew$y == "y")

levels(bacterianew$trt) <- c("placebo", "drug", "drugplus")

library(glmmBUGS)

bacrag <- glmmBUGS(formula = yInt ~ trt + week, data=bacterianew,
effects = "ID",
modelFile="model.bug", family="bernoulli")

names(bacrag$ragged)

source("getInits.R")

startingValues = bacrag$startingValues
------------------------

With WinBUGS, it runs well:

------------------------
library(R2WinBUGS)

bacResult1 = bugs(bacrag$ragged, getInits, model.file="model.bug", n.chain=3,
n.iter=2000, n.burnin=100, parameters=names(getInits()), n.thin=10)
------------------------

But with JAGS, I got error message "Error in FUN(50L[[1L]], ...) :
invalid first argument"

------------------------
library(R2jags)

jags.parms=names(getInits())

bacResult = jags(data=bacrag$ragged, n.chain=3, n.iter=2000,
model.file=model.bug)
------------------------

I hope somebody can help me figuring out how to make this work. Many thanks.

Best,
Shige
---------------------------------------------------------------------------------------

Here is the answer from Jens Åströ
---------------------------------------------------------------------------------------
Hi!

Not an expert by any means, but I got it to run by doing this:

1) change inprod2 to inprod in the model file (JAGS does not have
inprod2 function)

2) change ~dflat() to other uninformative prior, e.g.
~dnorm(0.0,1.0E-6) (Jags does not have dflat distribution)

3) specify/compile the model with e.g.
bac.jags<-jags.model("model.bug",data=bacrag$ragged,n.chains=4)

4) update it, update(bac.jags,1000)

5) Collect coda samples,
bacResult<-coda.samples(temp,names(getInits()),n.iter=10000,thin=10)

This seemed to converge well, (gelman.diag(bacResult))

Good luck
----------------------------------------------------------------------------------------

It works.

To leave a comment for the author, please follow the link and comment on his blog: Shige's Research Blog.

R-bloggers.com offers daily e-mail updates about R news and tutorials on topics such as: visualization (ggplot2, Boxplots, maps, animation), programming (RStudio, Sweave, LaTeX, SQL, Eclipse, git, hadoop, Web Scraping) statistics (regression, PCA, time series, trading) and more...



If you got this far, why not subscribe for updates from the site? Choose your flavor: e-mail, twitter, RSS, or facebook...

Tags: , ,

Comments are closed.