rcpp - Rcpp 在不同类型的向量列表上嵌套循环：高效的模板解决方案？

Question

另一个有趣的 Rcpp 模板编程问题：我想计算 R data.frame 或矩阵中所有列的成对完整观测值。对于矩阵，由于RCPP_RETURN宏，代码很容易：

// [[Rcpp::plugins(cpp11)]]
#include <Rcpp.h>
using namespace Rcpp;

template <int RTYPE>
IntegerMatrix pwNobsmCppImpl(const Matrix<RTYPE>& x) {
  int l = x.nrow(), col = x.ncol();
  auto isnnanT = (RTYPE == REALSXP) ? [](typename Rcpp::traits::storage_type<RTYPE>::type x) { return x == x; } : 
       [](typename Rcpp::traits::storage_type<RTYPE>::type x) { return x != Vector<RTYPE>::get_na(); };
  IntegerMatrix out = no_init_matrix(col, col);
  for(int j = 0; j != col; ++j) {
    ConstMatrixColumn<RTYPE> colj = x( _ , j);
    int nj = std::count_if(colj.begin(), colj.end(), isnnanT);
    out(j, j) = nj;
    for(int k = j+1; k != col; ++k) {
      ConstMatrixColumn<RTYPE> colk = x( _ , k);
      int njk = 0;
      for(int i = l; i--; ) if(isnnanT(colj[i]) && isnnanT(colk[i])) ++njk; // fastest? or save logical vector with colj Non-missing?
      out(j, k) = out(k, j) = njk;
    }
  }
  out.attr("dimnames") = List::create(colnames(x), colnames(x));
  return out;
}

template <>
IntegerMatrix pwNobsmCppImpl(const Matrix<CPLXSXP>& x) {
  stop("Not supported SEXP type!");
}

template <>
IntegerMatrix pwNobsmCppImpl(const Matrix<VECSXP>& x) {
  stop("Not supported SEXP type!");
}

template <>
IntegerMatrix pwNobsmCppImpl(const Matrix<RAWSXP>& x) {
  stop("Not supported SEXP type!");
}

template <>
IntegerMatrix pwNobsmCppImpl(const Matrix<EXPRSXP>& x) {
  stop("Not supported SEXP type!");
}

// [[Rcpp::export]]
IntegerMatrix pwNobsmCpp(SEXP x){
  RCPP_RETURN_MATRIX(pwNobsmCppImpl, x);
}

现在棘手的是 data.frames / 列表，其列可以是不同的类型。在我看来，原生 Rcpp 解决方案是一个switch(TYPEOF(x[j]), ...)语句，它根据具体情况处理列表中不同类型的列。但是，在这种情况下，这将需要 2 个嵌套开关，从而导致大量代码重复。我想知道是否可以通过一些智能模板编程来避免这种情况。列表的基本代码是：

// [[Rcpp::export]]
IntegerMatrix pwNobslCpp(const List& x) {
  int l = x.size();
  IntegerMatrix out = no_init_matrix(l, l);
  for(int j = 0; j != l; ++j) {
    int RTYPEj = TYPEOF(x[j]); // The code fails because the SEXPTYPE of x[j] is unknown beforehand!
    auto isnnanTj = (RTYPEj == REALSXP) ? [](typename Rcpp::traits::storage_type<RTYPEj>::type x) { return x == x; } : 
      [](typename Rcpp::traits::storage_type<RTYPEj>::type x) { return x != Vector<RTYPEj>::get_na(); };
    Vector<RTYPEj> colj = x[j];
    int nj = std::count_if(colj.begin(), colj.end(), isnnanTj);
    int rowj = colj.size();
    out(j, j) = nj;
    for(int k = j+1; k != l; ++k) {
      int RTYPEk = TYPEOF(x[k]); // The code fails because the SEXPTYPE of x[k] is unknown beforehand!
      auto isnnanTk = (RTYPEk == REALSXP) ? [](typename Rcpp::traits::storage_type<RTYPEk>::type x) { return x == x; } : 
        [](typename Rcpp::traits::storage_type<RTYPEk>::type x) { return x != Vector<RTYPEk>::get_na(); };
      Vector<RTYPEk> colk = x[k];
      if(colk.size() != rowj) stop("All columns need to have the same length!");
      int njk = 0;
      for(int i = rowj; i--; ) if(isnnanTj(colj[i]) && isnnanTk(colk[i])) ++njk; // fastest? or save logical vector with 
j Non-missing?
      out(j, k) = out(k, j) = njk;
    }
  }
  out.attr("dimnames") = List::create(names(x), names(x));
  return out;
}

此代码当然会失败，因为其中包含的 Rcpp 模板未使用已知值激活。但是由于计算的成对性质，模板化整个代码也没有意义。没有可怕的嵌套开关，有什么办法可以做到这一点？

标签： rcpp