model {
for (i in 1:2){
count[i,1:5] ~ dmulti(q[i,1:5], M[i])
for (r in 1:5) {
q[i,r] <- phi[i, r] / sum(phi[i, ])
log(phi[i, r]) <- a[r] + b.treat[r] * treat[i]
}
}
for (r in 2:5){
a[r] ~ dnorm(0, 0.00001)
}
a[1] <- 0
b.treat[1] <- 0
b.treat[2] ~ dnorm(0, 0.00001)
or.treat <- exp(b.treat[2])
for (r in 3:5) {
b.treat[r] <- 0
}
treat[1] <- 0
treat[2] <- 1
}
Data:
list(count=structure(.Data=c(210,60,0,1,1,
66,32,0,0,2),.Dim=c(2,5)),
M=c(272, 100))
Inits:
list(a = c(NA, 0, 0, 0, 0), b.treat = c(NA, 0, NA, NA, NA))
node mean sd MC error 2.5% median 97.5% start sample
or.treat 1.718 0.4566 0.01769 0.9899 1.672 2.798 4001 6000
q[1,1] 0.7684 0.02534 9.308E-4 0.7169 0.7687 0.8142 4001 6000
q[1,2] 0.2206 0.02507 9.07E-4 0.1751 0.2204 0.2711 4001 6000
q[1,3] 2.679E-6 7.726E-5 1.35E-6 0.0 0.0 1.16E-7 4001 6000
q[1,4] 0.002745 0.002862 8.859E-5 7.079E-5 0.001846 0.01036 4001 6000
q[1,5] 0.008323 0.004679 1.075E-4 0.001878 0.007449 0.01956 4001 6000
q[2,1] 0.6702 0.04553 0.001336 0.5786 0.6709 0.7583 4001 6000
q[2,2] 0.3202 0.04607 0.001356 0.2312 0.319 0.4126 4001 6000
q[2,3] 2.18E-6 6.024E-5 1.051E-6 0.0 0.0 1.024E-7 4001 6000
q[2,4] 0.002394 0.002505 7.774E-5 6.522E-5 0.001599 0.008891 4001 6000
q[2,5] 0.00727 0.004142 9.662E-5 0.001613 0.006541 0.01733 4001 6000
Using the Poisson trick:...
model {
for(i in 1:2){
for (r in 1:5) {
count[i,r] ~ dpois(mu[i,r])
log(mu[i,r]) <- lambda[i] + a[r] + b.treat[r] * treat[i]
}
lambda[i] ~ dflat()
}
a[1] <- 0
for (r in 2:5){
a[r] ~ dnorm(0, 0.00001)
}
b.treat[1] <- 0
b.treat[2] ~ dnorm(0, 0.00001)
or.treat <- exp(b.treat[2])
for (r in 3:5) {
b.treat[r] <- 0
}
treat[1] <- 0
treat[2] <- 1
}
Data:
list(count=structure(.Data=c(210,60,0,1,1,
66,32,0,0,2),.Dim=c(2,5)))
Inits:
list(a = c(NA, 0, 0, 0, 0), b.treat = c(NA, 0, NA, NA, NA), lambda=c(0,0))