Не удается скопировать код RStan ESS из бумаги Vehtari - PullRequest
0 голосов
/ 01 мая 2020

Я пытаюсь повторить расчет ESS (эффективный размер выборки), используя метод Vehtari et al. в: Нормализация рангов, свертывание и локализация: улучшенная Rhat для оценки сходимости MCM C

Я работаю с кодом здесь: https://github.com/avehtari/rhat_ess/blob/master/code/monitornew.R

  # Geyer's initial positive sequence
  rho_hat_t <- rep.int(0, n_samples)
  t <- 0
  rho_hat_even <- 1
  rho_hat_t[t + 1] <- rho_hat_even
  rho_hat_odd <- 1 - (mean_var - mean(acov[t + 2, ])) / var_plus # 251
  rho_hat_t[t + 2] <- rho_hat_odd
  while (t < nrow(acov) - 5 && !is.nan(rho_hat_even + rho_hat_odd) &&
         (rho_hat_even + rho_hat_odd > 0)) {
    t <- t + 2
    rho_hat_even = 1 - (mean_var - mean(acov[t + 1, ])) / var_plus # 256
    rho_hat_odd = 1 - (mean_var - mean(acov[t + 2, ])) / var_plus # 257
    if ((rho_hat_even + rho_hat_odd) >= 0) {
      rho_hat_t[t + 1] <- rho_hat_even
      rho_hat_t[t + 2] <- rho_hat_odd
    }
  }

Я могу следовать коду из статьи, кроме случаев, когда мы дойдем до уравнения 10 в статье (вычисление перекрестной автокорреляции). Код (строки 251, 256 и 257) отображается в виде:

1 - (mean_var - mean(acov[t + 1, ])) / var_plus

, который близок к уравнению 10, за исключением пропущенных в уравнении 10 членов 's':

\sqrt{foo}

Я не вижу нигде в коде, что это каким-то образом учитывается в другом месте в способе вычисления. Я попытался вставить термины «s» обратно в эти строки кода, и это сильно влияет на окончательное значение ESS.

Может ли кто-нибудь помочь мне понять несоответствие между бумагой и кодом?

Спасибо.

1 Ответ

1 голос
/ 01 мая 2020

В формуле, приведенной в статье, s^2 - это оценка дисперсии, а rho - оценка автокорреляции. Таким образом, s^2 * rho является оценкой автоковариации, которую вы видите в коде.

...