在大型数组中使用 apply 的更快方法（可能涉及并行化）？

Question

我在 R 中有一个大数组，比如说，

A = array(rnorm(100*100*30*30*100),dim=c(100,100,30,30,100))

我想找到第 4 维中小于 1 的第一个条目：

first_entry = apply(A,c(1,2,3,5),function(x) min(which(x < 1)))

因为数组太大，所以这需要一段时间，在我的计算机上大约需要 20 秒，并且我将来会为很多数组执行此操作。 for 循环的速度大致相同。我也尝试过使用 parApply，但它需要大约相同的时间，甚至更长。可能我的函数不够“复杂”，无法通过并行化来实现速度增益。有没有更快的方法来做到这一点？实际上，我希望能够理想地将另一个数组 B（与 A 具有相同维度）中的值设置为 0。所以，类似，

B[first_entry] = 0

请注意，鉴于“apply”（上面）的当前输出，这不起作用，因为first_entry的尺寸为100x100x30x100。

Answer 1

我稍微简化了你的例子（我的机器可能比你的慢！），并添加了一些时间测量。当我使用 parApplym 而不是 applym 时，代码会快一些（parApply 为 2.4 秒，apply 为 3.8 秒）。

如果 parApply 确实在您的计算机上产生类似的结果，您是否考虑以多线程方式运行更广泛的任务？您提到您将对很多数组执行此操作。因此，您是否可以同时运行多个查找，而不是加快查找过程？

# Create data
A = array(rnorm(100*100*30),dim=c(100,100,30))

# Using apply
start.time <- Sys.time()
first_entry = apply(A,c(1,2,3),function(x) min(which(x < 1)))
print(Sys.time() - start.time)

# Load packages
library(parallel)

# Create cluster and export
nrCores <- detectCores()
cl <- makeCluster(nrCores)
clusterExport(cl=cl, varlist=c("A"))

# Using parApply
start.time <- Sys.time()
first_entry = parApply(cl=cl,A,c(1,2,3),function(x) min(which(x < 1)))
print(Sys.time() - start.time)

Answer 2

在第四维上循环将允许提前打破。如果有很大比例的值符合您正在检查的条件（如示例所示），这将提供显着的加速。

旁注：要获得第一种情况，

which.max(x < 1)

会比

which(x < 1)[1]

更快；但是，如果没有一个

x

小于 1，则两者都会导致问题。下面的方法可以避免这些问题。

A <- array(rnorm(100*100*30*30*100),dim=c(100,100,30,30,100))

system.time(first_entry <- apply(A, c(1:3,5), \(x) which.max(x < 1)))
#>    user  system elapsed 
#>   83.48    2.89   86.43
system.time({
  d <- dim(A)
  i <- collapse::alloc(1L, prod(d[-4])) # or `rep_len(1L, prod(d[-4]))`
  j <- which(A[,,,1,] >= 1)
  for (k in 2:d[4]) {
    b <- A[,,,k,][j] < 1
    i[j[b]] <- k
    j <- j[!b]
    if (!length(j)) break
  }
  i[j] <- d[4] + 1L # in case the condition is not met
  dim(i) <- d[-4]
})
#>    user  system elapsed 
#>    2.25    0.63    2.89

检查这两种方法是否等效。

identical(first_entry, i)
#> [1] TRUE

根据我们的索引更新数组

B

B <- array(rnorm(100*100*30*30*100),dim=c(100,100,30,30,100))

system.time({
  s <- split(1:length(i), i)
  for (k in 1:min(length(s), d[4])) B[,,,k,][s[[k]]] <- 0
})
#>    user  system elapsed 
#>    6.39    3.46    9.86

或者，我们可以直接更新

B

而不计算

i

。

system.time({
  d <- dim(A)
  b <- A[,,,1,] < 1
  B[,,,1,][b] <- 0
  j <- which(!b)
  for (k in 2:d[4]) {
    b <- A[,,,k,][j] < 1
    B[,,,k,][j[b]] <- 0
    j <- j[!b]
    if (!length(j)) break
  }
})
#>    user  system elapsed 
#>    6.63    0.93    7.57