向量的稀疏矩阵

Question

我有一个带有值 ( val) 的向量和一个指示组成员身份的向量 ( group)：

vec   <- 1:9
group <- rep(1:3, c(2,4,3))

假设我们有K组和总值N，因此两个向量都有 length N。目标是有效构建稀疏“块对角”矩阵，其中第一列保存组 1 的值，第二列保存组 2 的值，依此类推。但是，这些值不应该“重叠”，因为每行应该只有一个值，请参阅下面的解决方案。我需要用非常大的尺寸来执行此操作数千K次N。因此，以下基于循环的解决方案不够高效：

K     <- length(unique(group))
N     <- length(group)
M     <- matrix(0, N, K)

for(k in 1:K){
  
 M[group == k, k] <- vec[group == k]
        
}

Matrix::Matrix(M, sparse = T)

9 x 3 sparse Matrix of class "dgCMatrix"
           
 [1,] 1 . .
 [2,] 2 . .
 [3,] . 3 .
 [4,] . 4 .
 [5,] . 5 .
 [6,] . 6 .
 [7,] . . 7
 [8,] . . 8
 [9,] . . 9

出于内存原因，理想的情况是直接基于稠密N时间K矩阵构造稀疏矩阵，而无需中间步骤。

---

编辑

对于上面给出的小例子，事实证明基于循环的解决方案非常有效：

Unit: microseconds
     expr     min       lq     mean   median       uq      max neval cld
      ben 734.280 771.7000 826.8372 787.5230 805.2710 3185.158   100   b
      CJR 711.187 745.1855 813.9948 766.9960 781.7495 4832.476   100   b
 original 199.714 221.9520 235.4320 227.9395 236.7065  379.757   100  a 

然而，当转向高维示例（N = 10,000 且 K = 1,000）时，CJR 的解决方案在速度方面是获胜者：

Unit: milliseconds
     expr        min         lq       mean     median         uq        max neval cld
      ben 128.529311 133.308972 140.032070 135.921289 139.272589 289.668852   100  b 
      CJR   1.841474   2.055513   2.261732   2.201557   2.395925   6.330544   100 a  
 original  93.387806 118.348398 171.380301 125.884493 244.421699 365.871433   100   c

Answer 1

Matrix::.bdiag()可以让您直接从矩阵列表构造块对角（稀疏）矩阵：

\n

mm <- lapply(split(vec, group), matrix)\nMatrix::.bdiag(mm)\n

\n

.bdiag(mm)大约相当于do.call(Matrix::bdiag, mm)：?bdiag说

\n

\n

\xe2\x80\x98bdiag()\xe2\x80\x99 的值继承自类 \xe2\x80\x98CsparseMatrix\xe2\x80\x99，而 \xe2\x80\x98.bdiag()\xe2\x80\x99返回 \xe2\x80\x98TsparseMatrix\xe2\x80\x99。

\n

\n

（前者是排序压缩的面向列的形式，后者是三元组形式：?"TsparseMatrix-class"表示“一旦创建了[面向三元组的矩阵]，然而，该矩阵通常被强制为\xe2\x80\x98CsparseMatrix\xe2” \x80\x99 用于进一步操作。\')

\n

?bdiag还有一个注意事项：

\n

\n

该函数已被编写并且对于相对较少的块矩阵（通常本身是稀疏的）的情况是有效的。

\n

\n

因此，这个解决方案肯定会比您现有的解决方案更好，但可能还可以进一步改进。

\n