WinBUGS
で単純なJolly-Seberモデルを実行できますが、Jags
では実行できません。私はJags
で線形回帰を実行することができ、R
はJags
を見つけて実行することができます。したがって、私は、Jags
がモデルコード内の1つ(またはそれ以上)の行を解釈できない可能性があると考えています。以下のコードを調べて、Jags
で実行するように修正する方法をお勧めします。当初私はおそらくprod
関数がJags
で利用できないと思われました。しかし、Jags
マニュアルを検索すると、Jags
にはprod
機能が含まれています。Jolly-Seberはジャグ付きではなくWinBUGSで動作します
これは私が考えることができる最も簡単な実例ですが、可能であればさらに簡素化します。
例のデータセットが下部にあります。モデルコードは、KeryとSchaub(2012)から変更されています。ここで
# BUGS code
sink("C:/Users/mmiller/Documents/simple R programs/my.model.txt")
cat("
model {
for (i in 1:M) {
for (t in 1:(n.occasions-1)) {
phi[i,t] <- mean.phi
}
for (t in 1:n.occasions) {
p[i,t] <- mean.p
}
}
mean.phi ~ dunif(0, 1)
mean.p ~ dunif(0, 1)
for (t in 1:n.occasions) {
gamma[t] ~ dunif(0, 1)
}
for (i in 1:M) {
z[i,1] ~ dbern(gamma[1])
mu1[i] <- z[i,1] * p[i,1]
y[i,1] ~ dbern(mu1[i])
for (t in 2:n.occasions) {
q[i,t-1] <- 1-z[i,t-1]
mu2[i,t] <- phi[i,t-1] * z[i,t-1] + gamma[t] * prod(q[i,1:(t-1)])
z[i,t] ~ dbern(mu2[i,t])
mu3[i,t] <- z[i,t] * p[i,t]
y[i,t] ~ dbern(mu3[i,t])
}
}
}
",fill=TRUE)
sink()
# run R2WinBUGS
setwd('C:/Users/mmiller/Documents/simple R programs')
library(R2WinBUGS)
data <- list(y = my.data, n.occasions = dim(my.data)[2], M = dim(my.data)[1])
inits <- function() {list(mean.phi = runif(1, 0, 1), mean.p = runif(1, 0, 1), z = my.data)}
parameters <- c("mean.p", "mean.phi")
bugs.out <- bugs(data, inits, parameters,
"C:/Users/mmiller/Documents/simple R programs/my.model.txt",
n.thin=1, n.chains=2, n.burnin=500, n.iter=1000, debug=FALSE, working.dir=getwd())
print(bugs.out, digits=2)
# run R2jags
library('R2jags')
data <- list(y = my.data, n.occasions = dim(my.data)[2], M = dim(my.data)[1])
inits <- function() {list(mean.phi = runif(1, 0, 1), mean.p = runif(1, 0, 1), z = my.data)}
parameters <- c("mean.p", "mean.phi")
jags.out2 <- jags(data = data, inits = inits, parameters,
model.file = "C:/Users/mmiller/Documents/simple R programs/my.model.txt",
n.thin=1, n.chains=2, n.burnin=500, n.iter=1000, working.dir=getwd())
print(jags.out2, digits=2)
例データセットされています
my.data <- read.table(text = '
1 1 1 0 0 0 0
0 1 1 1 1 0 0
0 1 0 0 0 0 0
0 0 1 1 0 1 0
0 1 1 0 0 0 0
1 0 0 0 0 0 0
1 0 0 0 0 0 0
1 0 0 0 0 0 0
1 0 0 0 0 0 0
1 0 1 0 0 0 0
1 0 0 0 0 0 0
1 0 1 1 1 0 0
1 0 0 0 0 0 0
1 0 0 0 0 0 0
1 1 0 0 0 0 0
0 0 0 0 1 0 0
1 0 0 0 0 0 0
1 1 0 0 0 0 0
1 0 1 0 0 0 0
0 1 0 0 0 0 0
0 1 1 0 0 0 0
1 0 0 0 0 0 0
0 1 1 0 0 0 0
0 0 1 0 0 0 0
0 0 0 1 0 0 0
0 1 0 0 0 0 0
0 0 0 0 1 1 0
0 1 1 0 0 0 0
0 1 1 0 0 0 0
0 0 1 1 0 0 0
0 0 0 1 0 0 0
0 0 1 1 0 0 0
0 0 1 1 1 0 0
0 0 1 1 0 0 0
0 0 1 0 0 0 0
0 0 0 1 0 1 1
0 0 0 0 1 1 0
0 0 0 1 0 0 0
0 0 0 1 0 0 0
0 0 0 0 1 0 0
0 0 0 1 0 0 0
0 0 0 0 1 1 1
0 0 0 1 0 0 0
0 0 0 1 0 1 1
0 0 0 1 0 0 0
0 0 0 0 1 0 0
0 0 0 1 0 0 0
0 0 0 0 1 0 0
0 0 0 0 0 1 1
0 0 0 0 0 0 1
0 0 0 0 1 1 0
0 0 0 0 0 1 1
0 0 0 0 0 1 1
0 0 0 0 0 1 1
0 0 0 0 0 1 0
0 0 0 0 0 1 0
0 0 0 0 0 1 1
0 0 0 0 0 0 1
0 0 0 0 0 0 1
0 0 0 0 0 0 1
0 0 0 0 0 1 1
0 0 0 0 0 1 0
0 0 0 0 0 0 1
0 0 0 0 0 0 1
0 0 0 0 0 0 1
0 0 0 0 0 0 1
0 0 0 0 0 0 1
0 0 0 0 0 0 1', header = FALSE)
nz <- 300
my.data <- rbind(my.data, matrix(0, ncol = ncol(my.data), nrow = nz))
dim(my.data)
head(my.data)
my.data <- as.matrix(my.data)
ここではJags
で実行するん線形回帰である:
# Linear regression in JAGS using R2jags
library('R2jags')
x <- rnorm(10)
mu <- -3.2 + 1.5 * x
y <- rnorm(10, mu, sd = 4)
cat ("
model {
for (i in 1:10) {
y[i] ~ dnorm(mu[i], tau)
mu[i] <- beta0 + beta1*x[i]
}
beta0 ~ dnorm(0, .01)
beta1 ~ dnorm(0, .01)
sigma ~ dunif(0,100)
tau <- 1/(sigma * sigma)
}
", file = 'C:/Users/mmiller/Documents/simple R programs/normal.txt')
data <- list(y=y, x=x)
inits <- function()
list(beta1 = rnorm(1), beta0 = rnorm(1), sigma = runif(1,0,2))
parameters <- c("beta0", "beta1", "sigma", "tau")
out <- jags(data = data, inits = inits, parameters,
model.file = 'C:/Users/mmiller/Documents/simple R programs/normal.txt',
n.thin=1, n.chains=2, n.burnin=2000, n.iter=6000, working.dir=getwd())
print(out, digits=2)
通常、ジャグのエラーメッセージはかなり有益です。これを実行しようとしたときに何を伝えましたか? –