2

我正在尝试使用 R2OpenBUGS 在 OpenBUGS 中运行微分方程求解器。我已经用 OpenBUGS/Diff/Examples 文件夹中的指数衰减示例 (Example01.odc) 进行了尝试。

这个例子的第一部分是一个模拟,第二部分是我尝试应用的推理版本。两者都可以直接在 OpenBUGS 中运行。

当我尝试从 R 运行它时,日志中出现错误:

expected the collection operator c error pos 1443
variable ngrid is not defined 

我认为第二个是第一个的后续。

下面是我编写的从 R 运行它的代码:

# example of adjusting a differential equation using OpenBUGS   
library(R2OpenBUGS)
library(coda)
workingdir <- "M:/data/R/OpenBUGSDiff"
setwd(workingdir)
# time points
tgrid <- c(0.00000, 0.20000, 0.40000, 0.60000, 0.80000,
       1.00000, 1.20000, 1.40000, 1.60000, 1.80000,
       2.00000, 2.20000, 2.40000, 2.60000, 2.80000,
       3.00000, 3.20000, 3.40000, 3.60000, 3.80000,
       4.00000, 4.20000, 4.40000, 4.60000, 4.80000,
       5.00000, 5.20000, 5.40000, 5.60000, 5.80000,
       6.00000, 6.20000, 6.40000, 6.60000, 6.80000,
       7.00000, 7.20000, 7.40000, 7.60000, 7.80000,
       8.00000, 8.20000, 8.40000, 8.60000, 8.80000,
       9.00000, 9.20000, 9.40000, 9.60000, 9.80000,
       10.00000)
obs <- c(  0.7887,1.241,0.7051,0.7388,0.3903,
    0.2632,0.1174,0.549,-0.1767,0.02938,
    0.154,0.1465,0.07878,-0.5569,0.01293,
    0.2905,-0.2665,-0.3881,0.02517,-0.138,
    0.4004,0.2859,-0.1217,0.3961,0.3813,
    0.1846,-0.3581,0.3293,0.04089,0.01972,
    0.3203,0.5294,-0.1389,-0.3732,0.1341,
    -0.02432,0.2261,-0.3612,0.3131,-0.258,
    0.02948,-0.0208,0.1066,0.3796,-0.2645,
    0.1035,0.1001,-0.2415,0.06739,-0.1554,
    -0.2388)
# inference model 
Modele <- function() 
{
  solution[1:ngrid, 1] <- 
    ode(init[1], 
        tgrid[1:ngrid], 
        D(C[1], t), 
        origin, 
        tol)
  D(C[1], t) <- -lambda * C[1] 
  log.lambda ~ dnorm(0.0, tau.lambda); 
  lambda <- exp(log.lambda)
  for (i in 1:ngrid) 
  {
    obs[i] ~ dnorm(solution[i, 1], tau) 
  }
  tau ~ dgamma(a, b)
}
write.model(Modele,"Diff.txt")
# data definition
origin = 0.0
tol = 1.0E-3
ngrid = 51 
init = c(1.0)
a = 0.001
b = 0.001
tau.lambda = 0.01
data  <- list(obs=obs,
          tgrid=tgrid,
          origin=origin,
          tol=tol,
          ngrid=ngrid,
          init=init,
          a=a,
          b=b,
          tau.lambda=tau.lambda)
inits <- function(){
  list(log.lambda = 1.0,tau = 0.01)
}
diff.inf <- bugs(data,inits,model.file = "Diff.txt",
              parameters = c("lambda","tau"),
              n.chains = 1, n.iter = 3000,n.burnin=500,n.thin=10,
              working.directory=workingdir,
              debug = TRUE,
              codaPkg=T)

这是来自 BUGS 的完整错误消息:

model is syntactically correct
expected the collection operator c error pos 1443
variable ngrid is not defined
model must have been compiled but not updated to be able to change RN generator
BugsCmds:NoCompileInits
model must be compiled before generating initial values
model must be initialized before updating
model must be initialized before monitors used
model must be initialized before monitors used
model must be initialized before monitors used
model must be initialized before DIC can be monitored
model must be initialized before updating
model must be initialized before monitors used
DIC monitor not set

我试图用“obs”替换模型中的“C[]”,因为我看不到“C”的定义位置,但这并没有改变任何东西。

我检查了科学数字符号中的“e”而不是“E”,这在其他地方已被建议在 BUGS 中创建“变量未定义”错误,但这似乎不是问题。

欢迎任何想法,或任何其他显示如何在 R 中使用“解决方案”或“颂歌”的代码。

提前致谢

奥利弗

4

1 回答 1

2

存在与init节点相关的 BUGS 错误。在上面的模型中,您有一个与之关联的索引 ( init[1]),但在数据中,您只提供一个值init = c(1.0)。我认为 BUGS 在这里与索引混淆(即它需要数据中有两个或多个初始值,而您只提供一个)。当您在模型中init[1]替换为时,问题就消失了。init

#BUGS model
model{
  solution[1:ngrid, 1] <-  ode(init,  tgrid[1:ngrid],    D(C[1], t),    origin,     tol)
  D(C[1], t) <- -lambda * C[1] 
  log.lambda ~ dnorm(0.0, tau.lambda); 
  lambda <- exp(log.lambda)
  for (i in 1:ngrid) {
    obs[i] ~ dnorm(solution[i, 1], tau) 
  }
  tau ~ dgamma(a, b)
}

#data
list(obs = c(0.7887, 1.241, 0.7051, 0.7388, 0.3903, 0.2632, 0.1174, 0.549, -0.1767, 0.02938, 0.154, 0.1465, 0.07878, -0.5569, 0.01293, 0.2905, -0.2665, -0.3881, 0.02517, 0.138, 0.4004, 0.2859, -0.1217, 0.3961, 0.3813, 0.1846, -0.3581, 0.3293, 0.04089, 0.01972, 0.3203, 0.5294, -0.1389, -0.3732, 0.1341, -0.02432, 0.2261, -0.3612, 0.3131, -0.258, 0.02948, -0.0208, 0.1066, 0.3796, -0.2645, 0.1035, 0.1001, -0.2415, 0.06739, -0.1554, -0.2388), 
tgrid = c(0, 0.2, 0.4, 0.6, 0.8, 1, 1.2, 1.4, 1.6, 1.8, 2,  2.2, 2.4, 2.6, 2.8, 3, 3.2, 3.4, 3.6, 3.8, 4, 4.2, 4.4, 4.6,  4.8, 5, 5.2, 5.4, 5.6, 5.8, 6, 6.2, 6.4, 6.6, 6.8, 7, 7.2,  7.4, 7.6, 7.8, 8, 8.2, 8.4, 8.6, 8.8, 9, 9.2, 9.4, 9.6, 9.8,  10),
origin = 0, tol = 0.001, ngrid = 51, init = 1, a = 0.001,  b = 0.001, tau.lambda = 0.01))

#inits
list(log.lambda = 1.0,tau = 0.01)

如果您将上面的模型替换为以前的模型 ( Modele),一切都应该从 R 顺利运行。

一般观点:我发现直接在 OpenBUGS 中调试(使用上面的脚本和 GUI 来编译、加载数据和初始化)比通过 R2OpenBUGS 更容易。

于 2014-05-23T09:29:55.087 回答