r - 将Lambda向量传递给Rcpp的rpois

在R中，rpois可以传递描述多个Poisson分布的lambda向量，例如

rpois(5, (1:5)*1000)

# [1] 1043 1974 3002 3930 4992

在上面，输出矢量的每个元素都是从不同的泊松分布中提取的，分别使用1000、2000、3000、4000和5000。

如果我有一个arma::mat（因为在其他地方我正在使用多维数据集而使用这些元素）包含Poisson分布的lambda，那么将这些（一次一行）传递给Rcpp中的rpois的最佳方法是什么？

这是一个玩具示例，以及随后出现的错误消息的摘录：

library(inline)
library(RcppArmadillo)
code <- "
  using namespace Rcpp;
  using namespace arma;
  arma_rng::set_seed(42); // Dirk's seed of choice

  mat lam = randu(5, 5); // ignore the fact these are all 0-1
  mat out(5, 5);

  for (int i = 0; i < 5; i++) {
    out.row(i) = rpois(5, lam.row(i));
  }

  return(wrap(out));
"

f <- cxxfunction(body=code, plugin="RcppArmadillo")

# cannot convert 'arma::subview_row<double>' to 'double' for argument '2'
#   to 'Rcpp::NumericVector Rcpp::rpois(int, double)'

我必须承认，我对C ++中类型转换的理解很差。我正在尝试做的事情是否可能（我的猜测是否定的，因为似乎rpois期望是两倍），还是我需要遍历矩阵的每个单元格，每次生成一个偏差？

最佳答案

在C / C ++中，您可以访问至少2个泊松偏差生成例程（在此online manual中搜索rpois）。

它们的声明如下：

double R::rpois(double mu);
NumericVector Rcpp::rpois(int n, double mu);

它们均不允许在mu中传递> 1个值（也称为lambda）。第一个功能是R的原始例程，用于实现rpois，正如我们从R stats包中所知道的那样（该包已向量化所有参数）。给定单个mu，它将返回单个（伪）随机偏差。

第二个是所谓的Rcpp糖功能。它允许一次计算n偏差并将其作为NumericVector返回（再次使用R::rpois）。

换句话说，您应该通过调用for使用两个嵌套的R::rpois循环来填充矩阵。不要害怕这种方法，这就是C ++。 :)