scipy.stats.

permutation_test

scipy.stats.permutation_test(data, statistic, *, permutation_type='independent', vectorized=None, n_resamples=9999, batch=None, alternative='two-sided', axis=0, random_state=None)

对给定的统计量在提供的数据上执行置换检验。

对于独立样本统计，零假设是数据是从同一分布中随机抽取的。对于配对样本统计，可以测试两个零假设：数据是随机配对的，或者数据是随机分配给样本的。

参数:

数据类似数组的迭代对象

包含样本，每个样本是一个观测值数组。样本数组的维度必须沿着 axis 之外的维度兼容广播。

统计可调用

要计算假设检验p值的统计量。statistic 必须是一个可调用的函数，它接受样本作为单独的参数（例如 statistic(*data)）并返回计算出的统计量。如果 vectorized 设置为 True，statistic 还必须接受一个关键字参数 axis，并且能够向量化地沿着样本数组的提供的 axis 计算统计量。

permutation_type{‘independent’, ‘samples’, ‘pairings’}, 可选

要执行的排列类型，根据零假设。前两种排列类型用于配对样本统计，其中所有样本包含相同数量的观测值，并且沿 axis 具有相应索引的观测值被认为是配对的；第三种用于独立样本统计。

• 'samples' : 观测值被分配到不同的样本中，但仍与来自其他样本的相同观测值配对。这种排列类型适用于配对样本假设检验，如Wilcoxon符号秩检验和配对t检验。

• 'pairings' : 观察值与不同的观察值配对，但它们仍保持在同一个样本中。这种排列类型适用于关联/相关性检验，统计量如斯皮尔曼的 \(\rho\) 、肯德尔的 \(\tau\) 和皮尔逊的 \(r\) 。

• 'independent' (默认) : 观测值被分配到不同的样本中。样本可能包含不同数量的观测值。这种排列类型适用于独立样本假设检验，如Mann-Whitney \(U\) 检验和独立样本t检验。

  请参阅下面的“注释”部分，以获取有关排列类型的更详细描述。

矢量化bool, 可选

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

n_resamplesint 或 np.inf, 默认值: 9999

用于近似零分布的随机排列（重采样）的数量。如果大于或等于不同排列的数量，将计算精确的零分布。请注意，不同排列的数量随着样本大小的增加而迅速增长，因此只有非常小的数据集才能进行精确测试。

批处理int, 可选

每次调用 statistic 时处理的排列数量。内存使用量为 O( batch * n )，其中 n 是所有样本的总大小，无论 vectorized 的值如何。默认值为 None，在这种情况下，batch 是排列的数量。

替代方案{‘双侧’, ‘小于’, ‘大于’}, 可选

计算p值的备择假设。对于每个备择假设，p值在精确检验中定义如下。

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

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

• 'two-sided' (默认) : 两次较小的p值。

请注意，随机化测试的p值是根据[R5641c5b1ce56-2]_和[R5641c5b1ce56-3]_中建议的保守（高估）近似值计算的，而不是根据[R5641c5b1ce56-4]_中建议的无偏估计量计算的。也就是说，在计算随机化零分布中与观察到的检验统计量同样极端的比例时，分子和分母的值都增加了一。这种调整的解释是，观察到的检验统计量的值总是被包含为随机化零分布的一个元素。用于双侧p值的约定并不普遍；如果偏好不同的定义，则返回观察到的检验统计量和零分布。

轴int, 默认值: 0

计算统计量时所沿的（广播后的）样本轴。如果样本的维度数不同，则在考虑 axis 之前，会在维度较少的样本前添加单例维度。

random_state{None, int,}

numpy.random.RandomState}, 可选

用于生成排列的伪随机数生成器状态。

如果 random_state 是 None``（默认），则使用 `numpy.random.RandomState` 单例。如果 `random_state` 是整数，则使用一个新的 ``RandomState 实例，并以 random_state 为种子。如果 random_state 已经是 Generator 或 RandomState 实例，则使用该实例。

返回:

resPermutationTestResult

一个带有属性的对象：

统计浮点数或ndarray

数据的观测测试统计量。

p值浮点数或ndarray

给定替代假设的p值。

null_distributionndarray

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

注释

此函数支持的三种排列测试类型如下所述。

未配对统计 (permutation_type='independent'):

与此排列类型相关的零假设是，所有观测值都是从相同的潜在分布中抽取的，并且它们已被随机分配到其中一个样本中。

假设 data 包含两个样本；例如 a, b = data。当 1 < n_resamples < binom(n, k) 时，其中

• k 是 a 中的观察次数。

• n 是 a 和 b 中的观测总数，并且

• binom(n, k) 是二项式系数（n 选 k），

数据被合并（连接），随机分配到第一个或第二个样本，并计算统计量。这个过程重复进行，permutation 次，生成在零假设下的统计量分布。原始数据的统计量与这个分布进行比较，以确定 p 值。

当 n_resamples >= binom(n, k) 时，执行精确检验：数据在每个不同的方式中*分区*到样本之间，并且形成精确的零分布。请注意，对于给定的数据在样本之间的分区，只考虑每个样本内数据的一个排序/排列。对于不依赖于样本内数据顺序的统计量，这大大减少了计算成本，而不影响零分布的形状（因为每个值的频率/计数受到相同因素的影响）。

对于 a = [a1, a2, a3, a4] 和 b = [b1, b2, b3]，这种排列类型的一个例子是 x = [b3, a1, a2, b2] 和 y = [a4, b1, a3]。因为在精确测试中，每个样本内的数据只有一种排序/排列被考虑，像 x = [b3, a1, b2, a2] 和 y = [a4, a3, b1] 这样的重采样将 不 被视为与上述例子不同。

permutation_type='independent' 不支持单样本统计，但它可以应用于两个以上样本的统计。在这种情况下，如果 n 是每个样本中观测数量的数组，不同分区的数量为:

np.prod([binom(sum(n[i:]), sum(n[i+1:])) for i in range(len(n)-1)])

配对统计，置换配对 (permutation_type='pairings'):

与此排列类型相关的零假设是，每个样本内的观察值是从相同的潜在分布中抽取的，并且与其他样本元素的配对是随机分配的。

假设 data 只包含一个样本；例如 a, = data，我们希望考虑 a 的所有元素与第二个样本 b 的元素的所有可能配对。设 n 为 a 中的观察数，这也必须等于 b 中的观察数。

当 1 < n_resamples < factorial(n) 时，a 的元素会被随机排列。用户提供的统计量接受一个数据参数，例如 a_perm，并计算统计量，考虑 a_perm 和 b。这个过程会重复执行 permutation 次，生成在零假设下的统计量分布。原始数据的统计量会与这个分布进行比较，以确定p值。

当 n_resamples >= factorial(n) 时，将执行精确检验：a 在每种不同的方式中恰好排列一次。因此，统计量 对于 a 和 b 之间的每种唯一样本配对恰好计算一次。

对于 a = [a1, a2, a3] 和 b = [b1, b2, b3]，这种排列类型的一个例子是 a_perm = [a3, a1, a2]，而 b 保持其原始顺序。

permutation_type='pairings' 支持包含任意数量样本的 data，每个样本必须包含相同数量的观测值。data 中提供的所有样本都是 独立 排列的。因此，如果 m 是样本数量，n 是每个样本中的观测值数量，那么在精确检验中的排列数量为:

factorial(n)**m

请注意，如果一个两样本统计量，例如，本质上不依赖于观察值提供的顺序 - 只依赖于观察值的*配对* - 那么应该只在 data 中提供两个样本中的一个。这大大减少了计算成本，而不会影响零分布的形状（因为每个值的频率/计数受到相同因素的影响）。

配对统计，置换样本 (permutation_type='samples'):

与此排列类型相关的零假设是，每对观察值来自相同的潜在分布，并且它们被分配到的样本是随机的。

假设 data 包含两个样本；例如 a, b = data。设 n 为 a 中的观测数，该数目必须也等于 b 中的观测数。

当 1 < n_resamples < 2**n 时，a 和 b 的元素在样本之间随机交换（保持它们的配对），并计算统计量。这个过程重复进行，permutation 次，生成在零假设下统计量的分布。原始数据的统计量与这个分布进行比较，以确定p值。

当 n_resamples >= 2**n 时，将执行精确检验：观测值在保持配对的情况下，被分配到两个样本中的每一种不同方式中，恰好一次。

对于 a = [a1, a2, a3] 和 b = [b1, b2, b3]，这种排列类型的一个例子是 x = [b1, a2, b3] 和 y = [a1, b2, a3]。

permutation_type='samples' 支持包含任意数量样本的 data，每个样本必须包含相同数量的观测值。如果 data 包含多个样本，则 data 内的配对观测值在样本之间独立交换。因此，如果 m 是样本数量，n 是每个样本内的观测值数量，那么精确测试中的排列数量为:

factorial(m)**n

几种配对样本的统计检验，如Wilcoxon符号秩检验和配对样本t检验，可以在仅考虑两个配对元素之间的*差异*的情况下进行。因此，如果``data``只包含一个样本，那么零分布是通过独立地改变每个观测值的*符号*来形成的。

警告

p 值是通过计算零分布中与统计量的观测值同样极端或更极端的元素数量来计算的。由于使用了有限精度算术，一些统计函数在理论值完全相等时返回数值上不同的值。在某些情况下，这可能导致计算的 p 值出现较大误差。permutation_test 通过将零分布中与检验统计量的观测值“接近”（在相对容差为 100 倍不精确数据类型的浮点 epsilon 内）的元素视为等于检验统计量的观测值来防止这种情况。然而，建议用户检查零分布以评估这种比较方法是否合适，如果不合适，则手动计算 p 值。请参见下例。

参考文献

[1]

18. 1. 费舍尔。《实验设计》，第六版（1951年）。

[2]

B. Phipson and G. K. Smyth. “Permutation P-values Should Never Be Zero: Calculating Exact P-values When Permutations Are Randomly Drawn.” Statistical Applications in Genetics and Molecular Biology 9.1 (2010).

[3]

M. D. Ernst. “Permutation Methods: A Basis for Exact Inference”. Statistical Science (2004).

[4]

B. Efron and R. J. Tibshirani. An Introduction to the Bootstrap (1993).

示例

Try it in your browser!

假设我们希望测试两个样本是否来自同一分布。假设我们对潜在的分布一无所知，并且在观察数据之前，我们假设第一个样本的均值将小于第二个样本的均值。我们决定使用样本均值之间的差异作为检验统计量，并且我们认为p值为0.05在统计上是显著的。

为了提高效率，我们将定义检验统计量的函数以矢量化方式编写：样本 x 和 y 可以是 ND 数组，并且统计量将沿着 axis 计算每个轴切片。

>>> import numpy as np
>>> def statistic(x, y, axis):
...     return np.mean(x, axis=axis) - np.mean(y, axis=axis)

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

>>> from scipy.stats import norm
>>> rng = np.random.default_rng()
>>> x = norm.rvs(size=5, random_state=rng)
>>> y = norm.rvs(size=6, loc = 3, random_state=rng)
>>> statistic(x, y, 0)
-3.5411688580987266

确实，检验统计量为负，表明 x 所基于的分布的真实均值小于 y 所基于的分布的真实均值。为了确定如果两个样本来自同一分布，这种情况偶然发生的概率，我们进行置换检验。

>>> from scipy.stats import permutation_test
>>> # because our statistic is vectorized, we pass `vectorized=True`
>>> # `n_resamples=np.inf` indicates that an exact test is to be performed
>>> res = permutation_test((x, y), statistic, vectorized=True,
...                        n_resamples=np.inf, alternative='less')
>>> print(res.statistic)
-3.5411688580987266
>>> print(res.pvalue)
0.004329004329004329

在零假设下，获得小于或等于观测值的检验统计量的概率为0.4329%。这低于我们选择的5%的阈值，因此我们认为这是反对零假设而支持备择假设的有力证据。

由于上述样本的规模较小，permutation_test 可以执行精确检验。对于较大的样本，我们采用随机排列检验。

>>> x = norm.rvs(size=100, random_state=rng)
>>> y = norm.rvs(size=120, loc=0.2, random_state=rng)
>>> res = permutation_test((x, y), statistic, n_resamples=9999,
...                        vectorized=True, alternative='less',
...                        random_state=rng)
>>> print(res.statistic)
-0.4230459671240913
>>> print(res.pvalue)
0.0015

在零假设下，获得小于或等于观测值的检验统计量的近似概率为0.0225%。这再次低于我们选择的5%的阈值，因此我们再次有显著证据拒绝零假设，支持备择假设。

对于大样本和大量的排列，结果与相应的渐近检验（独立样本t检验）相当。

>>> from scipy.stats import ttest_ind
>>> res_asymptotic = ttest_ind(x, y, alternative='less')
>>> print(res_asymptotic.pvalue)
0.0014669545224902675

提供了测试统计量的排列分布，以便进一步研究。

>>> import matplotlib.pyplot as plt
>>> plt.hist(res.null_distribution, bins=50)
>>> plt.title("Permutation distribution of test statistic")
>>> plt.xlabel("Value of Statistic")
>>> plt.ylabel("Frequency")
>>> plt.show()

如果统计量由于机器精度有限而存在不准确性，检查零分布是至关重要的。考虑以下情况：

>>> from scipy.stats import pearsonr
>>> x = [1, 2, 4, 3]
>>> y = [2, 4, 6, 8]
>>> def statistic(x, y, axis=-1):
...     return pearsonr(x, y, axis=axis).statistic
>>> res = permutation_test((x, y), statistic, vectorized=True,
...                        permutation_type='pairings',
...                        alternative='greater')
>>> r, pvalue, null = res.statistic, res.pvalue, res.null_distribution

在这种情况下，由于数值噪声，零分布的某些元素与相关系数 r 的观测值不同。我们手动检查零分布中与检验统计量的观测值几乎相同的元素。

>>> r
0.7999999999999999
>>> unique = np.unique(null)
>>> unique
array([-1. , -1. , -0.8, -0.8, -0.8, -0.6, -0.4, -0.4, -0.2, -0.2, -0.2,
    0. ,  0.2,  0.2,  0.2,  0.4,  0.4,  0.6,  0.8,  0.8,  0.8,  1. ,
    1. ])  # may vary
>>> unique[np.isclose(r, unique)].tolist()
[0.7999999999999998, 0.7999999999999999, 0.8]  # may vary

如果 permutation_test 天真地进行比较，那么空分布中值为 0.7999999999999998 的元素将不会被视为极端或比统计量的观测值更极端，因此计算出的 p 值会太小。

>>> incorrect_pvalue = np.count_nonzero(null >= r) / len(null)
>>> incorrect_pvalue
0.14583333333333334  # may vary

相反，permutation_test 将空分布中在 max(1e-14, abs(r)*1e-14) 范围内的元素视为等于统计量 r 的观测值。

>>> correct_pvalue = np.count_nonzero(null >= r - 1e-14) / len(null)
>>> correct_pvalue
0.16666666666666666
>>> res.pvalue == correct_pvalue
True

这种比较方法在大多数实际情况下预计是准确的，但建议用户通过检查接近统计量观测值的零分布元素来评估这一点。此外，考虑使用可以使用精确算术计算的统计量（例如整数统计量）。

Go BackOpen In Tab