Rcpp::table handles NAs inconsistently

[kevinushey]

E.g.

#include <Rcpp.h>
using namespace Rcpp;

template <typename T>
inline IntegerVector do_counts(const T& x) {
  return table(x);
}

// [[Rcpp::export]]
IntegerVector counts(SEXP x) {
  switch (TYPEOF(x)) {
  case INTSXP: return do_counts<IntegerVector>(x);
  case REALSXP: return do_counts<NumericVector>(x);
  case STRSXP: return do_counts<CharacterVector>(x);
  case LGLSXP: return do_counts<LogicalVector>(x);
  default: {
    stop("Unrecognized SEXP type");
    return R_NilValue;
  }
  }
}

/*** R
set.seed(123)
x_num <- rnorm(5)
x_num[ sample(1:length(x_num), 2) ] <- NA
print(x_num)
x_int <- as.integer( round(x_num) )
x_char <- as.character(x_int)
x_lgl <- x_num > 0

counts(x_num)
counts(x_int)
counts(x_char)
counts(x_lgl)
*/

gives me

> set.seed(123)

> x_num <- rnorm(5)

> x_num[ sample(1:length(x_num), 2) ] <- NA

> print(x_num)
[1] -0.56047565          NA  1.55870831  0.07050839          NA

> x_int <- as.integer( round(x_num) )

> x_char <- as.character(x_int)

> x_lgl <- x_num > 0

> counts(x_num)
-0.560476  0.070508  1.558708 
        2         2         1 

> counts(x_int)
<NA>   -1    0    2 
   2    1    1    1 

> counts(x_char)
  -1    0    2 <NA> 
   1    1    1    2 

> counts(x_lgl)
 <NA> FALSE  TRUE 
    2     1     2 

Two things:

1. The output for numeric case is wrong,
2. The location of the NA element is different between integer and character. (not such a big problem, but might be worth standardizing)

As far as I can see, table works fine in all the cases where there are no 'special' values (NA, NaN, Inf, -Inf) for numerics.

[kevinushey]

An option: we could just use std::map, and have nicely designed comparators for each case. A useful example (inspired by this SO answer):

There are some deficiencies for very small / very large doubles (no use of exponential notation, doesn't work for small numbers). EDIT: fixed by using R's 'coerce to character' to get names as R would construct them from a double.

#include <Rcpp.h>
using namespace Rcpp;

template <typename T>
struct NACompare;

template <>
struct NACompare<int> {
  inline bool operator()(int left, int right) const {
    if (left == NA_INTEGER) return false;
    if (right == NA_INTEGER) return true;
    return left < right;
  }
};

template <>
struct NACompare<double> {  
  inline bool operator()(double left, double right) const {
    bool leftNaN = (left != left);
    bool rightNaN = (right != right);
    if (leftNaN != rightNaN) {
      return leftNaN < rightNaN;
    } else {
      return left < right;
    }
  }
};

template <>
struct NACompare<String> {
  inline bool operator()(const String& left, const String& right) const {
    return left < right;
  }
};

template <typename T, typename U>
inline IntegerVector do_counts(const T& x) {
  std::map< U, int, NACompare<U> > counts;
  int n = x.size();
  for (int i=0; i < n; ++i) {
    ++counts[ x[i] ];
  }
  return wrap(counts);
}

template <>
inline IntegerVector do_counts<NumericVector, double>(const NumericVector& x) {
  std::map< double, int, NACompare<double> > counts;
  int n = x.size();
  for (int i=0; i < n; ++i) {
    ++counts[ x[i] ];
  }
  IntegerVector output = wrap(counts);

  // explicitly use R's coercion mechanism to get good names
  int m = counts.size();
  NumericVector keys = no_init(m);
  typedef std::map< double, int, NACompare<double> >::iterator MapItr;
  int i = 0;
  for (MapItr it = counts.begin(); it != counts.end(); ++it) {
    keys[i] = it->first;
    ++i;
  }
  CharacterVector names = Rf_coerceVector(keys, STRSXP);
  output.attr("names") = names;
  return output;
}

// [[Rcpp::export]]
IntegerVector counts(SEXP x) {
  switch (TYPEOF(x)) {
  case REALSXP: {
    IntegerVector output = do_counts<NumericVector, double>(x);
    // fix names
    for (int i=0; i < output.size(); ++i) {
      if (CHAR(STRING_ELT(output.attr("names"), i)) == "-0") {
        SET_STRING_ELT(output.attr("names"), i, Rf_mkChar("-0"));
      }
    }
    return output;
  }
  case STRSXP: return do_counts<CharacterVector, String>(x);
  case INTSXP: return do_counts<IntegerVector, int>(x);
  case LGLSXP: {
    IntegerVector output = do_counts<LogicalVector, int>(x);
    SET_STRING_ELT( output.attr("names"), 0, Rf_mkChar("FALSE") );
    SET_STRING_ELT( output.attr("names"), 1, Rf_mkChar("TRUE") );
    return output;
  }
  default: {
    stop("unrecognized SEXP type");
    return R_NilValue;
  }
  }
}

/*** R
set.seed(123)
x <- round( rnorm(1E5), 1 )
x[sample(length(x), 100)] <- NA
x_num <- x
x_int <- as.integer(x)
x_char <- as.character(x)
x_lgl <- x > 0
all.equal( counts(x), c(table(x, useNA='ifany') ) )
all.equal( counts(x_int), c(table(x_int, useNA="ifany")))
all.equal(counts(x_char), c(table(x_char, useNA="ifany")))
all.equal(counts(x_lgl), c(table(x_lgl, useNA="ifany")))
if (require(microbenchmark)) microbenchmark( times=5, counts(x), c(table(x, useNA='ifany')) )
*/

[kevinushey]

Ok. I am finishing up a solution and will send a pull request soon.

The crux of the problem: if Rcpp is using std::map we have to use a comparator explicitly or else insertions fail on NA, NaN.