# set the necessary libraries
library(lattice)
library(Matrix)
library(lme4)
library(coda)

library(R2WinBUGS)
library(BRugs)
library(MASS)
library(foreign)
library(car)
library(arm)

modelSetSeed(1237) # set the seed
# Step 1
# Write a WinBUGS model to a temporary file
BRugsmodel <- function()
{
   for( i in 1 : N ) {
      mu[i] <- alpha + beta * (x[i] - xbar)
   }
   for( i in 1 : N ) {
      Y[i] ~ dnorm(mu[i],tau)
   }
   alpha ~ dnorm( 0.0,1.0E-6)
   beta ~ dnorm( 0.0,1.0E-6)
   tau ~ dgamma(0.001,0.001)
   sigma <- 1 / tau
}

model.file <- file.path(tempdir(), "regrmodel.txt")
write.model(BRugsmodel, model.file)
file.show(model.file)
modelCheck(model.file) # check model file

# Step 2
# Set data
model.data <- list(N=5,xbar=3,
    x = c(1,2,3,4,5),
    Y = c(1,3,3,3,5))
modelData(bugsData(model.data,fileName=tempfile())) # read data file

# Step 3
# Compile data
modelCompile(numChains=1)       # compile model with 1 chain

# Step 4
# Set inits
model.inits <- function() list(alpha = 0, beta = 0, tau = 1)
modelInits(bugsInits(model.inits,fileName = tempfile()))

# Step 5
# Burnin
modelUpdate(2000) # 2000 iterations

# Step 6
# Set Parameters
model.parameters <- c("alpha", "beta", "tau")
samplesSet(model.parameters) # set parameters to be monitored
modelUpdate(2000) # 2000 more iterations ....

# Step 7
# Look at Sample Stats
samplesStats("*") # the summarized results

# Step 8
# Look at graphical info
samplesHistory("*", mfrow = c(4, 2)) # histories,
samplesDensity("*") # densities
samplesBgr("*") # bgr statistics
samplesAutoC("*", 1) # autocorrelations (1st chain)

# Step 9
# Clear monitor and free memory
samplesClear("*") # clear monitor

# do it all in one hit
fit.model = BRugsFit(data = model.data,
 inits = model.inits,
 para = model.parameters,
 nBurnin = 10000, nIter = 10000,
 modelFile = model.file, numChains = 1)