时间戳

首页 > 解决方案 > 一些来自矩阵的输入和其他来自向量的输入，如何替换二维循环？

r - 一些来自矩阵的输入和其他来自向量的输入，如何替换二维循环？

问题描述

我有一个我不想编辑的功能。一些输入是时间相关的（这里显示为向量），一些是时间相关的，还取决于另一个变量 Nj。

我目前正在循环时间 (Ni) 并循环 Nj 并单独计算每个值。据我所知，apply 系列函数仅适用于所有输入具有相同维度的场景。还有另一种方法可以做到这一点吗？

Ni <- 10
Nj <- 10

a <- matrix(1:100/100, Ni, Nj)
b <- matrix(runif(100)*500, Ni, Nj)
c <- runif(Ni)
d <- c + runif(Ni)
e <- runif(1)*100
f <- c(0.3, 0.7)

funky <- function(a, b, c, d, e, f) {

  firstLine <- a / b
  secondLine <- firstLine * c
  thirdLine <- (secondLine + 45) / d
  fourthLine <- thirdLine + e
  result <- c(f[1] * fourthLine, f[2] * fourthLine)

  result

}

resultMatrix1 <- matrix(numeric(), Ni, Nj)
resultMatrix2 <- matrix(numeric(), Ni, Nj)

for (i in 1:Ni) {
  for (j in 1:Nj) {

    result <- funky(a[i, j],
                    b[i, j],
                    c[i],
                    d[i],
                    e, 
                    f
    )

    resultMatrix1[i, j] <- result[1]
    resultMatrix2[i, j] <- result[2]
  }
}

这是我刚刚拼凑的一些代码，它显示了我对输入维度的意思。问题是我使用的实际功能不是很快，我填充的实际结果网格大约是 100*150，运行大约需要半小时。

标签： r

解决方案

您可以矢量化其中一个循环。我将for在i.

set.seed(1234)

Ni <- 10
Nj <- 10

a <- matrix(1:100/100, Ni, Nj)
b <- matrix(runif(100)*500, Ni, Nj)
c <- runif(Ni)
d <- c + runif(Ni)
e <- runif(1)*100
f <- c(0.3, 0.7)

funky <- function(a, b, c, d, e, f) {

  firstLine <- a / b
  secondLine <- firstLine * c
  thirdLine <- (secondLine + 45) / d
  fourthLine <- thirdLine + e
  result <- c(f[1] * fourthLine, f[2] * fourthLine)

  result

}

resultMatrix1 <- matrix(numeric(), Ni, Nj)
resultMatrix2 <- matrix(numeric(), Ni, Nj)

for (i in 1:Ni) {
  for (j in 1:Nj) {

    result <- funky(a[i, j],
                    b[i, j],
                    c[i],
                    d[i],
                    e, 
                    f
    )

    resultMatrix1[i, j] <- result[1]
    resultMatrix2[i, j] <- result[2]
  }
}

resultMatrix1
resultMatrix2


funky2 <- function(a, b, c, d, e, f) {

  firstLine <- a / b
  secondLine <- firstLine * c
  thirdLine <- (secondLine + 45) / d
  fourthLine <- thirdLine + e
  result <- matrix(c(f[1] * fourthLine, f[2] * fourthLine), ncol = 2)

  result
}

rmat <- array(NA, dim = c(2, Ni, Nj))
for(j in 1:Nj) {
  result <- funky2(a[, j], b[, j], c, d, e, f)
  rmat[1, , j] <- result[, 1]
  rmat[2, , j] <- result[, 2]
}


identical(resultMatrix1, rmat[1, , ])
#[1] TRUE
identical(resultMatrix2, rmat[2, , ])
#[1] TRUE


JoePye <- function(a, b, c, d, e, f, Ni, Nj){
  resultMatrix1 <- matrix(numeric(), Ni, Nj)
  resultMatrix2 <- matrix(numeric(), Ni, Nj)

  for (i in 1:Ni) {
    for (j in 1:Nj) {

      result <- funky(a[i, j],
                      b[i, j],
                      c[i],
                      d[i],
                      e, 
                      f
      )

      resultMatrix1[i, j] <- result[1]
      resultMatrix2[i, j] <- result[2]
    }
  }

  list(Mat1 = resultMatrix1,
       Mat2 = resultMatrix2)
}
RuiB <- function(a, b, c, d, e, f, Ni, Nj){
  rmat <- array(NA, dim = c(2, Ni, Nj))
  for(j in 1:Nj) {
    result <- funky2(a[, j], b[, j], c, d, e, f)
    rmat[1, , j] <- result[, 1]
    rmat[2, , j] <- result[, 2]
  }
  rmat
}

library(microbenchmark)

mb10 <- microbenchmark(
  JoePye = JoePye(a, b, c, d, e, f, Ni, Nj),
  RuiB = RuiB(a, b, c, d, e, f, Ni, Nj)
)
mb10
#Unit: microseconds
#   expr     min       lq     mean   median       uq     max neval cld
# JoePye 473.950 479.2840 496.5165 483.6325 495.7015 758.162   100   b
#   RuiB 150.502 157.8205 174.8730 160.1100 165.2400 647.653   100  a

这是一个 3 倍的加速。现在用更大的输入进行测试。

set.seed(1234)

Ni <- 1e2
Nj <- 2e1

a <- matrix(seq.int(Ni*Nj)/100, Ni, Nj)
b <- matrix(runif(Ni*Nj)*500, Ni, Nj)
c <- runif(Ni)
d <- c + runif(Ni)
e <- runif(1)*100
f <- c(0.3, 0.7)

res1 <- JoePye(a, b, c, d, e, f, Ni, Nj)
res2 <- RuiB(a, b, c, d, e, f, Ni, Nj)

identical(res1[[1]], res2[1, , ])
#[1] TRUE
identical(res1[[2]], res2[2, , ])
#[1] TRUE


mb100 <- microbenchmark(
  JoePye = JoePye(a, b, c, d, e, f, Ni, Nj),
  RuiB = RuiB(a, b, c, d, e, f, Ni, Nj),
  times = 10
)
mb100
#Unit: microseconds
#   expr      min       lq      mean   median       uq       max neval cld
# JoePye 9198.846 9248.114 9421.2359 9352.244 9426.161 10147.642    10   b
#   RuiB  478.564  490.404  533.8198  522.573  594.841   602.938    10  a

随着更大的输入，加速因子上升到 18 倍。

推荐阅读