将 arma::cube 子视图转换为 NumericVector 以使用糖

Transform arma::cube subview into NumericVector to use sugar

我将 3D 数组从 R 传递到 C++ 并 运行 传递到类型转换问题。我们如何将 RcppArmadillo 中的 arma::cube subviews 转换为 NumericVectors 以使用 Rcpp 中的糖函数对它们进行操作,例如 which_min?

假设您有一个带有一些数字条目的 3D 立方体 Q。我的目标是获取每行 i 和每个第三维 k 的列条目的最小值的索引。在 R 语法中,这是 which.min(Q[i,,k]).

例如 i = 1k = 1

cube Q = randu<cube>(3,3,3);
which_min(Q.slice(1).row(1)); // this fails

我认为转换为 NumericVector 可以解决问题,但此转换失败了

which_min(as<NumericVector>(Q.slice(1).row(1))); // conversion failed

我怎样才能让它工作?谢谢你的帮助。

这里有几个选项:

  1. 您可以为此使用 Armadillo 函数,成员函数 .index_min()(请参阅 Armadillo 文档 here)。
  2. 您可以使用 Rcpp::wrap(),其中 "transforms an arbitrary object into a SEXP"arma::cube subviews 变成 Rcpp::NumericVector 并使用糖函数 Rcpp::which_min().

最初我只是将第一个选项作为答案,因为它似乎是完成 objective 的更直接的方法,但我添加了第二个选项(在对答案的更新中),因为我现在考虑任意转换可能是您好奇的一部分。

我将以下 C++ 代码放入文件 so-answer.cpp:

// [[Rcpp::depends(RcppArmadillo)]]
#include <RcppArmadillo.h>

// [[Rcpp::export]]
Rcpp::List index_min_test() {
    arma::cube Q = arma::randu<arma::cube>(3, 3, 3);
    int whichmin = Q.slice(1).row(1).index_min();
    Rcpp::List result = Rcpp::List::create(Rcpp::Named("Q") = Q,
                                           Rcpp::Named("whichmin") = whichmin);
    return result;
}

// [[Rcpp::export]]
Rcpp::List which_min_test() {
    arma::cube Q = arma::randu<arma::cube>(3, 3, 3);
    Rcpp::NumericVector x = Rcpp::wrap(Q.slice(1).row(1));
    int whichmin = Rcpp::which_min(x);
    Rcpp::List result = Rcpp::List::create(Rcpp::Named("Q") = Q,
                                           Rcpp::Named("whichmin") = whichmin);
    return result;
}

我们有一个函数使用 Armadillo 的 .index_min() 和一个使用 Rcpp::wrap() 的函数来启用 Rcpp::which_min().

然后我使用 Rcpp::sourceCpp() 编译它,使函数可用于 R,并演示如何使用几个不同的种子调用它们:

Rcpp::sourceCpp("so-answer.cpp")
set.seed(1)
arma <- index_min_test()
set.seed(1)
wrap <- which_min_test()
arma$Q[2, , 2]
#> [1] 0.2059746 0.3841037 0.7176185
wrap$Q[2, , 2]
#> [1] 0.2059746 0.3841037 0.7176185
arma$whichmin
#> [1] 0
wrap$whichmin
#> [1] 0
set.seed(2)
arma <- index_min_test()
set.seed(2)
wrap <- which_min_test()
arma$Q[2, , 2]
#> [1] 0.5526741 0.1808201 0.9763985
wrap$Q[2, , 2]
#> [1] 0.5526741 0.1808201 0.9763985
arma$whichmin
#> [1] 1
wrap$whichmin
#> [1] 1
library(microbenchmark)
microbenchmark(arma = index_min_test(), wrap = which_min_test())
#> Unit: microseconds
#>  expr    min      lq     mean  median      uq    max neval cld
#>  arma 12.981 13.7105 15.09386 14.1970 14.9920 62.907   100   a
#>  wrap 13.636 14.3490 15.66753 14.7405 15.5415 64.189   100   a

reprex package (v0.2.1)

创建于 2018-12-21