#The model:
model{

   for(i in 1:n){
 y[i]~dnorm(mu[i],taue)
      mu[i] <- beta0 + beta1*x[i,1] + beta2*x[i,2] + beta3*x[i,3]
#      mu[i] <- beta0 + beta1*x[i,1] + beta2*x[i,2] 

   }
   
    beta0~dnorm(betamean0,tau0)
	beta1~dnorm(betamean1,tau1)
	beta2~dnorm(betamean2,tau2)
	beta3~dnorm(betamean3,tau3)
 
betamean0 ~ dnorm(0,0.0001)
betamean1 ~ dnorm(0,0.0001)
betamean2 ~ dnorm(0,0.0001)
betamean3 ~ dnorm(0,0.0001)

tau0 ~ dgamma(1,0.001)
tau1 ~ dgamma(1,0.001)
tau2 ~ dgamma(1,0.001)
tau3 ~ dgamma(1,0.001)

taue ~ dgamma(1,0.001)
 
}