welcome:
location: RatsExample

# 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 Yij is the weight of the ith rat measured at age xj .

 Weight Yij of rat i on day xj xj = 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

Yij ~ Normal( αi + βi (xj - xbar ), τc )

αi ~ Normal( αc , τa )

βi ~ Normal( βc, τb )

where xbar = 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 xj's around their mean to reduce dependence between αi and &betai 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 xbar.

## 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.

## 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.

RatsExample (last edited 2009-10-26 15:59:15 by BobOH)