0

我正在用 JAGS 中的 BUGS 代码在 R 中开发贝叶斯层次模型。

在我的模型中,我有两个矩阵,它们在相同的矩阵位置包含彼此相关的信息。我的信息按行排列。我按行对第一个矩阵 Distmat 应用数学运算:

diffmat[i,j] <- abs(Distmat[birthterr[i],j] - Dist[i]) 

我有兴趣在新向量中记录每行 diffmat 中每个最小值的列位置,然后将此向量应用于第二个矩阵。这在使用函数“which”或“which.min”的常规 R 代码中相对容易:

a <- numeric()
for (i in 1:dim(diffmat)[1])
  for (j in 1:dim(diffmat)[2])
a[i] <- which.min(diffmat[i,])

然后将向量“a”应用于第二个矩阵(terrmat)以获得与 Distmat 位置相关的值:

b <- numeric(0)
for (i in 1:dim(diffmat)[1])
  for (j in 1:dim(diffmat)[2])
b[i] <- terrmat[i, a[i]]

但是,显然 BUGS 代码无法识别 which 或 which.min(),我正在努力寻找一种将这些矩阵行位置存储在向量中的方法。也许有一个非常简单的解决方案,但我真的被困在那里。希望我的信息足够清楚。

任何建议将不胜感激。谢谢你的时间!

4

1 回答 1

0

这是一个最小的工作示例。这里的类似物就是x你的diffmat. 我在这里随机绘制它,但如果您以其他方式定义它,它应该仍然有效。本质上,您对x每行中的值进行排名并创建一个新矩阵e,该矩阵是一个虚拟矩阵,如果x[i,j]排名为 1,则编码为 1,否则为 0。然后,您从1:ncol(terrmat)假设中获取该值的内积和一个值向量,terrmat并且diffmat具有相同的维度。然后,这将为您提供观察的第一个排名值的列索引i。下面ymat的例子是你terrmat要去的地方。我认为任何实际大小的问题都会很慢,但它似乎可以从下面的输出中工作。

dl <- list(
ymat = matrix(1:3, ncol=3, nrow=5, byrow=TRUE),
yinds = 1:3
)
mods <- "model{
  for(i in 1:5){
    for(j in 1:3){
      x[i,j] ~ dnorm(0,1)
      e[i,j] <- equals(rx[i,j], 1)
    }
    rx[i,1:3] <- rank(x[i,1:3])
    ind[i] <- inprod(e[i,], yinds)
    yval[i] <- ymat[i,ind[i]]
  }
  
}"

library(runjags)

out <- run.jags(mods, data=dl, monitor="yval")

out
# 
# JAGS model summary statistics from 20000 samples (chains = 2; adapt+burnin = 5000):
#   
#         Lower95 Median Upper95   Mean      SD Mode     MCerr MC%ofSD SSeff       AC.10    psrf
# yval[1]       1      2       3 1.9973  0.8139    2 0.0058421     0.7 19409  -0.0077146 0.99996
# yval[2]       1      2       3 2.0067 0.81605    3 0.0057704     0.7 20000  0.00049096  1.0003
# yval[3]       1      2       3 1.9895  0.8142    2 0.0057573     0.7 20000  0.00066309       1
# yval[4]       1      2       3 1.9973 0.81638    1 0.0057727     0.7 20000 -0.00040016 0.99998
# yval[5]       1      2       3  1.993 0.81611    1 0.0057708     0.7 20000  -0.0027988 0.99996
# 
# Total time taken: 0.7 seconds


于 2021-03-16T18:02:26.757 回答