时间戳

首页 > 解决方案 > 具有不同值的 Rcpp 函数填充矩阵

c++ - 具有不同值的 Rcpp 函数填充矩阵

问题描述

我正在构建一个过程，它将实例化一个 NumericMatrix 并用 Sorenson-Dice 相似系数（一个相似矩阵）填充它。矩阵本身具有可变维度，并且取决于正在处理的元素的数量。通常有超过 100 个单独的元素可以随时进行比较（因此矩阵维度通常为 100+ x 100+）。到目前为止我构建的将创建矩阵，计算系数，然后用这些计算值填充矩阵。但是，当我重复运行该函数时，我注意到矩阵中的值在每次运行之间发生变化，这不是预期的行为，因为被比较的数据在每次运行之间都没有变化或重新排序。我也得到大于 1 的相似性，这绝对不应该发生。我有四个功能，

这是C++代码：

// function to calculate the denominator of the dice coefficient
int diceDenomcpp(NumericVector val1, NumericVector val2){
  
  
  int val1Len = na_omit(val1).size();
  int val2Len = na_omit(val2).size();
  int bands = 0;
  
  
  bands = val1Len + val2Len;
  // return the computed total data points within both arrays
  
  
  return bands;
}

//######################################################################
//######################################################################
//######################################################################

// function to calculate the numerator for the dice coefficient
int diceNumcpp(NumericVector iso1, NumericVector iso2){
  
  // declare and initialize vectors with the element band data
  // remove any NA values within each vector
  NumericVector is1 = na_omit(iso1);
  NumericVector is2 = na_omit(iso2);
  
  // declare and initialize some counter variables
  int n = 0;
  int m = 0;
  int match = 0;
  
  // loop through the first element's first datum and check for matching datum
  // with the second element then continue to loop through each datum within each element 
  while (n<=is1.size()){
    if (m>=is2.size()){
      n++;
      m=0;
    }
    // if a suitable match is found, increment the match variable
    if((fabs(is1[n]-is2[m])/is1[n])<0.01 && (fabs(is1[n]-is2[m])/is2[m])<0.01){
      match++;
      
    }
    m++;
  }
  return match;
}

//########################################################################
//########################################################################
//########################################################################

// function to put the coefficient together
double diceCoefcpp(NumericVector val1, NumericVector val2){
  
  NumericVector is1 = clone(val1);
  NumericVector is2 = clone(val2);
  double dVal;
  double num = 2*diceNumcpp(is1, is2);
  double denom = diceDenomcpp(is1, is2);
  
  dVal = num/denom;
  
  return dVal;
  
}

//#######################################################################
//#######################################################################
//#######################################################################


// function to build the similarity matrix with the coefficients

NumericMatrix simMatGencpp(NumericMatrix df){
  
  // clone the input data frame
  NumericMatrix rapdDat = clone(df);

  // create a data frame for the output 
  NumericMatrix simMat(rapdDat.nrow(),rapdDat.nrow());
    std::fill(simMat.begin(), simMat.end(), NumericVector::get_na());
  
  // declare and initialize the iterator
  int i = 0;

  // declare and initialize the column counter
  int col = 0;  
  
  // declare an initialize the isolate counter
  int iso = 0;
  
  //simMat(_,0)=rapdDat(_,0);
  
  while (iso < rapdDat.nrow()){
    if (iso+i > rapdDat.nrow()){
      col++;
      i=0;
      iso++;
    }
    if (iso+i < rapdDat.nrow()){
      simMat(iso+i, col) = diceCoefcpp(rapdDat(iso,_), rapdDat(iso+i,_));
      
    }
    i++;
  }
  
  
  //Rcout << "SimMatrix:" << simMat << "\n";
  
  return simMat;
}

这是输入数据应该是什么样子的示例。. .

sampleData

    band1  band2  band3  band4  band5  band6
1   593.05 578.04 439.01     NA     NA     NA
2   589.07 567.03     NA     NA     NA     NA
3   591.04 575.10 438.12     NA     NA     NA
4   591.04     NA     NA     NA     NA     NA
5   588.08 573.18     NA     NA     NA     NA
6   591.04 576.09 552.10     NA     NA     NA
7  1805.00 949.00 639.19 589.07 576.09 440.06
8   952.00 588.08 574.14 550.04     NA     NA
9  1718.00 576.09 425.01     NA     NA     NA
10 1708.00 577.05 425.01     NA     NA     NA

对于足够小的数据集，输出 simMatGencpp() 函数每次都会产生相同的结果，但是当数据集变大时，值将开始在运行之间发生变化。

我已经尝试在各个元素上独立运行 diceNumcpp()、diceDenomcpp() 和 diceCoefcpp() 函数，并且每次都能获得一致的预期输出。但是，一旦我使用 simMatGencpp() ，输出就会再次变得混乱。所以我尝试循环每个单独的函数，如下所示。

例子：

for(i in 1:100){
  print(diceNumcpp(sampleData[7,], sampleData[3,]))
}

上面的预期输出应该是 3，但有时是 4。每次我运行这个循环时，无论哪个时间 4 输出都会发生变化，有时是第二次迭代，有时是第 14 次，或者根本没有，或者连续三次。

我的第一个想法是，也许由于垃圾收集并没有完全发生在 c++ 中，可能之前运行的函数调用将旧向量留在内存中，因为输出对象的名称在运行之间没有变化。但是这篇文章说，当函数退出时，在函数调用范围内创建的任何对象也会被销毁。

当我仅在 R 代码中编写相同的解决方案时，运行时很糟糕，但它每次都会始终返回具有相同值的矩阵或示例向量。

我不知所措。任何人可以在这个主题上提供的任何帮助或启发都将受到极大的欢迎！

谢谢你的帮助。

2020-08-19 更新

我希望这将有助于为更精通 c++ 的人提供一些见解，以便您对可能发生的事情有一些额外的想法。我有一些示例数据，类似于上面显示的内容，即 187 行长，这意味着这些数据的相似矩阵将有 17578 个元素。我一直在比较这个解决方案的 R 版本和这个解决方案的 c++ 版本，使用这样的代码和示例数据：

# create the similarity matrix with the R-solution to compare iteratively
# with another R-solution similarity matrix
simMat1 <- simMatGen(isoMat)
resultsR <- c()
for(i in 1:100){
  
  simMat2 <- simMatGen(isoMat)

  # check for any mis-matched elements in each matrix
  resultsR[[i]]<-length(which(simMat1 == simMat2)==TRUE)

  #######################################################################
  # everytime this runs I get the expected number of true values 17578
  # and check this by subtracting the mean(resultsR) from the expected 
  # number of true values of 17578 
}

mean(resultsR)

现在，当我对 c++ 版本执行相同的过程时，情况会发生巨大而迅速的变化。我用 64 位和 32 位 R-3.6.0 都试过这个，只是因为。

simMat1 <- simMatGen(isoMat)
isoMat <- as.matrix(isoMat)
resultscpp <- c()
for(i in 1:10000){
  
  simMat2 <- simMatGencpp(isoMat)
  resultscpp[[i]]<-length(which(simMat1 == simMat2)==TRUE)

  ############  64 bit R  ##############
  # first iteration length(which(simMat1 == simMat2)==TRUE)-17578 equals 2
  # second iteration 740 elements differ: length(which(simMat1 == simMat2)==TRUE)-17578 equals 740 
  # third iteration 1142 elements differ
  # after 100 iterations the average difference is 2487.7 elements
  # after 10000 iterations the average difference is 2625.91 elements
  
  ############  32 bit R  ##############
  # first iteration difference = 1
  # second iteration difference = 694
  # 100 iterations difference = 2520.94
  # 10000 iterations difference = 2665.04
}

mean(resultscpp)

这是 sessionInfo()

R version 3.6.0 (2019-04-26)
Platform: i386-w64-mingw32/i386 (32-bit)
Running under: Windows 10 x64 (build 17763)

Matrix products: default

locale:
[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United States.1252    LC_MONETARY=English_United States.1252
[4] LC_NUMERIC=C                           LC_TIME=English_United States.1252    

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

loaded via a namespace (and not attached):
 [1] Rcpp_1.0.5        rstudioapi_0.10   magrittr_1.5      usethis_1.5.0     devtools_2.1.0    pkgload_1.0.2     R6_2.4.0          rlang_0.4.4      
 [9] tools_3.6.0       pkgbuild_1.0.3    sessioninfo_1.1.1 cli_1.1.0         withr_2.1.2       remotes_2.1.0     assertthat_0.2.1  digest_0.6.20    
[17] rprojroot_1.3-2   crayon_1.3.4      processx_3.3.1    callr_3.2.0       fs_1.3.1          ps_1.3.0          testthat_2.3.1    memoise_1.1.0    
[25] glue_1.3.1        compiler_3.6.0    desc_1.2.0        backports_1.1.5   prettyunits_1.0.2

标签： c++rrcpp

解决方案

在这里犯了一个菜鸟c++错误。

在 diceNumcpp() 中，我没有进行任何检查，以免意外引用数组中的越界元素。

// if a suitable match is found, increment the match variable
    if((fabs(is1[n]-is2[m])/is1[n])<0.01 && (fabs(is1[n]-is2[m])/is2[m])<0.01){
      match++;
}

改为：

// if a suitable match is found, increment the match variable
    if(n<=(is1.size()-1) && (m<=is2.size()-1)){ // <- here need to make sure it stays inbounds 
     if((fabs(is1[n]-is2[m])/is1[n])<0.01 && (fabs(is1[n]-is2[m])/is2[m])<0.01){
       match++;
     }
    }

运行 1000 次后，每次都能得到正确的结果。

每天学些新东西。

干杯。

推荐阅读