Слишком низкие значения psrf в сводке runjags?

Runjags сообщает об очень низком уровне psrf = 1.0047 для цепочки, у которой явно есть проблемы с конвергенцией:

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

> print(o, vars = "q_date2")

JAGS model summary statistics from 3000 samples (chains = 3; adapt+burnin = 1000):

         Lower95   Median   Upper95      Mean       SD    MCerr MC%ofSD SSeff    AC.10   psrf
q_date2 -0.17611 -0.10467 0.0053844 -0.087376 0.063296 0.023726    37.5     7 0.022495 1.0047

Когда я пытаюсь вычислить psrf с помощью кода, я получаю результат, который выглядит гораздо более разумным:

> gelman.diag(as.mcmc.list(o)[,'q_date2'], transform=FALSE, autoburnin=FALSE)
Potential scale reduction factors:

     Point est. Upper C.I.
[1,]       3.54       7.94

Так почему же psrf, о котором сообщают беглецы, так низки? Это проблема бегунов или я что-то не так делаю?

Я использую текущую версию runjags (1.2.1-0) в R 3.1.0.

РЕДАКТИРОВАТЬ: во время создания резюме я получил предупреждения - извините, что не упомянул их раньше:

Warning messages:
1: In autocorrs[x$stochastic] <- x$autocorr[4, ] :
  number of items to replace is not a multiple of replacement length
2: In psrfs[x$stochastic] <- x$psrf$psrf[, 1] :
  number of items to replace is not a multiple of replacement length

r jags runjags
person Tomas    schedule 20.11.2014    source источник

Ответы (1)


arrow_upward
1
arrow_downward

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

Тем временем вы можете либо (а) игнорировать psrf, указанный в выводе для summary(), и вместо этого использовать RJout$psrf (или свой собственный код), либо (б) удалить отслеживаемые переменные (в данном случае M), которые вызывая проблему. Разрабатываемая версия runjags имеет еще лучшие решения (сводная статистика и графики (пере)рассчитываются после возврата модели) и должна появиться в CRAN в ближайшие пару месяцев.

Это также хорошее напоминание о том, что ручная проверка графиков трассировки является неотъемлемой частью анализа MCMC :)

person Matt Denwood    schedule 20.11.2014
comment
Спасибо! И на самом деле, если я запущу gelman.diag на MCMC с переменной M, то есть gelman.diag(as.mcmc.list(outRJ)), он скажет Ошибка в chol.default(W): старший минор порядка 1 не является положительно определенным. Так что это, вероятно, причина, почему эта проблема возникла! - person Tomas; 20.11.2014