= Name of Contributed Example = == Description of problem == This example is taken from section 6 of Gelfand et al (1990), and concerns 30 young rats whose weights were measured weekly for five weeks. Part of the data is shown below, where Y,,ij,, is the weight of the ith rat measured at age x,,j,, . || |||||||||| Weight Y,,ij,, of rat i on day x,,j,, || || || x,,j,, = 8 || 15 || 22 || 29 || 36 || || Rat 1 || 151 || 199 || 246 || 283 || 320 || || Rat 2 || 145 || 199 || 249 || 293 || 354 || || ... || || || || || || || Rat 1 || 153 || 200 || 244 || 286 || 324 || A plot of the 30 growth curves suggests some evidence of downward curvature. The model is essentially a random effects linear growth curve Y,,ij,, ~ Normal( α,,i,, + β,,i,, (x,,j,, - x,,bar,, ), τ,,c,, ) α,,i,, ~ Normal( α,,c,, , τ,,a,, ) β,,i,, ~ Normal( β,,c,,, τ,,b,, ) where x,,bar,, = 22, and τ represents the precision (1/variance) of a normal distribution. We note the absence of a parameter representing correlation between α,,i,, and β,,i,, unlike in Gelfand ''et al.'' 1990. However, see the Birats example in Volume 2 which does explicitly model the covariance between α,,i,, and β,,i,,. For now, we standardise the x,,j,,'s around their mean to reduce dependence between α,,i,, and &beta,,i,, in their likelihood: in fact for the full balanced data, complete independence is achieved. (Note that, in general, prior independence does not force the posterior distributions to be independent). α,,c,, τ,,a,, β,,c,,, τ,,b,, and τ,,c,, are given independent "noninformative" priors. Interest particularly focuses on the intercept at zero time (birth), denoted α,,0,, = α,,c,, - β,,c,, x,,bar,,. == BUGS Code == {{{ model { for( i in 1 : N ) { for( j in 1 : T ) { Y[i , j] ~ dnorm(mu[i , j],tau.c) mu[i , j] <- alpha[i] + beta[i] * (x[j] - xbar) culmative.Y[i , j] <- culmative(Y[i , j], Y[i , j]) post.pv.Y[i , j] <- post.p.value(Y[i , j]) prior.pv.Y[i , j] <- prior.p.value(Y[i , j]) replicate.post.Y[i , j] <- replicate.post(Y[i , j]) pv.post.Y[i , j] <- step(Y[i , j] - replicate.post.Y[i , j]) replicate.prior.Y[i , j] <- replicate.prior(Y[i , j]) pv.prior.Y[i , j] <- step(Y[i , j] - replicate.prior.Y[i , j]) } alpha[i] ~ dnorm(alpha.c,alpha.tau) beta[i] ~ dnorm(beta.c,beta.tau) } tau.c ~ dgamma(0.001,0.001) sigma <- 1 / sqrt(tau.c) alpha.c ~ dnorm(0.0,1.0E-6) alpha.tau ~ dgamma(0.001,0.001) beta.c ~ dnorm(0.0,1.0E-6) beta.tau ~ dgamma(0.001,0.001) alpha0 <- alpha.c - xbar * beta.c } }}} Note the use of a very flat but conjugate prior for the population effects: a locally uniform prior could also have been used. == Data == [[/RatsData|Data]] == Inits == [[/RatsInits|Initial values]] == Results == A 1000 update burn in followed by a further 10000 updates gave the parameter estimates: |||| mean|| sd || MC_error || val2.5pc || median || val97.5pc || start || sample || || alpha0 || 106.6 || 3.655 || 0.04079 || 99.44 || 106.5 || 113.8 || 1001 || 10000 || || beta.c || 6.185 || 0.1061 || 0.00132 || 5.975 || 6.185 || 6.394 || 1001 || 10000 || || sigma || 6.086 || 0.4606 || 0.007398 || 5.255 || 6.061 || 7.049 || 1001 || 10000 || These results may be compared with Figure 5 of Gelfand et al 1990 --- we note that the mean gradient of independent fitted straight lines is 6.19.