1
$\begingroup$

I want to analyze a multilevel multidimensional model in WinBUGS. the model is as below (N=2362 students responding to K=45 items of a test, students are nested within J=116 schools):

model{
#responses
for(i in 1:N){
    for(j in 1:K){
        logit(p[i,j])<- a1[j]*th[i,1]+a2[j]*th[i,2]-b[j]
        y[i,j]~dbern(p[i,j] )
    }
    th[i,1:2]~dmnorm(mu[sc[i],1:2],tau.p[1:2,1:2])
}
#school level
for(j in 1:J){  
    mu[j,1:2]~dmnorm(m[j,1:2],tau.s[1:2,1:2])
}    

#priors
for(j in 1:J){
    m[j,1:2]~dmnorm(m0[1:2],cov[1:2,1:2])
}

tau.p[1:2,1:2]~dwish(cov[1:2,1:2],2)
tau.s[1:2,1:2]~dwish(cov[1:2,1:2],2)
sigma.p[1:2,1:2]<-inverse(tau.p[,])
sigma.s[1:2,1:2]<-inverse(tau.s[,])
s2p<-sum(sigma.p[,])
s2s<-sum(sigma.s[,])
rho<-(s2s)/(s2s+s2p)

a1[1]~dlnorm(0,4)
a2[1]<-0
b[1]~dnorm(0,1)
for(s in 2:K) {
    a1[s]~dlnorm(0,4)
    a2[s]~dlnorm(0,4)
    b[s]~dnorm(0,1)
}    
}

I've set these functions as initial values:

ini<-function(){
list(tau.p=matrix(rgamma(4,100,100),2,2),
tau.s=matrix(rgamma(4,100,100),2,2),
th=rmvnorm(N,mean=c(0,0),sigma=diag(2)),
m=rmvnorm(J,mean=c(0,0),sigma=diag(2)),
mu=rmvnorm(J,mean=c(0,0),sigma=diag(2)),
a1=rlnorm(K,0, 0.4),
a2=c(NA,rlnorm(K-1,0, 0.4)),
b=rnorm(45,0,0.5))
}

I use rube package in R to check and run my analysis and everything looks fine. When I run the model I receive "Trap 66 (postcondition violated)" or "undefined real result". I think the problem is from the initials but I have no idea how to solve it.

Any idea?

|cite|improve this question
$\endgroup$
  • $\begingroup$ any chance you can provide a slimmed down data set, otherwise it is difficult for one to replicate your error and help you debug? $\endgroup$ – gjabel Oct 10 '13 at 16:07
  • 1
    $\begingroup$ Just as an FYI you shouldn't use "t" as a variable name since it's the name of an internal R function ?t. $\endgroup$ – emhart Oct 10 '13 at 22:58
  • $\begingroup$ @gjabel it's a huge data set.I'm not sure how I can provide a slim data.Do you the problem due to data? $\endgroup$ – Amin Oct 11 '13 at 16:25
  • 1
    $\begingroup$ I just posted a blog about this problem on aheblog.com/2014/01/08/…. Hope this helps. Pepijn $\endgroup$ – user37573 Jan 20 '14 at 20:06
  • 1
    $\begingroup$ I'd suggest trying JAGS. The error messages are far more transparent and you can use essentially the same code (maybe a small tweak here and there). Even if you are wedded to WinBUGS for some reason, you might get some clue as to what is wrong. $\endgroup$ – guy Jan 20 '14 at 20:51
2
$\begingroup$

I've run into this problem before. Often times it's due to a problem with passing negative values to a distribution that doesn't allow them. So perhaps your log normal or wishart distributions are getting negative values somehow.

Another issue is that you are giving priors on your variance parameters that is quite big. This can cause unexpected values to be passed to your log-normald and wishart variables. I'd start by severely constraining the variance priors and slowly expand the parameter space. In my experience fitting BUGS models is both and art and a science.

|cite|improve this answer
$\endgroup$
  • $\begingroup$ How can I narrow down priors of variances while I'm using Wishart Dist.? The cov[1:2,1:2] is an identity matrix and I set df of Wishart Dist. as 2 because of the dimensions of th. I completely agree with you about trickiness of BUGS models. $\endgroup$ – Amin Oct 11 '13 at 16:32
  • $\begingroup$ I'd try modifying the inits you draw from the rgamma in your init function as a place to start. $\endgroup$ – emhart Oct 11 '13 at 18:56

Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

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