Sobol#
- class scipy.stats.qmc.Sobol(d, *, scramble=True, bits=None, rng=None, optimization=None)[源代码]#
用于生成(扰乱的)Sobol’序列的引擎。
Sobol’ 序列是低差异、准随机数。可以使用两种方法绘制点
random_base2
: 安全地绘制 \(n=2^m\) 个点。此方法保证序列的平衡特性。random
: 从序列中绘制任意数量的点。请参阅下面的警告。
- 参数:
- dint
序列的维度。最大维度为 21201。
- scramblebool, 可选
如果为 True,则使用 LMS+shift 扰乱。否则,不进行扰乱。默认值为 True。
- bitsint, 可选
生成器的位数。控制可以生成的最大点数,即
2**bits
。最大值为 64。它不对应于返回类型,返回类型始终为np.float64
,以防止点重复自身。默认值为 None,为了向后兼容,对应于 30。1.9.0 版本中新增。
- optimization{None, “random-cd”, “lloyd”}, 可选
是否使用优化方案来提高采样后的质量。请注意,这是一个后处理步骤,不能保证样本的所有属性都会被保留。默认值为 None。
random-cd
:坐标的随机排列,以降低中心差异。基于中心差异的最佳样本会不断更新。与使用其他差异度量相比,基于中心差异的采样在 2D 和 3D 子投影中显示出更好的空间填充鲁棒性。lloyd
:使用修改后的 Lloyd-Max 算法扰动样本。该过程收敛到等距样本。
1.10.0 版本中新增。
- rng
numpy.random.Generator
, 可选 伪随机数生成器状态。当 rng 为 None 时,会使用来自操作系统的熵创建一个新的
numpy.random.Generator
。除了numpy.random.Generator
之外的类型会传递给numpy.random.default_rng
以实例化Generator
。在 1.15.0 版本中更改: 作为从使用
numpy.random.RandomState
过渡到numpy.random.Generator
的 SPEC-007 过渡的一部分,此关键字已从 seed 更改为 rng。在过渡期间,这两个关键字将继续工作,但一次只能指定一个。在过渡期后,使用 seed 关键字的函数调用将发出警告。在弃用期之后,seed 关键字将被删除。
说明
Sobol’ 序列 [1] 在 \([0,1)^{d}\) 中提供 \(n=2^m\) 个低差异点。扰乱它们 [3] 使它们适用于奇异被积函数,提供了一种误差估计方法,并且可以提高它们的收敛速度。实现的扰乱策略是(左)线性矩阵扰乱 (LMS),然后是数字随机偏移 (LMS+shift) [2]。
Sobol’ 序列有许多版本,具体取决于它们的“方向数”。此代码使用 [4] 中的方向数。因此,最大维度数为 21201。方向数已通过搜索标准 6 预先计算,可以在 https://web.maths.unsw.edu.au/~fkuo/sobol/ 获取。
警告
Sobol’ 序列是一个求积规则,如果使用的样本大小不是 2 的幂,或者跳过第一个点,或者稀疏化序列 [5],则它们会失去其平衡特性。
如果 \(n=2^m\) 个点不足,则应取 \(2^M\) 个点,其中 \(M>m\)。当扰乱时,独立副本 R 的数量不必是 2 的幂。
Sobol’ 序列生成到一定数量 \(B\) 的位。生成 \(2^B\) 个点后,序列将重复。因此,会引发错误。可以使用参数 bits 控制位数。
参考资料
[1]I. M. Sobol’, “The distribution of points in a cube and the accurate evaluation of integrals.” Zh. Vychisl. Mat. i Mat. Phys., 7:784-802, 1967.
[2]J. Matousek, “On the L2-discrepancy for anchored boxes.” J. of Complexity 14, 527-556, 1998.
[3]Art B. Owen, “Scrambling Sobol and Niederreiter-Xing points.” Journal of Complexity, 14(4):466-489, December 1998.
[4]S. Joe and F. Y. Kuo, “Constructing sobol sequences with better two-dimensional projections.” SIAM Journal on Scientific Computing, 30(5):2635-2654, 2008.
[5]Art B. Owen, “On dropping the first Sobol’ point.” arXiv:2008.08051, 2020.
示例
从 Sobol’ 的低差异序列生成样本。
>>> from scipy.stats import qmc >>> sampler = qmc.Sobol(d=2, scramble=False) >>> sample = sampler.random_base2(m=3) >>> sample array([[0. , 0. ], [0.5 , 0.5 ], [0.75 , 0.25 ], [0.25 , 0.75 ], [0.375, 0.375], [0.875, 0.875], [0.625, 0.125], [0.125, 0.625]])
使用差异准则计算样本的质量。
>>> qmc.discrepancy(sample) 0.013882107204860938
要继续现有的设计,可以通过再次调用
random_base2
来获得额外的点。或者,您可以跳过一些点,如>>> _ = sampler.reset() >>> _ = sampler.fast_forward(4) >>> sample_continued = sampler.random_base2(m=2) >>> sample_continued array([[0.375, 0.375], [0.875, 0.875], [0.625, 0.125], [0.125, 0.625]])
最后,可以将样本缩放到边界。
>>> l_bounds = [0, 2] >>> u_bounds = [10, 5] >>> qmc.scale(sample_continued, l_bounds, u_bounds) array([[3.75 , 3.125], [8.75 , 4.625], [6.25 , 2.375], [1.25 , 3.875]])
方法
fast_forward
(n)将序列快进 n 个位置。
integers
(l_bounds, *[, u_bounds, n, ...])从 l_bounds(包括)到 u_bounds(不包括)绘制 n 个整数,或者如果 endpoint=True,则从 l_bounds(包括)到 u_bounds(包括)。
random
([n, workers])在半开区间
[0, 1)
中绘制 n 个点。random_base2
(m)从 Sobol' 序列中绘制点。
reset
()将引擎重置为基本状态。