8

我们可以使用set.seed()在 R 中设置随机种子,这具有全局效果。这是一个最小的例子来说明我的目标:

set.seed(0)
runif(1)
# [1] 0.8966972

set.seed(0)
f <- function() {
  # I do not want this random number to be affected by the global seed
  runif(1)
}
f()
# [1] 0.8966972

基本上我希望能够避免全局随机种子(即.Random.seed)在本地环境中的影响,例如 R 函数,这样我就可以实现某种用户无法控制的随机性。例如,即使用户有set.seed(),他每次调用这个函数时仍然会得到不同的输出。

现在有两种实现方式。第一个依赖于set.seed(NULL)每次我想获得一些随机数时让 R 重新初始化随机种子:

createUniqueId <- function(bytes) {
  withPrivateSeed(
    paste(as.hexmode(sample(256, bytes, replace = TRUE) - 1), collapse = "")
  )
}
withPrivateSeed <- function(expr, seed = NULL) {
  oldSeed <- if (exists('.Random.seed', envir = .GlobalEnv, inherits = FALSE)) {
    get('.Random.seed', envir = .GlobalEnv, inherits = FALSE)
  }
  if (!is.null(oldSeed)) {
    on.exit(assign('.Random.seed', oldSeed, envir = .GlobalEnv), add = TRUE)
  }
  set.seed(seed)
  expr
}

即使我将种子设置为 0,您也可以看到我得到了不同的 id 字符串,并且全局随机数流仍然是可重现的:

> set.seed(0)
> runif(3)
[1] 0.8966972 0.2655087 0.3721239
> createUniqueId(4)
[1] "83a18600"
> runif(3)
[1] 0.5728534 0.9082078 0.2016819

> set.seed(0)
> runif(3)  # same
[1] 0.8966972 0.2655087 0.3721239
> createUniqueId(4)  # different
[1] "77cb3d91"
> runif(3)
[1] 0.5728534 0.9082078 0.2016819

> set.seed(0)
> runif(3)
[1] 0.8966972 0.2655087 0.3721239
> createUniqueId(4)
[1] "c41d61d8"
> runif(3)
[1] 0.5728534 0.9082078 0.2016819

第二个实现可以Github 上找到。比较复杂,基本思路是:

  • set.seed(NULL)使用(in .onLoad())在包启动期间初始化随机种子
  • 将随机种子存储在单独的环境中 ( .globals$ownSeed)
  • 每次我们想要生成随机数时:
    1. 将本地种子分配给全局随机种子
    2. 生成随机数
    3. 将新的全局种子(由于步骤 2 已更改)分配给本地种子
    4. 将全局种子恢复到其原始值

现在我的问题是这两种方法在理论上是否等效。第一种方法的随机性依赖于createUniqueId()调用时的当前时间和进程ID,而第二种方法似乎依赖于加载包时的时间和进程ID。对于第一种方法,是否有可能createUniqueId()在同一个 R 进程中同时发生两个调用,以便它们返回相同的 id 字符串?

更新

在下面的答案中,Robert Krzyzanowski 提供了一些set.seed(NULL)可能导致严重 ID 冲突的经验证据。我为它做了一个简单的可视化

createGlobalUniqueId <- function(bytes) {
  paste(as.hexmode(sample(256, bytes, replace = TRUE) - 1), collapse = "")
}
n <- 10000
length(unique(replicate(n, createGlobalUniqueId(5))))
length(unique(x <- replicate(n, createUniqueId(5))))
# denote duplicated values by 1, and unique ones by 0
png('rng-time.png', width = 4000, height = 400)
par(mar = c(4, 4, .1, .1), xaxs = 'i')
plot(1:n, duplicated(x), type = 'l')
dev.off()

来自 set.seed(NULL) 的随机数

当线到达图的顶部时,这意味着生成了一个重复值。但是,请注意这些重复不是连续出现的,即any(x[-1] == x[-n])通常是FALSE。可能存在与系统时间相关的重复模式。由于我对基于时间的随机种子的工作原理缺乏了解,我无法进一步调查,但您可以在此处此处查看相关的 C 源代码。

4

2 回答 2

7

我认为在你的函数中只有一个独立的 RNG 会很好,它不受全局种子的影响,但会有自己的种子。事实证明,randtoolbox提供了这个功能:

library(randtoolbox)
replicate(3, {
  set.seed(1)
  c(runif(1), WELL(3), runif(1))
})   
#            [,1]      [,2]      [,3]
#[1,] 0.265508663 0.2655087 0.2655087
#[2,] 0.481195594 0.3999952 0.9474923
#[3,] 0.003865934 0.6596869 0.4684255
#[4,] 0.484556709 0.9923884 0.1153879
#[5,] 0.372123900 0.3721239 0.3721239

顶部和底部行受种子影响,而中间行是“真正随机的”。

基于此,这是您的功能的实现:

sample_WELL <- function(n, size=n) {
  findInterval(WELL(size), 0:n/n)
}

createUniqueId_WELL <- function(bytes) {
  paste(as.hexmode(sample_WELL(256, bytes) - 1), collapse = "")
}

length(unique(replicate(10000, createUniqueId_WELL(5))))
#[1] 10000

# independency on the seed: 
set.seed(1)
x <- replicate(100, createGlobalUniqueId(5))
x_WELL <- replicate(100, createUniqueId_WELL(5))
set.seed(1)
y <- replicate(100, createGlobalUniqueId(5))
y_WELL <- replicate(100, createUniqueId_WELL(5))
sum(x==y)
#[1] 100
sum(x_WELL==y_WELL)
#[1] 0

编辑

要理解为什么我们会得到重复的键,我们应该看看当我们调用set.seed(NULL). 所有与RNG相关的代码都是用C编写的,所以我们应该直接去svn.r-project.org/R/trunk/src/main/RNG.c参考函数do_setseed。如果seed = NULLthen 清楚地TimeToSeed被调用。有一条评论指出它应该位于 datetime.c 中,但是,它可以在svn.r-project.org/R/trunk/src/main/times.c中找到。

导航 R 源代码可能很困难,所以我将函数粘贴在这里:

/* For RNG.c, main.c, mkdtemp.c */
attribute_hidden
unsigned int TimeToSeed(void)
{
    unsigned int seed, pid = getpid();
#if defined(HAVE_CLOCK_GETTIME) && defined(CLOCK_REALTIME)
    {
    struct timespec tp;
    clock_gettime(CLOCK_REALTIME, &tp);
    seed = (unsigned int)(((uint_least64_t) tp.tv_nsec << 16) ^ tp.tv_sec);
    }
#elif defined(HAVE_GETTIMEOFDAY)
    {
    struct timeval tv;
    gettimeofday (&tv, NULL);
    seed = (unsigned int)(((uint_least64_t) tv.tv_usec << 16) ^ tv.tv_sec);
    }
#else
    /* C89, so must work */
    seed = (Int32) time(NULL);
#endif
    seed ^= (pid <<16);
    return seed;
}

所以每次我们调用 时set.seed(NULL),R 都会执行以下步骤:

  1. 以秒和纳秒为单位占用当前时间(如果可能,以#if defined块为单位的平台依赖性)
  2. 将位移位应用于纳秒,位'xor'es 结果以秒为单位
  3. 将位移位应用于 pid 并与之前的结果进行位“异或”
  4. 将结果设置为新种子

好吧,现在几乎很明显,当生成的种子发生碰撞时,我们会得到重复的值。我的猜测是当两个呼叫在 1 秒内发生时会发生这种情况,因此 tv_sec 是恒定的。为了确认这一点,我引入了一个滞后:

createUniqueIdWithLag <- function(bytes, lag) {
  Sys.sleep(lag)
  createUniqueId(bytes)
}
lags <- 1 / 10 ^ (1:5)
sapply(lags, function(x) length(unique(replicate(n, createUniqueIdWithLag(5, x)))))
[1] 1000 1000  996  992  990

令人困惑的是,即使与纳秒相比,延迟也很大,我们仍然会遇到碰撞!让我们进一步挖掘它,我为种子编写了一个“调试模拟器”:

emulate_seed <- function() {
  tv <- as.numeric(system('echo $(($(date +%s%N)))', intern = TRUE))
  pid <- Sys.getpid()
  tv_nsec <- tv %% 1e9
  tv_sec <- tv %/% 1e9
  seed <- bitwXor(bitwShiftL(tv_nsec, 16), tv_sec)
  seed <- bitwXor(bitwShiftL(pid, 16), seed)
  c(seed, tv_nsec, tv_sec, pid)
}

z <- replicate(1000, emulate_seed())
sapply(1:4, function(i) length(unique(z[i, ])))
# unique seeds, nanosecs, secs, pids:
#[1]  941 1000   36    1

这真的很令人困惑:纳秒都是唯一的,这并不能保证最终种子的唯一性。为了证明这种效果,这里有一个副本:

#            [,1]        [,2] 
#[1,] -1654969360 -1654969360
#[2,]   135644672   962643456
#[3,]  1397894128  1397894128 
#[4,]        2057        2057
bitwShiftL(135644672, 16)
#[1] -973078528
bitwShiftL(962643456, 16)
#[1] -973078528

最后一点:这两个数字的二进制表示和移位是

00001000000101011100011000000000 << 16 => 1100011000000000 + 16 zeroes
00111001011000001100011000000000 << 16 => 1100011000000000 + 16 zeroes

所以是的,这确实是一次不受欢迎的碰撞。

好吧,毕竟,最后的结论是:set.seed(NULL)容易受到高负载的影响,并且不保证在处理多个连续调用时不会发生冲突!

于 2014-04-18T07:39:52.133 回答
4

对于第一种方法,似乎确实有可能createUniqueId()在同一个 R 进程中同时发生两个调用并返回相同的 ID 字符串。

length(unique(sapply(seq_len(100000), function(.) createUniqueId(5))))
# [1] 93906
createGlobalUniqueId <- function(bytes) paste(as.hexmode(sample(256, bytes, replace = TRUE) - 1), collapse = "")
length(unique(sapply(seq_len(100000), function(.) createGlobalUniqueId(5))))
# [1] 100000

因此,如果您不想要 ID 冲突,我会采用第二种方法。

于 2014-04-16T00:41:50.963 回答