密集与稀疏矩阵代数的速度

Question

我将在 R 中使用一个非常大 (7 e6 x 4.5 e3) 但非常稀疏的矩阵。所以我试图了解如何有效地处理稀疏矩阵。我有两个相关的问题。

首先：我了解到Matrix 包会自动链接到 LAPACK 和 SuiteSparse 编译的 dll。 （我在 Windows 中工作。）我的印象是，与使用 LAPACK 套件使用密集矩阵相比，使用 SuiteSparse 例程会缩短执行时间。但下面的测试运行表明，矩阵的稀疏版本的运行时间比密集版本慢很多。

> library(Matrix)
> sparse <- sparseMatrix(1:4, 1:4, x=rnorm(4))
> dense <- as.matrix(sparse)
> x <- 1:4
> system.time(for (i in 1:10000) sparse %*% x)
   user  system elapsed 
   0.23    0.00    0.23 
> system.time(for (i in 1:10000) dense %*% x)
   user  system elapsed 
      0       0       0 
> system.time(for (i in 1:1000) solve(sparse))
   user  system elapsed 
   3.94    0.00    3.94 
> system.time(for (i in 1:1000) solve(dense))
   user  system elapsed 
   0.05    0.00    0.05

a) Matrix 自动与上述两个编译库连接是否正确？如果没有，我如何链接到这些 DLL？ b) 使用稀疏矩阵代数实际上是否比使用密集矩阵代数慢很多？

第二：我已经安装了RcppEigen 和RcppArmadillo 包。我已经能够使用RcppArmadillo 编译一个测试程序（使用 Dirk Eddelbuettel 和 Conrad Sanderson 的论文）。但在我的一生中，我还没有找到与RcppEigen 类似的介绍，这将为我提供一些可用于入门的模型代码。你们中的任何人都可以指出一个类似于 Eddelbuettel 和 Sanderson 论文的文件，它可以帮助我开始使用 RcppEigen 吗？

Answer 1

（评论有点太长了。）我将首先对更大的矩阵进行分析；我可以想象，当矩阵很小且不是 非常 稀疏（例如，在这种情况下，25% 的单元格不为零）时，稀疏算法处于劣势。在下面的示例（1000x1000 矩阵）中，稀疏求解器比密集求解器快 26 倍。您可能会发现 Matrix 例程对于您的目的来说足够快，而无需承担学习 (Rcpp)Eigen/(Rcpp)Armadillo 的额外认知开销 ...

library(rbenchmark)
library(Matrix)
set.seed(101)
sparse <- sparseMatrix(1:1000,1:1000,x=rnorm(1000))
dense <- as.matrix(sparse)
benchmark(solve(sparse),solve(dense),replications=20,
          columns = c(
       "test", "replications", "elapsed", "relative", "user.self"))
##            test replications elapsed relative user.self
## 2  solve(dense)           20   6.932   26.868     6.692
## 1 solve(sparse)           20   0.258    1.000     0.256