Dear all,
I have problems understanding how to do indexing of random effects in
the multivariate normal distribution. I have used the Jaws data as an
example, and added a variable that is assumed to mimic the levels of a
random effect. I cannot see why the code below is wrong, any help is
appreciated:
model
{
for (i in 1:N) {
Y[i, 1:T] ~ dmnorm(mu[i,1:T], Omega.err[ , ])
mu[i,1:T] <- my[1:T]+popeff[popid[i],1:T]
}
for(j in 1:T) {
my[j] ~ dnorm(0.0, 0.001)
my2[j] ~ dnorm(0.0, 0.001)
}
for(k in 1:5) {
popeff[k,1:T] ~ dmnorm(mu.pop[k,1:T], Omega.pop[ , ])
mu.pop[k,1:T] <- my2[1:T]
}
Omega.err[1 : T , 1 : T] ~ dwish(R.err[ , ], 4)
Cov.err[1 : T , 1 : T] <- inverse(Omega.err[ , ])
Omega.pop[1 : T , 1 : T] ~ dwish(R.pop[ , ], 4)
Cov.pop[1 : T , 1 : T] <- inverse(Omega.pop[ , ])
}
# Data
list(T = 4, N = 20, Y = structure(
.Data = c(47.8, 48.8, 49.0, 49.7,
46.4, 47.3, 47.7, 48.4,
46.3, 46.8, 47.8, 48.5,
45.1, 45.3, 46.1, 47.2,
47.6, 48.5, 48.9, 49.3,
52.5, 53.2, 53.3, 53.7,
51.2, 53.0, 54.3, 54.5,
49.8, 50.0, 50.3, 52.7,
48.1, 50.8, 52.3, 54.4,
45.0, 47.0, 47.3, 48.3,
51.2, 51.4, 51.6, 51.9,
48.5, 49.2, 53.0, 55.5,
52.1, 52.8, 53.7, 55.0,
48.2, 48.9, 49.3, 49.8,
49.6, 50.4, 51.2, 51.8,
50.7, 51.7, 52.7, 53.3,
47.2, 47.7, 48.4, 49.5,
53.3, 54.6, 55.1, 55.3,
46.2, 47.5, 48.1, 48.4,
46.3, 47.6, 51.3, 51.8),
.Dim = c(20, 4)),
R.err = structure(
.Data = c(1, 0, 0, 0,
0, 1, 0, 0,
0, 0, 1, 0,
0, 0, 0, 1), .
Dim = c(4, 4)),
R.pop = structure(
.Data = c(1, 0, 0, 0,
0, 1, 0, 0,
0, 0, 1, 0,
0, 0, 0, 1), .
Dim = c(4, 4)),
popid = c(1,
1,
1,
1,
2,
2,
2,
2,
3,
3,
3,
3,
4,
4,
4,
4,
5,
5,
5,
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
|