example-MixSIR.R

library(runjags)
library(R2WinBUGS)
library(R2jags)
library(BRugs)
library(gtools)
library(gdata)

num.iso = 2   # number of isotopes
num.prey = 3  # number of prey species
num.pop = 3   # number of populations

# prey * isotopes
u = array(0,dim=c(num.prey,num.iso,num.pop))
sigma2 = array(0,dim=c(num.prey,num.iso,num.pop))


   X = as.matrix(read.table("test_fake_data.txt",header=F))
   Region = sort(rep(seq(1,3),30))
   N = dim(X)[1]
   MU = read.table("test_fake_means.txt",header=F)
   SIG = read.table("test_fake_sds.txt",header=F)

   u = array(0,dim=c(num.prey,num.iso,num.pop))
   sigma2 = array(0,dim=c(num.prey,num.iso,num.pop))
   for(i in 1:3) {
      for(j in 1:2) {
         u[,j,i] = MU[,j]
         sigma2[,j,i] = SIG[,j]^2
      }
   }
   
# These are the parameters that need to be set by the user
  mcmc.chainLength <- as.integer(10000)  # post-burn
  mcmc.burn <- as.integer(5000)
  mcmc.thin = 1;
  mcmc.chains = 3       # needs to be at least 2 for DIC
  dic.samples = 5000
  alpha=rep(1,num.prey)
  alphaMix = c(1,1)
  R = diag(num.prey)
  dic = 0
  calc.DIC = TRUE   
        
  jags.data = list("u", "sigma2", "N", "num.prey", "num.iso", "num.pop","X","alpha","Region")  
  jags.params = c("p")
  jags.inits <- function(){
     list("p"=c(0.33,0.33,0.34))	
  }
  model.loc = "MixSIR.txt"
  jags.1 = jags(jags.data, inits = jags.inits, parameters.to.save= jags.params, model.file=model.loc, n.chains = mcmc.chains, n.burnin = mcmc.burn, n.thin = mcmc.thin, n.iter = mcmc.chainLength, DIC = TRUE)