Linear regression in OpenBUGS

By | August 19, 2012

I always wondered why is it so difficult to find an OpenBUGS example of simple linear regression on the Web. Curiously, such example is even missing in the OpenBUGS help. The only nice example so far is in the book by Marc Kéry. I have simplified the code. You need to have OpenBUGS (or WinBUGS) installed, together with the R2OpenBUGS package. Then just copy-paste the code into R.
Good luck!

```# generation of artificial data
x <- sort(rnorm(100))
y <- 5 + 3*x + rnorm(100)
plot(x,y)

# putting the data in the OpenBugs-friendly format
my.data <- list(x=x, y=y, N=100)

library(R2OpenBUGS)

# sinking the model into a file
sink("regression.txt")
cat("
model
{
# priors
beta0 ~ dnorm(0,0.001)
beta1 ~ dnorm(0,0.001)
sigma ~ dunif(0,100)

tau <- 1/(sigma*sigma)

# likelihood
for(i in 1:N)
{
y[i] ~ dnorm(mu[i], tau)
mu[i] <- beta0 + beta1*x[i]
}
}
")
sink()

# setting the parameters to be monitored
params <- c("beta0", "beta1")

# specifying the initial MCMC values
inits <- function()
{
list (beta0=rnorm(1,0,10000),
beta1=rnorm(1,0,10000))
}

# calling OpenBugs to sample from the posterior distributions
# (you may need to change the OpenBUGS.pgm directory)
res <- bugs(data=my.data, inits=inits, parameters.to.save=params,
n.iter=2000, n.chains=3, n.burnin=1000,
model.file="regression.txt", debug=FALSE, codaPkg=TRUE,
OpenBUGS.pgm="C:/Program Files (x86)/OpenBUGS/OpenBUGS321/OpenBUGS.exe",
working.directory=getwd())

# and, finally, getting the results
1. Abdul Awal Opu