问题描述
我正在研究一个 Rcpp 稀疏矩阵类,它同时使用 Rcpp::IntegerVector
(行/列指针)和模板化的 std::vector<T>
。基本原理是,可以通过简单地将它们作为指向 R 对象的指针来避免在极大稀疏矩阵中深度复制整数指针向量(@i
,@p
)的开销,并且始终如一地,微基准显示与转换为 Eigen::SparseMatrix<T>
和 arma::SpMat<T>
相比,这种方法几乎只需要一半的时间,同时使用更少的内存。
Bare-bones Rcpp 稀疏矩阵类
namespace SpRcpp {
template<typename T>
class SpMatrix {
public:
Rcpp::IntegerVector i,p;
std::vector<T> x;
unsigned int n_row,n_col;
// constructor for the class from an Rcpp::S4 object
SpMatrix(Rcpp::S4& mat) {
Rcpp::IntegerVector dims = mat.slot("Dim");
n_row = (unsigned int)dims[0];
i = mat.slot("i");
p = mat.slot("p");
x = Rcpp::as<std::vector<T>>(mat.slot("x"));
};
// other constructors,class methods,iterators,etc.
};
}
示例用法:
//[[Rcpp::export]]
std::vector<float> SpRcpp_SpMatrix(Rcpp::S4& mat) {
SpRcpp::SpMatrix<float> A(mat);
return A.x;
}
这有效!
但是,我希望 Rcpp 将 S4 dgCMatrix 对象(例如)隐式转换为 SpRcpp::SpMatrix
对象以启用如下功能:
隐式 Rcpp 包装
//[[Rcpp::export]]
std::vector<float> SpRcpp_SpMatrix2(SpRcpp::SpMatrix<float>& mat) {
return mat.x;
}
这就是使用 Rcpp::as
的情况。
我尝试了以下方法:
namespace Rcpp {
namespace traits {
template <typename T>
class Exporter< SpRcpp::SpMatrix<T> > {
public:
Exporter(SEXP x) { Rcpp::S4 mat = x; }
SpRcpp::SpMatrix<T> get() {
return SpRcpp::SpMatrix<T>(mat);
}
private: Rcpp::S4 mat;
};
}
}
这可以编译,我知道 SpRcpp::SpMatrix<T>(Rcpp::S4& x)
构造函数可以工作,但是当我尝试将 dgCMatrix 输入 SpRcpp_SpMatrix()
时,我收到错误:
Error in SpRcpp_SpMatrix2(A) : Not an S4 object.
我认为这是因为我在所有类声明之前都声明了以下内容:
#include <RcppCommon.h>
#include <Rcpp.h>
根据文档 here、RcppGallery example 和 RcppArmadillo 实现,#include <Rcpp.h>
不得位于 Rcpp::as
和 Rcpp::wrap
函数之前,但就我而言,我不能这样做,因为我的类定义需要 Rcpp.h
。
问题:当 Rcpp::as
类依赖于 {{1} 时,如何从 SpRcpp::SpMatrix
dgCMatrix 创建 Rcpp::S4
的 SpMatrix
}?
解决方法
创建一个 Rcpp SparseMatrix 类其实很简单!我多虑了。
#include <rcpp.h>
// Rcpp for sparse matrices (spRcpp)
namespace Rcpp {
class SparseMatrix {
public:
Rcpp::IntegerVector i,p;
Rcpp::NumericVector x;
int n_rows,n_cols;
// constructor
SparseMatrix(Rcpp::S4 mat) {
Rcpp::IntegerVector dim = mat.slot("Dim");
i = mat.slot("i");
p = mat.slot("p");
x = mat.slot("x");
n_rows = (int)dim[0];
n_cols = (int)dim[1];
};
};
}
namespace Rcpp {
template <> Rcpp::SparseMatrix as(SEXP mat) {
return Rcpp::SparseMatrix(mat);
}
}
//[[Rcpp::export]]
Rcpp::NumericVector toRcppSparseMatrix(Rcpp::SparseMatrix& A) {
return A.x;
}
给定 Matrix::dgCMatrix
、mat
,调用 toRcppSparseMatrix(mat)
会在 1-2 微秒内返回 2500 万个值的非零值。这与 RcppArmadillo 或 RcppEigen 稀疏矩阵转换形成对比,后者需要大约 250 毫秒才能对同一矩阵进行深度复制并在内存中运行。
正如 Dirk 所建议的,使用 RcppArmadillo ivec
和 dvec
非常有效,但仍会创建浅拷贝,这会导致大约 100 毫秒的运行时间并消耗一些内存。
显然,上述方法仅限于双精度类型,因此没有深拷贝就无法进行浮点运算。