Solved – Undefined real result error at WinBUGS

distributionserrormultilevel-analysispriorwinbugs

I am currently working on my thesis and interested in estimating a multilevel differential item functioning model and I using at WinBUGS. Until I had done model check-up, there are no errors. However, when I tried to "update" the sample, suddenly the trap message said "undefined real result" popped up. What am I doing wrong? Also, I have tried different prior and initial values but I could not solve the problem!

   # Model
   Model
   {
   for (l in 1:10){
   y[l] ~ dbern(p[l])
   logit(p[l])<- u2[stu[l]] - beta[x[l]] + gamma[tea[l], x[l]]*grp[l] + alpha1[x[l]] *geo[l] +
   alpha2[x[l]]*conf[l] + alpha3[x[l]]*ses[l]
  }
  for (t in 1:5){
  for (i in 1:5){
  gamma[t,i] ~ dnorm(gamma.hat[t,i], tau.gamma[i])
  gamma.hat[t,i]<-pi1[i] + pi2[i]*inq[t]
  }
  }
  # fixed effect prior
   for (i in 1:5){
   beta[i] ~ dnorm(0, .0001)
   alpha1[i] ~ dnorm(0, .0001)
   alpha2[i] ~ dnorm(0, .0001)
   alpha3[i] ~ dnorm(0, .0001)
   pi1[i] ~ dnorm(0, .0001)
   pi2[i] ~ dnorm(0, .0001)
    }
   # Random effect prior
   for (s in 1:10){
    u2[s] ~ dnorm(0,tau.u2)
    }
   tau.u2 <- pow(sigma.u2, -2)
   sigma.u2 ~ dunif (0, 100)
   for (i in 1:5){
   tau.gamma[i] <- pow(sigma.gamma[i],-2)
   sigma.gamma[i] ~ dunif(0, 100)
   }
   }

   # Data
   list(y=c(0, 1, 1, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 1, 1, 0, 0, 0, 0, 1, 1, 0,  1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 0, 1, 1, 0), ses=c(2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1), conf=c(1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3), geo=c(3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 3, 3, 3, 3, 3), grp=c(1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0), inq=c(1, 2, 2, 1, 2), stu=c(1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 6, 6, 6, 6, 6, 7, 7, 7, 7, 7, 8, 8, 8, 8, 8, 9, 9, 9, 9, 9, 10, 10, 10, 10, 10), x=c(1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5), tea=c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5))

   #Initital values
   list(beta=c(0, 0, 0, 0, 0), alpha1=c(0, 0, 0, 0, 0),  alpha2=c(0, 0, 0, 0, 0),  alpha3=c(0, 0, 0, 0, 0), sigma.gamma=c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1), u2=c(0, 0, 0, 0, 0), pi1=c(0, 0, 0, 0, 0), pi2=c(0, 0, 0, 0, 0), sigma.u2=1, gamma=structure(
.Data=c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), .Dim=c(5, 5)))

Best Answer

There's a lot of unknown parameters with very weak priors, and only 10 observations. Did you mean to only use the first 10 out of the 50 y's? Even with 50 observations it still looks like a complicated model to try to identify. It's probably failed numerically because the initial values are implausible given the data, or the priors are too diffuse for the data to give any information. I'd suggest simplifying it until you have some idea what the parameter values should be, and then building it up step by step.

Related Question