我在Rcpp
中编写了一个函数qSelectMbycol
,它在O(n)时间内返回每列的k
最大元素。该功能运行正常。如果我尝试执行相同的操作,但对行而不是列进行操作(函数 qSelectMbyrow
),则会返回错误 “错误:Mat::init():请求的大小与列 vector 不兼容布局”
。有人认为我做错了什么吗?我将此文件另存为“qselect.cpp
”:
// [[Rcpp::depends(RcppArmadillo)]]
#define RCPP_ARMADILLO_RETURN_COLVEC_AS_VECTOR
#include <RcppArmadillo.h>
using namespace arma;
// [[Rcpp::export]]
arma::vec qSelectMbycol(arma::mat& M, const int k) {
// ARGUMENTS
// M: matrix for which we want to find the k-th largest elements of each column
// k: k-th statistic to look up
arma::mat Y(M.memptr(), M.n_rows, M.n_cols);
// we apply over columns
int c = M.n_cols;
arma::vec out(c);
int i;
for (i = 0; i < c; i++) {
arma::vec y = Y.col(i);
std::nth_element(y.begin(), y.begin() + k - 1, y.end());
out(i) = y(k-1); // the k-th largest value of each column
}
return out;
}
// [[Rcpp::export]]
arma::vec qSelectMbyrow(arma::mat& M, const int k) {
// ARGUMENTS
// M: matrix for which we want to find the k-th largest elements of each row
// k: k-th statistic to look up
arma::mat Y(M.memptr(), M.n_rows, M.n_cols);
// we apply over rows
int r = M.n_rows;
arma::vec out(r);
int i;
for (i = 0; i < r; i++) {
arma::vec y = Y.row(i); // this line throws the error "error: Mat::init(): requested size is not compatible with column vector layout"
std::nth_element(y.begin(), y.begin() + k - 1, y.end());
out(i) = y(k-1); // should give k-th largest value of each row
}
return out;
}
示例:
n=500
p=100
set.seed(1)
M=matrix(rnorm(n, mean = 100, sd = 1),n,1)
library(Rcpp)
library(RcppArmadillo)
Rcpp::sourceCpp('qselect.cpp')
qSelectMbycol(M,5) # works OK
qSelectMbyrow(M,5) # throws error "error: Mat::init(): requested size is not compatible with column vector layout"
我也尝试过插入
typedef std::vector<double> stdvec;
并将线设置 vector y
替换为
arma::vec y = arma::conv_to<stdvec>::from(Y.row(i));
在我的 qSelectMbyrow
函数中,虽然该函数随后运行,但与在列上应用相比,它运行缓慢,并且如果运行 100 次,还会使我的 R session 崩溃。
最佳答案
问题是 arma::vec
实际上是 arma::colvec
(参见 the docs )。所以,我们可以通过改变来解决这个问题
arma::vec y = Y.row(i);
(这是不兼容的,因为它认为您想要一个具有一列的矩阵,但您试图给它一个具有一行的矩阵)
arma::rowvec y = Y.row(i);
关于c++ - 使用 RcppArmadillo 在矩阵列上应用函数可以工作,但在行上应用时会返回错误,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/54852368/