model
{
for( i in 1 : N ) {
    for( j in 1 : T ) {
        Y[i , j] ~ dnorm(mu[i , j],eps.tau)
        mu[i , j] <- beta[i,1] + beta[i,2] * (x[j]-11)
    }
    beta[i,1:2] ~ dmnorm(beta.mu[1:2],iSigma[1:2,1:2])
}
beta.mu[1:2] ~ dmnorm(mean[1:2], prec[1:2, 1:2])
iSigma[1:2, 1:2] ~ dwish(R[1:2, 1:2], r)
Sigma[1:2, 1:2] <- inverse(iSigma[1:2, 1:2]) 
eps.tau <- exp(logtau)
logtau ~ dflat()
sigma <- 1 / sqrt(eps.tau)
}
list(x = c(8,10,12,14), N = 11, T = 4,  
Y = structure(
.Data =   c(21,20,21.5,23,
21,21.5,24,25.5,
20.5,24,24.5,26,
23.5,24.5,25,26.5,
21.5,23,22.5,23.5,
20,21,21,22.5,
21.5,22.5,23,25,
23,23,23.5,24,
20,21,22,21.5,
16.5,19,19,19.5,
24.5,25,28,28),
.Dim = c(11,4)),mean = c(0, 0),r=4,
R = structure(.Data = c(1, 0, 0, 0.1), 
.Dim = c(2, 2)),
prec = structure(.Data = c(1.0E-6, 0,0,1.0E-6), 
.Dim = c(2, 2))))
list(beta = structure(.Data = c(18,18,18,18,18,18,18,18,18,18,18,
     .5,.5,.5,.5,.5,.5,.5,.5,.5,.5,.5), .Dim=c(11,2)), beta.mu = c(18,.5), 
     iSigma = structure(.Data = c(1, 0, 0, 0.1), .Dim = c(2, 2)), logtau = 0)