受講される方は事前に環境の準備をしていただき,当日実際に手元で計算を進めていただくと,より理解が深まるかと思います。
Windowsユーザの方は以下の点にご注意ください
Rのサイトはこちらです。Rの最新バージョンは3.5.2ですので,ご準備ください。
RStudioのサイトはこちらです。RStudioの最新バージョンは1.1.463ですので,ご準備ください。OSに応じて最新版のインストーラをダウンロードし,実行してください。
RおよびRStudioのインストールが終われば,RStudioを起動し,次のコードを実行してください。
関連パッケージも含め,多くのパッケージが導入されますのでしばらくお待ちください。
install.packages('tidyverse')
install.packages('rstan')
install.packages('bayesplot')
install.packages('brms')
すでにRやRStudio,パッケージの導入をしたことがある人も,次のコードでパッケージを最新版にアップデートすることができますので,ご確認ください。
update.packages(ask=F)
環境の準備ができれば,次のコードをコピー&ペーストで実行してください。
library(rstan)
rstan_options(auto_write = TRUE)
options(mc.cores = parallel::detectCores())
stancode <- "
data {
int<lower=0> J; // number of schools
real y[J]; // estimated treatment effects
real<lower=0> sigma[J]; // s.e. of effect estimates
}
parameters {
real mu;
real<lower=0> tau;
real eta[J];
}
transformed parameters {
real theta[J];
for (j in 1:J)
theta[j] = mu + tau * eta[j];
}
model {
target += normal_lpdf(eta | 0, 1);
target += normal_lpdf(y | theta, sigma);
}
"
schools_dat <- list(J = 8, y = c(28, 8, -3, 7, -1, 1, 18, 12),
sigma = c(15, 10, 16, 11, 9, 11, 10, 18))
model.ex <- stan_model(model_code = stancode, model_name = "school")
fit.samp <- sampling(model.ex, data = schools_dat, iter = 1000,
chains = 4)
fit.samp
## Inference for Stan model: school.
## 4 chains, each with iter=1000; warmup=500; thin=1;
## post-warmup draws per chain=500, total post-warmup draws=2000.
##
## mean se_mean sd 2.5% 25% 50% 75% 97.5% n_eff Rhat
## mu 7.90 0.15 5.06 -1.77 4.54 7.83 11.22 17.64 1195 1.00
## tau 6.42 0.19 5.24 0.26 2.39 5.12 9.02 19.28 734 1.00
## eta[1] 0.41 0.02 0.91 -1.44 -0.16 0.44 1.02 2.11 2063 1.00
## eta[2] 0.01 0.02 0.89 -1.68 -0.58 0.00 0.60 1.82 1886 1.00
## eta[3] -0.17 0.02 0.91 -1.90 -0.78 -0.19 0.44 1.62 1962 1.00
## eta[4] 0.00 0.02 0.88 -1.73 -0.58 0.00 0.58 1.77 1743 1.00
## eta[5] -0.38 0.02 0.88 -2.11 -0.96 -0.37 0.20 1.44 1961 1.00
## eta[6] -0.17 0.02 0.90 -1.93 -0.75 -0.18 0.42 1.61 2361 1.00
## eta[7] 0.35 0.02 0.92 -1.52 -0.24 0.39 0.97 2.11 1617 1.00
## eta[8] 0.05 0.02 0.91 -1.74 -0.59 0.03 0.68 1.90 2157 1.00
## theta[1] 11.37 0.20 7.96 -1.36 6.14 10.42 15.48 30.10 1569 1.00
## theta[2] 7.91 0.13 6.39 -5.00 3.78 7.68 11.94 21.28 2328 1.00
## theta[3] 6.21 0.18 7.77 -13.05 2.18 6.79 11.11 20.26 1907 1.00
## theta[4] 7.81 0.14 6.56 -5.63 3.74 7.96 12.06 20.68 2236 1.00
## theta[5] 5.20 0.13 6.35 -7.56 1.20 5.53 9.41 16.56 2558 1.00
## theta[6] 6.36 0.14 6.76 -8.12 2.23 6.75 10.75 18.99 2183 1.00
## theta[7] 10.78 0.17 6.83 -0.80 6.12 10.17 15.03 25.92 1577 1.00
## theta[8] 8.39 0.17 7.81 -7.47 3.83 8.05 12.77 24.73 2171 1.00
## lp__ -39.54 0.10 2.56 -45.18 -41.16 -39.38 -37.74 -35.11 593 1.01
##
## Samples were drawn using NUTS(diag_e) at Tue Mar 5 18:21:17 2019.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at
## convergence, Rhat=1).
なにやら上記のような文字列が出ていたら準備完了です。
Rのコードについては,本サイトからコピー&ペーストで使ってもらって結構ですが,Rの基本的な使い方,RStudioの基本的な使い方,ggplotによる描画の方法などについて,「RユーザのためのRStudio入門」 を参考に予習しておいていただけると,理解が進むと思います。