【2回目】100日後にベイズ統計ができるゆーみん【12日目】

こんばんは！今日は早めに帰れたので今日中にアップできそうです。

4.4.12 transformed parameters ブロックとgenerated quantitiesブロック

処理はR側でするのではなく、
できるだけStan側で実行するのが望ましいそうな。
特にtransformed parameters ブロックとgenerated quantitiesブロックが利用できるそうです。

data {
 int N;
 real X[N];
 real Y[N];
}

parameters {
 real a;
 real b;
 real<lower=0> sigma;
}

model {
 for (n in 1:N) {
   Y[n] ~ normal(a + b*X[n], sigma);
 }
}

↑これを拡張して

data {
 int N;
 real X[N];
 real Y[N];
 int N_new;
 real X_new[N_new];
}

parameters {
 real a;
 real b;
 real<lower=0> sigma;
}

transformed parameters {
 real y_base[N];
 for (n in 1:N)
   y_base[n] = a + b*X[n];
}

model {
 for (n in 1:N)
   Y[n] ~ normal(y_base[n], sigma);
}

generated quantities {
 real y_base_new[N_new];
 real y_new[N_new];
 for (n in 1:N_new) {
   y_base_new[n] = a + b*X_new[n];
   y_new[n] = normal_rng(y_base_new[n], sigma);
 }
}

細かく見ていきましょう！

・transformed parameters ブロック

transformed parameters {
 real y_base[N];
 for (n in 1:N)
   y_base[n] = a + b*X[n];
}

dataブロック、parameterブロックで宣言されたパラメータを四則演算やlogなどの関数を使って新たなサンプリングする変数を作ることができるのである！

・for文の{}の省略

 for (n in 1:N)
 
 for (n in 1:N)

for文の中身が一文だけの場合は{}が省略できるぞ！（いるのか）

・generated quantitiesブロックの使用

generated quantities {
 real y_base_new[N_new];
 real y_new[N_new];
 for (n in 1:N_new) {
   y_base_new[n] = a + b*X_new[n];
   y_new[n] = normal_rng(y_base_new[n], sigma);
 }
}

ここではdataブロック、parameterブロック、transformed parameters ブロックで宣言されたパラメータを四則演算やlogなどの関数を使って新たなサンプリングする変数を作ることができるのである！

は？transformed parameters ブロックの存在意義なんやねん？と思うやないですか、

あれなんですって、このブロックは事後確率と完全に切り離されていて計算が速いんだって。

このブロック内では〜が使えないのである。
「〜分布名（）」の代わりに「分布名_rng()」と言う関数を用いて乱数を生成して＝で代入するんだって。へぇ（なんで「〜分布名（）」じゃだめなんかと思ってる顔をしている）。

はい、今日はここまでです〜

実は分からなくて飛ばしてるところがあるので明日時間があれば取り組もうかなと思います。

では、また明日お会いしましょう〜