How can I set a random seed using the jags() function?

Each time I run my JAGS model using the jags() function, I get very different values of fitted parameters. However, I want other people to reproduce my results.

I tried to add set.seed(123), but it didn't help. This link describes how to achieve my goal using the run.jags() function. I wonder how I can do similar things using jags(). Thank you!

Below is my model in R:

##------------- read data -------------##
m <- 6 
l <- 3 
node <- read.csv("answer.csv", header = F)
n <- nrow(node)

# values of nodes
## IG
IG <- c(c(0.0, 1.0, 0.0), c(0.0, 0.0, 1.0), c(1.0, 0.0, 0.0), c(1.0, 0.0, 0.0), c(0.0, 1.0, 0.0), c(0.0, 0.0, 1.0))
IG <- matrix(IG, nrow=6, ncol=3, byrow=T)
V_IG <- array(0, dim=c(n, m, l)) 
for (i in 1:n){
  for (j in 1:m){
    for (k in 1:l)
    {
      V_IG[i,j,k] <- IG[j,k] # alternatively, V[i,j,k] <- PTS[j,k]
    }
  }
} 

## PTS
PTS <- c(c(1.0, 0.5, 0.0), c(1.0, 0.0, 0.5), c(1.0, 1.0, 0.0), c(1.0, 0.0, 1.0), c(0.0, 0.5, 1.0), c(0.0, 1.0, 0.5))
PTS <- matrix(PTS, nrow=m, ncol=3, byrow=T)
V_PTS <- array(0, dim=c(n, m, l))
for (i in 1:n){
  for (j in 1:m){
    for (k in 1:l)
    {
      V_PTS[i,j,k] <- PTS[j,k] 
    }
  }
} 

##------------- fit model -------------##
set.seed(123)
data <- list("n", "m", "V_IG", "V_PTS", "node") 
myinits <- list(list(tau = rep(1,n), theta = rep(0.5,n))) 
parameters <- c("tau", "theta")

samples <- jags(data, inits=myinits, parameters,
                model.file ="model.txt", n.chains=1, n.iter=10000, 
                n.burnin=1, n.thin=1, DIC=T)

And my model file model.txt:

model{
    # data: which node (1, 2, 3) was chosen by each child in each puzzle
    for(i in 1:n) # for each child
    {
        for (j in 1:m) # for each problem
        {
            # node chosen
            node[i,j] ~ dcat(mu[i,j,1:3])
            mu[i,j,1:3] <- exp_v[i,j,1:3] / sum(exp_v[i,j,1:3])
              for (k in 1:3) {
                exp_v[i,j,k] <- exp((V_IG[i,j,k]*theta[i] + V_PTS[i,j,k]*(1-theta[i]))/tau[i])
             }
    }
}
    # priors on tau and theta
    for (i in 1:n)
    {
        tau[i] ~ dgamma(0.001,0.001)
        theta[i] ~ dbeta(1,1)
    }
}

Here is a toy example for linear regression. First the model:

model{

  a0 ~ dnorm(0, 0.0001)
  a1 ~ dnorm(0, 0.0001)
  tau ~ dgamma(0.001,0.001)

  for (i in 1:100) {

    y[i] ~ dnorm(mu[i], tau)
    mu[i] <- a0 + a1 * x[i]
  }
}

Now we generate some data and you the set.seed function to generate identical results from multiple calls to the jags function.

# make the data and prepare what we need to fit the model
x <- rnorm(100)
y <- 1 + 1.2 * x + rnorm(100)

data <- list("x", "y")
parameters <- c("a0", "a1", "tau")
inits = list(list(a0 = 1, a1=0.5, tau = 1))

# First fit
set.seed(121)
samples <- jags(data, inits, 
  parameters,model.file = "./sov/lin_reg.R",
  n.chains = 1, n.iter = 5000, n.burnin = 1, n.thin = 1)

# second fit
set.seed(121) # with set.seed at same value
samples2 <- jags(data, inits, 
  parameters,model.file = "./sov/lin_reg.R",
  n.chains = 1, n.iter = 5000, n.burnin = 1, n.thin = 1)

If we pull out the draws for one of the parameters from samples and samples2 we can see that they have generated the same values.

a0_1 <- samples$BUGSoutput$sims.list$a0

a0_2 <- samples2$BUGSoutput$sims.list$a0

head(cbind(a0_1, a0_2))
          [,1]      [,2]
[1,] 1.0392019 1.0392019
[2,] 0.9155636 0.9155636
[3,] 0.9497509 0.9497509
[4,] 1.0706620 1.0706620
[5,] 0.9901852 0.9901852
[6,] 0.9307072 0.9307072