rstanを利用してベイズ推定を行いたい。

Question

前提・実現したいこと

ベイズ推定を行いたいです。
Rおよびstanに関してはほぼほぼ初心者です。
根底からの理解は求めておらず、推定が完了できれば良いです。

環境導入、サンプルコードによる推定の実施は完了しています。

発生している問題・エラーメッセージ

自分で書いたコード、モデルで推定を実施しようとすると、以下のエラーメッセージが出ます。

Stan model "モデルの名称" does not contain samples.

知りたいこと

サンプルコードと見比べたところ、データをstanに渡せていない訳ではないと考えられる。
・このエラーメッセージが出るということは、どのような現象が起きているのか。
・stanのmodelの設定自体が間違っているのか。

ということを知りたいです。

参考にしたコード

model.stan

data {
    int t;
    int Dat[t]; 
}

parameters {
    vector <lower = -5, upper = 5>[t] mu;
    real<lower=0> sigma;
}

transformed parameters {
    vector<lower=0, upper=1>[t] theta;
    theta = inv_logit(mu);
}

model {
    sigma ~ normal(0, 1);
    for (i in 2:t){
        mu[i] ~ normal(mu[i-1], sigma);
    }
    Dat ~ bernoulli(theta);
}

Rコード

#data to be passed on stan
t=length(Dat)
datastan=list(t=t,Dat=Dat)

#
# stan model
#
#model<-stan_model("GaussianLogitDynamic.stan")
#saveRDS(model,"GaussianLogitDynamic.rds")
GLD<-readRDS("GaussianLogitDynamic.rds")

#
# stan fit
#
fit = sampling(GLD, data = datastan, chains = 4, iter = 4000, warmup = 2000, thin=4, seed=1234)

#Gelman-Rubin 
launch_shinystan(fit)

URL

https://mrunadon.github.io/GaussianLogitDynamic/

自分の書いたコード

pre_model.stan

data{
    int t;
    int Dat[t];
}

parameters{
    real v1;
    real v2;
    real<lower=0> sigma;
}

transformed parameters{
    vector<lower=0,upper=1>[t]  x;
    vector<lower=0,upper=1> [t]  y;
    vector<lower=0,upper=1>[t] theta;
    x[1]=0;
    for(i in 2:t){
      if(Dat[i]==1){
        x[i]=x[i-1]+v2;
      }  
      if(Dat[i]==0){
          x[i]=x[i-1]+v1;
      }
    }
    
    y[1]=0.5;
    theta[t]=(1-y[t])*x[t];
}

model{
  v1 ~ normal(0,1);
  v2 ~ normal(0,1);
  sigma ~ normal(0,1);
  for(i in 2:t){
      y[i] ~ normal(y[i-1],sigma);
  }
  Dat ~ bernoulli(theta);
}

Rコード
データは事前に指定しています。

#data to be passed on stan
t=length(Dat)
datastan=list(t=t,Dat=Dat)

#
# stan model
#
#model<-stan_model("pre_model.stan")
#saveRDS(model,"pre_model.rds")
GLD<-readRDS("pre_model.rds")

#
# stan fit
#
fit = sampling(GLD, data = datastan, chains = 4, iter = 4000, warmup = 2000, thin=4, seed=1234)

#Gelman-Rubin 
launch_shinystan(fit)

Comments

[mather]

参考にしたサンプルコードのURL、自分の書いたコード、をすべて質問内容に提示しないと回答できません。

[kzzzzz]

追加で記載いたしました。
お力を貸していただけると幸いです。

[mather]

そもそも、pre_model.stan の読み込みで警告が発生します。
stanファイルのコンパイルがコメントアウトされているようですが、ちゃんとこのstanファイルで動作確認しましたか？

[kzzzzz]

お返事ありがとうございます。
動作確認、未実施でした。
```
DIAGNOSTIC(S) FROM PARSER:
Info (non-fatal):
Left-hand side of sampling statement (~) may contain a non-linear transform of a parameter or local variable.
If it does, you need to include a target += statement with the log absolute determinant of the Jacobian of the transform.
Left-hand-side of sampling statement:
y[i] ~ normal(...)

recompiling to avoid crashing R session
```
というエラーメッセージが表示されるのですが、これはそもそもモデルが成り立っていないという理解で正しいでしょうか。
その成り立っていない点というのもエラーメッセージから読み取れるようでしたら、ご教示いただけると幸いです。

[mather]

transformed parametersはparametersから導出される値のことなので、 y[i] に関してサンプリングするのは不自然です。
v1, v2, sigma から y が導出されるように指定すべきでしょう。

[kzzzzz]

回答ありがとうございます。
基本の理解ができておらず、恥ずかしい限りです。
いただいた回答の理解としては、
推定したい値←（確率的）←transformed parameters←（確定的）←parameters
という関係性が成り立っていないため、モデルとして成り立っていない、という内容で間違いないでしょうか？

また、今回の目的としては、xおよびyを推定することなのですが、
どのようなモデルを記述すれば、良いと思われるでしょうか。
以下、想定してる関係性です。
二値データ←発生確率←（確定的）←xおよびy「推定したい値」←（確率的）←v1,v2,sigma

何かアイデア等ございましたら、お力添えいただけると幸いです。

[mather]

その理解でいいと思います。
「今回の目的としては、xおよびyを推定すること」なのであれば、x, y をパラメータにするべきだと思います。
アイデアと言うか、そもそも「現象として現れるデータをどんなモデルの結果値として推定しようとしているか」という説明を書いたほうがいいと思います。それがわからないとアドバイスできません。

[kzzzzz]

回答ありがとうございます。
ご指摘ありがとうございます。

表現したいモデルの要件としては以下の通りです。
・現象として現れるdataは確率θのベルヌーイ分布から現れる
・θ[t]はx[t]とy[t]によって確定的に決定される
・x[t]:
data[t-1]==1ならば x[t]=x[t-1]+v1
data[t-1]==0ならば x[t]=x[t-1]+v2
・y[t]は平均値y[t-1]、分散σの正規分布で推移
・v1,v2,σはそれぞれ平均値0、分散1の事前分布を持つ

[kzzzzz]

いただいたアドバイスをもとにやや修正を加え動作確認を行ったところ、以下のようなメッセージが表示されます。
：動作確認後のメッセージ
hash mismatch so recompiling; make sure Stan code ends with a blank line
：モデルのコード
data{
int t;
int Dat[t];
}

parameters{
real v1;
real v2;
real<lower=0> sigma;
vector<lower=0,upper=1> [t] y;

}

transformed parameters{
vector<lower=0,upper=1>[t] x;
vector<lower=0,upper=1>[t] theta;
for(i in 2:t){
if(Dat[i-1]==1){
x[i]=x[i-1]+v2;
}
if(Dat[i-1]==0){
x[i]=x[i-1]+v1;
}
}
theta[t]=(1-y[t])*x[t];
}

model{
v1 ~ normal(0,1);
v2 ~ normal(0,1);
sigma ~ normal(0,1);
for(i in 2:t){
y[i] ~ normal(y[i-1],sigma);
}
Dat ~ bernoulli(theta);
}

ーーーーーーーーーーーーーーーーーーーーーー
最終行は空白で終えているのですが、どのようなエラーなのでしょうか。

[mather]

機械翻訳でもいいので、内容をまず読んでみてください。
> hash mismatch so recompiling; make sure Stan code ends with a blank line
「ハッシュが一致しないためリコンパイルします。Stanコードは空行で終わっていることを確認してください。」
これはただのログです。コンパイルエラーではありません。

[kzzzzz]

回答ありがとうございます。
先ほどのモデルコードで推定を実施したところ、
>Stan model 'model' does not contain samples.
というエラーメッセージが出ます。
コンソールを遡ると
>Chain 4: Rejecting initial value:
>Chain 4: Error evaluating the log probability at the initial value.
>Chain 4: Exception: modelea826ff04640_model_namespace::write_array: x[1] is nan, but must be greater than or equal to 0 (in 'modelea826ff04640_model' at line 15)
というメッセージが表示されています。

xの初期値が見つからないということだと考え、
inits=list(x=0,y=0.5)
fit = sampling(GLD, data = datastan, chains = 4, iter = 1000, warmup = 100, thin=2, seed=1234,init=inits)
として実行したところ、Rが落ちてしまいます。
何が原因でしょうか。

Answer 1

stanを触ったことないので的外れかもしれませんが。

transformed parameters{
    vector<lower=0,upper=1>[t]  x;
    vector<lower=0,upper=1> [t]  y;

ここのVector<lower=0,upper=1>[t]って、>[ の間にスペース入れても動くの？