R 中使用 set.seed 并行 foreach 共享内存

Question

正如这个问题中所见，在Windows中不可能在R中使用共享内存运行并行进程。因此，我设计了以下方法，使用一系列

set.seed()

来模拟之间的“共享内存”我的核心。

我想模拟一系列随机变量矩阵，其中：

1. 线条代表贷款的固定数量。
2. 各列代表 Monte Calro 模拟的可变数量。

nruns <- 4
nsim.runs <- c(5000, 10000)
nport.total <- 10
n.sectors <- 2
K <- n.sectors

ISCorrMatrix <- matrix(data=0.5, nrow=n.sectors, ncol=n.sectors, byrow=TRUE)
diag(ISCorrMatrix) <- 1

Alpha.Matrix <- t(chol(ISCorrMatrix))

Alpha.Matrix.extended <- NULL

for (i in 1:n.sectors) {
  new <- matrix(data=Alpha.Matrix[i, ], nrow=nport.total/2, ncol=n.sectors, 
                byrow=TRUE)
  Alpha.Matrix.extended <- rbind(Alpha.Matrix.extended, new)
}

beta.star <- runif(n.sectors, min=0.1, max=0.9)  ## so just draw them at random.
beta.star <- beta.star/sqrt(sum(beta.star^2))  ## to normalise the draw.

sum(beta.star^2)

## Initialize parallel backend
num_cores <- length(nsim.runs)
library(doParallel)
cl <- makeCluster(num_cores)
registerDoParallel(cl)

## Set up a reproducible RNG stream
clusterSetRNGStream(cl)

parallel_results <- foreach(core=1:length(nsim.runs), .combine="c") %dopar% {
  nsim <- nsim.runs[core]
  results <- list()
  set.seed(123)
  for (run in 1:nruns) { 
    Z.systematic <- matrix(data=rnorm(nsim*K), nrow=K, ncol=nsim, byrow=FALSE)
    Y.systematic <- Alpha.Matrix %*% Z.systematic
    Y.systematic.star <- beta.star %*% Z.systematic
    Y.systematic.extended <- c()
    for (i in 1:n.sectors) {
      new <- matrix(data=rep(Y.systematic[i, ], nport.total/2), 
                    nrow=nport.total/2, ncol=nsim, byrow=TRUE)
      Y.systematic.extended <- rbind(Y.systematic.extended, new)
    }
    Y.systematic.extended.star <- 
      matrix(data=Y.systematic.star, nrow=nport.total, ncol=nsim, byrow=TRUE)
    set.seed(123 + run*100000) 
    Z.Default <- matrix(data=rnorm(nsim*nport.total), nrow=nport.total,
                        ncol=nsim, byrow=FALSE) 
    set.seed(123 + run*100000)
    LGD.Basel <- matrix(data=rbeta(n=nsim*nport.total, shape1=2, shape2=3),
                        nrow=nport.total, ncol=nsim, byrow=FALSE)
    results[[paste0("run_", run, "_Nsim_", nsim)]] <- 
      list(Y.systematic=Y.systematic.extended.star, 
           Z.Default=Z.Default, 
           LGD=LGD.Basel)
  }
  return(results)  ## The return call of the parallelisation.
}

## Stop the parallel backend
stopCluster(cl)

每个核心分配了不同数量的蒙特卡洛模拟（我知道这不是最佳的，但这是不爆炸 RAM 的唯一方法），从较小的数量到较大的数量，例如 {5000, 10000} .

为了避免模拟噪声，我希望第二个核心（使用 10000 nsim 的那个）使用与前一个核心中模拟的相同的“第一个”5000 RV，依此类推。

同时，如果我重新运行相同的 5000 次模拟，我希望在同一个核心内实现这一点，以确保模拟值之间的可变性。

我的代码似乎实现了这一点（这里我只显示 Y 和 Z，但也适用于 LGD）：

跨 nsim 的相同运行的稳定性：Y

> parallel_results$run_1_Nsim_5000$Y.systematic[1:2,4998:5000]
           [,1]     [,2]      [,3]
[1,] -0.2397531 1.548914 0.2781637
[2,] -0.2397531 1.548914 0.2781637
> parallel_results$run_1_Nsim_10000$Y.systematic[1:2,4998:5002]
           [,1]     [,2]      [,3]    [,4]       [,5]
[1,] -0.2397531 1.548914 0.2781637 2.23994 -0.6209843
[2,] -0.2397531 1.548914 0.2781637 2.23994 -0.6209843

同一 nsim 上不同运行的差异：Y

> parallel_results$run_1_Nsim_5000$Y.systematic[1:2,1:2]
           [,1]      [,2]
[1,] -0.3996827 0.2019272
[2,] -0.3996827 0.2019272
> parallel_results$run_2_Nsim_5000$Y.systematic[1:2,1:2]
           [,1]        [,2]
[1,] -0.4580132 -0.08322681
[2,] -0.4580132 -0.08322681

跨 nsim 的相同运行的稳定性：Z

> parallel_results$run_1_Nsim_5000$Z.Default[1:2,4998:5000]
          [,1]         [,2]       [,3]
[1,] -1.077092 -0.006225626 -0.6152041
[2,]  2.033250 -1.350659878 -0.7040274
> parallel_results$run_1_Nsim_10000$Z.Default[1:2,4998:5002]
          [,1]         [,2]       [,3]        [,4]      [,5]
[1,] -1.077092 -0.006225626 -0.6152041  0.57472772  2.593232
[2,]  2.033250 -1.350659878 -0.7040274 -0.08696669 -1.456689

同一 nsim 上不同运行的差异：Z

> parallel_results$run_1_Nsim_5000$Z.Default[1:2,1:2]
              [,1]       [,2]
[1,]  0.0008058532 -1.0539683
[2,] -1.3029004815 -0.2031874
> parallel_results$run_2_Nsim_5000$Z.Default[1:2,1:2]
            [,1]      [,2]
[1,]  0.50913790 -1.120769
[2,] -0.03775884  0.894895

这些抽签是否真正独立（我已经检查了汇总统计数据和更正，它们似乎是......）？我是否遗漏了与 R 中的 RNG（及其与并行包的交互）相关的注意事项？我在堆栈中进行了广泛的搜索，但无法找到符合我的应用程序的东西。任何帮助或评论将不胜感激。

Answer 1

不确定我是否理解正确，但你正在做这样的事情，

> library(doParallel)
> cl <- makeCluster(num_cores)
> registerDoParallel(cl)
> 
> foreach(core=1:3, .combine="list") %dopar% {
+   rn <- array(, c(2, 2))
+   for (run in seq_len(nrow(rn))) {
+     set.seed(41 + run)
+     rn[, run] <- rnorm(2)
+   }
+   rn
+ }
[[1]]
[[1]][[1]]
           [,1]        [,2]
[1,]  1.3709584 -0.03751376
[2,] -0.5646982 -1.57460441

[[1]][[2]]
           [,1]        [,2]
[1,]  1.3709584 -0.03751376
[2,] -0.5646982 -1.57460441


[[2]]
           [,1]        [,2]
[1,]  1.3709584 -0.03751376
[2,] -0.5646982 -1.57460441

产生相似的值。如果这不是您的意图，您可能想要添加迭代

i

，尽管这并不理想，因为随机种子通常实际上是伪随机的。

> 
> stopCluster(cl)
> 
> library(parallel)
> cl <- makeCluster(num_cores)
> 
> parLapply(cl, 1:3, \(i) {
+   rn <- array(, c(2, 2))
+   for (run in seq_len(nrow(rn))) {
+     set.seed(41 + run + i)
+     rn[, run] <- rnorm(2)
+   }
+   rn}
+ )
[[1]]
            [,1]       [,2]
[1,] -0.03751376 0.65391826
[2,] -1.57460441 0.01905227

[[2]]
           [,1]       [,2]
[1,] 0.65391826  0.3407997
[2,] 0.01905227 -0.7033403

[[3]]
           [,1]       [,2]
[1,]  0.3407997 -0.8989165
[2,] -0.7033403  0.2121311

> 
> stopCluster(cl)