在R中对多个矩阵进行乘法运算

Question

在R中对多个矩阵进行乘法运算

4

我想要将同样大小的多个矩阵与一个初始向量相乘。在下面的例子中，p.state 是一个包含 m 个元素的向量，tran.mat 是一个列表，其中每个成员都是一个 m x m 的矩阵。

for (i in 1:length(tran.mat)){
  p.state <- p.state %*% tran.mat[[i]]
}

以上代码给出了正确答案，但在length(tran.mat)较大时可能会很慢。我想知道是否有更高效的方法来实现这个问题？

以下是一个示例，其中m=3和length(mat)=10可以生成此结果：

p.state <- c(1,0,0)
tran.mat<-lapply(1:10,function(y){apply(matrix(runif(9),3,3),1,function(x){x/sum(x)})})

for (i in 1:length(tran.mat)){
  p.state <- p.state %*% tran.mat[[i]]
}

print(p.state)

注意： tran.mat 不一定是一个列表，它只是目前被写成一个列表。

在一些评论之后进行编辑：

当 m 很小的时候，可以使用 Reduce。然而，当 m=6 时，循环表现优于上述两种解决方案。请注意，保留了 html 标签。

p.state1 <- p.state <- c(1,0,0,0,0,0)
tran.mat<-lapply(1:10000,function(y){t(apply(matrix(runif(36),6,6),1,function(x){x/sum(x)}))})

tst<-do.call(c, list(list(p.state), tran.mat))

benchmark(
  'loop' = {
    for (i in 1:length(tran.mat)){
      p.state <- p.state %*% tran.mat[[i]]
    }
  },
  'reduce' = {
    p.state1 %*% Reduce('%*%', tran.mat)
   },
  'reorder' = {
    Reduce(`%*%`,tran.mat,p.state1)
  }

这将导致

        test replications elapsed relative user.self sys.self user.child     sys.child
   1    loop          100    0.87    1.000      0.87        0         NA        NA
   2  reduce          100    1.41    1.621      1.39        0         NA        NA
   3 reorder          100    1.00    1.149      1.00        0         NA        NA

- wearetheboro

3个回答

1

我将使用来自Rfast包的函数来减少乘法的执行时间。不幸的是，循环的时间无法缩短。

名为Rfast::eachcol.apply的函数是您目的的绝佳解决方案。您的乘法也是crossprod函数，但对于我们的目的而言速度较慢。

以下是一些辅助函数：

mult.list<-function(x,y){
    for (xm in x){
        y <- y %*% xm
    }
    y
}

mult.list2<-function(x,y){
    for (xm in x){
        y <- Rfast::eachcol.apply(xm,y,oper="*",apply="sum")
    }
    y
}

Here is an example:

x<-list()
y<-rnomr(1000)
for(i in 1:100){
    x[[i]]<-Rfast::matrnorm(1000,1000)
}


microbenchmark::microbenchmark(R=a<-mult.list(x,y),Rfast=b<-mult.list2(x,y),times = 10)
 Unit: milliseconds
     expr        min         lq        mean     median         uq        max neval
        R 410.067525 532.176979 633.3700627 649.155826 699.721086 916.542414    10
    Rfast 239.987159 251.266488 352.1951486 276.382339 458.089342 741.340268    10

all.equal(as.numeric(a),as.numeric(b))
[1] TRUE

参数oper用于每个元素的操作，而apply用于每列的操作。在大矩阵中应该很快。我无法在我的笔记本电脑上测试更大的矩阵。

- Manos Papadakis

1

我不确定这样做会更快，但它更短:

prod <- Reduce("%*%", L)

all.equal(prod, L[[1]] %*% L[[2]] %*% L[[3]] %*% L[[4]])
## [1] TRUE

注意

我们使用了以下测试输入：

m <- matrix(1:9, 3)
L <- list(m^0, m, m^2, m^3)

- G. Grothendieck

你也在字符串中引用函数名吗？:( - Konrad Rudolph

我强烈反对。这会破坏类型系统并产生语义后果。语法上的差异（或者说相似性）纯粹是偶然的，实际上非常误导人。你不可能“调用一个字符串”，这个概念本身就是无意义的。而且 R 通过 match.fun 和语法灵活性在幕后进行魔法操作，这对初学者和中级 R 用户造成了极大的伤害，因为他们误解了正在发生的事情。按照你的论点，写 "a" = 1 是“可以”的。你肯定同意那段代码是荒谬的吧。 - Konrad Rudolph

那是推卸责任。我的“观点”是有理有据的，不是凭空臆断的。这就是为什么我很惊讶看到像你这样经验丰富的R用户持不同意见。当您将函数名作为参数传递时，您是否引用所有函数名？还是只有运算符？ - Konrad Rudolph

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- qdread · Accepted Answer

更快的方法是使用Reduce()在矩阵列表上进行顺序矩阵乘法。

这样可以获得大约4倍的加速。以下是您的代码示例，测试了1000个元素而不是10个，以便更容易地看到性能提升。

代码

library(rbenchmark)

p.state <- c(1,0,0)
tran.mat<-lapply(1:1000,function(y){apply(matrix(runif(9),3,3),1,function(x){x/sum(x)})})

benchmark(
  'loop' = {
    for (i in 1:length(tran.mat)){
      p.state <- p.state %*% tran.mat[[i]]
    }
  },
  'reduce' = {
    p.state %*% Reduce('%*%', tran.mat)
  }
)

输出

    test replications elapsed relative user.self sys.self user.child sys.child
1   loop          100    0.23    3.833      0.23        0         NA        NA
2 reduce          100    0.06    1.000      0.07        0         NA        NA

你可以看到，reduce方法的速度大约快了3.8倍。