scipy.stats.sampling.

离散指南表#

class scipy.stats.sampling.DiscreteGuideTable(dist, *, domain=None, guide_factor=1, random_state=None)#

离散引导表方法。

离散引导表方法从任意但有限的概率向量中采样。它使用大小为 \(N\) 的概率向量或具有有限支持的概率质量函数来生成随机数。离散引导表的设置非常慢（与向量长度成线性关系），但提供了非常快的采样。

参数:
dist类数组或对象，可选

分布的概率向量（PV）。如果PV不可用，则期望一个具有 pmf 方法的类实例。PMF的签名期望为：def pmf(self, k: int) -> float。即它应接受一个Python整数并返回一个Python浮点数。

领域int, 可选

支持PMF。如果概率向量 (pv) 不可用，则必须给定一个有限域。即PMF必须有一个有限的支持。默认是 None。当 None 时：

  • 如果分布对象 dist 提供了 support 方法，它将用于设置分布的域。

  • 否则，支持被假定为 (0, len(pv)) 。当此参数与概率向量一起传递时，domain[0] 用于将分布从 (0, len(pv)) 重新定位到 (domain[0], domain[0]+len(pv)) ，而 domain[1] 被忽略。有关更详细的解释，请参见注释和教程。

guide_factor: int, 可选

引导表相对于PV长度的尺寸。较大的引导表会导致生成时间更快，但需要更昂贵的设置。不推荐使用大于3的尺寸。如果相对尺寸设置为0，则使用顺序搜索。默认值为1。

random_state{None, int,}

用于生成均匀随机数流的底层 NumPy 随机数生成器的 NumPy 随机数生成器或种子。如果 random_state 是 None（或 np.random），则使用 numpy.random.RandomState 单例。如果 random_state 是整数，则使用一个新的 RandomState 实例，并以 random_state 为种子。如果 random_state 已经是 Generator 或 RandomState 实例，则使用该实例。

方法

ppf(u)

给定分布的PPF。

rvs([size, random_state])

来自发行版的示例。

set_random_state([random_state])

设置底层均匀随机数生成器。

注释

当有限概率向量可用或分布的概率质量函数（PMF）可用时，此方法有效。如果仅提供PMF，则还必须给出PMF的*有限*支持（域）。建议首先通过在支持的每个点上评估PMF来获得概率向量，然后使用它。

DGT 从任意但有限概率向量中采样。随机数通过反演方法生成，即。

  1. 生成一个随机数 U ~ U(0,1)。

  2. 找到最小的整数 I 使得 F(I) = P(X<=I) >= U。

步骤 (2) 是关键步骤。使用顺序搜索需要 O(E(X)) 次比较，其中 E(X) 是分布的期望值。然而，索引搜索使用一个引导表来跳转到某个 I’ <= I 接近 I 的位置，以在常数时间内找到 X。实际上，当引导表的大小与概率向量相同时（这是默认设置），预期的比较次数减少到 2。对于更大的引导表，这个数字变得更小（但总是大于 1），对于更小的表，这个数字变得更大。对于表大小为 1 的极限情况，算法简单地执行顺序搜索。

另一方面，引导表的设置时间为 O(N)，其中 N 表示概率向量的长度（对于长度为 1 的情况，不需要预处理）。此外，对于非常大的引导表，内存效应可能会降低算法的速度。因此，我们不建议使用比给定概率向量大三倍以上的引导表。如果只需要生成少量随机数，（远）小的表大小会更好。引导表相对于给定概率向量长度的尺寸可以通过 guide_factor 参数设置。

如果给定一个概率向量，它必须是一个一维的非负浮点数数组，且不包含任何 inf 或 nan 值。此外，必须至少有一个非零项，否则会引发异常。

默认情况下，概率向量从0开始索引。然而，这可以通过传递一个 domain 参数来改变。当 domain 与PV一起给出时，它会将分布从 (0, len(pv)) 重新定位到 (domain[0], domain[0] + len(pv))。在这种情况下，domain[1] 被忽略。

参考文献

[1]

UNU.RAN 参考手册，第5.8.4节，“DGT - （离散）引导表方法（索引搜索）” https://statmath.wu.ac.at/unuran/doc/unuran.html#DGT

[2]

H.C. Chen 和 Y. Asau (1974). 关于从经验分布生成随机变量，AIIE Trans. 6, pp. 163-166.

示例

>>> from scipy.stats.sampling import DiscreteGuideTable
>>> import numpy as np

要使用概率向量创建随机数生成器，请使用：

>>> pv = [0.1, 0.3, 0.6]
>>> urng = np.random.default_rng()
>>> rng = DiscreteGuideTable(pv, random_state=urng)

RNG 已经设置好。现在，我们可以使用 rvs 方法从分布中生成样本：

>>> rvs = rng.rvs(size=1000)

为了验证随机变量是否遵循给定的分布，我们可以使用卡方检验（作为拟合优度的度量）：

>>> from scipy.stats import chisquare
>>> _, freqs = np.unique(rvs, return_counts=True)
>>> freqs = freqs / np.sum(freqs)
>>> freqs
array([0.092, 0.355, 0.553])
>>> chisquare(freqs, pv).pvalue
0.9987382966178464

由于p值非常高，我们未能拒绝原假设，即观察到的频率与预期频率相同。因此，我们可以安全地假设变量是从给定的分布中生成的。请注意，这仅证明了算法的正确性，而不是样本的质量。

如果无法提供 PV，也可以传递一个具有 PMF 方法和有限域的类实例。

>>> urng = np.random.default_rng()
>>> from scipy.stats import binom
>>> n, p = 10, 0.2
>>> dist = binom(n, p)
>>> rng = DiscreteGuideTable(dist, random_state=urng)

现在，我们可以使用 rvs 方法从分布中采样，并测量样本的拟合优度：

>>> rvs = rng.rvs(1000)
>>> _, freqs = np.unique(rvs, return_counts=True)
>>> freqs = freqs / np.sum(freqs)
>>> obs_freqs = np.zeros(11)  # some frequencies may be zero.
>>> obs_freqs[:freqs.size] = freqs
>>> pv = [dist.pmf(i) for i in range(0, 11)]
>>> pv = np.asarray(pv) / np.sum(pv)
>>> chisquare(obs_freqs, pv).pvalue
0.9999999999999989

为了检查样本是否从正确的分布中抽取，我们可以可视化样本的直方图：

>>> import matplotlib.pyplot as plt
>>> rvs = rng.rvs(1000)
>>> fig = plt.figure()
>>> ax = fig.add_subplot(111)
>>> x = np.arange(0, n+1)
>>> fx = dist.pmf(x)
>>> fx = fx / fx.sum()
>>> ax.plot(x, fx, 'bo', label='true distribution')
>>> ax.vlines(x, 0, fx, lw=2)
>>> ax.hist(rvs, bins=np.r_[x, n+1]-0.5, density=True, alpha=0.5,
...         color='r', label='samples')
>>> ax.set_xlabel('x')
>>> ax.set_ylabel('PMF(x)')
>>> ax.set_title('Discrete Guide Table Samples')
>>> plt.legend()
>>> plt.show()
../../_images/scipy-stats-sampling-DiscreteGuideTable-1_00_00.png

要设置引导表的大小，请使用 guide_factor 关键字参数。这将设置引导表相对于概率向量的大小。

>>> rng = DiscreteGuideTable(pv, guide_factor=1, random_state=urng)

要计算参数为 \(n=4\) 和 \(p=0.1\) 的二项分布的PPF：我们可以设置如下引导表：

>>> n, p = 4, 0.1
>>> dist = binom(n, p)
>>> rng = DiscreteGuideTable(dist, random_state=42)
>>> rng.ppf(0.5)
0.0
On this page