model kidney;

const

N = 38,  # of patients
M = 2;   # number of observations per patient

var

t[N,M],                       # failture time
t.cen[N,M],                    # censoring time
mu[N,M], r,                   # Weibull parameters
b[N],                         # random effects for patients
tau,                          # precision for random effects 
sigma,                        # 1/sqrt(tau)
age[N,M], sex[N], disease[N],  # covariates
beta.age, beta.sex,            # regression coefficients
beta.disease[4], alpha;        # regression coefficients



data t, t.cen, age, sex, disease in "kidney.dat";
inits in "kidney.in";
{
for (i in 1:N) {
 for (j in 1:M) {

     # survival times bounded below by censoring times:
     t[i,j] ~ dweib(r,mu[i,j]) I(t.cen[i,j],);

     log(mu[i,j]) <- alpha + beta.age*age[i,j]
                     + beta.sex*sex[i]
                     + beta.disease[disease[i]] + b[i];

     }
     # Random effects:
     b[i] ~ dnorm(0.0,tau);
}

# Priors:

alpha ~ dnorm(0.0, 0.0001);
beta.age ~ dnorm(0.0, 0.0001);
beta.sex ~ dnorm(0.0, 0.0001);
beta.disease[1] <- 0;  # corner-point constraint
for (k in 2:4) {
  beta.disease[k] ~ dnorm(0.0,0.0001);
}
tau ~ dgamma(1.0E-3,1.0E-3);
r ~ dgamma(1.0, 1.0E-3);
sigma <- 1/sqrt(tau);  # s.d. of random effects
}