#setwd("c:/Users/gustaf/Documents/Workspace")
library("arm")    ### includes the R -> BUGS stuff

### useful functions
logit <- function(p) { log(p)-log(1-p) }
expit <- function(x) {1/(1+exp(-x)) }
trpl <- function(x) {c(quantile(x, .025), mean(x), quantile(x, .975))}

### three pre-cursors to calling BUGS
### first, data to pass from R to BUGS
dattobugs <- list(
  wtot0=390, n0=1250, wtot1=430, n1=1250,
  a.sn=50.4, b.sn=10.2,
  a.sp=50.4, b.sp=10.2)

### second, scheme to generate initial values for MCMC
  inittobugs <- function() {
    list(SN=0.8, SP=0.8, r0=0.1, r1=0.1)
  }

### third, parameters to monitor
  paramstobugs <- c("SN", "SP", "r0", "r1", "lor")

### now call BUGS, model in file "adjust.txt"
bugsfit1 <- bugs(dattobugs, inittobugs, paramstobugs,
            "adjust.txt",
            n.chains=1, n.thin=10,
            n.iter=100000, bugs.directory="c:/WinBUGS14/",
            debug=F)

### can add debug=T above to see more output as it goes

### make output `live'
attach.bugs(bugsfit1)

### summary of inference on log-OR
print(trpl(lor))