0

我最近开始使用 Rcpp 包将我的一些 R 代码段写入 C++。

给定一个数据矩阵，我有以下 Rcpp 函数，它计算一些观察的协方差的内核重新加权估计。

cppFunction('
        NumericVector get_cov_1obs(NumericMatrix cdata, int ID, float radius){

        int nrow = cdata.nrow(), ncol = cdata.ncol();
        float norm_ = 0;
        float w;
        NumericMatrix out(ncol, ncol);

        NumericMatrix outer_prod(ncol, ncol);

        for (int i=0; i<ncol;i++){
        for (int j=0;j<ncol;j++){
        out(i,j) = 0;
        outer_prod(i,j) = 0;
        }
        }

        for (int i=0; i<nrow;i++){
        w =  exp( -(i-ID)*(i-ID)/(2*radius));
        norm_ += w;
        for (int j=0; j<ncol;j++){
        for (int k=0;k<ncol;k++){
        outer_prod(j,k) = cdata(i,j) * cdata(i,k);
        }
        }

        for (int j=0; j<ncol;j++){
        for (int k=0;k<ncol;k++){
        out(j,k) += outer_prod(j,k)*w;
        }
        }
        }

        for (int i=0; i<ncol;i++){
        for (int j=0;j<ncol;j++){
        out(i,j) /= norm_;
        }
        }

        return out;
        }')

我想快速估计数据集中所有观察值的内核重新加权协方差矩阵，并将它们存储为数组。由于 Rcpp 不处理数组，我编写了以下 R 函数：

get_kern_cov_C = function(data, radius){
  # data is data for which we wish to estimate covariances
  # radius is the radius of the gaussian kernel

  # calculate covariances:
  kern_cov = array(0, c(ncol(data),ncol(data),nrow(data))) 
  for (i in 1:nrow(data)){
    kern_cov[,,i] = get_cov_1obs(cdata=data, ID = i-1, radius=radius)
  }
  return(kern_cov)
}

这似乎工作正常（而且比 R 快得多）但是问题是我时不时地（似乎是随机的）我收到以下形式的错误：

Error in kern_cov[, , i] = get_cov_1obs(cdata = data, ID = i - 1, radius = radius) : 
  incompatible types (from X to Y)

其中X是内置的或 NULL 并且Y是双精度的。

我大致了解为什么会发生这种情况（我试图将内置/NULL 变量放入双精度），但我不确定代码中是否存在错误。我怀疑这可能与内存管理有关，因为它只会时不时地发生。

4

1 回答 1

0

您也可以在 C(++) 级别测试 NULL，在这种情况下可能应该这样做。

至于它发生的原因：恐怕你需要调试它。

于 2013-05-27T14:08:20.070 回答