/* This is the original form of the inhaler example in which the ordered
   logistic distribution is constructed in bugs code. This is superseded
   by the dordered distribution in the glm module from JAGS 3.4.0.
*/
var 
  pattern[Npattern,T], # response pattern e.g. (1,1) or (2,4) etc. 
  Ncum[Npattern,T],    # cumulative total
  response[N,T],       # response for patient i in period t
  p[N,T,(Ncut+1)],     # prob of response j for patient i in period t 
  Q[N,T,Ncut],         # cumulative prob of response worse than j
                       # for patient i in period t
  group[N],            # treatment group (1=AB; 2=BA)
  mu[G,T],             # logistic mean for group g & period t
  treat[2,T], beta,    # treatment effect
  period[2,T], pi,     # period effect
  carry[2,T], kappa,   # carryover effect
  a[Ncut],             # cut points for latent response variable
  b[N],                # subject random effect
  tau,                 # precision of subject effects
  sigma, log.sigma;
data {
   # Construct individual response data from contingency table
   for (i in 1:Ncum[1,1]) { 
      group[i] <- 1; for (t in 1:T) { response[i,t] <- pattern[1,t] }
   }
   for (i in (Ncum[1,1]+1):Ncum[1,2]) { 
      group[i] <- 2; for (t in 1:T) { response[i,t] <- pattern[1,t] }
   }

   for (k in 2:Npattern) {
      for(i in (Ncum[k-1,2]+1):Ncum[k,1]) {
         group[i] <- 1; for (t in 1:T) { response[i,t] <- pattern[k,t] }
      }
      for(i in (Ncum[k,1]+1):Ncum[k,2]) {
         group[i] <- 2; for (t in 1:T) { response[i,t] <- pattern[k,t] }
      }
   }
}
model {
  for (i in 1:N) {
     for (t in 1:T) {
        for (j in 1:Ncut) {
#  
# Cumulative probability of worse response than j
#
           logit(Q[i,t,j]) <- -(a[j] + mu[group[i],t] + b[i]); 
        }
#
# Probability of response = j
#
        p[i,t,1] <- 1 - Q[i,t,1];
        for (j in 2:Ncut) { p[i,t,j] <- Q[i,t,j-1] - Q[i,t,j] }
        p[i,t,(Ncut+1)] <- Q[i,t,Ncut];
      
        response[i,t] ~ dcat(p[i,t,]);
     }
#
# Subject (random) effects
#
     b[i] ~ dnorm(0.0, tau);
  }

#
# Fixed effects
#
  for (g in 1:G) {
     for(t in 1:T) { 
        # logistic mean for group i in period t
        mu[g,t] <- beta*treat[g,t]/2 + pi*period[g,t]/2 + kappa*carry[g,t]; 
     }
  }                                                             
  beta ~ dnorm(0, 1.0E-06);
  pi ~ dnorm(0, 1.0E-06);
  kappa ~ dnorm(0, 1.0E-06);

# ordered cut points for underlying continuous latent variable  
  a <- sort(a0)
  for(i in 1:3) {
     a0[i] ~ dnorm(0, 1.0E-6);
  }
 
  tau ~ dgamma(0.001, 0.001);
  sigma <- sqrt(1/tau);
  log.sigma <- log(sigma);

}