Оценка неизвестной переменной отклика в JAGS — обучение без учителя

Я пытаюсь оценить значения отклика процентного покрытия (COV) по известным параметрам распределения. Я могу сделать это, указав данные ответа как NA в OpenBUGS (например, код ниже), но JAGS этого не допустит. Кто-нибудь знает, как я могу добиться этого в JAGS?

Я думаю, что это относится к категории «неконтролируемого статистического обучения».

model {
  for (i in 1:n.sites) {  # loop around sites
    # specify prior distribution forms for effectively unknown percentage  cover
    COV[i] ~ dbeta(a[i], b[i])T(r1[i], r2[i]) 
  }  
}

# DATA
list(n.sites=5, COV=c(NA, NA, NA, NA, NA), a=c(7.1,2.2,13,10,13),
     b=c(25,11,44,27,44), r1=c(0.05,0.1,0.2,0.1,0.2),
     r2=c(0.15,0.3,0.6,0.3,0.6) )

# INITS
list(COV=c(0.1, 0.2, 0.4, 0.2, 0.4))

parameters jags winbugs
person Chris Jones    schedule 16.02.2015    source источник

Ответы (1)


arrow_upward
0
arrow_downward

Если вы просто хотите смоделировать значения для COV, которые соответствуют вашему усеченному бета-распределению, вы можете опустить COV из своего списка данных. Например:

library(R2jags)

cat('
model {
  for (i in 1:n.sites) {
    COV[i] ~ dbeta(a[i], b[i])T(r1[i], r2[i]) 
  }  
}', file={M <- tempfile()})


dat <- list(n.sites=5, a=c(7.1, 2.2, 13, 10, 13), b=c(25, 11, 44, 27, 44), 
            r1=c(0.05, 0.1, 0.2, 0.1, 0.2), r2=c(0.15, 0.3, 0.6, 0.3, 0.6))

j <- jags(dat, NULL, 'COV', M, 1, 10000, DIC=FALSE, n.burnin=0, n.thin=1)

Вот несколько верхних строк матрицы из 10000 смоделированных COV векторов...

head(j$BUGSoutput$sims.list$COV)

##           [,1]      [,2]      [,3]      [,4]      [,5]
## [1,] 0.1169165 0.2889155 0.2543063 0.2083161 0.2426788
## [2,] 0.1494647 0.1430956 0.2867575 0.2410594 0.2795923
## [3,] 0.1200414 0.2093230 0.2736719 0.2189734 0.2469634
## [4,] 0.1472082 0.1442609 0.2911482 0.2625216 0.2714883
## [5,] 0.1403574 0.1100977 0.2556352 0.1918480 0.2353231
## [6,] 0.1310404 0.1677148 0.3011752 0.1974136 0.2131811

ИЗМЕНИТЬ

Поскольку вы просто берете образцы из известного дистрибутива, вы также можете сделать это в чистом R (пакет distr предоставляет некоторые функции, которые помогают с выборками из усеченных дистрибутивов).

library(distr)
n <- 10000 # How many samples?
COV <- mapply(function(shape1, shape2, min, max) {
  r(Truncate(Beta(shape1, shape2), min, max))(n)
}, shape1=c(7.1, 2.2, 13, 10, 13), shape2=c(25, 11, 44, 27, 44), 
   min=c(0.05, 0.1, 0.2, 0.1, 0.2), max=c(0.15, 0.3, 0.6, 0.3, 0.6))

Выше вы должны передать векторы равной длины shape1, shape2, min и max в mapply, которые будут генерировать n случайных бета-распределенных переменных для каждого shape1 и соответствующих shape2, min и max по очереди.

Мы можем сравнить плотность ядра столбцов COV (наши образцы чистого R) с плотностью столбцов j$BUGSoutput$sims.list$COV (наши образцы JAGS).

par(mfrow=c(3, 2), mar=c(3, 0.5, 0.5, 0.5))
sapply(1:5, function(i) {
  djags <- density(j$BUGSoutput$sims.list$COV[, i])
  dr <- density(COV[, i])
  plot(djags, lwd=4, col='gray80', main='', ylab='', xlab='', yaxt='n', 
       ylim=c(0, max(djags$y, dr$y)))
  lines(dr)
})
plot.new()
legend('topleft', c('JAGS', 'R'), col=c('gray80', 'black'), lwd=c(4, 1), bty='n')

введите здесь описание изображения

person jbaums    schedule 16.02.2015
comment
Это отлично, спасибо JB! Оказывается, аргумент DIC=F здесь имеет решающее значение, можете ли вы сказать мне, почему он не будет работать, если DIC=T? - person Chris Jones; 17.02.2015
comment
@Chris - DIC по умолчанию равно TRUE. Когда TRUE, это говорит JAGS рассчитать отклонение модели, которое представляет собой сумму отклонений по наблюдаемым стохастическим узлам. В этом случае у вас нет наблюдаемых узлов — вы не подгоняете модель к данным, а скорее выбираете COV из определенного распределения. Для JAGS нет истинных COV, с которыми можно было бы сравнить смоделированные данные, и поэтому невозможно рассчитать отклонение (т. е. отсутствует отклонение). - person jbaums; 17.02.2015