通过替换双层循环实现更快的矩阵乘法

7

我有一个数据框,它看起来有点像以下代码产生的结果(但要大得多)

set.seed(10)    
mat <- matrix(rbinom(200, size=1, prob = .5), ncol = 10)

在这些列中,数字1表示一个观察值对特定问题感兴趣。我想生成一个网络,比较所有的观察值,并统计每个二元组共同感兴趣的问题数量。

我已经编写了下面的代码,看起来工作得很好:

mat2 <- matrix(NA,20,20)

for(i in 1:nrow(mat)){
    for(j in 1:nrow(mat)){
       mat2[i,j] <- sum(as.numeric(mat[i,]==1) + as.numeric(mat[j,]==1) == 2)
    }
 }

我将每一条记录与其他所有记录进行比较,只有当两个记录都有一个输入(即它们对此感兴趣),这时才会计为该主题的共同兴趣。

我的问题是我的数据集非常大,而循环现在已经运行了数小时。

是否有人有想法可以避免使用循环来完成这个任务?


3
try mat %*% t(mat) - Tom
@Flow请查看复制嵌套循环的问题。https://dev59.com/Ybfna4cB1Zd3GeqPnSVz - simar
3
在这种情况下,应该使用tcrossprod(mat==1)。它非常快速。请将其发布为答案! - jogo
感谢大家的评论和答案,非常有帮助!我甚至没有考虑过矩阵乘法,尽管这些解决方案非常简单和好。再次感谢你们在这里的帮助。 - Flow
2个回答

5

这应该会更快:

tmat <- t(mat==1)
mat4 <- apply(tmat, 2, function(x) colSums(tmat & x))

5

我将继续支持@jogo的评论,因为这是目前最快捷的方法(感谢提示,我也会在生产中使用它)。

set.seed(10)    
mat <- matrix(rbinom(200, size=1, prob = .5), ncol = 10)
mat2 <- matrix(NA,20,20)
binary_mat <- mat == 1
tmat <- t(mat==1)

microbenchmark::microbenchmark(
  "loop" = for(i in 1:nrow(mat)){
             for(j in 1:nrow(mat)){
               mat2[i,j] <- sum(as.numeric(mat[i,]==1) + as.numeric(mat[j,]==1) == 2)
             }
           }, 
  "apply" = mat4 <- apply(tmat, 2, function(x) colSums(tmat & x)), 
  "matrix multiplication" = mat5 <- mat %*% t(mat),
  "tcrossprod" = tcrossprod(mat),
  "tcrossprod binary" = tcrossprod(binary_mat)
)

在我的机器上,这个基准测试的结果是
Unit: microseconds
                  expr       min        lq        mean    median         uq       max neval cld
                  loop 16699.634 16972.271 17931.82535 17180.397 17546.1545 31502.706   100   b
                 apply   322.942   330.046   395.69045   357.886   368.8300  4299.228   100  a 
 matrix multiplication    21.889    28.801    36.76869    39.360    43.9685    50.689   100  a 
            tcrossprod     7.297     8.449    11.20218     9.984    14.4005    18.433   100  a 
     tcrossprod binary     7.680     8.833    11.08316     9.601    12.0970    35.713   100  a 

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接