Rcpp功能比相同的R功能低(Rcpp function to be SLOWER than same R function)

 我一直在编写一个R函数来计算关于某些分布的积分，参见下面的代码。 
 
EVofPsi = function(psi, probabilityMeasure, eps=0.01, ...){

distFun = function(u){
 probabilityMeasure(u, ...)
}
xx = yy = seq(0,1,length=1/eps+1)
summand=0

for(i in 1:(length(xx)-1)){
  for(j in 1:(length(yy)-1)){
    signPlus = distFun(c(xx[i+1],yy[j+1]))+distFun(c(xx[i],yy[j]))
    signMinus = distFun(c(xx[i+1],yy[j]))+distFun(c(xx[i],yy[j+1]))
    summand = c(summand, psi(c(xx[i],yy[j]))*(signPlus-signMinus))
  }
}
sum(summand)
}
 
 它工作正常，但它很慢。 通常听说用C ++等编译语言对函数进行重新编程会加快速度，特别是因为上面的R代码涉及双循环。 我也是这样，使用Rcpp： 
 
#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
double EVofPsiCPP(Function distFun, Function psi, int n, double eps) {

  NumericVector xx(n+1);
  NumericVector yy(n+1);
  xx[0] = 0;
  yy[0] = 0;

  // discretize [0,1]^2
  for(int i = 1; i < n+1; i++) {
      xx[i] = xx[i-1] + eps;
      yy[i] = yy[i-1] + eps;
  }

  Function psiCPP(psi);
  Function distFunCPP(distFun);
  double signPlus;
  double signMinus;
  double summand = 0;

  NumericVector topRight(2); 
  NumericVector bottomLeft(2);
  NumericVector bottomRight(2);
  NumericVector topLeft(2);

  // compute the integral
  for(int i=0; i<n; i++){
    //printf("i:%d \n",i);
    for(int j=0; j<n; j++){
      //printf("j:%d \n",j);
      topRight[0] = xx[i+1];
      topRight[1] = yy[j+1];
      bottomLeft[0] = xx[i];
      bottomLeft[1] = yy[j];
      bottomRight[0] = xx[i+1];
      bottomRight[1] = yy[j];
      topLeft[0] = xx[i];
      topLeft[1] = yy[j+1];
      signPlus = NumericVector(distFunCPP(topRight))[0] +  NumericVector(distFunCPP(bottomLeft))[0];
      signMinus = NumericVector(distFunCPP(bottomRight))[0] + NumericVector(distFunCPP(topLeft))[0];
      summand = summand + NumericVector(psiCPP(bottomLeft))[0]*(signPlus-signMinus);
      //printf("summand:%f \n",summand);
    }
  }
  return summand;
}
 
 我很开心，因为这个C ++函数可以正常工作。 但是，当我测试了这两个函数时，C ++运行了SLOWER： 
 
sourceCpp("EVofPsiCPP.cpp")
pFGM = function(u,theta){
  u[1]*u[2] + theta*u[1]*u[2]*(1-u[1])*(1-u[2])
}
psi = function(u){
  u[1]*u[2]
}
print(system.time(
for(i in 1:10){
  test = EVofPsi(psi, pFGM, 1/100, 0.2)  
}
))
test

print(system.time(
  for(i in 1:10){
    test = EVofPsiCPP(psi, function(u){pFGM(u,0.2)}, 100, 1/100)  
  }
))
 
 那么，是否有一位善意的专家愿意向我解释这一点？ 我是否像猴子一样编码，并且有办法加快这一功能？ 此外，我还有第二个问题。 事实上，我可以用SEXP替换输出类型double，并且参数类型Function by SEXP也是如此，它似乎没有改变任何东西。 那么区别是什么呢？ 
 
 非常感谢你，吉尔达斯 

I have been coding a R function to compute an integral with respect to certain distributions, see code below.
 
EVofPsi = function(psi, probabilityMeasure, eps=0.01, ...){

distFun = function(u){
 probabilityMeasure(u, ...)
}
xx = yy = seq(0,1,length=1/eps+1)
summand=0

for(i in 1:(length(xx)-1)){
  for(j in 1:(length(yy)-1)){
    signPlus = distFun(c(xx[i+1],yy[j+1]))+distFun(c(xx[i],yy[j]))
    signMinus = distFun(c(xx[i+1],yy[j]))+distFun(c(xx[i],yy[j+1]))
    summand = c(summand, psi(c(xx[i],yy[j]))*(signPlus-signMinus))
  }
}
sum(summand)
}
 
It works fine, but it is pretty slow. It is common to hear that re-programming the function in a compiled language such as C++ would speed it up, especially because the R code above involves a double loop. So did I, using Rcpp:
 
#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
double EVofPsiCPP(Function distFun, Function psi, int n, double eps) {

  NumericVector xx(n+1);
  NumericVector yy(n+1);
  xx[0] = 0;
  yy[0] = 0;

  // discretize [0,1]^2
  for(int i = 1; i < n+1; i++) {
      xx[i] = xx[i-1] + eps;
      yy[i] = yy[i-1] + eps;
  }

  Function psiCPP(psi);
  Function distFunCPP(distFun);
  double signPlus;
  double signMinus;
  double summand = 0;

  NumericVector topRight(2); 
  NumericVector bottomLeft(2);
  NumericVector bottomRight(2);
  NumericVector topLeft(2);

  // compute the integral
  for(int i=0; i<n; i++){
    //printf("i:%d \n",i);
    for(int j=0; j<n; j++){
      //printf("j:%d \n",j);
      topRight[0] = xx[i+1];
      topRight[1] = yy[j+1];
      bottomLeft[0] = xx[i];
      bottomLeft[1] = yy[j];
      bottomRight[0] = xx[i+1];
      bottomRight[1] = yy[j];
      topLeft[0] = xx[i];
      topLeft[1] = yy[j+1];
      signPlus = NumericVector(distFunCPP(topRight))[0] +  NumericVector(distFunCPP(bottomLeft))[0];
      signMinus = NumericVector(distFunCPP(bottomRight))[0] + NumericVector(distFunCPP(topLeft))[0];
      summand = summand + NumericVector(psiCPP(bottomLeft))[0]*(signPlus-signMinus);
      //printf("summand:%f \n",summand);
    }
  }
  return summand;
}
 
I'm pretty happy since this C++ function works fine. However, when I tested both functions, the C++ one ran SLOWER:
 
sourceCpp("EVofPsiCPP.cpp")
pFGM = function(u,theta){
  u[1]*u[2] + theta*u[1]*u[2]*(1-u[1])*(1-u[2])
}
psi = function(u){
  u[1]*u[2]
}
print(system.time(
for(i in 1:10){
  test = EVofPsi(psi, pFGM, 1/100, 0.2)  
}
))
test

print(system.time(
  for(i in 1:10){
    test = EVofPsiCPP(psi, function(u){pFGM(u,0.2)}, 100, 1/100)  
  }
))
 
So, is there some kind expert around willing to explain me this? Did I code like a monkey and is there a way to speed up that function? Moreover, I would have a second question. Indeed, I could have replaced the output type double by SEXP, and the argument types Function by SEXP as well, it doesn't seem to change anything. So what is the difference?
 
Thank you very much in advance, Gildas 

原文：https://stackoverflow.com/questions/17958168