model { for (i in 1:8) { y[i] ~ dbin(p[i], n[i]) logit(p[i]) <- alpha + beta*(x[i] - mean(x[])) phat[i] <- y[i]/n[i] yhat[i] <- n[i]*p[i] } alpha ~ dnorm(0, 0.0001) beta ~ dnorm(0, 0.0001) ED95 <- (logit(0.95) - alpha)/beta + mean(x[]) } list(x = c(1.6907, 1.7242, 1.7552, 1.7842, 1.8113, 1.8369, 1.8610, 1.8839), n = c(59, 60, 62, 56, 63, 59, 62, 60), y = c(6, 13, 18, 28, 52, 53, 61, 60)) list(alpha = 50, beta = 0)