使用 foreach 循环和并行处理生成矩阵

Question

我正在尝试转换我当前用于跨大型矩阵运行进程的 for 循环。当前的 for 循环在 30 x 30 的部分中找到最大值并创建一个具有最大值的新矩阵。

for 循环的当前代码如下所示：

mat <- as.matrix(CHM) # CHM is the original raster image
maxm <- matrix(nrow=nrow(mat)/30, ncol=ncol(mat)/30) # create new matrix with new dimensions

for(i in 1:dim(maxm)[1]) {
  for(j in 1:dim(maxm)[2]) {
    row <- 30 * (i - 1) + 1
    col <- 30 * (j - 1) + 1
    maxm[i,j] <- max(CHM[row:(row + 29), col:(col + 29)])
   }
 }

我想将其转换为 foreach 循环以使用并行处理。我已经制作了以下代码，但是这很有效。我不确定如何在 foreach 循环中生成新矩阵：

ro<-nrow(mat)/30
co<-ncol(mat)/30
maxm <- matrix(nrow=nrow(mat)/30, ncol=ncol(mat)/30)

foreach(i=ro, .combine='cbind') %:%
  foreach(j=co, .combine='c') %dopar% {
    row <- 30 * (i - 1) + 1
    col <- 30 * (j - 1) + 1
    maxm[i,j]<-(max(CHM[row:(row + 29), col:(col + 29)]))
  }

请给点建议！

Answer 1

让我试着在这里得到答案。

据我所知，R 使用集群系统进行并行计算，每个节点都使用自己的环境。因此，foreach-%dopar%，首先，将所有当前的 .globalEnv 复制到每个集群节点，然后尝试运行循环体中编写的代码。代码执行后没有备份。您只会得到result = foreach(...) { } 的结果。因此，每个节点中的代码 maxm[i,j]<-(max(CHM[row:(row + 29), col:(col + 29)])) 仅更改矩阵的本地副本，仅此而已。 因此，“正确”的代码可能是这样的：

mat <- as.matrix(CHM);
ro<-nrow(mat)/30;
co<-ncol(mat)/30;

maxm = foreach(i=1:ro, .combine='cbind') %:% 
{ 
   result = foreach(j = 1:co, .combine='c') %dopar% 
            { 
                row <- 30 * (i - 1) + 1; 
                col <- 30 * (j - 1) + 1; 
                max(CHM[row:(row + 29), col:(col + 29)]); 
            } 
   result; 
} 

也许还需要使用as.matrix 来表示maxm。

Answer 2

在并行执行任何操作之前，应该尝试查看是否可以进行任何矢量化。一旦完成，问题“并行化是否合理？”

在这个特定示例中，并行化不太可能像您期望的那样快，因为在每次迭代中您都将输出保存到一个公共对象中。 R 在并行化中通常不支持这一点，相反，人们应该在所谓的“令人尴尬的可并行化”问题中寻求并行化，直到人们更好地理解并行问题的工作原理。简而言之：除非您知道自己在做什么，否则不要对 R 中的数据执行并行更改。它不太可能更快。

也就是说，在您的情况下，它实际上变得非常棘手。您似乎正在执行“滚动最大窗口”，并且输出应保存在组合矩阵中。将数据直接保存到其他矩阵的另一种方法是返回一个包含 3 列的矩阵 x、i、j，其中后两个是指示 x 值的索引应该放在里面。

正如 Dmitriy 在他的回答中指出的那样，为了使其正常工作，需要将数据导出到每个 cluster（并行会话），以便我们可以使用它。然后下面的例子展示了如何执行并行化

首先：创建集群并导出数据集

set.seed(1)
#Generate test example
n <- 3000
dat <- matrix(runif(n^2), ncol = n)
library(foreach)
library(doParallel)
#Create cluster
cl <- parallel::makeCluster(parallel::detectCores())
#Register it for the foreach loop
doParallel::registerDoParallel(cl)
#Export the dataset (could be done directly in the foreach, but this is more explicit)
parallel::clusterExport(cl, "dat")

接下来我们来到foreach 循环。请注意，根据文档，嵌套的foreach 循环应使用%:% 标签分隔，如下面的示例所示：

output <- foreach(i = 1:(nrow(dat)/30), .combine = rbind, .inorder = FALSE) %:% 
    foreach(j = 1:(ncol(dat)/30), .combine = rbind, .inorder = FALSE) %dopar%{
        row <- 30 * (i - 1) + 1
        col <- 30 * (j - 1) + 1
        c(x = max(dat[row:(row + 29), col:(col + 29)]), i = i, j = j)
    }

注意.inorder = FALSE。当我返回索引时，我不关心顺序，只关心速度。 最后但同样重要的是，我们需要创建矩阵。 Matrix 包函数Matrix::SparseMatrix 允许指定值和索引。

output <- Matrix::sparseMatrix(output[,"i"], output[,"j"], x = output[,"x"])

这仍然很慢。对于n = 3000，执行计算大约需要 6 秒 + 导出数据的开销并不小。但它可能比使用顺序循环的相同方法更快。