7
$\begingroup$

I'm trying to specify a model in JAGS/rjags with one between subjects factor (a, with two levels - Y, N) interacting with one repeated measures continuous variable x plus subject varying slopes and intercepts that correlate. I can specify this model simply enough with the lmer function:

lmer(y ~ a + x + a:x + (1 + a | id))

My JAGS/rjags is very rusty (or very fresh). The below seems to me to be fitting a model with subject varying intercepts and subject varying slopes while estimating the slope for both levels of a, but I'm not sure I'm doing what I think I'm doing. There's also no correlation specified between the two.

modelstring = "
model {
  for ( i in 1:Ntotal ) {
    y[i] ~ dnorm( mu[i] , tau )
    mu[i] <- a1[aLvl[i]] + s1[sLvl[i]] + a2[aLvl[i]] * x[i] + s2[sLvl[i]] * x[i]
  }
  # Prior:
  tau <- pow( sigma , -2 )
  sigma ~ dunif(0,1000)
  for ( j in 1:2 ) { 
      a1[j] ~ dnorm( 0.0 , aTau ) 
      a2[j] ~ dnorm( 0.0 , aTau ) 
  }
  aTau <- 1 / pow( aSD , 2 )
  aSD <- abs( aSDunabs ) + .1
  aSDunabs ~ dt( 0 , 1.0E-7 , 2 )
  #
  for ( j in 1:NsLvl ) { 
    s1[j] ~ dnorm( 0.0 , sTau ) 
    s2[j] ~ dnorm( 0.0 , sTau )
  }
  sTau <- 1 / pow( sSD , 2 )
  sSD <- abs( sSDunabs ) + .1
  sSDunabs ~ dt( 0 , 1.0E-7 , 2 )
}
"

The framework for this comes from Kruschke and this has been of some help too. I would appreciate some pointers or links to examples of similar analyses.

$\endgroup$

1 Answer 1

7
$\begingroup$

I eventually figured this one out with much help from Doing Bayesian Data Analysis (Kruschke) and Data Analysis Using Regression and Multilevel/Hierarchical Models (Gelman). This model gives varying intercepts and slopes and the correlation between them.

y = dependent variable
sLvl = participant id at each data point
aLvlx = between subjects factor for each id
NaLvl = Number of levels for the between subject's factor
Ntotal = total length of data in long form

modelstring = "
model {
for( r in 1 : Ntotal ) {
    y[r] ~ dnorm( mu[r] , tau )
    mu[r] <- b0[ sLvl[r] ] + b1[ sLvl[r] ] * x[r]
}
#General priors
tau ~ dgamma( sG , rG )
sG <- pow(m,2)/pow(d,2)
rG <- m/pow(d,2)
m ~ dgamma(1, 0.001)
d ~ dgamma(1, 0.001)
#Subject level priors
for ( s in 1 : NsLvl ) {
    b0[s] <- B[s,1]
    b1[s] <- B[s,2]
    B[s, 1:2] ~ dmnorm( B.hat[s, ], Tau.B[ , ] )
    B.hat[s,1] <- hix1[aLvlx[s]]
    B.hat[s,2] <- hix2[aLvlx[s]]        
}
Tau.B[1:2 , 1:2] <- inverse(Sigma.B[,])
Sigma.B[1,1] <- pow(tau0G, 2)
Sigma.B[2,2] <- pow(tau1G, 2)
Sigma.B[1,2] <- rho * tau0G * tau1G
Sigma.B[2,1] <- Sigma.B[1,2]

tau0G ~ dunif(0.001,100)
tau1G ~ dunif(0.001,100)
rho ~ dunif(-1,1)
#Between subjects level priors
for ( k in 1:NaLvl ) { 
    hix1[k] ~ dnorm( 0 , 0.0000001 )
    hix2[k] ~ dnorm( 0 , 0.0000001 ) 
}
}"
writeLines(modelstring,con="model.txt")
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.