事後点推定
ベイズ回帰モデルにおける、身長 \(X\) による体重 \(Y\) の尤度を振り返りましょう。$Y \sim N(m, s^2)$(ここで $m = a + b X$)。事後分布の100,000回のRJAGSシミュレーション weight_sim_big がワークスペースに用意されており、マルコフ連鎖の出力を格納したデータフレームも含まれています。
> head(weight_chains, 2)
a b s iter
1 -113.9029 1.072505 8.772007 1
2 -115.0644 1.077914 8.986393 2
切片と傾きのパラメータ \(a\) と \(b\) の事後平均は、体重と身長の関係における事後平均トレンドを表します。一方、\(a\) と \(b\) の完全な事後分布は、妥当なパラメータの範囲、すなわちトレンドに関する事後不確実性を表します。以下では、このトレンドとその不確実性を確認しましょう。bdims データもワークスペースに用意されています。
この演習はコースの一部です
RJAGS によるベイズモデリング
演習の手順
weight_sim_bigの連鎖に対してsummary()統計量を求めましょう。- \(b\) の事後平均は
summary()の表1に示されています。生のweight_chainsを使ってこの計算を確認しましょう。 bdimsのwgtとhgtのデータを散布図で描きましょう。geom_abline()を使って事後平均トレンドを重ねて表示しましょう。wgtとhgtの散布図をもう一枚作成しましょう。weight_chainsの最初の20セットの \(a\) と \(b\) パラメータ値で定義される20本の回帰直線を重ねて表示しましょう。
実践的なインタラクティブ演習
このサンプルコードを完成させて、この演習に挑戦してみましょう。
# Summarize the posterior Markov chains
# Calculate the estimated posterior mean of b
mean(___)
# Plot the posterior mean regression model
ggplot(bdims, aes(x = ___, y = ___)) +
geom_point() +
geom_abline(intercept = mean(___), slope = mean(___), color = "red")
# Visualize the range of 20 posterior regression models
ggplot(bdims, aes(x = ___, y = ___)) +
geom_point() +
geom_abline(intercept = ___[1:20], slope = ___[1:20], color = "gray", size = 0.25)