Ошибка JAGS - возможный направленный цикл с участием некоторых или всех следующих узлов - PullRequest
Новая
       109

Ошибка JAGS - возможный направленный цикл с участием некоторых или всех следующих узлов

0 голосов
/ 29 июня 2018

Полный набор данных содержит ~ 11 000 строк. Я выполняю код с K = 400, проверяя, что код выполняется.

Все строки относятся к конкретной ячейке на карте и содержат информацию, извлеченную из изображений Sentinel-2 и цифровой карты высот.

Подгруппа из 117 ячеек также содержит ковариаты среды обитания, зарегистрированные во время экскурсии. Таким образом, некоторые из столбцов, включая переменные ответа (S1 и S2) и tussac, характеризуются высокой долей NA.

Код: 

add_c4 <- "model{
for(i in 1:K) {
S1[i]~dpois(lambda1[i])
lambda1[i]<-exp(a0+a1*DEM_slope[i]+a2*DEM_elevation[i]+a3*tussac[i]+a4*S2[i])

S2[i]~dpois(lambda2[i])
lambda2[i]<-exp(c0+c1*DEM_slope[i]+c2*DEM_elevation[i]+c3*tussac[i]+c4*S1[i])

muLogit_tussac[i]<-b0 + sentinel1[i] + sentinel3[i] + sentinel7[i] + sentinel8[i] + sentinel9[i] + DEM_slope[i]

Logit_tussac[i]~dnorm(muLogit_tussac[i], tau)
logit(tussac[i])<-Logit_tussac[i]
}

# Priors

a0~dnorm(0, 10)
a1~dnorm(0, 10)
a2~dnorm(0, 10)
a3~dnorm(0, 10)
a4~dnorm(0, 10)

b0~dnorm(0, 10)
b1~dnorm(0, 10)
b2~dnorm(0, 10)
b3~dnorm(0, 10)

c0~dnorm(0, 10)
c1~dnorm(0, 10)
c2~dnorm(0, 10)
c3~dnorm(0, 10)
c4~dnorm(0, 10)

tau~dgamma(0.001, 0.001)

#data# S1, S2, K, sentinel1, sentinel3, sentinel7, sentinel8, sentinel9, DEM_slope, DEM_elevation
#inits# a0, a2, a3, a4, b0, b1, b2, b3, c0, c2, c3, c4
#monitor# a0, a1, a2, a3, a4, b0, b1, b2, b3, tau, ped, dic, c0, c1, c2, c3, c4
}"

Когда я включаю 'c4 * S1 [i]', я получаю следующую ошибку: 

Possible directed cycle involving some or all of the following nodes

Затем происходит перечисление всех значений S1, S2, lambda1 и lambda2.

Удаление 'c4 * S1 [i]' приводит к выполнению кода.

Я просмотрел следующие темы:

Возможная направленная ошибка цикла в JAGS

https://stats.stackexchange.com/questions/220312/coding-a-jags-error-model-for-a-dependent-variable-that-has-increasing-variance

Проблемы в них, похоже, были вызваны тем, что плакат использовал «y» с обеих сторон уравнения. Я думаю, что моя проблема вызвана тем фактом, что a4 отправляет код в секцию S2, а c4 отправляет его обратно в секцию S1, что немного похоже на направленный цикл. Есть идеи как обойти это?

Я включил верхние строки набора данных на случай, если он пригодится: 

S1 S2 Logit_tussac moisture DEM_slope DEM_aspect DEM_elevation sentinel1 sentinel2 sentinel3 sentinel4 sentinel5 sentinel6 sentinel7 sentinel8 sentinel9 sentinel10
NA NA     NA            NA  2.434239   168.5011   0.588606366    0.0413    0.0499    0.0531    0.1035    0.1862    0.1968    0.1808    0.1318    0.0400     0.0199
NA NA     NA            NA  3.705001   178.1289   1.007037127    0.0966    0.1108    0.1212    0.0855    0.0917    0.1063    0.0937    0.1842    0.0341     0.0161
NA NA     NA            NA  5.006181   180.0000   1.883010797    0.1309    0.1472    0.1361    0.0855    0.0917    0.1063    0.0937    0.1572    0.0341     0.0161
NA NA     NA            NA  5.006181   180.0000   2.758984468    0.0542    0.0512    0.0472    0.0145    0.0127    0.0092    0.0166    0.0510    0.0148     0.0080

Подмножество набора данных, чтобы только 117 строк содержали данные удаленного и локального зондирования: 

S1 S2 Logit_tussac moisture DEM_slope DEM_aspect DEM_elevation sentinel1 sentinel2 sentinel3 sentinel4 sentinel5 sentinel6 sentinel7 sentinel8 sentinel9 sentinel10
NA NA        NA        NA   14.917334   256.1612      12.24432    0.0513    0.0588    0.0541    0.1145    0.1676    0.1988    0.1977    0.1658    0.1566     0.0770
0  0  -9.210240         1   23.803741   225.1231      16.88028    0.1058    0.1370    0.2139    0.2387    0.2654    0.2933    0.3235    0.2928    0.3093     0.1601
NA NA        NA        NA   20.789165   306.0945      18.52480    0.0287    0.0279    0.0271    0.0276    0.0290    0.0321    0.0346    0.0452    0.0475     0.0219
NA NA -9.210240         1    6.689442   287.9641      36.08975    0.0462    0.0679    0.1274    0.1535    0.1797    0.2201    0.2982    0.2545    0.4170     0.2252
0  0  -9.210240         1   25.476444   203.0659      23.59964    0.0758    0.1041    0.1326    0.1571    0.2143    0.2486    0.2939    0.2536    0.3336     0.1937
1  0  -1.385919         3    1.672511   270.0000      39.55215    0.0466    0.0716    0.1227    0.1482    0.2215    0.2715    0.3334    0.2903    0.3577     0.1957

1 Ответ

0 голосов
/ 07 сентября 2018

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

Например, возьмем следующую модель, содержащую направленный цикл: 

model <- 'model{

    for(i in 1:N){
        a[i] ~ dnorm(b[i], tau)
        b[i] ~ dnorm(a[i], tau)
    }

    tau ~ dgamma(0.01,0.01)

    #monitor# tau
    #data# N

}'

N <- 10
runjags::run.jags(model)

Нет разумного способа оценить эту модель, и JAGS скажет вам об этом. Но теоретически можно оценить эту модель: 

model <- 'model{

    for(i in 1:N){
        a[i] ~ dnorm(b[i], tau)
        b[i] ~ dnorm(a[i], tau)
    }

    tau ~ dgamma(0.01,0.01)

    #monitor# tau
    #data# N, a

}'

N <- 10
a <- rnorm(N)
runjags::run.jags(model)

Что изменилось, так это то, что все a [] теперь фиксированы (наблюдаются), поэтому мы действительно можем оценить эту модель. Но JAGS все еще обнаруживает направленный цикл, поэтому требуется обходной путь: 

model <- 'model{

    for(i in 1:N){
        a[i] ~ dnorm(b[i], tau)
        b[i] ~ dnorm(aa[i], tau)
    }

    tau ~ dgamma(0.01,0.01)

    #monitor# tau
    #data# N, a, aa

}'

N <- 10
a <- rnorm(N)
aa <- a
runjags::run.jags(model)

Это скрывает направленный цикл, заставляя ЯГОВ думать, что a [] и aa [] не связаны. Но это работает только тогда, когда все a [] наблюдаются / фиксируются, в противном случае отсутствующие aa [] не оцениваются или не определяются в модели. В вашем случае S1 [] и S2 [] кажутся частично наблюдаемыми, поэтому этот трюк не сработает, если вы просто не пропустите строки / наблюдения с отсутствующими S1 или S2 (что может оказаться невозможным, поскольку вы говорите, что у них высокий доля НК).

В противном случае вам придется как-то переформулировать свою модель, чтобы разорвать направленный цикл. Это потребует размышлений о биологическом процессе, лежащем в основе вашей системы, и о том, как представить требуемые отношения, не создавая направленного цикла, поэтому мы не можем помочь в этом.

Надеюсь, это поможет,

Мэтт

...