var
 Y[N,M], age[M],               # jaw bone length in mm and age
 mu[M],Omega[M,M],Sigma2[M,M], # mean, precision & covariance matrix for Y
 beta0,                        # regression coefficients for location models
 R[M,M],                       # prior guess at magnitude of Sigma2 
 resid[N,M],resid2[N,M],RSS,   # residuals and residual sum of squares
 M;    # log-likelihood terms and deviance
model {
  for (i in 1:N) {
     Y[i,] ~ dmnorm(mu[], Omega[,]);  # The 4 measurements for each  
  }                                   # boy are multivariate normal

  for(j in 1:M) {     # location model for mean bone length at each age
     mu[j] <- beta0;                                         # constant
  }

  beta0 ~ dnorm(0.0, 0.001); 
  Omega[,] ~ dwish(R[,], 4);	# between-child variance in length at each age	
  Sigma2[,] <- inverse(Omega[,]);

  for (i in 1:N) {
     for  (j in 1:M) {
        resid[i,j] <- Y[i,j] - mu[j];         # residuals
        resid2[i,j] <- resid[i,j]^2;     # squared residuals
     }
  }
  RSS <- sum(resid2[,]);                    # Residual Sum of Squares
}