使用 NumPy 随机数生成器的最佳实践

鉴于实现真正的随机性所面临的实际挑战，科学界采用名为伪随机数生成器 (RNG) 的确定性算法来创建模拟随机性的序列。这些生成器用于模拟、实验和分析，在这些场景中，拥有看似不可预测的数字至关重要。我想在这里分享我关于伪 RNG，尤其是 NumPy 中可用的那些，的最佳实践的经验。

伪 RNG 通过确定性算法更新其内部状态来工作。此内部状态使用称为种子的值进行初始化，并且每次更新都会产生一个看起来是随机生成的数字。这里的关键是该过程是确定性的，这意味着如果您从同一个种子开始并应用相同的算法，您将获得相同的内部状态序列（以及数字）。尽管有这种确定性，但生成的数字表现出随机性的属性，看起来不可预测且均匀分布。用户可以手动指定种子，从而在一定程度上控制生成的序列，或者他们可以选择让 RNG 对象从系统熵自动推导出种子。后一种方法通过将外部因素纳入种子来增强不可预测性。

我假设您对 NumPy 有一定的了解，并且使用的是 NumPy 1.17 或更高版本。之所以如此，是因为在版本 1.17 的 random 模块中引入了非常棒的新功能。由于 numpy 通常被导入为 np，因此有时我会使用 np 代替 numpy。最后，在本文的其余部分，RNG 将始终表示伪 RNG。

主要信息#

1. 避免使用全局 NumPy RNG。这意味着您应该避免使用 np.random.seed 和 np.random.* 函数（例如 np.random.random）来生成随机值。
2. 使用 np.random.default_rng 函数创建一个新的 RNG 并将其传递。
3. 注意并行随机数生成，并使用 NumPy 提供的策略。

请注意，在旧版本的 NumPy（<1.17）中，创建新 RNG 的方法是使用 np.random.RandomState，它基于流行的 Mersenne Twister 19937 算法。此函数在较新版本的 NumPy 中仍然可用，但现在建议使用 default_rng 代替，后者返回统计上更好的 PCG64 RNG 的实例。您可能仍然会在测试中看到使用 np.random.RandomState，因为它在不同 NumPy 版本之间具有很强的稳定性保证。

使用 NumPy 生成随机数#

当您在 Python 脚本中导入 numpy 时，后台会创建一个 RNG。此 RNG 用于您使用 np.random.random 等函数生成新的随机值时。在这里，我将此 RNG 称为全局 NumPy RNG。

虽然不推荐，但使用 np.random.seed 函数在脚本开头重置此全局 RNG 的种子是一种常见做法。在开头修复种子可以确保脚本可重复：每次运行脚本时都会产生相同的值和结果。但是，虽然有时很方便，但使用全局 NumPy RNG 是一种不好的做法。一个简单的原因是使用全局变量会导致不希望发生的副作用。例如，有人可能会使用 np.random.random，而不知道全局 RNG 的种子在代码库中的其他地方被设置了。引用 Robert Kern 的 NumPy 增强提案 (NEP) 19

np.random.* 方便函数背后的隐式全局 RandomState 会导致问题，尤其是在涉及线程或其他并发形式时。全局状态始终有问题。我们明确建议在涉及可重复性时避免使用这些方便函数。[…] 获取可重复的伪随机数的最佳实践是使用种子实例化一个生成器对象并将其传递。

简而言之

• 不要使用 np.random.seed（它会重新播种已经创建的全局 NumPy RNG），然后使用 np.random.* 函数，而应该创建一个新的 RNG。
• 您应该在脚本开头创建一个 RNG（如果您希望可重复性，请使用您自己的种子），并在脚本的其余部分使用此 RNG。

要创建一个新的 RNG，您可以使用 default_rng 函数，如 random 模块文档的介绍 中所述

import numpy as np

rng = np.random.default_rng()
rng.random()  # generate a floating point number between 0 and 1

如果您想使用种子来确保可重复性，NumPy 文档 建议使用一个大的随机数，其中“大”意味着至少 128 位。使用大的随机数的第一个原因是，这增加了拥有与其他人不同的种子，从而获得独立结果的概率。第二个原因是，仅依赖于种子的小数会导致偏差，因为它们没有完全探索 RNG 的状态空间。此限制意味着由 RNG 生成的第一个数字可能看起来不像预期的那样随机，因为无法访问第一个内部状态。例如，某些数字永远不会作为第一个输出产生。一种可能性是在 RNG 的状态空间中随机选择种子，但是 根据 Robert Kern 的说法，128 位随机数足够大¹。要为您的种子生成一个 128 位随机数，您可以依赖于 secrets 模块

import secrets

secrets.randbits(128)

运行此代码时，我会得到 65647437836358831880808032086803839626，用作种子的数字。此数字是随机生成的，因此您需要复制粘贴由 secrets.randbits(128) 返回的值，否则您每次运行代码时都会有不同的种子，从而破坏可重复性。

import numpy as np

seed = 65647437836358831880808032086803839626
rng = np.random.default_rng(seed)
rng.random()

只播种 RNG 一次（并将该 RNG 传递）的原因是，使用像 default_rng 返回的那样的优质 RNG，您可以确保生成数字的良好随机性和独立性。但是，如果操作不当，使用多个 RNG（每个 RNG 使用自己的种子创建）可能会导致随机数流的独立性低于使用相同种子创建的随机数流²。也就是说，如 Robert Kern 所解释，使用 NumPy 1.17 中引入的 RNG 和播种策略，创建使用系统熵的 RNG，即多次使用 default_rng(None)，被认为是相当安全的。但是，如后文所述，在并行运行作业并依赖 default_rng(None) 时要小心。只播种 RNG 一次的原因是，获取优质种子可能很耗时。一旦您拥有用于实例化生成器的优质种子，您就可以继续使用它。

传递 NumPy RNG#

当您编写在单独使用以及在更复杂的脚本中使用的函数时，能够传递种子或您已创建的 RNG 会很方便。default_rng 函数可以让您轻松实现这一点。如上所述，此函数可用于根据您选择的种子创建新的 RNG（如果您将种子传递给它），或者在传递 None 时从系统熵创建新的 RNG，但您也可以传递已经创建的 RNG。在这种情况下，返回的 RNG 是您传递的 RNG。

import numpy as np


def stochastic_function(high=10, rng=None):
    rng = np.random.default_rng(rng)
    return rng.integers(high, size=5)

您可以将 int 种子或您已创建的 RNG 传递给 stochastic_function。确切地说，对于某些类型的 RNG（例如使用 default_rng 本身创建的 RNG），default_rng 函数会返回与传递给它的 RNG 完全相同的 RNG。您可以参考 default_rng 文档，了解有关您可以传递给此函数的参数的更多详细信息³。

并行处理#

在将 RNG 与并行处理结合使用时，您必须小心。让我们考虑蒙特卡罗模拟的背景：您有一个返回随机输出的随机函数，并且您想要生成这些随机输出很多次，例如为了计算经验平均值。如果函数的计算成本很高，则加速计算时间的简单解决方案是使用并行处理。根据您使用的并行处理库或后端，可能会观察到不同的行为。例如，如果您没有自己设置种子，则分叉的 Python 进程可能会使用相同的随机种子（例如从系统熵生成的），从而产生完全相同的输出，这对计算资源是一种浪费。使用 Joblib 并行处理库时，一个很好的示例说明了这一点 这里。

如果您在主脚本开头修复种子以确保可重复性，然后将您播种的 RNG 传递给每个要并行运行的进程，大多数情况下这不会得到您想要的结果，因为此 RNG 将被深度复制。因此，每个进程都会产生相同的结果。一种解决方案是创建与并行进程一样多的 RNG，每个 RNG 使用不同的种子。现在的问题是，您无法像您想象的那样轻松地选择种子。当您选择两个不同的种子来实例化两个不同的 RNG 时，您如何知道这些 RNG 生成的数字看起来是统计独立的？² 为并行进程设计独立 RNG 一直是一个重要的研究问题。例如，参见 L’Ecuyer 等人 (2017) 的 并行计算机的随机数：需求和方法，重点关注 GPU，了解不同方法的良好概述。

从 NumPy 1.17 开始，实例化独立的 RNG 变得非常容易。根据您使用的 RNG 类型，不同的策略可用，如 NumPy 文档的并行随机数生成部分所述。其中一种策略是使用SeedSequence，这是一种算法，可以确保将不良输入种子转换为良好的初始 RNG 状态。更准确地说，这确保了您的 RNG 不会出现退化行为，并且后续的数字将看起来是随机且独立的。此外，它确保将接近的种子映射到截然不同的初始状态，从而导致 RNG 非常有可能相互独立。您可以参考SeedSequence 生成的文档，了解如何从SeedSequence或现有 RNG 生成独立 RNG 的示例。这里我将展示如何将它应用于上面提到的joblib 示例。

import numpy as np
from joblib import Parallel, delayed


def stochastic_function(high=10, rng=None):
    rng = np.random.default_rng(rng)
    return rng.integers(high, size=5)


seed = 319929794527176038403653493598663843656
# creating the RNG that is passed around.
rng = np.random.default_rng(seed)
# create 5 independent RNGs
child_rngs = rng.spawn(5)

# use 2 processes to run the stochastic_function 5 times with joblib
random_vector = Parallel(n_jobs=2)(
    delayed(stochastic_function)(rng=child_rng) for child_rng in child_rngs
)
print(random_vector)

通过使用固定种子，每次运行此代码时始终获得相同的结果，并且通过使用rng.spawn，您在每次调用stochastic_function时都有一个独立的 RNG。请注意，这里您也可以从使用种子创建的SeedSequence生成，而不是创建 RNG。但是，一般来说您传递的是 RNG，因此我假设只能访问 RNG。另请注意，从版本 1.25 的 NumPy⁴开始，才能从 RNG 生成。

希望这篇博文能帮助您了解使用 NumPy RNG 的最佳方法。新的 Numpy API 为您提供了所需的所有工具。以下资源可供进一步阅读。最后，我要感谢 Pamphile Roy、Stefan van der Walt 和 Jarrod Millman 对他们的宝贵反馈和评论，这些评论极大地改进了这篇博文的原始版本。

资源#

NumPy RNG#

• NumPy 随机模块的文档是查找信息和我在本文中分享的大部分信息的最佳地点。
• 关于随机数生成器策略的 NumPy 增强提案 (NEP) 19，它导致了 NumPy 1.17 中引入的更改。
• 关于check_random_state函数和 RNG 良好实践的NumPy 问题，尤其是 Robert Kern 的这条评论。
• 还可以查看 Robert Kern 对我关于SeedSequence能做什么和不能做什么的问题的这个回答。这也解释了为什么建议使用非常大的随机数作为种子。
• 如何为整个执行设置随机状态？来自 scikit-learn FAQ。
• 关于不同库用来为 RNG 种子的 API 的持续讨论。

一般意义上的 RNG#

• 并行计算机的随机数：需求和方法，重点关注 GPU，作者 L’Ecuyer 等人 (2017 年)。
• 要了解更多关于 NumPy 中使用的默认 RNG (名为 PCG) 的信息，我推荐PCG 论文，其中还包含大量关于一般意义上的 RNG 的有用信息。pcg-random.org 网站也充满了关于 RNG 的有趣信息。