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) 个单独的样本（例如，如果 rvs 包含两个可调用对象且 data 包含两个样本，则为 statistic(samples1, sample2)）并返回结果统计量的可调用对象。如果将 vectorized 设置为 True，则 statistic 还必须接受关键字参数 axis，并且可以向量化以计算 data 中样本的给定 axis 的统计量。

vectorizedbool，可选

如果将 vectorized 设置为 False，则不会向 statistic 传递关键字参数 axis，并且希望仅针对 1D 样本计算统计量。如果为 True，则会向 statistic 传递关键字参数 axis，并且希望在传递 ND 样本数组时沿 axis 计算统计量。如果为 None（默认），则如果 axis 是 statistic 的参数，则将 vectorized 设置为 True。使用向量化统计量通常会减少计算时间。

n_resamplesint，默认值：9999

从 rvs 的每个可调用对象中抽取的样本数。等效地，在零假设下用作蒙特卡洛零分布的统计量值数量。

batchint，可选

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

alternative{‘two-sided’，‘less’，‘greater’}

计算 p 值的备择假设。对于每个备择假设，p 值定义如下。

'greater' ：大于或等于检验统计量的观察值的零分布百分比。
'less' ：小于或等于检验统计量的观察值的零分布百分比。
'two-sided' ：上述 p 值中较小值的两倍。

axisint，默认值：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)。

示例

假设我们希望测试一个小的样本是否是从正态分布中抽取的。我们决定将样本的偏度用作检验统计量，并且我们将考虑 0.05 的 p 值为统计学上的显着性。

>>> 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 的 p 值基本匹配，后者依赖于基于样本偏度的检验统计量的渐近分布。

>>> 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