1

我试图在 Stan 中运行这个模型。我有一个正在运行的 JAGS 版本(返回高度自相关的参数),并且我知道如何将其公式化为双指数的 CDF(具有两个速率),这可能会毫无问题地运行。但是,我想将此版本用作类似但更复杂模型的起点。

到目前为止,我怀疑这样的模型在 Stan 中是不可能的。可能由于采用布尔值之和引入的离散性,Stan 可能无法计算梯度。

有谁知道是这种情况,还是我在这个模型中以错误的方式做了其他事情?我将得到的错误粘贴到模型代码下方。

非常感谢提前一月

Model:

data {
    int y[11]; 
    int reps[11];
    real soas[11]; 

}
parameters {
    real<lower=0.001,upper=0.200> v1;
    real<lower=0.001,upper=0.200> v2;

}


model {
    real dif[11,96];
    real cf[11];

    real p[11];

    real t1[11,96];
    real t2[11,96];

    for (i in 1:11){
        for (r in 1:reps[i]){     
            t1[i,r]  ~ exponential(v1);
            t2[i,r]  ~ exponential(v2);
            dif[i,r] <-  (t1[i,r]+soas[i]<=(t2[i,r]));

            }
        cf[i] <- sum(dif[i]);
        p[i]  <-cf[i]/reps[i];
        y[i] ~ binomial(reps[i],p[i]); 
    }

}

这是一些虚拟数据:

psy_dat = { 
         'soas' :  numpy.array(range(-100,101,20)),
            'y' :  [47, 46, 62, 50, 59, 47, 36, 13, 7, 2, 1],
         'reps' :  [48, 48, 64, 64, 92, 92, 92, 64, 64, 48, 48]
          }

以下是错误:

DIAGNOSTIC(S) FROM PARSER:
Warning (non-fatal): Left-hand side of sampling statement (~) contains a non-linear transform of a parameter or local variable.
You must call increment_log_prob() with the log absolute determinant of the Jacobian of the transform.
Sampling Statement left-hand-side expression:
get_base1(get_base1(t1,i,"t1",1),r,"t1",2) ~ exponential_log(...)
Warning (non-fatal): Left-hand side of sampling statement (~) contains a non-linear transform of a parameter or local variable.
You must call increment_log_prob() with the log absolute determinant of the Jacobian of the transform.
Sampling Statement left-hand-side expression:
get_base1(get_base1(t2,i,"t2",1),r,"t2",2) ~ exponential_log(...)

在运行时:

Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
 stan::prob::exponential_log(N4stan5agrad3varE): Random variable is nan:0, but must not be nan!
 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
 Rejecting proposed initial value with zero density.


Initialization between (-2, 2) failed after 100 attempts. 
 Try specifying initial values, reducing ranges of constrained values, or   reparameterizing the model

这是此模型的工作 JAGS 版本:

   model {
   for ( n in 1 : N  ) { 
     for (r in 1 : reps[n]){
       t1[r,n] ~ dexp(v1)
       t2[r,n] ~ dexp(v2)
       c[r,n] <- (1.0*((t1[r,n]+durs[n])<=t2[r,n]))
     } 
     p[n] <- max((min(sum(c[,n]) /  (reps[n]),0.99999999999999)),   1-0.99999999999999)) 
     y[n] ~ dbin(p[n],reps[n])
   }

   v1 ~ dunif(0.0001,0.2)
   v2 ~ dunif(0.0001,0.2)
   }

关于 min() 和 max():见这篇文章https://stats.stackexchange.com/questions/130978/observed-node-inconsistent-when-binomial-success-rate-exactly-one?noredirect=1 #comment250046_130978

4

1 回答 1

1

我仍然不确定您要估算的模型是什么(最好发布 JAGS 代码),但是您上面的内容在 Stan 中无法使用。Stan 在您必须声明然后定义对象的意义上更接近 C++。在您的 Stan 程序中,您有两个声明 real t1[11,96]; real t2[11,96]; ,但没有t1or的定义t2。因此,它们被初始化为NaN,当你这样做时 t1[i,r] ~ exponential(v1); ,它会被解析为类似于 for(i in 1:11) for(r in 1:reps[i]) lp__ += log(v1) - v1 * NaN wherelp__是一个内部符号,它保存对数后验的值,它变成 NaN,它不能对参数进行 Metropolis 风格的更新。

也许您的意思是 fort1t2to be unknown 参数,在这种情况下,它们应该在parameters块中声明。以下 [EDITED] Stan 程序是有效的并且应该可以工作,但它可能不是您想到的程序(它对我来说没有多大意义,并且不连续性dif可能会阻止 Stan 有效地从这个后验分布中采样)。 data { int<lower=1> N; int y[N]; int reps[N]; real soas[N]; } parameters { real<lower=0.001,upper=0.200> v1; real<lower=0.001,upper=0.200> v2; real t1[N,max(reps)]; real t2[N,max(reps)]; } model { for (i in 1:N) { real dif[reps[i]]; for (r in 1:reps[i]) { dif[r] <- (t1[i,r]+soas[i]) <= t2[i,r]; } y[i] ~ binomial(reps[i], (1.0 + sum(dif)) / (1.0 + reps[i])); } to_array_1d(t1) ~ exponential(v1); to_array_1d(t2) ~ exponential(v2); }

于 2015-01-16T05:42:56.320 回答