将lambda的向量传递给Rcpp的rpois

jba*_*ums 5 random r rcpp

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

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

# [1] 1043 1974 3002 3930 4992
Run Code Online (Sandbox Code Playgroud)

在上面,输出矢量的每个元素是从不同的泊松分布中绘制的,分别具有1000,2000,3000,4000和5000的平均值.

如果我有一个arma::mat(使用这些,因为在其他地方,我使用立方体)包含泊松分布的lambda,那么将这些(一次一行)传递rpois到Rcpp内的最佳方法是什么?

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

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)'
Run Code Online (Sandbox Code Playgroud)

我必须承认我对c ++中类型转换的理解很差.我正在尝试做什么(我的猜测是否定的,因为它似乎rpois期望加倍),或者我是否需要迭代矩阵的每个单元格,每次产生一个偏差?

gag*_*ews 6

从C/C++开始,您可以访问至少2个Poisson偏差生成例程(rpois在此在线手册中搜索).

他们的声明如下:

double R::rpois(double mu);
NumericVector Rcpp::rpois(int n, double mu);
Run Code Online (Sandbox Code Playgroud)

它们都不允许在mu(aka lambda)中传递> 1个值.第一个函数是R的原始例程,用于实现rpois我们从R stats包中知道的(与其所有参数一起向量化的包).给定单个mu,它返回单个(伪)随机偏差.

第二个是所谓的Rcpp糖功能.它允许n一次计算偏差并将它们作为a NumericVector(再次,通过使用R::rpois)返回.

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