gpt4 book ai didi

r - 列表示 3d 矩阵(立方体)Rcpp

转载 作者:行者123 更新时间:2023-12-04 09:39:09 30 4
gpt4 key购买 nike

我有一个程序,我需要在 Rcpp 中重复计算立方体 X(nRow, nCol, nSlice) 的每个切片的列均值,结果均值形成矩阵 M(nCol, nSlice)。以下代码产生错误:

#include <RcppArmadillo.h>

// [[Rcpp::depends(RcppArmadillo)]]
using namespace Rcpp;
using namespace arma;

// [[Rcpp::export]]

mat cubeMeans(arma::cube X){
int nSlice = X.n_slices;
int nCol = X.n_cols;
int nRow = X.n_rows;
arma::vec Vtmp(nCol);
arma::mat Mtmp(nRow, nCol);
arma::mat Means(nCol, nSlice);
for (int i = 0; i < nSlice; i++){
Mtmp = X.slice(i);
for(int j = 0; j < nCol; j++){
Vtmp(j) = sum(Mtmp.col(j))/nRow;
}
Means.col(i) = Vtmp;
}
return(wrap(Means));
}

'/Rcpp/internal/Exporter.h:31:31: error: no matching function for call to 'arma::Cube::Cube(SEXPREC*&)'

我想不通。当函数的输入是一个矩阵(并返回一个向量)时,我没有得到错误。但是,我将上述功能作为我的主程序的一部分,即

#include <RcppArmadillo.h>

// [[Rcpp::depends(RcppArmadillo)]]
using namespace Rcpp;
using namespace arma;

mat cubeMeans(arma::cube X){
int nSlice = X.n_slices;
...
return(Means);
}

// [[Rcpp::export]]

main part of program

程序编译成功,但速度非常慢(几乎与使用 colMeans 的 R 版本程序一样慢)。有没有更好的方法来计算多维数据集上的列均值,为什么会出现编译错误?

如果有任何帮助,我将不胜感激。

问候,

最佳答案

我在尝试使用 arma::cube 时也收到此错误作为 Rcpp 函数参数。 根据编译器错误,我相信这是因为没有 Rcpp::wrap<arma::cube>当前已定义(这是处理传递给函数的 R 对象所必需的)。 † 在网上阅读了几个相关示例后,看起来典型的解决方法是在 R 中阅读 array作为 NumericVector , 因为它保留了它的 dims属性,使用这些设置您的 arma::cube方面。尽管需要额外的一两步来解决 缺失 wrap 的问题特化 †,我放在一起的 Armadillo 版本似乎比我的 R 解决方案快很多:

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

// [[Rcpp::export]]
arma::mat cube_means(Rcpp::NumericVector vx) {

Rcpp::IntegerVector x_dims = vx.attr("dim");
arma::cube x(vx.begin(), x_dims[0], x_dims[1], x_dims[2], false);

arma::mat result(x.n_cols, x.n_slices);
for (unsigned int i = 0; i < x.n_slices; i++) {
result.col(i) = arma::conv_to<arma::colvec>::from(arma::mean(x.slice(i)));
}

return result;
}

/*** R

rcube_means <- function(x) t(apply(x, 2, colMeans))

xl <- array(1:10e4, c(100, 100 ,10))
all.equal(rcube_means(xl), cube_means(xl))
#[1] TRUE

R> microbenchmark::microbenchmark(
"R Cube Means" = rcube_means(xl),
"Arma Cube Means" = cube_means(xl),
times = 200L)
Unit: microseconds
expr min lq mean median uq max neval
R Cube Means 6856.691 8204.334 9843.7455 8886.408 9859.385 97857.999 200
Arma Cube Means 325.499 380.540 643.7565 416.863 459.800 3068.367 200

*/

我在这里利用了 arma::mean 的事实arma::mat 的函数重载s 将默认计算列均值(arma::mean(x.slice(i), 1) 将为您提供该切片的行均值)。


编辑: † 再想一想,我不确定这是否与 Rcpp::wrap 有关或不是 - 但问题似乎与缺失 Exporter<> 有关特化arma::cube - Rcpp 的 Exporter.h 第 31 行:

template <typename T>
class Exporter{
public:
Exporter( SEXP x ) : t(x){}
inline T get(){ return t ; }

private:
T t ;
} ;

无论如何,NumericVector/我使用的设置尺寸方法目前看来是实用的解决方案。


根据您在问题中描述的输出维度,我假设您希望结果矩阵的每一列都是相应数组切片的列均值向量(第 1 列 = 切片 1 的列均值,等等... ), 即

R> x <- array(1:27, c(3, 3, 3))
R> rcube_means(x)
[,1] [,2] [,3]
[1,] 2 11 20
[2,] 5 14 23
[3,] 8 17 26
R> cube_means(x)
[,1] [,2] [,3]
[1,] 2 11 20
[2,] 5 14 23
[3,] 8 17 26

但如果需要,您可以轻松更改它。

关于r - 列表示 3d 矩阵(立方体)Rcpp,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/30969799/

30 4 0
Copyright 2021 - 2024 cfsdn All Rights Reserved 蜀ICP备2022000587号
广告合作:1813099741@qq.com 6ren.com