Rcpp 中的 na.locf 和 inverse.rle答案

【问题标题】：na.locf and inverse.rle in RcppRcpp 中的 na.locf 和 inverse.rle
【发布时间】：2014-06-02 22:05:51
【问题描述】：

我想检查na.locf（来自zoo 包）、rle 和inverse.rle 在RCpp 中是否存在任何预先存在的技巧？

我写了一个循环来实现，例如我对na.locf(x, na.rm=FALSE, fromLast=FALSE)的实现如下：

#include <Rcpp.h>
using namespace Rcpp;

//[[Rcpp::export]]
NumericVector naLocf(NumericVector x) {
  int n=x.size();
  for (int i=1;i<n;i++) {
    if (R_IsNA(x[i]) & !R_IsNA(x[i-1])) {
      x[i]=x[i-1];
    }
  }
  return x;
}

我只是想知道，由于这些是非常基本的功能，有人可能已经在RCpp 中以更好的方式（可能是避免循环）或更快的方式实现了它们？

【问题讨论】：

在 C++ 代码中循环会很快。请注意，当有人在 R 中编写矢量化代码时，C/C++/FORTRAN 代码中仍然存在您看不到的循环。

标签： r rcpp

【解决方案1】：

我要说的唯一一件事是，您只需要为每个值测试两次NA，而您只需要执行一次。测试NA 不是免费操作。也许是这样的：

//[[Rcpp::export]]
NumericVector naLocf(NumericVector x) {
    int n = x.size() ;
    double v = x[0]
    for( int i=1; i<n; i++){
        if( NumericVector::is_na(x[i]) ) {
            x[i] = v ;
        } else {
            v = x[i] ;    
        }
    }

    return x;
}

然而，这仍然会做一些不必要的事情，比如每次我们只能在最后一次没有看到NA 时设置v。我们可以试试这样：

//[[Rcpp::export]]
NumericVector naLocf3(NumericVector x) {
    double *p=x.begin(), *end = x.end() ;
    double v = *p ; p++ ;

    while( p < end ){
        while( p<end && !NumericVector::is_na(*p) ) p++ ;
        v = *(p-1) ;
        while( p<end && NumericVector::is_na(*p) ) {
            *p = v ;
            p++ ;
        }
    }

    return x;
}

现在，我们可以尝试一些基准测试：

x <- rnorm(1e6)
x[sample(1:1e6, 1000)] <- NA 
require(microbenchmark)
microbenchmark( naLocf1(x), naLocf2(x), naLocf3(x) )
#  Unit: milliseconds
#       expr      min       lq   median       uq      max neval
# naLocf1(x) 6.296135 6.323142 6.339132 6.354798 6.749864   100
# naLocf2(x) 4.097829 4.123418 4.139589 4.151527 4.266292   100
# naLocf3(x) 3.467858 3.486582 3.507802 3.521673 3.569041   100

【讨论】：