#From Thompson SG, Nixon RM and Grieve R, Addressing the issues that arise in analysing multicentre cost
 data, with application to a multinational study, Health Economics, 2005

model{
for(i in 1:N){
#The Generalized linear mixed model regression
    totcos[i]~dgamma(shape[centre[i]],rate[i]) 
    rate[i]<-shape[centre[i]]/ztheta[i]
    log(ztheta[i])<-lmu[centre[i]]+lIncont[centre[i]]*(incont[i]-incont.bar)
				  +lParaly*(paraly[i]-paraly.bar)
                       +lStype2*(stype2[i]-stype2.bar)+lStype3*(stype3[i]-stype3.bar)
				   +lLiving2*(living2[i]-living2.bar)+lLiving3*(living3[i]-living3.bar)

}

for(j in 1:13){
    Incont[j]<-exp(lIncont[j])
}

Paraly<-exp(lParaly)
Stype2<-exp(lStype2)
Living2<-exp(lLiving2)
Stype3<-exp(lStype3)
Living3<-exp(lLiving3)
Gdp<-exp(lGdp)

#Variable means
incont.bar<-mean(incont[])
paraly.bar<-mean(paraly[])
stype2.bar<-mean(stype2[])
stype3.bar<-mean(stype3[])
living2.bar<-mean(living2[])
living3.bar<-mean(living3[])
gdp.bar<-mean(gdp[]) 

#Random effect for the mean - the model mixes better if the Gdp effect is put at this level
for(j in 1:13){
    lmean.mu[j]<-lmu.mu+lGdp*(gdp[j]-gdp.bar)
    lmu[j]~dnorm(lmean.mu[j], tau.mu)
}
tau.mu<-1/ss.mu
ss.mu<-s.mu*s.mu

#Random effect for incontinence
for(j in 1:13){
    lIncont[j]~dnorm(mu.lincont, tau.lincont)
}
tau.lincont<-1/ss.lincont
ss.lincont<-s.lincont*s.lincont

#Likelihood
for(i in 1:N){
zloglik[i]<- -loggam(shape[centre[i]]) + shape[centre[i]]*log(rate[i]) + (shape[centre[i]]-1)*log(totcos[i]) - rate[i]*totcos[i]
}
loglik<--2*sum(zloglik[])

#Priors
for(j in 1:13){
    shape[j]~dunif(0,10)
}

lParaly~dunif(-0.5,0.5)
lStype2~dunif(-0.5,0.5)
lLiving2~dunif(-0.5,0.5)
lStype3~dunif(-0.5,0.5)
lLiving3~dunif(-0.5,0.5)
mu.lincont~dunif(0,1)
s.lincont~dunif(0.01,1)
lGdp~dnorm(0,1.0E-6)
lmu.mu~dnorm(0,1.0E-6)
s.mu~dunif(0.01,10)


}