与 R 并行置换矩阵列

Question

我正在尝试使用 R 按列排列矩阵。 但是，它需要很长时间（矩阵是 68k x 32k 整数）。

我想并行执行（因为每一列都是独立排列的）。我怎样才能实现它？它应该与 R 中的 embarrassingly parallel for 有关，但我没有找到解决方案。

目前，我的功能如下：

permMTX <- function(x) {
    nr <- nrow(x)
    nc <- ncol(x)
    # I'd like to parallelize this for, since each
    # column can be permuted independently
    for (i in 1:nc) {
        x[,i] <- x[sample(nr),i]
    }
    x
} 

Answer 1

免责声明：我是 bigstatsr 包的作者。

您可以使用共享内存（存储在磁盘上的矩阵）并执行以下操作：

# devtools::install_github("privefl/bigstatsr")
library(bigstatsr)

# matrix on disk
mat <- FBM(68e3, 32e2, backingfile = "test")
# inialize with 1:nrow(mat) for each column
system.time(
  big_apply(mat, a.FUN = function(X, ind) {
    print(min(ind))
    X[, ind] <- rep(rows_along(X), length(ind))
    NULL
  }, a.combine = 'c')
) # 15 sec

# permute each column, in parallel
system.time(
  big_apply(mat, a.FUN = function(X, ind) {
    print(min(ind))
    X[, ind] <- apply(X[, ind], 2, sample)
    NULL
  }, a.combine = 'c', ncores = nb_cores())
) # 27 sec

十分之一的数据需要 27 秒，整个数据集需要 378 秒（在只有 2 个物理内核和 8GB RAM 的笔记本电脑上）。

Answer 2

解决方案

首先，我会使用矢量化，这应该会提高效率。

permMTX = function(x) apply(x, 2L, sample)

然后我们可以使用库parallel 来并行化该函数：

library(parallel)

parPermMTX = function(x, cluster) parApply(cl = cluster, X = x, MARGIN = 2L, FUN = sample)

---

用法

使用parallel，您必须在使用前注册一个集群。这是一个例子：

cl = makeCluster(detectCores(logical = FALSE))
parPermMTX(diag(10), cl)
#     [,1] [,2] [,3] [,4] [,5]
#[1,]    0    1    0    0    0
#[2,]    0    0    0    0    0
#[3,]    0    0    0    0    0
#[4,]    1    0    0    1    1
#[5,]    0    0    1    0    0

parallel 的工作方式（产生多个 R 进程）您必须确保您有足够的内存来容纳数据的多个副本。

我认为建议也将数据导出到进程中，您只需调用即可做到这一点

clusterExport(cl, varlist = "exampleData")

虽然它确实在我端并行运行，但它根本不比简单地使用 apply 快，但我无法使用与你相同规模的数据进行测试，所以我不能确定它会起作用的。

这是因为sample 已经进行了高度优化，因此生成进程的开销比简单地调用sample 更大。见Why is the parallel package slower than just using apply?

在我的系统上，对 68E3 个整数进行 32E3 次采样大约需要 40 秒：

microbenchmark(sample(68E3), times = 32E3)
#Unit: milliseconds
#          expr      min       lq     mean   median       uq      max neval
# sample(68000) 1.132273 1.192923 1.290838 1.227912 1.286229 7.880191 32000

也许你的内存用完了，并且使用了硬盘缓存，这真的很慢。

---

第二个解决方案

那么，如果我们尝试将尽可能多的对sample 的调用顺序分配给单个进程呢？这是我在这里尝试过的：

parPermMTX2 = function(x, cluster) do.call(cbind, parLapply(cl = cluster, X = list(x[,seq(floor(ncol(x)/2))], x[,floor(ncol(x)/2)+seq(ceiling(ncol(x)/2))]), fun = permMTX))

我们将x分成两半，然后分别调用permMTX，然后与cbind重组。

遗憾的是，这样我也无法获得更好的性能。所以，虽然我回答了你的问题，但我不确定它是否有任何帮助。