All geek questions in one place

Fast check for missing values ​​in Rcpp

This question is related to the NA values ​​in the Rcpp conditional expression .

Basically, I have Rcpp code that iterates over several (double) elements. And I need to check if there are no values ​​for each element (and I cannot use the vector). Let me count the number of missing values ​​in a vector, as well as the minimum reproducible example:

#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
int nb_na(const NumericVector& x) {
  int n = x.size();
  int c = 0;
  for (int i = 0; i < n; i++) if (R_IsNA(x[i])) c++;
  return c;
}

// [[Rcpp::export]]
int nb_na3(const NumericVector& x) {
  int n = x.size();
  int c = 0;
  for (int i = 0; i < n; i++) if (x[i] == 3) c++;
  return c;
}

// [[Rcpp::export]]
LogicalVector na_real(NumericVector x) {
  return x == NA_REAL;
}

Then in R we get:

> x <- rep(c(1, 2, NA), 1e4)

> x2 <- replace(x, is.na(x), 3)

> microbenchmark::microbenchmark(
+   nb_na(x),
+   nb_na3(x2)
+ )
Unit: microseconds
       expr     min      lq      mean  median       uq      max neval
   nb_na(x) 135.633 135.982 153.08586 139.753 140.3115 1294.928   100
 nb_na3(x2)  22.490  22.908  30.14005  23.188  23.5025  684.026   100

> all.equal(nb_na(x), nb_na3(x2))
[1] TRUE

> na_real(x[1:3])
[1] NA NA NA

As noted in a related question, you cannot just check x[i] == NA_REALbecause it always returns a missing value. However, use is R_IsNA(x[i])much slower than checking for equality with a numerical value (e.g., 3).

, , , . , .

+4
2

NaN , . .

. NumericVector::is_na R_IsNA, .

, .. if (R_IsNA(x[i])) c++; c += NumericVector::is_na(x[i]). :

// [[Rcpp::export]]
int nb_na4(const NumericVector& x) {
  int n = x.size();
  int c = 0;
  for (int i = 0; i < n; i++) c += NumericVector::is_na(x[i]) ;
  return c;
}

int x[i] , std::count_if. raison d'être. :

// [[Rcpp::export]]
int nb_na5(const NumericVector& x) {
  return std::count_if(x.begin(), x.end(), NumericVector::is_na ) ;
}

, - , , tbb RcppParallel.

// [[Rcpp::export]]
int nb_na6(const NumericVector& x) {
  return tbb::parallel_reduce( 
    tbb::blocked_range<const double*>(x.begin(), x.end()),
    0, 
    [](const tbb::blocked_range<const double*>& r, int init) -> int {
      return init + std::count_if( r.begin(), r.end(), NumericVector::is_na );
    }, 
    []( int x, int y){ return x+y; }
  ) ;
}

:

library(microbenchmark)

bench <- function(n){
  x <- rep(c(1, 2, NA), n)
  microbenchmark(
    nb_na = nb_na(x), 
    nb_na4 = nb_na4(x), 
    nb_na5 = nb_na5(x), 
    nb_na6 = nb_na6(x)
  )
}
bench(1e5)

:

> bench(1e4)
Unit: microseconds
expr    min      lq      mean  median       uq     max neval  cld
nb_na  84.358 94.6500 107.41957 110.482 118.9580 137.393   100    d
nb_na4 59.984 69.4925  79.42195  82.442  85.9175 106.567   100  b  
nb_na5 65.047 75.2625  85.17134  87.501  93.0315 116.993   100   c 
nb_na6 39.205 51.0785  59.20582  54.457  68.9625  97.225   100 a   

> bench(1e5)
Unit: microseconds
expr     min       lq     mean   median       uq      max neval  cld
nb_na  730.416 732.2660 829.8440 797.4350 872.3335 1410.467   100    d
nb_na4 520.800 521.6215 598.8783 562.7200 657.1755 1059.991   100  b  
nb_na5 578.527 579.3805 664.8795 626.5530 710.5925 1166.365   100   c 
nb_na6 294.486 345.2050 368.6664 353.6945 372.6205  897.552   100 a

- , long long, 64- NA_REAL:

  > devtools::install_github( "ThinkR-open/seven31" )
  > seven31::reveal(NA, NaN, +Inf, -Inf )
  0 11111111111 ( NaN ) 0000000000000000000000000000000000000000011110100010 : NA
  0 11111111111 ( NaN ) 1000000000000000000000000000000000000000000000000000 : NaN
  0 11111111111 ( NaN ) 0000000000000000000000000000000000000000000000000000 : +Inf
  1 11111111111 ( NaN ) 0000000000000000000000000000000000000000000000000000 : -Inf

:

// [[Rcpp::export]]
int nb_na7( const NumericVector& x){
  const long long* p = reinterpret_cast<const long long*>(x.begin()) ;
  long long na = *reinterpret_cast<long long*>(&NA_REAL) ;

  return std::count(p, p + x.size(), na ) ;

}

:

// [[Rcpp::export]]
int nb_na8( const NumericVector& x){
  const long long* p = reinterpret_cast<const long long*>(x.begin()) ;
  long long na = *reinterpret_cast<long long*>(&NA_REAL) ;

  auto count_chunk = [=](const tbb::blocked_range<const long long*>& r, int init) -> int {
    return init + std::count( r.begin(), r.end(), na);
  } ;

  return tbb::parallel_reduce( 
    tbb::blocked_range<const long long*>(p, p + x.size()),
    0, 
    count_chunk, 
    []( int x, int y){ return x+y; }
  ) ;

}

  > bench(1e5)
  Unit: microseconds
     expr     min       lq     mean   median       uq      max neval    cld
    nb_na 730.346 762.5720 839.9479 857.5865 881.8635 1045.048   100      f
   nb_na4 520.946 521.6850 589.0911 578.2825 653.4950  832.449   100    d  
   nb_na5 578.621 579.3245 640.9772 616.8645 701.8125  890.736   100     e 
   nb_na6 291.115 307.4300 340.1626 344.7955 360.7030  484.261   100   c   
   nb_na7 122.156 123.4990 141.1954 132.6385 149.7895  253.988   100  b    
   nb_na8  69.356  86.9980 109.6427 115.2865 126.2775  182.184   100 a     

  > bench(1e6)
  Unit: microseconds
     expr      min        lq      mean    median        uq      max neval  cld
    nb_na 7342.984 7956.3375 10261.583 9227.7450 10869.605 79757.09   100    d
   nb_na4 5286.970 5721.9150  7659.009 6660.2390  9234.646 31141.47   100   c 
   nb_na5 5840.946 6272.7050  7307.055 6883.2430  8205.117 10420.48   100   c 
   nb_na6 2833.378 2895.7160  3891.745 3049.4160  4054.022 18242.26   100  b  
   nb_na7 1661.421 1791.1085  2708.992 1916.6055  2232.720 60827.63   100 ab  
   nb_na8  650.639  869.6685  1289.373  939.0045  1291.025 10223.29   100 a

NA.

:

#include <Rcpp.h>
#include <RcppParallel.h>

// [[Rcpp::depends(RcppParallel)]]
// [[Rcpp::plugins(cpp11)]]
using namespace Rcpp;

// [[Rcpp::export]]
int nb_na(const NumericVector& x) {
  int n = x.size();
  int c = 0;
  for (int i = 0; i < n; i++) if (R_IsNA(x[i])) c++;
  return c;
}

// [[Rcpp::export]]
int nb_na4(const NumericVector& x) {
  int n = x.size();
  int c = 0;
  for (int i = 0; i < n; i++) c += NumericVector::is_na(x[i]) ;
  return c;
}

// [[Rcpp::export]]
int nb_na5(const NumericVector& x) {
  return std::count_if(x.begin(), x.end(), NumericVector::is_na ) ;
}

// [[Rcpp::export]]
int nb_na6(const NumericVector& x) {
  return tbb::parallel_reduce( 
    tbb::blocked_range<const double*>(x.begin(), x.end()),
    0, 
    [](const tbb::blocked_range<const double*>& r, int init) -> int {
      return init + std::count_if( r.begin(), r.end(), NumericVector::is_na );
    }, 
    []( int x, int y){ return x+y; }
  ) ;
}

// [[Rcpp::export]]
int nb_na7( const NumericVector& x){
  const long long* p = reinterpret_cast<const long long*>(x.begin()) ;
  long long na = *reinterpret_cast<long long*>(&NA_REAL) ;

  return std::count(p, p + x.size(), na ) ;

}

// [[Rcpp::export]]
int nb_na8( const NumericVector& x){
  const long long* p = reinterpret_cast<const long long*>(x.begin()) ;
  long long na = *reinterpret_cast<long long*>(&NA_REAL) ;

  auto count_chunk = [=](const tbb::blocked_range<const long long*>& r, int init) -> int {
    return init + std::count( r.begin(), r.end(), na);
  } ;

  return tbb::parallel_reduce( 
    tbb::blocked_range<const long long*>(p, p + x.size()),
    0, 
    count_chunk, 
    []( int x, int y){ return x+y; }
  ) ;

}

/*** R
library(microbenchmark)

bench <- function(n){
  x <- rep(c(1, 2, NA), n)
  microbenchmark(
    nb_na = nb_na(x), 
    nb_na4 = nb_na4(x), 
    nb_na5 = nb_na5(x), 
    nb_na6 = nb_na6(x), 
    nb_na7 = nb_na7(x), 
    nb_na8 = nb_na8(x)
  )
}
bench(1e5)
bench(1e6)
*/
+6

(IEEE) , . R.

, ALTREP R - , , , / - , . , ALTREP, https://github.com/HenrikBengtsson/Wishlist-for-R/issues/12

+4

Source: https://habr.com/ru/post/1688020/

More articles:


All Articles