Вывод процесса чистой бортом с JAGS
Я пытаюсь использовать JAGS для определения уровня рождаемости в (стохастическом) процессе чистого рождения.
На языке химии эта модель эквивалентна реакции: X->2X со скоростью альфа *X (также можно рассматривать как модель цепной реакции)
Это код R, который я использую для генерации процесса (в фиксированное время), и код зазубрин, чтобы сделать вывод параметра alpha.
library(rjags)
y <- 1; # Starting number of "individuals"
N <- 25 # number of time samplings
alpha <- 0.2 # per-capita birth rate
# Generate the time series
for(i in 2:N) {
y<- c(y,y[i-1]+rpois(1,alpha*y[i-1]))
};
# The jags code
model_string <- "model{
for(i in 2:N) {
New[i] ~ dpois(alpha*y[i-1])
y[i] <- y[i-1] + New[i]
}
alpha ~ dunif(0, 2)
}"
# Create and run the jags model
model <- jags.model(textConnection(model_string), data = list(y = y,N = N), n.chains = 3, n.adapt= 10000)
update(model, 5000); # Burnin for 10000 samples
mcmc_samples <- coda.samples(model, variable.names=c("alpha"), n.iter=5000)
Когда я запускаю код, я получаю следующую ошибку:
Error in jags.model(textConnection(model_string), data = list(y = y, N = N), :
RUNTIME ERROR:
Compilation error on line 4.
y[2] is a logical node and cannot be observed
Я пробовал разные вещи, такие как помещение alpha*y[i-1] в новую переменную (скажем, lambda[i]) или изменение вызовов New[i] на New[i-1], но ничего не получалось. Есть идеи, почему это не удается? Еще один умный способ сделать это?
Заранее спасибо.