稀疏矩阵上的 R 扫描答案

【问题标题】：R sweep on a Sparse Matrix稀疏矩阵上的 R 扫描
【发布时间】：2019-08-19 19:37:33
【问题描述】：

我正在尝试将 sweep 函数应用于稀疏矩阵 (dgCMatrix)。不幸的是，当我这样做时，我得到了一个内存错误。似乎扫描正在将我的稀疏矩阵扩展为一个完整的密集矩阵。

如果有一种简单的方法来执行此功能而不会破坏我的记忆？

这就是我想要做的。

sparse_matrix <- sweep(sparse_matrix, 1, vector_to_multiply, '*')

【问题讨论】：

假设vector_to_multiply的长度等于sparse_matrix的行数，你也许可以做到sparse_matrix* vector_to_multiply

标签： r sparse-matrix

【解决方案1】：

我赞同@user20650 的建议，即使用mat * vec 形式的直接乘法，它通过隐式回收vec 将矩阵mat 的每一列与向量vec 相乘。

处理时间分析

我知道您在这里的主要要求是内存，但是执行 microbenchmark 比较 sweep 和密集矩阵和稀疏矩阵的直接乘法方法很有趣：

# Sample data
library(Matrix)  
set.seed(2018)
mat <- matrix(sample(c(0, 1), 10^6, replace = T), nrow = 10^3)
mat_sparse <- Matrix(mat, sparse = T)
vec <- 1:dim(mat)[1]

library(microbenchmark)
res <- microbenchmark(
    sweep_dense = sweep(mat, 1, vec, '*'),
    sweep_sparse = sweep(mat_sparse, 1, vec, '*'),
    mult_dense = mat * vec,
    mult_sparse = mat_sparse * vec
)
res
Unit: milliseconds
         expr        min         lq       mean     median        uq      max
  sweep_dense   8.639459  10.038711  14.857274  13.064084  18.07434  32.2172
 sweep_sparse 116.649865 128.111162 162.736864 135.932811 155.63415 369.3997
   mult_dense   2.030882   3.193082   7.744076   4.033918   7.10471 184.9396
  mult_sparse  12.998628  15.020373  20.760181  16.894000  22.95510 201.5509

library(ggplot2)
autoplot(res)

平均而言，涉及稀疏矩阵的运算实际上比具有密集矩阵的运算要慢一些。但是请注意，直接乘法比sweep 快得多。

内存分析

我们可以使用memprof 来分析不同方法的内存使用情况。

library(profmem)
mem <- list(
    sweep_dense = profmem(sweep(mat, 1, vec, '*')),
    sweep_sparse = profmem(sweep(mat_sparse, 1, vec, '*')),
    mult_dense = profmem(sweep(mat * vec)),
    mult_sparse = profmem(sweep(mat_sparse * vec)))
lapply(mem, function(x) utils:::format.object_size(sum(x$bytes), units = "Mb"))
#$sweep_dense
#[1] "15.3 Mb"
#
#$sweep_sparse
#[1] "103.1 Mb"
#
#$mult_dense
#[1] "7.6 Mb"
#
#$mult_sparse
#[1] "13.4 Mb"

说实话，我很惊讶与稀疏矩阵的直接乘法的内存印记不小于涉及密集矩阵的内存印记。也许样本数据过于简单。可能值得用您的实际数据（或其代表性子集）来探索这一点。

【讨论】：

在你的例子中，稀疏矩阵并不是真的稀疏（50%）；如果你增加稀疏度，相对速度和内存将会提高
谢谢@user20650；如果/当我有空闲时间时，我会尝试提出一个更具代表性的例子。
很棒的答案。谢谢你们俩。我学到了很多东西。

【解决方案2】：

我正在处理 NLP 问题中的一个大且非常稀疏的 dgTMatrix 矩阵（200k 行和 10k 列）。经过几个小时思考一个好的解决方案后，我为稀疏矩阵创建了一个替代的sweep 函数。它非常快速且内存高效。将所有矩阵行乘以权重数组仅需 1 秒且不到 1G 的内存。对于margin = 1，它适用于dgCMatrix 和dgTMatrix。

如下：

sweep_sparse <- function(x, margin, stats, fun = "*") {
   f <- match.fun(fun)
   if (margin == 1) {
      idx <- x@i + 1
   } else {
      idx <- x@j + 1
   }
   x@x <- f(x@x, stats[idx])
   return(x)
}

【讨论】：