使设置随机种子的函数独立

Making functions that set the random seed independent

有时我想编写一个随机函数,该函数始终 return 为特定输入提供相同的输出。我总是通过在函数顶部设置随机种子然后继续执行来实现这一点。考虑以这种方式定义的两个函数:

sample.12 <- function(size) {
  set.seed(144)
  sample(1:2, size, replace=TRUE)
}
rand.prod <- function(x) {
  set.seed(144)
  runif(length(x)) * x
}

sample.12 returns从集合{1, 2}rand.prod中随机抽取的指定大小的向量将指定向量的每个元素乘以均匀选择的随机值来自 [0, 1]。通常我希望 x <- sample.12(10000) ; rand.prod(x) 有一个 "step" 分布,pdf 3/4 在 [0, 1] 范围内,1/4 在 [1, 2] 范围内,但由于我不幸的选择上面相同的随机种子我看到了不同的结果:

x <- sample.12(10000)
hist(rand.prod(x))

在这种情况下,我可以通过将其中一个函数中的随机种子更改为其他值来解决此问题。例如,在 rand.prod 中使用 set.seed(10000) 我得到预期的分布:

Previously on SO this solution of using different seeds has been accepted as the best approach to generate independent random number streams. However, I find the solution to be unsatisfying because streams with different seeds could be related to one another (possibly even highly related to one another);事实上,根据 ?set.seed:

,它们甚至可能产生相同的流

There is no guarantee that different values of seed will seed the RNG differently, although any exceptions would be extremely rare.

有没有办法在 R 中实现一对随机函数:

  1. 始终return特定输入的相同输出,并且
  2. 通过不仅仅是使用不同的随机种子来加强随机源之间的独立性?

我对此进行了深入研究,看起来 rlecuyer 包提供了独立的随机流:

Provides an interface to the C implementation of the random number generator with multiple independent streams developed by L'Ecuyer et al (2002). The main purpose of this package is to enable the use of this random number generator in parallel R applications.

第一步是独立流的全局初始化:

library(rlecuyer)
.lec.CreateStream(c("stream.12", "stream.prod"))

然后需要修改每个函数以将适当的流重置为其开始状态 (.lec.RestartStartStream),将 R 随机数生成器设置为适当的流 (.lec.CurrentStream),然后设置R 随机数生成器回到调用函数之前的状态 (.lec.CurrentStreamEnd).

sample.12 <- function(size) {
  .lec.ResetStartStream("stream.12")
  .lec.CurrentStream("stream.12")
  x <- sample(1:2, size, replace=TRUE)
  .lec.CurrentStreamEnd()
  x
}
rand.prod <- function(x) {
  .lec.ResetStartStream("stream.prod")
  .lec.CurrentStream("stream.prod")
  y <- runif(length(x)) * x
  .lec.CurrentStreamEnd()
  y
}

这满足"always returns the same output given the same input"要求:

all.equal(rand.prod(sample.12(10000)), rand.prod(sample.12(10000)))
# [1] TRUE

在我们的示例中,流似乎也是独立运行的:

x <- sample.12(10000)
hist(rand.prod(x))

请注意,这不会在我们的脚本运行中给出一致的值,因为每次调用 .lec.CreateStream 都会给出不同的初始状态。为了解决这个问题,我们可以记录每个流的初始状态:

.lec.GetState("stream.12")
# [1] 3161578179 1307260052 2724279262 1101690876 1009565594  836476762
.lec.GetState("stream.prod")
# [1]  596094074 2279636413 3050913596 1739649456 2368706608 3058697049

然后我们可以将脚本开头的流初始化更改为:

library(rlecuyer)
.lec.CreateStream(c("stream.12", "stream.prod"))
.lec.SetSeed("stream.12", c(3161578179, 1307260052, 2724279262, 1101690876, 1009565594, 836476762))
.lec.SetSeed("stream.prod", c(596094074, 2279636413, 3050913596, 1739649456, 2368706608, 3058697049))

现在对 sample.12rand.prod 的调用将匹配对脚本的调用。