【R Cpp Error系列】 Rcpp::sugar::Vectorized＜exp, true, Rcpp::sugar::Times_Vector_Vector＜14, true, Rcpp::

目录

一、代码遇到的错误

​编辑

 二、运行正确的类似代码

三、对照修改自己的代码

---

一、代码遇到的错误

运行下面代码：

 sourceCpp(code='
 #include <Rcpp.h>
 using namespace Rcpp;
 
 // [[Rcpp::export]]
 NumericMatrix test(NumericMatrix z, NumericVector b_hat){
 double exp_bz;
 exp_bz = exp(z.row(1)*b_hat);
    
 return exp_bz;
 }
 
')

得到如下错误：

exp()函数是一个向量化函数，也就是说，exp()中可以接收一个元素，也可以接受一个向量。R是向量化运算的，在Rcpp中也引进了向量化运算。这个概念是【Rcpp糖】。

错误原因：

上述代码错误的原因是什么呢？exp()是一个向量化函数，上述代码exp中接收的参数是一个向量，即z的第一行元素分别乘b_hat的每个元素，（其中z是一个n by p的矩阵，b_hat是一个长度为p的向量），然而exp函数的输入是一个向量，但是输出对象exp_bz被定义为一个double数，是一维的。所以产生了上述错误。

在R中可以通过exp.bz <- as.numeric(exp(z %*% b.hat)) # n vector直接获得一个列向量。

Rcpp糖的介绍具体参考：56 Rcpp糖 | R语言教程 (pku.edu.cn)

在C++中，向量和矩阵的运算通常需要逐个元素进行， 或者调用相应的函数。 Rcpp通过C++的表达式模板(expression template)功能和懒惰求值(lazy evaluation)功能， 可以在C++中写出像R中对向量和矩阵运算那样的向量化表达式。 这称为Rcpp糖(sugar)。

 二、运行正确的类似代码

sourceCpp(code='
#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
NumericVector distC(double x, NumericVector y) {
  int n = y.size();
  NumericVector out(n);
  for(int i = 0; i < n; ++i) {
    out[i] = sqrt(pow(y[i] - x, 2.0));
  }
  return out;
}

')

三、对照修改自己的代码

sourceCpp(code='
 #include <Rcpp.h>
 using namespace Rcpp;
 
 // [[Rcpp::export]]
 NumericVector test(NumericMatrix z, NumericVector b_hat){
 int n = z.nrow();
 NumericVector exp_bz(n);
 for(int i=0; i<n; ++i){
 exp_bz(i) = exp(z(i,1)*b_hat(1)+z(i,2)*b_hat(2));
 }
    return exp_bz;
 }
 
')

参考：

1. Advanced R and R Packages Workshop - C++ and Rcpp
2. 56 Rcpp糖 | R语言教程 (pku.edu.cn)