階層ベイズのモデル比較

階層ベイズのモデル比較をまとめてみました。
（備忘録として）

【DIC】
1. 個人のβの平均値を計算
このβとデータを使って、再度、対数尤度（mean_para.loglike）を計算
mean_para.Deviance <- -2 * mean_para.loglike

2. 各MCMCの対数尤度の平均（mean.loglike）を計算
mean.Deviance <- -2 * mean.loglike

3. pD & DIC
pD <- mean.Deviance - mean_para.Deviance
DIC <- mean.Deviance + pD

pDは、AICでいうところの罰則項みたいなもの。
ただ、パラメータ数が増えても、pDが必ず増えない場合もある。

【log marginal density】
bayesmには、log marginal densityを計算する関数が用意されているので、それを使えば一発で答えが出ます。

# LMD(log marginal density)
LMD <- llike[-(1:mcmc_start)]
logMargDenNR(LMD)

【WAIC】
渡辺先生が提唱されている指標。
検証用に2つ作ってみました。
（値は、同じです。）

# WAIC 1
library("matrixStats")
log_lik <- llike[-(1:mcmc_start)]
dim(log_lik) <- if (length(dim(log_lik))==1) c(length(log_lik),1) else
c(dim(log_lik)[1], prod(dim(log_lik)[2:length(dim(log_lik))]))

S <- nrow(log_lik)
n <- ncol(log_lik)

lpd <- colMeans(log_lik)
p_waic <- colVars(log_lik) # 正規化後のデータの各列の分散

elpd_waic <- lpd - p_waic
waic <- -2*elpd_waic

# WAIC 2
lppd <- sum(colMeans(log_lik))
p_waic <- sum(apply(temp,2,var))
waic2 <- -2 * lppd + 2 * p_waic

【ホールドイン、ホールドアウト】
機械学習などでおなじみの、モデルを作っていないデータを使ってどれくらいの予測精度になるかという指標。

まぁ、分かりやすいですね。