scipy.stats.

monte_carlo_test#

scipy.stats.monte_carlo_test(data, rvs, statistic, *, vectorized=None, n_resamples=9999, batch=None, alternative='two-sided', axis=0)[源代码]#

执行蒙特卡罗假设检验。

data 包含一个样本或一个或多个样本的序列。 rvs 指定 data 中样本在零假设下的分布。将给定 data 的 statistic 的值与蒙特卡罗零分布进行比较：使用 rvs 生成的 n_resamples 组样本的每个样本的统计量值。这给出了 p 值，即在零假设下观察到如此极端的检验统计量值的概率。

参数:
data类数组或类数组序列

观测值的数组或数组序列。

rvs可调用的对象或可调用的对象元组

在零假设下生成随机变量的可调用的对象或可调用的对象序列。rvs 的每个元素必须是可以接受关键词参数 size（例如 rvs(size=(m, n))）并且返回该形状的 N 维数组样本的可调用的对象。如果 rvs 是一个序列，则 rvs 中的可调用的对象数量必须匹配 data 中的样本数量，即 len(rvs) == len(data)。如果 rvs 是一个可调用的对象，则会将 data 视为一个样本。

statistic可调用的对象

要计算其假设检验的 p 值的统计数据。statistic 必须是可以接受一个样本（例如 statistic(sample)）或 len(rvs) 个独立样本（例如 statistic(samples1, sample2)（如果 rvs 包含两个可调用的对象且 data 包含两个样本））并返回值的统计数据。如果 vectorized 设置为 True，则 statistic 还必须接受关键词参数 axis，而且还需要进行矢量化处理才能在 data 中样本的提供的 axis 上计算统计数据。

vectorized布尔值，可选

如果 vectorized 设置为 False，则不会将关键词参数 axis 传递给 statistic，并且预计仅为一维样本计算统计数据。如果为 True，则会将关键词参数 axis 传递给 statistic，并且预计在将 ND 样本数组传递过来时，沿着 axis 计算统计数据。如果为 None（默认值），则在 statistic 的参数中有 axis 时，vectorized 会设置为 True。使用经过矢量化处理的统计数据通常可以减少计算时间。

n_resamples整数，默认值：9999

从 rvs 的每个可调用的对象中提取样本的数量。对应地，这是在蒙特卡罗零分布中用作零假设下的统计值数量。

批次int，可选

在对 statistic 进行每次调用时，蒙特卡洛样本处理的数量。内存使用是 O( batch * sample.size[axis] )。默认为 None，在这种情况下 batch 等于 n_resamples。

备选方案{‘双边’，‘较小’，‘较大’}

将为此计算 p 值的备选假设。对于每个备选方案，p 值定义如下。

  • 'greater'：大于或等于观测值检验统计量的空分布中的百分比。

  • 'less'：小于或等于观测值检验统计量的空分布中的百分比。

  • 'two-sided'：上述两个 p 值中较小者乘以 2。

轴int，默认值：0

用于计算统计量的data （或data 中的每个样本）轴。

返回值:
resMonteCarloTestResult

具有属性的对象

statisticfloat 或 ndarray

观测所得data 的检验统计量。

pvaluefloat 或 ndarray

给定备选方案的 p 值。

null_distributionndarray

在空假设下生成检验统计量的值。

警告

p 值的计算方法是计算与观测值统计量一样极端或比其更极端的空分布元素数量。由于使用有限精度算术，某些统计量函数会在理论值绝对相等的情况下返回数值上不同的值。在某些情况下，这会导致计算出的 p 值出现较大误差。 monte_carlo_test 通过将空分布中与检验统计量观测值 “接近”（在不准确类型的浮点 epsilon 的 100 倍的相对容差范围内）的元素视为等于检验统计量观测值来防止这种情况的发生。但是，建议用户检查空分布以评估此比较方法是否合适，如果不合适，请手动计算 p 值。

参考

[1]

B. Phipson 和 G. K. Smyth。“置换 p 值绝不应为零：在置换随机抽取时计算确切的 p 值”。统计学在遗传学和分子生物学中的应用程序 9.1 (2010)。

示例

假设我们希望测试小型样本是否来自于正态分布。我们决定使用样本的偏度作为检验统计量，并将 p 值 0.05 视为有统计学意义。

>>> import numpy as np
>>> from scipy import stats
>>> def statistic(x, axis):
...     return stats.skew(x, axis)

收集数据后，我们计算检验统计量的观测值。

>>> rng = np.random.default_rng()
>>> x = stats.skewnorm.rvs(a=1, size=50, random_state=rng)
>>> statistic(x, axis=0)
0.12457412450240658

为了确定，如果样本来自于正态分布，则在偶然情况下观测到如此极端的偏度值的概率，我们可以执行蒙特卡罗假设检验。该检验会从正态分布中随机抽取许多样本，计算每个样本的偏度，并将我们的原始偏度与此分布进行比较，以确定近似 p 值。

>>> from scipy.stats import monte_carlo_test
>>> # because our statistic is vectorized, we pass `vectorized=True`
>>> rvs = lambda size: stats.norm.rvs(size=size, random_state=rng)
>>> res = monte_carlo_test(x, rvs, statistic, vectorized=True)
>>> print(res.statistic)
0.12457412450240658
>>> print(res.pvalue)
0.7012

在零假设下，获得小于或等于观测值的检验统计量的概率约为 70%。这大于我们选择的 5% 的阈值，因此我们不能将其视为反对零假设的有力证据。

请注意，此 p 值与 scipy.stats.skewtest 基本相同，它依赖于基于样本偏度的检验统计量的渐近分布。

>>> stats.skewtest(x).pvalue
0.6892046027110614

这种渐近近似对小样本量无效，但 monte_carlo_test 可以用于任何大小的样本。

>>> x = stats.skewnorm.rvs(a=1, size=7, random_state=rng)
>>> # stats.skewtest(x) would produce an error due to small sample
>>> res = monte_carlo_test(x, rvs, statistic, vectorized=True)

提供检验统计量的蒙特卡罗分布以供进一步调查。

>>> import matplotlib.pyplot as plt
>>> fig, ax = plt.subplots()
>>> ax.hist(res.null_distribution, bins=50)
>>> ax.set_title("Monte Carlo distribution of test statistic")
>>> ax.set_xlabel("Value of Statistic")
>>> ax.set_ylabel("Frequency")
>>> plt.show()
../../_images/scipy-stats-monte_carlo_test-1.png
本页内容