在 R 的 JAGS 或 BUGS 中指定离散威布尔分布

我使用 R 中的 JAGS 将威布尔模型拟合到离散值。将威布尔模型拟合到连续数据没有问题，但当我切换到离散值时，我遇到了麻烦。

以下是在 JAGS 中拟合威布尔模型的一些数据和代码：

#draw data from a weibull distribution
y <- rweibull(200, shape = 1, scale = 0.9)
#y <- round(y)

#load jags, specify a jags model.
library(runjags)

j.model ="
model{
for (i in 1:N){
y[i] ~ dweib(shape[i], scale[i])

shape[i] <- b1
scale[i] <- b2
}

#priors
b1 ~ dnorm(0, .0001) I(0, )
b2 ~ dnorm(0, .0001) I(0, )
}
"

#load data as list
data <- list(y=y, N = length(y))

#run jags model.
jags.out <-     run.jags(j.model,
                         data=data,
                         n.chains=3,
                         monitor=c('b1','b2')
                        )
summary(jags.out)

这个模型很合适。但是，如果我使用将 y 值转换为离散值y <- round(y)，并运行相同的模型，它失败并出现错误Error in node y[7], Node inconsistent with parents。每次我尝试时，节点的具体数量都会改变，但它总是一个很小的数字。

我知道我可以通过在所有值中添加一个非常小的数字来完成此运行，但是，这并没有考虑到数据是离散的事实。我知道存在离散威布尔分布，但如何在 JAGS 中实现离散威布尔分布？

您可以使用“ones 技巧”在 JAGS 中实现离散威布尔分布。使用 PMFhere https://en.wikipedia.org/wiki/Discrete_Weibull_distribution我们可以创建一个函数来生成一些数据：

pmf_weib <- function(x, scale, shape){

  exp(-(x/scale)^shape) - exp(-((x+1)/scale)^shape)
}

# probability of getting 0 through 200 with scale = 7 and shape = 4
probs <- pmf_weib(seq(0,200), 7, 4) 

y <- sample(0:200, 100, TRUE, probs ) # sample from those probabilities

为了使“个数技巧”发挥作用，您通常必须将新的 pmf 除以某个大常数，以确保概率在 0 到 1 之间。虽然看起来离散威布尔的 pmf 已经确保了这一点，但我们仍然添加了一些无论如何，模型中的常数很大。所以，模型现在是这样的：

j.model ="
data{ 
C <- 10000
for(i in 1:N){
ones[i] <- 1
}
}
model{
for (i in 1:N){
discrete_weib[i] <- exp(-(y[i]/scale)^shape) - exp(-((y[i]+1)/scale)^shape)

ones[i] ~ dbern(discrete_weib[i]/C)
}

#priors
scale ~ dnorm(0, .0001) I(0, )
shape ~ dnorm(0, .0001) I(0, )
}
"

请注意，我们在数据参数中添加了 1) 一个由 1 组成的向量和一个大常数，2) 离散威布尔的 pmf，以及 3) 我们通过伯努利试验运行该概率。

您可以使用上面相同的代码来拟合模型，以下是显示模型成功恢复参数值（尺度 = 7 且形状 = 4）的摘要。

       Lower95   Median  Upper95     Mean        SD Mode       MCerr MC%ofSD SSeff
scale 6.968277 7.289216 7.629413 7.290810 0.1695400   NA 0.001364831     0.8 15431
shape 3.843055 4.599420 5.357713 4.611583 0.3842862   NA 0.003124576     0.8 15126