我需要估计一个过渡矩阵。由于我有很多数据,我尝试使用 foreach
并行运行它,并尝试了bigstatsr
的共享内存函数FBM()
。而且该函数似乎并不总是返回正确的结果。(有时确实如此。可能是该功能无法正常工作吗?
下面是代码正常工作时的示例:
x <- c(1,2,1,1,3,4,4,1,2,4,1,4,3,4,4,4,3,1,3,2,3,3,3,4,2,2,3)
n <- length(unique(x))
A <- matrix(nrow = n, ncol = n, 0)
for (t in 1:(length(x) - 1)) {A[x[t], x[t + 1]] <- A[x[t], x[t + 1]] + 1}
A
以下是并不总是返回正确结果的代码:
library(foreach)
library(doParallel)
library(bigstatsr)
cl <- makeCluster(8)
registerDoParallel(cl)
B <- FBM(n, n)
set.seed(3)
foreach (t = 1:(length(x) - 1)) %dopar% {B[x[t], x[t + 1]] <- B[x[t], x[t + 1]] + 1}
stopCluster(cl)
B[]
identical(A,B[])
使用库时也会发生同样的情况snow
library(snow)
library(bigstatsr)
cl <- makeCluster(8)
f.trans.m <- function(t) {
D[x[t], x[t + 1]] <<- D[x[t], x[t + 1]] + 1
}
D <- FBM(n, n)
clusterExport(cl, "f.trans.m")
clusterExport(cl, "D")
clusterExport(cl, "x")
parLapply(cl, seq(1,(length(x) - 1)), function(t) f.trans.m(t))
D[]
identical(A,D[])
我是否正确使用了该软件包,还是FBM()
中存在错误?
一个解决方案:
缺少包flock
提供的文件锁。
B <- FBM(n, n)
lock <- tempfile()
foreach (t = 1:(length(x) - 1)) %dopar% {
locked <- flock::lock(lock)
B[x[t], x[t + 1]] <- B[x[t], x[t + 1]] + 1
flock::unlock(locked)
}
对于这个特定示例,问题在于值的并行并发更新(参见 https://privefl.github.io/blog/a-guide-to-parallelism-in-r/#advanced-parallelism-synchronization(。
在这里,我根本不会使用并行性。我宁愿选择顺序(但矢量化(访问器。
我将首先重新组合索引以递增:
library(dplyr)
ind <- data.frame(i = x[-length(x)], j = x[-1]) %>%
group_by(i, j) %>%
count()
然后,我将使用两列矩阵访问器来更新相应的值,而无需使用 R 循环。
B <- FBM(n, n, init = 0)
ind2 <- as.matrix(ind[1:2])
B[ind2] <- B[ind2] + ind[[3]]