差分

このページの2つのバージョン間の差分を表示します。

この比較画面へのリンク

両方とも前のリビジョン前のリビジョン
次のリビジョン		前のリビジョン
stan:r [2018/04/18 19:49] – watalustan:r [2018/04/18 20:18] (現在) – watalu
行 1: 行 1:
 ==== StanをRから使う ==== ==== StanをRから使う ====
 +
 +Stanが生成するのは、事後分布に従う乱数列(数値のデータ)である。この設計思想は実に潔い。Stanのユーザは必ず、Stan以外にも何らかのソフトウェアを用いて、メディアンを計算したり、ヒストグラムを描いたりすることになる。
  
 Rを通じてStanを利用すると決めた場合は、Stanを単独でインストールするのではなく、次のような手順を踏んでインストールを進める。 Rを通じてStanを利用すると決めた場合は、Stanを単独でインストールするのではなく、次のような手順を踏んでインストールを進める。
行 19: 行 21:
 |Windows|[[https://cran.r-project.org/bin/windows/Rtools/|RTools]]| |Windows|[[https://cran.r-project.org/bin/windows/Rtools/|RTools]]|
 |macOS|[[https://itunes.apple.com/jp/app/xcode/id497799835|Xcode]]| |macOS|[[https://itunes.apple.com/jp/app/xcode/id497799835|Xcode]]|
-|Linux|[[https://cran.r-project.org/|R]]|GCCおよび必要なライブラリ一式|+|Linux|GCCおよび必要なライブラリ一式|
  
 === RとRStudio === === RとRStudio ===
行 59: 行 61:
 ここまで終われば、たくさんのパッケージがインストールされて、やっとStanをRから呼び出す準備が整う。 ここまで終われば、たくさんのパッケージがインストールされて、やっとStanをRから呼び出す準備が整う。
  
-=== RStanの使い方 ===+=== RStanの使い方(概要) === 
 + 
 +RStanというパッケージは、 
 +  - Stanモデリング言語で記述されたファイルから、Stanでモデルファイルを作る 
 +  - データとモデルを与える未知パラメータの事後分布からのサンプリングを行う 
 + 
 +という順序で動作する。 
 +まずRで 
 +<code> 
 +library(rstan) 
 +</code> 
 +を実行し、必要なパッケージを読み込む。 
 +R上でこの２つのステップを同時に実行するには 
 +<code> 
 +binomial.fit <- stan("binomial.stan",data=list(x=20,n=30)) 
 +</code> 
 +とすればよい。 
 +２段階を分けて実行するには 
 +<code> 
 +binomial.model <- stan_model("binomial.stan") 
 +binomial.fit <- sampling(binomial.model, data=list(x=20,n=30)) 
 +</code> 
 +とする。モデルファイルを作る際にコンパイルを行っていて、ここに少々の時間がかかる。一度コンパイルしたモデルファイルは、RDSという形式でRの外に保存することができる。 
 +<code> 
 +saveRDS(binomial.model,"binomial.rds") 
 +</code> 
 +読み込むには、次のreadRDS関数を用いる。 
 +<code> 
 +binomial.model <- readRDS("binomial.rds") 
 +</code> 
 + 
 +=== RStanの使い方(スケルトン) === 
 + 
 +以下では、Stanモデリング言語を記述したファイルをexample.stan、観測される確率変数をx、モデルの未知パラメータをthetaと置く。 
 + 
 +  - Stanモデリング言語で記述したコードを用意する。 
 +  - Rの中で観測値を入力する。(以下では、example.obsというオブジェクトとする) 
 +  - Rからモデルを作成するように指示を出す。<code>example.stan <- stan_model("example.stan")</code> 
 +  - 作成したモデルに基づいて、事後分布に従う乱数を生成するように指示を出す。<code>example.nuts <- sampling(example.stan, data=list(x=example.obs), chains=4)</code> 
 +  - NUTSが法則収束していることを確認する。<code>stan_trace(example.nuts, inc_warmup=TRUE) 
 +stan_dens(example.nuts, pars="theta", separate_chains=TRUE) 
 +stan_ac(example.nuts,pars="theta",separate_chains=TRUE) 
 +stan_rhat(example.nuts) 
 +stan_ess(example.nuts) 
 +</code> 
 +  - 事後分布を確認する。<code>print(example.nuts)</code> 
 +  -事後分布を目視で確認し、併せて先ほどの平均やメディアンが合致しているかどうかを確認する。<code>stan_hist(example.nuts, pars="theta")</code> 
 +  - 事後分布と事後信頼区間を得る。<code>stan_plot(example.nuts, point_est="median", show_density=TRUE, ci_level=0.95, outer_level=1.00,show_outer_line=TRUE)</code>