[R-sig-ME] mixture DIC model

[Mohammadreza Kordbagheri]

hello 
my name is alireza kordbagheri of shahid beheshti university tehran.

i have two question of you:

can you help me for my question:
one question: i read your article tittle:

i dont know, how draw heatmap ANOVA GAMMA.

two question:i writing code model with JAGS softwar.
i used of function GB2.

I have problem, what DIC estimation?
This article is open access at

https://projecteuclid.org/euclid.ba/1340370933.

library(rjags)

cat(file = "model.txt",
"
data {
  C <- 100
  for (i in 1:N) {
    zeros[i] <- 0
  }
}
model {
  for (i in 1:N) {
    zeros[i] ~ dpois(loglik[i] + C)
    loglik[i] <- -log(pi[1] * l1[i] + pi[2] * l2[i])
    l1[i] <- abs(a[1]) * pow(y[i], a[1] * p[1] - 1) / (pow(b1[i], a[1] * p[1]) * exp(loggam(p[1]) + loggam(q[1]) - loggam(p[1] + q[1])) * pow(1 + pow(y[i] / b1[i], a[1]), p[1] + q[1]))
    l2[i] <- abs(a[2]) * pow(y[i], a[2] * p[2] - 1) / (pow(b2[i], a[2] * p[2]) * exp(loggam(p[2]) + loggam(q[2]) - loggam(p[2] + q[2])) * pow(1 + pow(y[i] / b2[i], a[2]), p[2] + q[2]))
    b1[i]<- exp(mu[1] + alpha[1] * f[i] + beta1[1] * col1[i] + beta2[1] * col2[i] + beta3[1] * col3[i] + loggam(p[1]) + loggam(q[1]) - loggam(p[1] + 1 / a[1]) - loggam(q[1] - 1 / a[1]))
    b2[i]<- exp(mu[2] + alpha[2] * f[i] + beta1[2] * col1[i] + beta2[2] * col2[i] + beta3[2] * col3[i] + loggam(p[2]) + loggam(q[2]) - loggam(p[2] + 1 / a[2]) - loggam(q[2] - 1 / a[2]))
  }
  pi[1] ~ dunif(0, 1)
  pi[2] <- 1 - pi[1]
  for(k in 1:2) {
    mu[k] ~ dnorm(0, .01)
    alpha[k] ~ dnorm(0, .01)
    beta1[k] ~ dnorm(0, .01)
    beta2[k] ~ dnorm(0, .01)
    beta3[k] ~ dnorm(0, .01)
    a[k] ~ dnorm(0, .01)
    p[k] ~ dgamma(.01, .01)I(0, 90)
    q[k] ~ dgamma(.01, .001)I(1, 3.5)
  }
}
")

data <- matrix(c(
3323,  5, 1, 0, 0,
8332,  8, 0, 1, 1,
9572,  5, 1, 0, 0,
10172, 10, 1, 0, 1,
7631,  2, 1, 0, 0,
3855,  9, 0, 1, 1,
3252,  6, 1, 0, 1,
4433,  8, 0, 1, 1,
2188,  7, 1, 0, 0,
  333,  4, 1, 0, 0,
  199,  3, 0, 1, 1,
  692,  9, 0, 1, 0,
  311, 12, 1, 0, 0,
0.01,  5, 1, 0, 1,
  405,  6, 0, 1, 0,
  293,  6, 0, 1, 0,
  76,  7, 1, 0, 1,
  14,  9, 1, 0, 0,
3785, 21, 0, 1, 1,
10342, 11, 0, 1, 1), ncol = 5, byrow = TRUE)

bugs.data <- list(N = nrow(data),
                  y = data[, 1],
                  f = data[, 2],
                  col1 = data[, 3],
                  col2 = data[, 4],
                  col3 = data[, 5])
set.seed(123)
inits <- lapply(1:3, function(i) {
    list(mu = rnorm(2, 0, 1),
        alpha = rnorm(2, 0, 1),
        beta1 = rnorm(2, 0, 1),
        beta2 = rnorm(2, 0, 1),
        beta3 = rnorm(2, 0, 1),
        pi = c(runif(1, 0, 1), NA),
        a = runif(2, -4, -2),
        p = runif(2, 1, 2),
        q = runif(2, 1, 2),
        .RNG.name = "base::Mersenne-Twister",
        .RNG.seed = i)
    })
pars <- c("mu", "alpha", "beta1", "beta2", "beta3", "pi",
          "a", "p", "q")
model <- jags.model("model.txt",
                    data = bugs.data,
                    inits = inits,
                    n.chains = 3,
                    n.adapt = 1000)
update(model, 3000)
samp <- coda.samples(model,
                    variable.names = pars,
                    n.iter = 4000, thin = 4)
gelman.diag(samp, multivariate = FALSE)
traceplot(samp[, 15]) # pi[1]
traceplot(samp[, 16]) # pi[2]

can you help me.

Regards,
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: data.txt
URL: <https://stat.ethz.ch/pipermail/r-sig-mixed-models/attachments/20170414/b20eceb4/attachment-0001.txt>