R: RStanをためす：Taglibro de H：SSブログ

	ブログをはじめるログイン

2012年8月にみた鳥｜R: RStanをためす (2) ブログトップ

R: RStanをためす　[統計] [編集]

あたらしいBayes推定ソフトウェアであるStan 1.0.0と、そのRインターフェスのRStan 1.0.0が発表されている(Revolutions: RStan: Fast, multilevel Bayesian modeling in R)。

Stan → C++ → ネイティブバイナリ
と、コンパイルして高速に計算できということだ。

さっそくためしてみる。

library(rstan)

x <- rpois(30, 4)

model <- '
  data {
    int<lower=0> n;        // number of data
    int<lower=0> x[n];     // data
  } 
  parameters {
    real<lower=0> lambda;  // mean 
  } 
  model {
    for (i in 1:n) {
      x[i] ~ poisson(lambda);
    }
    lambda ~ gamma(0.001, 0.001);
  } 
'

data <- list(n = length(x), x = x)

fit <- stan(model_code = model, data = data, 
            warmup = 500, iter = 1500, chains = 4)

結果

> fit >- stan(model_code = model, data = data, 
+             warmup = 500, iter = 1500, chains = 4)

TRANSLATING MODEL 'anon_model' FROM Stan CODE TO C++ CODE NOW.
COMPILING THE C++ CODE FOR MODEL 'anon_model' NOW.
SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 1).
Iteration: 1500 / 1500 [100%]  (Sampling)

SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 2).
Iteration: 1500 / 1500 [100%]  (Sampling)

SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 3).
Iteration: 1500 / 1500 [100%]  (Sampling)

SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 4).
Iteration: 1500 / 1500 [100%]  (Sampling)

> 
> print(fit)
Inference for Stan model: anon_model.
4 chains: each with iter=1500; warmup=500; thin=1; 1500 iterations saved.

       mean se_mean  sd 2.5%  25%  50%  75% 97.5% n_eff Rhat
lambda  3.6       0 0.4  3.0  3.4  3.6  3.9   4.3   824    1
lp__   29.8       0 0.7 27.8 29.6 30.1 30.3  30.3   745    1

Sample were drawn using NUTS2 at Sat Sep  1 10:43:27 2012.
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).