Понимание оценки ошибок, используемых в пакете mcmc - PullRequest
Новая
       20

Понимание оценки ошибок, используемых в пакете mcmc

1 голос
/ 19 марта 2019

Я смотрю на Пример учебника MCMC из пакета mcmc в R .Ниже приведен код, где MCMC используется для вычисления параметров в примере логистической регрессии. 

library(mcmc) 

data(logit)
foo <- logit

out <- glm(y ~ x1 + x2 + x3 + x4, data = foo, family = binomial())

x <- foo
x$y <- NULL 
x <- as.matrix(x)
x <- cbind(1, x)
dimnames(x) <- NULL
y <- foo$y

lupost <- function(beta, x, y){
  eta <- x %*% beta
  p <- 1/(1 + exp(-eta))

  logl <- sum(log(p[y == 1])) + sum(log(1-p[y == 0]))

  return(logl + sum(dnorm(beta, 0, 2, log = TRUE)))
}

set.seed(42)
beta.init <- as.numeric(coefficients(out))
out <- metrop(lupost, beta.init, 1000, x=x, y = y)

Идея состоит в том, чтобы вычислить стандартную ошибку Монте-Карло.Поэтому используется 

> apply(out$batch, 2, mean)
[1] 0.6531950 0.7920342 1.1701075 0.5077331 0.7488265
[6] 0.5145751 0.7560775 1.4973807 0.3913837 0.7244162

.У меня вопрос, что означает последние 5 столбцов в выводе этой команды? В руководстве говорится, что это как-то E(X^2).Но в какой строке генерируется X?Что здесь X?

1 Ответ

1 голос
/ 19 марта 2019

Если я введу только тот код, который вы написали в своем вопросе выше, я получу только 5 цифр: 

apply(out$batch, 2, mean)
# [1] 0.5990174 0.8027482 1.1539329 0.4547702 0.7172724

Похоже, что учебник, которому вы следуете, - это mcmc::demo виньетка. Виньетка выполняет все шаги, которые вы опубликовали выше, , а затем вычисляет 

out <- metrop(out, nbatch = 100, blen = 100, outfun = function(z) c(z, z^2))

с последующим 

apply(out$batch, 2, mean)
# [1] 0.6531950 0.7920342 1.1701075 0.5077331 0.7488265
# [6] 0.5145751 0.7560775 1.4973807 0.3913837 0.7244162

Здесь вы можете видеть, что outfun вычисляет X и X ^ 2, и что вы получаете оценки E [X] и E [X ^ 2], когда берете среднее значение, используя apply().

Здесь X, похоже, взяты из апостериорного распределения параметров модели (бета). Обратите внимание, что байесовские апостериорные средние (первые 5 чисел) довольно близки к небайесовским точечным оценкам, рассчитанным с помощью glm в четвертой строке вашего кода.

...