使用 Rcpp 处理 NA 值

Dealing with NA values using Rcpp

我正在测试我的一段代码,如下所示:

#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
NumericMatrix testOutMat(const int& ncols, const int& nrows, const NumericVector& col_prob){

  //Store row and column positions
  NumericVector col_pos = no_init(nrows);
  NumericVector row_pos = no_init(nrows);

  int row_val;
  int nz_counter=0;
  for(int j=0; j<ncols; ++j){ 
    for(int i=0; i<nrows; ++i){
      row_val = R::rbinom(1,col_prob[j]);
      Rcout << "i,j: " << i << "," << j << std::endl;
      Rcout << "val: " << row_val << std::endl;
      if(row_val==1){ //if (i,j)th entry is a 1, save location
        row_pos[i] = i;
        col_pos[i] = j;
        nz_counter += 1;
      } else{ //assign as NA
        row_pos[i] = NA_REAL;
        col_pos[i] = NA_REAL;
      }
      Rcout << "row_pos[i]: " << row_pos[i] << std::endl;
      Rcout << "col_pos[i]: " << col_pos[i] << std::endl;
      Rcout << "num non-zeros: " << nz_counter << std::endl;
    }
  }

  NumericMatrix out = no_init(nz_counter,2);

  Rcout << "Printing output matrix" << std::endl;
  for(int i=0; i<nz_counter; ++i){
    if(!Rcpp::NumericVector::is_na(row_pos[i])){ 
      out(i,0) = row_pos[i];
      out(i,1) = col_pos[i];
    }
    Rcout << "row_pos[i]: " << row_pos[i] << std::endl;
    Rcout << "col_pos[i]: " << col_pos[i] << std::endl; 
  }

  return out;
}

/*** R
set.seed(1)
res <- testOutMat(ncols=5,nrows=5,col_prob = runif(20, 0.1, 0.2))
*/

从输出中,我知道条目 (i,j)={(0,0),(3,1)} 是非零的,因此 res 应该是第一行有 0 02x2 矩阵第二个是 3 1。但是,我得到了一些非常不同的东西:

     [,1] [,2]
[1,]   64 1024
[2,]    1    4

我怀疑这是由于我处理 NA 的方式所致。该函数的总体目标是为非零元素生成行和列索引(通过调用 rbinom 生成)。

我已经尝试调试了一段时间,但似乎无法修复。

这里的问题是你一遍又一遍地写 row_poscol_posncols 次)而没有任何形式的跟踪先前的结果。那,加上您的 no_init() 使用,是导致您看到的最终结果的原因。我们可以稍微更改您的代码以确保 row_poscol_pos 不会被覆盖:

#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
IntegerMatrix testOutMat(const int ncols, const int nrows,
                         const NumericVector& col_prob) {

    IntegerMatrix binomial_deviates(nrows, ncols);
    IntegerVector row_positions;
    IntegerVector col_positions;
    int nz_counter = 0;

    for ( int j = 0; j < ncols; ++j ) {
        binomial_deviates(_, j) = rbinom(nrows, 1, col_prob[j]);
        for ( int i = 0; i < nrows; ++i ) {
            if ( binomial_deviates(i, j) == 1 ) {
                row_positions.push_back(i);
                col_positions.push_back(j);
                nz_counter += 1;
            }
        }
    }

    IntegerMatrix out(nz_counter, 2);

    for ( int i = 0; i < nz_counter; ++i ) {
        out(i, 0) = row_positions[i];
        out(i, 1) = col_positions[i];
    }

    return out;
}

/*** R
set.seed(1)
res <- testOutMat(ncols=5,nrows=5,col_prob = runif(20, 0.1, 0.2))
*/

结果:

> set.seed(1)

> res <- testOutMat(ncols=5,nrows=5,col_prob = runif(20, 0.1, 0.2))
> res
     [,1] [,2]
[1,]    0    0
[2,]    3    1