代码之家 › 专栏 › 技术社区 › alberto

带有指向矩阵子视图的指针的rmultinom

armadillo rcpp

alberto · 技术社区 · 7 年前

要从Rcpp中的多项式分布中采样,我们可以做到:

int n = 100;
int k = 3;
arma::vec probs = {0.4,0.2,0.4}
arma::irowvec c(k);
Rcpp::rmultinom(n, probs.begin(), k, c.begin());

我想知道当C是一个矩阵时,我们是否也可以这样做。我试着

int n = 100;
int k = 3;
arma::vec probs = {0.4,0.2,0.4}
# C (arma::mat C) passed by reference to the function
Rcpp::rmultinom(n, probs.begin(), k, C.row(1).begin());

但它抛出了一个错误。有简单的方法吗?

我想尝试第二种方法,因为我有一个大的矩阵C,我通过引用传递给我的函数,然后我想在多项式之后更新它的行。

MWE:

#include <RcppArmadillo.h>

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

using namespace Rcpp;
using namespace arma;

void subfunction(const arma::imat& C) {

  int n = 100;
  int k = 3;
  arma::vec probs = {0.4,0.2,0.4};

  rmultinom(n, probs.begin(), k, C.col(1).begin());
}

// [[Rcpp::export]]
arma::imat myfunction(){

  arma::imat C = { {1, 2}, 
                  {3, 4},
                  {5, 6}};
  subfunction(C);
  Rcpp::Rcout << "C: " <<  C << std::endl;

}

错误是:

test_multim.cpp:18:44: error: âclass arma::subview_col<int>â 
has no member named âbeginâ
rmultinom(n, probs.begin(), k, Ct.col(1).begin());

1 回复 | 直到 7 年前

coatless 7 年前

这个问题有两个问题:

要求的类型 rmultinom 函数定义。
访问内存指针以获取 armadillo 矩阵

首先,请注意其中一个错误是:

mondayso.cpp:16:3: error: no matching function for call to 'Rf_rmultinom'
  rmultinom(1, probs.begin(), k, C.colptr(1));
  ^~~~~~~~~
/Library/Frameworks/R.framework/Resources/include/Rmath.h:468:6: note: candidate function not viable: no known conversion from 'double *' to 'int *' for 4th argument
void    rmultinom(int, double*, int, int*);
        ^
1 error generated.

本质上 rmultinorm 作用必须在第四个参数中传入一个整数。由于施工 arma::mat 存在 double 默认情况下,矩阵的类型是不合适的。在这种情况下 C 矩阵必须为 arma::imat 因为它使用 arma::sword 或签名 int 组成部分

接下来,数据 armadillo matrices is stored in a column-by-column order (参见 Wikipedia's entry 容易地 柱通过 .colptr . 这解决了出现的第二个错误:

error: no member named 'begin' in 'arma::subview_row<int>'
  rmultinom(n, probs.begin(), k, C.row(1).begin());
                                 ~~~~~~~~ ^
1 error generated.

说到这里,我已经构建了一个示例来促进转换。

#include <RcppArmadillo.h>

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

// [[Rcpp::plugins(cpp11)]]

// [[Rcpp::export]]
arma::imat test() {

  int n = 100;
  int k = 3;
  arma::vec probs = {0.4,0.2,0.4};

  arma::imat C = { {1, 3, 5},
                  {2, 4, 6} };

  arma::imat Ct = C.t();

  // C++ indices start at 0 (thus, this is the second column!)
  rmultinom(n, probs.begin(), k, Ct.colptr(1));

  return Ct;
}

测验

set.seed(111)
test()
#      [,1] [,2]
# [1,]    1   43
# [2,]    3   18
# [3,]    5   39