R -RCPP-如何在C 代码中使用分布功能

r - Rcpp - how to use a distribution function in the C++ code

本文关键字:分布 功能 代码 -RCPP-      更新时间:2023-10-16

我正在编写N(0,1)分布的Metropolis-Hastings算法:

#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
NumericVector metropolis(int R, double b, double x0){
    NumericVector y(R);
    y(1) = x0;
    for(int i=1; i<R; ++i){
       y(i) = y(i-1);
       double xs = y(i)+runif(1, -b,b)[0];
       double a = dnorm(xs, 0, 1)[0]/dnorm(y(i), 0, 1)[0];
       NumericVector A(2);
       A(0) = 1;
       A(1) = a;
       double random = runif(1)[0];
       if(random <= min(A)){
            y[i] = xs;
       }
   }
   return y;
}

但是,每当我尝试编译函数时,都会发生此错误:

第12行:无匹配函数呼叫'dnorm4'

我试图使用DNORS编写一个微不足道的函数,例如

NumericVector den(NumericVector y, double a, double b){
    NumericVector x = dnorm(y,a,b);
    return x;
}

它有效。有人知道为什么我在大都市代码中有这种错误吗?在r?

中,还有其他一些使用密度函数的方法

在RCPP中,有两组采样器-Scalar和Vector-由名称空间R::Rcpp::分开。他们被分割了:

  • 标量返回一个值(例如doubleint
  • 向量返回多个值(例如NumericVectorIntegerVector

在这种情况下,您想使用标量采样空间,而不是矢量采样空间。

这很明显,因为:

double a = dnorm(xs, 0, 1)[0]/dnorm(y(i), 0, 1)[0];

调用子集运算符[0]以获取向量中的唯一元素。


此问题的第二部分是

暗示的第四参数的缺失部分

第12行:无匹配函数呼叫'dnorm4'

如果您查看dnorm函数的R定义,则会看到:

dnorm(x, mean = 0, sd = 1, log = FALSE)

在这种情况下,除了第四个参数外,您都提供了所有内容。


结果,您可能需要以下内容:

// [[Rcpp::export]]
NumericVector metropolis(int R, double b, double x0){
    NumericVector y(R);
    y(1) = x0; // C++ indices start at 0, you should change this!
    for(int i=1; i<R; ++i){ // C++ indices start at 0!!!
        y(i) = y(i-1);
        double xs = y(i) + R::runif(-b, b);
        double a = R::dnorm(xs, 0.0, 1.0, false)/R::dnorm(y(i), 0.0, 1.0, false);
        NumericVector A(2);
        A(0) = 1;
        A(1) = a;
        double random = R::runif(0.0, 1.0);
        if(random <= min(A)){
            y[i] = xs;
        }
    }
    return y;
}

旁注:

c 索引从0 不是 1开始,因此,您上面的向量是稍有问题的,因为您通过在y(1)填充y矢量而不是y(0)。我将其作为练习供您纠正。