并行嵌套的 foreach 循环

parallel nested foreach loops

我正在尝试为 Metropolis-Hastings 算法编写 nested 并行 foreach 循环,但矩阵组合不正确。示例代码如下,最终矩阵 mtx2 的尺寸应与原始矩阵 mtx 的尺寸相同,但有些行随机更改。矩阵行应该如何组合?

我直接尝试了 foreach 包,但结果相同 - mtx2 将列组合了 5 次。

# library(doParallel)
library(foreach)

no_cores <- detectCores() - 2  
cl <- makeCluster(no_cores)  
registerDoParallel(cl)  

mtx <- matrix(data=rnorm(n=1e3*5,mean=0,sd=1),nrow=1e3,ncol=5)
mtx2 <- matrix(data=NA,nrow=1e3,ncol=5)

#basic for loop - slow for large number of rows
for(k in 1:nrow(mtx)){
  for(r in 1:5) {
    if(runif(n=1,min=0,max=1)>0.9){
      mtx2[k,] <- mtx[k,]*10
    }else{
      mtx2[k,] <- mtx[k,]
    }  
  }
}

#series version for de-bugging
mtx2 <-foreach(k=1:nrow(mtx),.combine="rbind") %do% {
  foreach(r=1:5,.combine="c") %do% {
    if(runif(n=1,min=0,max=1)>0.9){
      mtx[k,]*10
    }else{
      mtx[k,]
    }  
  }
}

#parallel version
mtx2 <-foreach(k=1:nrow(mtx),.combine="rbind") %:% {
  foreach(r=1:5,.combine="c") %dopar% {
    if(runif(n=1,min=0,max=1)>0.9){
      mtx[k,]*10
    }else{
      mtx[k,]
    }  
  }
}

mtx2 <- round(mtx2,2)

要扩展评论,您可以通过一次创建所有逻辑比较来跳过循环。在这里,我们创建 runif(nrow(mtx) * ncol(mtx)) 但只取每 5 个结果来匹配 for (r in 1:5) {...}

的 OP 内循环

关键在于,虽然 OP 问题无法找到在嵌套并行循环中更新矩阵的方法,但重构代码有时可以提供显着的性能提升。

nr = 1e4
nc = 5
mtx <- matrix(data=rnorm(n=nr*nc,mean=0,sd=1),nrow=nr,ncol=nc)

set.seed(123L)
lgl = matrix(runif(n = nr * nc), ncol = nc, byrow = TRUE)[, nc] > 0.9
mtx3 = sweep(mtx, 1L, 1 + 9 * lgl, FUN = '*')

all.equal(mtx2, mtx3) ##mtx2 was created with set.seed(123L)

# [1] TRUE

对于 100 万行,这要快得多:

system.time({
  lgl = matrix(runif(n = nr * nc), ncol = nc, byrow = TRUE)[, nc] > 0.9
  mtx3 = sweep(mtx, 1L, 1 + 9 * lgl, FUN = '*')
})

##    user  system elapsed 
##    0.27    0.00    0.27 

system.time({
  for(k in 1:nrow(mtx)){
    for(r in 1:5) {
      if(runif(n=1,min=0,max=1)>0.9){
        mtx2[k,] <- mtx[k,]*10
      }else{
        mtx2[k,] <- mtx[k,]
      }  
    }
  }
})

##    user  system elapsed 
##   14.09    0.03   14.12