RcppArmadillo 中的矢量化 log1p()

Posted 2023-02-17

【中文标题】RcppArmadillo 中的矢量化 log1p()

【英文标题】：Vectorized log1p() in RcppArmadillo

【发布时间】：2018-11-07 03:23:07

【问题描述】：

将log1p() 应用于整个arma::vec 的适当方法是什么？似乎有log() 和exp() 的矢量化版本，但没有log1p()。我发现NumericVector 有语法糖，所以我最终将arma::vec 转换为NumericVector，应用log1p()，然后再转换回来：

#include <RcppArmadillo.h>

using namespace Rcpp;

// [[Rcpp::depends(RcppArmadillo)]]

// [[Rcpp::export]]
arma::vec test_log1p( arma::vec v )

  // arma::vec res = log1p(v);         // results in a compilation error
  NumericVector v1 = log1p( wrap(v) );
  arma::vec res = as<arma::vec>(v1);
  return res;

有没有更优雅的方式来做到这一点？

【参考方案1】：

魔鬼在细节中。

对于初学者来说，RcppArmadillo 没有“糖”，所以你的推理是有缺陷的——你不能只看在 Rcpp::NumericVector 上工作的 Rcpp Sugar 函数。

再一次，一个人可以像你一样转换。但是您选择了昂贵的转换。查看显式重用内存的高级构造函数——不需要副本。

一个更简单和更直接（但本地）的方法就是添加一个小的本地内联函数。这就是我要做的 :) 几分钟后完成。

最后，我们有一些兄弟项目将 Rcpp Sugar 推广到任何可以使用迭代器的东西上。那是“大道”，它可以做一些新的发展。也许开始at this repo。

【讨论】：

谢谢，德克。感谢您的洞察力。

【参考方案2】：

使用可用于犰狳向量和矩阵的.transform() 或.for_each() 工具。示例：

v.transform( [](double val)  return log1p(val);  );

或

v.for_each( [](double& val)  val = log1p(val);  );  // note the & character

您可能需要使用std 前缀：std::log1p() 而不是log1p()。

【讨论】：

.for_each() 将成为我的建议 (+1)。请注意，OP 不需要使用std::log1p()，因为RcppArmadillo pulls in Rcpp 和Rcpp ensures log1p() is defined（参见第227-28 行）。

这就是我对“添加单行”的评论所想到的。我们可以几乎总是添加这样的东西，无论 Sugar 是否有预制的东西。