I am trying to fit a linear mixed model of the form
y_ij = X_ij \beta + delta_i + e_ij
where e_ij ~N(0,s^2_ij) with s_ij known
and delta_i~N(0,tau^2)
I have code in R (using method of moments or mle) to fit this.
I am trying to learn Bayesian statistics so I want to fit this
model in winbugs. I get horrible posterior point estmates. I could give
mle values as starting values, but that seems cheating a bit.
Am I doing anything incredibly moronic in my code below ?
Thanks,
Benn
model
{
for( i in 1 : N ) {
for( j in 1 : clustsize[i] ) {
Y[i , j] ~ dnorm(mu[i , j],varmatrix[i,j])
mu[i , j] <- delta[i]+beta0 + beta2*x[(i-1)*mx+j,2]
+beta3*x[(i-1)*mx+j,3]+beta4*x[(i-1)*mx+j,4]+beta5*x[(i-1)*mx+j,5]
}
delta[i]~dnorm(0.0,.000001)
}
beta0 ~ dnorm(0,.000001)
beta2 ~ dnorm(0,.000001)
beta3 ~ dnorm(0,.000001)
beta4~ dnorm(0,.000001)
beta5~ dnorm(0,.000001)
}
Data list(N= 11 ,
clustsize = c(1, 1, 12, 1, 2, 3, 3, 6, 5, 3, 1),
mx = 12 ,
Y = structure(.Data=c(0.08, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
0.19, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 0.8143174,
0.1528186, 0.1362733, 0.034752, 0.1494959, 0.1261007, 0.1255043,
0.0988546, 0.1320305, 0.1265011, 0.1297557, 0.117368, 0.57, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 0.15, 0.27, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, 0.0125672, 0.0159507, 0.0159507,
NA, NA, NA, NA, NA, NA, NA, NA, NA, 0.19, 0.21, 0.35, NA, NA,
NA, NA, NA, NA, NA, NA, NA, -0.01, 0.1, 0.53, -0.68, -0.16, -0.51,
NA, NA, NA, NA, NA, NA, -1.034374, -0.8834668, -1.075397, -0.9376762,
-1.038769, NA, NA, NA, NA, NA, NA, NA, 0.081, 0.138, 0.146, NA,
NA, NA, NA, NA, NA, NA, NA, NA, 0.49, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA), .Dim = c(11, 12))
,
varmatrix = structure(.Data=c(263.401498574366, NA, NA, NA, NA, NA, NA, NA,
NA,
NA, NA, NA, 76.1573363395064, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, 194.045721584835, 235.409569774326, 220.690108188553,
138.486761611400, 219.134195139914, 225.45857703789, 217.649383164342,
207.055159667532, 216.504935185474, 212.016249515029, 201.825762149932,
208.637968513130, 1111.11111111111, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, 172.660805895495, 53.2930775455279, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, 24975.7023429189, 15512.8997877784,
27245.5591609796, NA, NA, NA, NA, NA, NA, NA, NA, NA, 109.145949585032,
94.3224431426345, 68.5693865862892, NA, NA, NA, NA, NA, NA, NA,
NA, NA, 16667.5561648032, 166.675553040714, 5.93362597109009,
3.60457511722423, 65.1076381166684, 6.4081335336211, NA, NA,
NA, NA, NA, NA, 18.0124923451681, 17.7298649674369, 17.9013112392251,
20.8344366200993, 23.6343931776020, NA, NA, NA, NA, NA, NA, NA,
62.9881582262535, 44.4444444444444, 51.7571554267377, NA, NA,
NA, NA, NA, NA, NA, NA, NA, 28.1405535998313, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA), .Dim = c(11, 12))
,
x = structure(.Data=c(1, 0, 0, 1, 1, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1, 1,
1, 1, 1, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1, 1, 1, 1, 0, 1, 1,
1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1,
1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1,
0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, 0, 0, 1, 0, 1, 0, 0, 1, 0, 1, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 1, 1, 1,
0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, 0, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1,
0, 1, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 0,
0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 0,
1, 0, 1, 0, 0, 1, 0, 1, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0,
0, 0, 0, 1, 0, 0, 0, 0, 1, 0, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 0, 0, 1, 1, 1, 0,
0, 1, 1, 1, 0, 0, 1, 1, 1, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, 1, 1, 1, 1, 1, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA,
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA), .Dim = c(132,
5))
)
-------------------------------------------------------------------
This list is for discussion of modelling issues and the BUGS software.
For help with crashes and error messages, first mail [log in to unmask]
To mail the BUGS list, mail to [log in to unmask]
Before mailing, please check the archive at www.jiscmail.ac.uk/lists/bugs.html
Please do not mail attachments to the list.
To leave the BUGS list, send LEAVE BUGS to [log in to unmask]
If this fails, mail [log in to unmask], NOT the whole list
|