Monday, February 27, 2012

JAGS code for Exercise 12.1 (initializing vectors with nonstochastic components)

Exercise 12.1 provides code snippets for modifying the BUGS code, but not for modifying the JAGS code. Here are the corresponding JAGS modifications.The primary issue here is how to initialize vectors that have some stochastic and some nonstochastic components.

In the model specification, change the prior on the mu[j] as indicated in the exercise:

# Hyperpriors for mu and kappa:
## Four mu version:
#for ( j in 1:nCond ) {
# mu[j] ~ dbeta( a[j,mdlIdx] , b[j,mdlIdx] )
# Two mu version:
for ( j in 1:2 ) {
  mu[j] ~ dbeta( a[j,mdlIdx] , b[j,mdlIdx] )
mu[3] <- mu[2]
mu[4] <- mu[2]

From p. 315 of the book: "... when initializing chains in BUGS [and JAGS!], only stochastic variables have random chains and therefore only the stochastic variables are initialized. Any variable that has a value set by the assignment operator, “<-”, has a deterministic, not stochastic value. But sometimes, as in the present application, a vector has some stochastic components and some nonstochastic components. The way to initialize such vectors is to provide a specific value for the stochastic components, and an NA value to the nonstochastic components..."

Hence, in the initialization section, include the following code:

sqzData = .01+.98*dataList$nCorrOfSubj/dataList$nTrlOfSubj
mu = aggregate( sqzData , list(dataList$CondOfSubj) , "mean" )[,"x"]
sd = aggregate( sqzData , list(dataList$CondOfSubj) , "sd" )[,"x"]
kappa = mu*(1-mu)/sd^2 - 1
initsList = list( theta = sqzData ,

 # For four mu version: 
 # mu = mu ,
 # For two mu version:
 mu = c( mu[1] , mean( mu[2:4] ) , NA , NA ) , 
 mu0 = mean(mu) ,
 kappa = kappa ,
 mdlIdx = 1 )

And, in run-the-chains section, uncomment the explicit initialization:

jagsModel = jags.model( "model.txt" , data=dataList , inits=initsList ,

 n.chains=nChains , n.adapt=adaptSteps )

If you copy and paste from this page, it would be best not to directly paste into your R script because extraneous font formatting may bother the R interpreter. To avoid this problem, copy from here but then paste into a plain-text editor such as Notepad in Windows, which will strip away the font info. Then copy and paste from the plain-text editor into R.

No comments:

Post a Comment