`integrate_1d` 语法并以交互方式调试它

Question

我是第一次使用 rstan/Stan 的用户，尝试设置自定义可能性并使用 integrate_1d 集成它。基于documentation，我设置了我的被积函数如下：

// internal function to integrate over z 
// integrate_1d is super picky: the argument types need to be real, real, 
    real[], real[], int[]
// e.g., cannot use vectors in place of real[] or int[]

real integrand( real x,  // the observed Z-stat; value at which to evaluate integral
              real xc,
              real[] theta,
              real[] x_r,
              int[] x_i ){

    // separate the parameters
    real zeta = theta[1];
    real eta = theta[2];
    real vi = theta[3];

    // significance indicator
    int signif; 
    signif = fabs(x) > 1.96;

    return( (1/eta) * (1 - signif) + signif ) * exp( normal_lpdf( x | zeta, sqrt(vi) ) );
  }

我正在尝试测试和调试此功能，因为它最终会导致在我运行模型时评估对数似然性时出错。为此，我希望能够以交互方式运行该函数，传递我自己选择的参数。我正在尝试使用 rstan 的 expose_stan_functions 来执行此操作：

code = "functions{

  real integrand( real x,  // the observed Z-stat
                  real xc,
                  real[] theta,
                  real[] x_r,
                  int[] x_i ){

    // separate the parameters
    real zeta = theta[1];
    real eta = theta[2];
    real vi = theta[3];

    // significance indicator
    int signif; 
    signif = fabs(x) > 1.96;

    return( (1/eta) * (1 - signif) + signif ) * exp( normal_lpdf( x | zeta, sqrt(vi) ) );
  }
}"

expose_stan_functions(stanc(model_code = code))

integrand(1, {1,1,1}, real[], int[], 1e-8)

integrate_1d( integrand, negative_infinity(), positive_infinity(), {1, 1, 1}, real[], int[], 1e-8 )

但是对integrand 和integrate_1d 的调用都会引发错误

错误：意外的'，'

就好像我用根本不正确的语法调用函数一样。 （在对integrate_1d 的调用中，我不太确定如何处理xc。我省略了它，因为文档（第9.3.2.4 节）似乎建议它应该是NaN 以获得不定积分。）

如果expose_stan_functions 不是理想的方式，我也欢迎关于交互式调试integrand 的完全不同的建议。

Answer 1

integrate_1d 函数是一个 Stan / C++ 函数，所以如果你想从 R 调用它，你需要指定另一个调用它的函数。在你的情况下，你可以这样做

code = "functions{

  real integrand( real x,  // the observed Z-stat
                  real xc,
                  real[] theta,
                  real[] x_r,
                  int[] x_i ){

    // separate the parameters
    real zeta = theta[1];
    real eta = theta[2];
    real vi = theta[3];

    // significance indicator
    int signif; 
    signif = fabs(x) > 1.96;

    return( (1/eta) * (1 - signif) + signif ) * exp( normal_lpdf( x | zeta, sqrt(vi) ) );
  }

  real area(real[] theta, data real[] x_r) {
    int x_i[0];
    return integrate_1d(integrand, negative_infinity(), positive_infinity(),
                        theta, x_r, x_i, 1e-8);
  }
}"

library(rstan)
expose_stan_functions(stanc(model_code = code))
area(theta = rep(1, 3), x_r = double()) # yields 1