
pausa <- function() { cat ("\n >> press Enter to continue\n"); scan() }

mu.true    = 3
sigma.true = 2
sample.n   = 10
x = rnorm(sample.n, mu.true, sigma.true) 
print(round(x,3))

data = list(x=x)
inits = list(mu=mean(x), tau=1/var(x))

model = "tmp_model.bug"
write("
model{
    for (i in 1:length(x)) {
      x[i] ~ dnorm(mu, tau);
     }
    mu ~ dnorm(0.0, 1.0E-6);
    tau ~ dgamma(1.0, 1.0E-6);
    sigma <- 1.0/sqrt(tau);
    # predictions
    xf ~ dnorm(mu, tau);
}
", model)

library(rjags)

n.iter = 100000
jm <- jags.model(model, data, inits)
update(jm, 100)
chain <- coda.samples(jm, c("mu","sigma","xf"), n.iter=n.iter)

plot(chain)
print(summary(chain))

