[hepatitis0]    Hepatitis: a normal hierarchical
         model with measurement error

This example is taken from Spiegelhalter et al (1996) (chapter in Markov Chain Monte Carlo in Practice ) and concerns 106 children whose post-vaccination anti Hb titre was measured 2 or 3 times. Both measurements and times have been transformed to a log scale. One covariate y0 = log titre at baseline, is available.

The model is essentially a random effects linear growth curve

   Y ij ~ Normal( a i + b i ( t ij -t bar ), t )

    a i ~ Normal( a c , t a )

    b i ~ Normal( b c , t b )

where t represents the precision (1/variance) of a normal distribution. We note the absence of a parameter representing correlation between a i and b i unlike in Gelfand et al 1990. However, see the Birats example in Volume 2 which does explicitly model the covariance between a i and b i .

a c , t a , b c , t b , t are given independent ``noninformative'' priors.

Graphical model for hep example:
[hepatitis1]


BUGS language for hep example:
       model
       {
          for( i in 1 : N ) {
             for( j in 1 : T ) {
                Y[i , j] ~ dnorm(mu[i , j],tau)
                mu[i , j] <- alpha[i] + beta[i] * (t[i,j] - 6.5) +
                           gamma * (y0[i] - mean(y0[]))
             }
             alpha[i] ~ dnorm(alpha0,tau.alpha)
             beta[i] ~ dnorm(beta0,tau.beta)
          }
          tau ~ dgamma(0.001,0.001)
          sigma <- 1 / sqrt(tau)
          alpha0 ~ dnorm(0.0,1.0E-6)   
          tau.alpha ~ dgamma(0.001,0.001)
          beta0 ~ dnorm(0.0,1.0E-6)
          tau.beta ~ dgamma(0.001,0.001)
          gamma ~ dnorm(0.0,1.0E-6)
       }


Note the use of a very flat but conjugate prior for the population effects: a locally uniform prior could also have been used.

 Data ( click to open )

Inits for chain 1        Inits for chain 2    ( click to open )



Results

[hepatitis2]

With measurement error

[hepatitis3]

       model
      {
         tau.alpha ~ dgamma(0.001,0.001)
         alpha0 ~ dnorm( 0.0,1.0E-6)
         beta0 ~ dnorm( 0.0,1.0E-6)
         tau.beta ~ dgamma(0.001,0.001)
         for( i in 1 : N ) {
            alpha[i] ~ dnorm(alpha0,tau.alpha)
            beta[i] ~ dnorm(beta0,tau.beta)
            y0[i] ~ dnorm(mu0[i],tau)
            mu0[i] ~ dnorm(theta,psi)
         }
         for( j in 1 : T ) {
            for( i in 1 : N ) {
               Y[i , j] ~ dnorm(mu[i , j],tau)
               mu[i , j] <- alpha[i] + beta[i] * (t[i , j] - 6.5) +
                  gamma * (mu0[i] - mean(y0[]))
            }
         }
         tau ~ dgamma(0.001,0.001)
         sigma <- 1 / sqrt(tau)
         gamma ~ dnorm( 0.0,1.0E-6)
         theta ~ dnorm( 0.0,1.0E-6)
         psi ~ dgamma(0.001,0.001)
      }



Data ( click to open )

Inits for chain 1     Inits for chain 2    ( click to open )


Results

[hepatitis4]