制作独立设置随机种子的函数

13
有时候我需要编写一个随机函数,它对于特定的输入总是返回相同的输出。我一直通过在函数顶部设置随机种子然后继续实现这个功能。考虑以下两个以此方式定义的函数:
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会随机从集合{1, 2}中选取指定数量的元素并返回一个向量,而rand.prod会将给定向量中的每个元素乘以一个在区间[0, 1]内均匀分布的随机数。通常我期望x <- sample.12(10000) ; rand.prod(x)的结果呈现出“阶梯”分布,其中区间[0, 1]的概率密度函数为3/4,区间[1, 2]的概率密度函数为1/4。但由于上述不幸的相同随机种子选择,我看到了不同的结果:
x <- sample.12(10000)
hist(rand.prod(x))

输入图像描述

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

输入图像描述

之前在SO上,使用不同的种子被接受为生成独立随机数流的最佳方法。然而,我发现这个解决方案并不令人满意,因为具有不同种子的流可能与彼此相关(甚至可能是高度相关的);事实上,根据?set.seed,它们甚至可能产生相同的流:

没有保证不同的种子值将不同地播种RNG,尽管任何例外都极为罕见。

是否有一种方法在R中实现一对随机化函数:

  1. 始终为特定输入返回相同的输出,并且
  2. 通过不仅仅使用不同的随机种子来强制执行它们的随机性源之间的独立性?
1个回答

9
我已经进一步研究了这个问题,看起来rlecuyer包提供了独立的随机流:

提供了一个接口,用于访问L'Ecuyer等人(2002)开发的具有多个独立流的随机数生成器的C实现。此软件包的主要目的是使该随机数生成器能够在并行R应用程序中使用。

第一步是全局初始化独立流:
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
}

这符合“相同的输入始终返回相同的输出”的要求:
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的调用将在脚本的不同调用之间匹配。

1
很好的发现。为了完整起见,rlecuyer使用MRG32k3a此论文第1.1章),因此它可能也有其局限性(就像Mersenne-Twister一样)。不过在99%的情况下应该不会有太大问题。 - tonytonov

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接