前回（16）では、”Bayesian Linear Regression(BLR)” と通常の ”Linear Regression(LR)” の回帰係数（推定値）が一致せず、その理由がよく分かりませんでした（図１）。
図１ BLR と LR の回帰係数（推定値）
上段が LRの、下段が BLRの Coefficientsであり、Interceptが大きく違っています。
そこで、
今回はデータ解析環境「R」で BLRをやってみました。
「R」には、
MCMC（マルコフ連鎖モンテカルロシミュレーション）による推定方法が用意されています。
これは乱数を使ってシミュレーション計算する方法で、「R」の環境があれば簡単に実行できます。
「RStan」での方法もありますが、少々厄介です。
そこで、
今回はデータ解析環境「R」で BLRをやってみました。
「R」には、
MCMC（マルコフ連鎖モンテカルロシミュレーション）による推定方法が用意されています。
これは乱数を使ってシミュレーション計算する方法で、「R」の環境があれば簡単に実行できます。
「RStan」での方法もありますが、少々厄介です。
「R」の実行：
事前にpackage（MCMApack)をインストールしておいて下さい。
事前にpackage（MCMApack)をインストールしておいて下さい。
データは前回の「ZTT-Protein.csv」を Excelで開き「ZTT～gGlb」のすべてを選択しグリップボードにコピーします。
そして、下記により読込み実行します。
そして、下記により読込み実行します。
***
dat<- read.delim("clipboard", header=T)
head(dat)
attach(dat)
dat<- read.delim("clipboard", header=T)
head(dat)
attach(dat)
library(MCMCpack)
fit <- MCMCregress(ZTT ~ gGlb , data=dat, burnin = 1000, mcmc = 10000, thin = 2, b0 = 0, B0 = 0,
fit <- MCMCregress(ZTT ~ gGlb , data=dat, burnin = 1000, mcmc = 10000, thin = 2, b0 = 0, B0 = 0,
c0 = 0.001 d0= 0.001)
summary(fit)
summary(fit)
出力結果（図１）
***
ここで、
事前情報が一様分布を当てはめているなら、
「事前分布の平均：b0=0、事前分布の分散：B0=0、逆ガンマ分布のshapeパラメータ：c0=0.001、逆ガンマ分布のscaleパラメータ：d0=0.001」としておきます。
事前情報が一様分布を当てはめているなら、
「事前分布の平均：b0=0、事前分布の分散：B0=0、逆ガンマ分布のshapeパラメータ：c0=0.001、逆ガンマ分布のscaleパラメータ：d0=0.001」としておきます。
注釈：
burnin：バーイン期間、mcmc:くり返し数、thin：サンプリング間隔
burnin：バーイン期間、mcmc:くり返し数、thin：サンプリング間隔
引数の詳細は、「Functions in MCMCpack」を参考にして下さい。
図１（出力結果）のとおり、MCMCregressによる BLRは、ほぼ LRの結果に近いものになりました。
JASPとの回帰係数（Intercept）の違いは、良くわかりません（ご教示頂ければ幸いです）。
JASPとの回帰係数（Intercept）の違いは、良くわかりません（ご教示頂ければ幸いです）。
次回に続く！
情報統計研究所をご気軽にご利用くださいませ。