model{ theta[1] ~ dunif(theta[2], theta.max) theta[2] ~ dunif(theta[3], theta[1]) theta[3] ~ dunif(theta[9], theta[2]) theta[4] ~ dunif(theta[9], theta.max) theta[5] ~ dunif(theta[7], theta.max) theta[6] ~ dunif(theta[7], theta.max) theta[7] ~ dunif(theta[9], theta7max) theta7max <- min(theta[5], theta[6]) theta[8] ~ dunif(theta[9], theta.max) theta[9] ~ dunif(theta[10], theta9max) theta9max <-min(min(theta[3], theta[4]), min(theta[7], theta[8])) theta[10] ~ dunif(theta[11], theta[9]) theta[11] ~ dunif(0 ,theta[10]) bound[1] <- ranked(theta[1:8], 8) bound[2] <- ranked(theta[1:8], 1) bound[3] <- ranked(theta[9:11], 3) bound[4] <- ranked(theta[9:11], 1) for (j in 1 : 5){ theta[j + 11] ~ dunif(0, theta.max) within[j, 1] <- 1 - step(bound[1] - theta[j + 11]) for (k in 2 : 4){ within[j, k] <- step(bound[k - 1] - theta[j + 11]) - step(bound[k] - theta[j + 11]) } within[j, 5] <- step(bound[4] - theta[j + 11]) } for (i in 1:nDate){ X[i] ~ dnorm(mu[i], tau[i]) tau[i] <- 1/pow(sigma[i],2) mu[i] <- interp.lin(theta[i], calBP[], C14BP[]) # monitor the following variable to smooth density of theta theta.smooth[i] <- 10 * round(theta[i] / 10) } }