Реконструкция переменных из объектов mcmc

Я использую rjags в качестве семплера. В модели определены 3 матрицы. Функция coda.samples возвращает список образцов. Если я возьму первый образец списка, имена столбцов будут выглядеть примерно так:

> colnames(output[[1]])
"A[1,1]"  "A[2,1]"  "A[1,2]"  "A[2,2]" ... 
"B[1,1]"  "B[2,1]"  "B[3,1]"  "B[4,1]" ... 
"C[1,1]"  "C[2,1]"

Очевидно, что A, B и C являются матрицами в моей модели. Я хочу реконструировать их на основе среднего значения этих выборок. Я могу легко получить средние значения с помощью colMeans(output[[1]]), но я понятия не имею, как легко восстановить матрицы из этого вектора.

Хорошим способом реконструкции будет функция relist(). Итак, если бы у меня были матрицы A, B и C в списке L = list(A=A,B=B,C=C), я мог бы преобразовать этот список в вектор с помощью unlist() и преобразовать обратно с помощью relist(). Я ищу что-то похожее/готовое для объектов mcmc, но пока безрезультатно - не могу поверить, что я первый, кому это нужно. Очевидно, что relist(colMeans(output[[1]])) не работает.

Кто-нибудь может помочь мне с реконструкцией?

Редактировать: также обратите внимание, что функции relist() нужен только скелет, поэтому извлечение скелета из colnames(output[[1]]) также поможет. Или я усложняю?


person Davor Josipovic    schedule 11.04.2017    source источник


Ответы (1)


Не думаю, что relist() поможет...

Я предполагаю, что ваш объект output является объектом класса mcmc.list, как определено в пакете R coda, а output[[1]] — объект класса mcmc, представляющий образец для первой цепочки MCMC.

Я почти уверен, что coda не понимает, что например "A[1,1]" — это матрица JAGS, она просто обрабатывает ее как имя переменной. Поэтому вам придется перебирать соответствующие переменные и накладывать структуру самостоятельно.

В идеале вы должны обернуть это в функцию, подобную следующей:

getMatrix <- function(output, varname, rows, cols) {
  unname(
    sapply(1:cols, function(j)
      sapply(1:rows, function(i)
        summary(output[,sprintf("%s[%s,%s]", varname, i, j)])[[1]][1]
      )
    )
  )
}

Таким образом, например, если ваша матрица B, хранящаяся в output[[1]], имеет 3 строки и 4 столбца, вы должны написать:

getMatrix(output[[1]], "B", 3, 4)

чтобы получить средства как матричный объект в R.

person Theodore Lytras    schedule 11.04.2017