Stan による正規分布のパラメーター推定

正規分布

正規分布は、平均 μ と標準偏差 σ によって表される分布である。多くのデータが正規分布に従うために、よく利用される。ここで、平均 12.2、標準偏差 1.5 の正規分布から乱数を 100 個生成して、Stan でこれらの乱数が生成された正規分布のパラメーターを推定する例を示す。

x <- rnorm(n = 100, mean = 12.2, sd = 1.5)

Stan を利用して x を生成した正規分布の平均と標準偏差を推定する。Stan コードの data ブロックには、これから入力する値である x を定義する。Stan 自身は、x に何個のデータが入っているのかを計算できないので、外部からその個数を与える必要がある。このため、data ブロックには x の他に、x の要素数も合わせて定義する。次に、parameters ブロックには、これから推定したい正規分布のパラメーターである平均 μ と標準偏差 σ を定義する。最後に、model ブロックで、入力値変数とパラメーター変数を使って正規分布モデルを構築する。

data {
    int n;
    real x[n];
}
parameters {
    real mu;
    real sigma;
}
model {
    x ~ normal(mu, sigma);
}

R で Stan コードを呼び出して実行する。

d <- list(x = x, n = length(x))
fit <- stan(file = 'normdisp.stan', data = d)
fit
## Inference for Stan model: tmp.
## 4 chains, each with iter=2000; warmup=1000; thin=1;
## post-warmup draws per chain=1000, total post-warmup draws=4000.
## 
##         mean se_mean   sd   2.5%    25%    50%    75%  97.5% n_eff Rhat
## mu     12.17    0.00 0.16  11.84  12.06  12.17  12.27  12.49  3210    1
## sigma   1.59    0.00 0.11   1.39   1.51   1.59   1.67   1.84  2966    1
## lp__  -96.00    0.03 1.00 -98.64 -96.41 -95.70 -95.28 -95.01  1396    1
## 
## Samples were drawn using NUTS(diag_e) at Wed Jan 10 15:08:55 2018.
## 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).

各回 MCMC サンプリングの結果をプロットする場合は、stan_trace 関数を利用し、パラメーターの事後確率のヒストグラムをプロットする場合は stan_hist 関数を利用する。

stan_trace(fit)
stan_hist(fit)
stan による正規分布のパラメーター推定結果(サンプリング軌跡)
stan による正規分布のパラメーター推定結果(事後確率の分布)