Укажите дискретное распределение Вейбулла в JAGS или BUGS для R.

Я подгоняю модель Вейбулла к дискретным значениям, используя JAGS в R. У меня нет проблем с подгонкой Вейбулла к непрерывным данным, но у меня возникают проблемы, когда я переключаюсь на дискретные значения.

Вот некоторые данные и код для модели вейбулла в 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?


r jags winbugs
person colin    schedule 25.08.2017    source источник

Ответы (1)


arrow_upward
2
arrow_downward

Вы можете использовать «уловку с единицами» для реализации дискретного распределения Вейбулла в JAGS. Используя pmf здесь, мы можем создать функцию для генерации некоторых данных:

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

Чтобы «уловка с единицами» сработала, вам, как правило, нужно разделить вашу новую ПДС на некоторую большую константу, чтобы гарантировать, что вероятность находится между 0 и 1. Хотя кажется, что ПДС дискретного Вейбулла уже обеспечивает это, мы все же добавили некоторые в любом случае большая константа в модели. Итак, вот как модель выглядит сейчас:

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) вектор из единиц и большую константу в аргументе данных, 2) плотность массы дискретного Вейбулла и 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
person mfidino    schedule 28.08.2017