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)[source]#

执行蒙特卡洛假设检验。

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
在此页面上