回帰の事前分布を可視化する
前の演習では、体重 \(Y\) を身長 \(X\) で予測するベイズ回帰モデル($Y \sim N(m, s^2)$、平均 $m = a + bX$)の各パラメータ(\(a\), \(b\), $s$)について、10,000 個の samples をシミュレートしました。samples の各行に含まれる \(a\), \(b\), \(s\) の値は、事前分布として妥当な回帰シナリオの 1 つを表しています。これらの事前シナリオの広がりを探るために、最初の 12 組の事前パラメータ \(a\), \(b\), \(s\) のそれぞれから、身長と体重のペアを 50 組シミュレートします。
この演習はコースの一部です
RJAGS によるベイズモデリング
演習の手順
samplesの最初の 12 組の事前パラメータをn = 50回複製したデータフレームprior_simulationを作成します(合計 600 行になります)。prior_simulationの 600 行それぞれについて、以下を実行します。- \(N(170, 10^2)\) モデルから
heightの値をシミュレートします。 - \(N(a + bX, s^2)\) モデルから
weightの値をシミュレートします。ここで \(X\) は身長、\((a, b, s)\) は事前パラメータの組です。
- \(N(170, 10^2)\) モデルから
- これで 12 組のパラメータセットそれぞれについて、
heightとweightのペアが 50 組シミュレートされました。ggplot()を使って、各パラメータのsetごとに 50 組の散布図を作成しましょう。weightを y 軸に配置することを忘れずに!
実践的なインタラクティブ演習
このサンプルコードを完成させて、この演習に挑戦してみましょう。
# Replicate the first 12 parameter sets 50 times each
prior_scenarios_rep <- bind_rows(replicate(n = ___, expr = samples[1:12, ], simplify = FALSE))
# Simulate 50 height & weight data points for each parameter set
prior_simulation <- prior_scenarios_rep %>%
mutate(height = rnorm(n = 600, mean = ___, sd = ___)) %>%
mutate(weight = rnorm(n = 600, mean = ___, sd = ___))
# Plot the simulated data & regression model for each parameter set
ggplot(prior_simulation, aes(x = ___, y = ___)) +
geom_point() +
geom_smooth(method = "lm", se = FALSE, size = 0.75) +
facet_wrap(~ set)