SciPy 1.12 中文文档(六十三)
scipy.stats.weightedtau
scipy.stats.weightedtau(x, y, rank=True, weigher=None, additive=True)
计算 Kendall 的加权版本 (\tau)。
加权 (\tau) 是 Kendall (\tau) 的加权版本,在此版本中,高权重的交换比低权重的交换更具影响力。默认参数计算指数加法版本的指数,(\tau_\mathrm h),已被证明在重要和不重要元素之间提供了最佳平衡[1]。
加权是通过一个等级数组和一个称重函数定义的,该函数为每个元素分配基于等级的权重(较重要的等级与较小的值相关联,例如,0 是最高可能的等级),然后交换的权重是交换元素等级的权重的和或乘积。默认参数计算 (\tau_\mathrm h):在等级为 (r) 和 (s)(从零开始)的元素之间的交换的权重为 (1/(r+1) + 1/(s+1))。
指定等级数组只有在您有一个外部重要性标准的情况下才有意义。如果像通常发生的那样,您没有一个具体的等级标准在脑海中,那么加权 (\tau) 就是通过使用 (x, y) 和 (y, x) 递减字典序排名得到的值的平均值来定义的。这是默认参数的行为。请注意,这里用于排名的约定(较低的值意味着更高的重要性)与其他 SciPy 统计函数使用的约定相反。
参数:
x, y数组样本
得分数组,形状相同。如果数组不是 1-D,则将其展平为 1-D。
rank整数数组或布尔值的数组,可选
给每个元素分配一个非负的等级。如果为 None,则将使用递减字典序排名 (x, y):更高等级的元素将是具有更大 x 值的元素,使用 y 值来打破并列(特别地,交换 x 和 y 将产生不同的结果)。如果为 False,则将直接使用元素索引作为等级。默认为 True,此时该函数返回使用 (x, y) 和 (y, x) 递减字典序排名得到的值的平均值。
weigher可调用对象,可选
该称重函数必须将非负整数(零表示最重要的元素)映射到非负权重。默认情况下,None 提供双曲线加权,即,排名 (r) 被映射到权重 (1/(r+1))。
additive布尔值,可选
如果为 True,则交换的权重通过添加交换元素的等级的权重来计算;否则,权重将相乘。默认为 True。
返回:
res: SignificanceResult
包含属性的对象:
statisticfloat
加权的τ相关指数。
pvaluefloat
目前为np.nan,因为统计量的空分布未知(即使在加性双曲线情况下也是如此)。
另见
kendalltau
计算 Kendall's tau。
spearmanr
计算 Spearman 等级相关系数。
theilslopes
计算一组点(x,y)的 Theil-Sen 估计器。
注意
此函数使用基于(O(n \log n))的归并排序算法[1],这是肯德尔τ的 Knight 算法的加权扩展[2]。它可以通过将additive和rank设置为 False 来计算 Shieh 的加权τ[3],用于排名之间无并列(即排列)的情况,因为[1]中给出的定义是 Shieh 的一般化。
NaNs 被认为是最小可能的分数。
0.19.0 版中的新功能。
参考文献
[1] (1,2,3)
Sebastiano Vigna,《带有并列的排名的加权相关指数》,《第 24 届国际万维网会议论文集》,第 1166-1176 页,ACM,2015 年。
[2]
W.R. Knight,《一种计算 Kendall's Tau 的计算机方法,适用于非分组数据》,《美国统计协会杂志》,第 61 卷,第 314 号,第一部分,第 436-439 页,1966 年。
[3]
Grace S. Shieh,《加权的肯德尔τ统计量》,《统计与概率信函》,第 39 卷,第 1 期,第 17-24 页,1998 年。
示例
>>> import numpy as np
>>> from scipy import stats
>>> x = [12, 2, 1, 12, 2]
>>> y = [1, 4, 7, 1, 0]
>>> res = stats.weightedtau(x, y)
>>> res.statistic
-0.56694968153682723
>>> res.pvalue
nan
>>> res = stats.weightedtau(x, y, additive=False)
>>> res.statistic
-0.62205716951801038
NaNs 被认为是最小可能的分数:
>>> x = [12, 2, 1, 12, 2]
>>> y = [1, 4, 7, 1, np.nan]
>>> res = stats.weightedtau(x, y)
>>> res.statistic
-0.56694968153682723
这恰好是 Kendall's tau:
>>> x = [12, 2, 1, 12, 2]
>>> y = [1, 4, 7, 1, 0]
>>> res = stats.weightedtau(x, y, weigher=lambda x: 1)
>>> res.statistic
-0.47140452079103173
>>> x = [12, 2, 1, 12, 2]
>>> y = [1, 4, 7, 1, 0]
>>> stats.weightedtau(x, y, rank=None)
SignificanceResult(statistic=-0.4157652301037516, pvalue=nan)
>>> stats.weightedtau(y, x, rank=None)
SignificanceResult(statistic=-0.7181341329699028, pvalue=nan)
scipy.stats.somersd
原文:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.stats.somersd.html#scipy.stats.somersd
scipy.stats.somersd(x, y=None, alternative='two-sided')
计算 Somers' D,一种有序关联的非对称度量。
像 Kendall's (\tau) 一样,Somers' (D) 是两个排名之间对应的一种度量。这两个统计量都考虑了两个排名 (X) 和 (Y) 中协调和不协调对的差异,并且都被归一化,使得接近 1 的值表示强烈一致,接近-1 的值表示强烈不一致。它们在归一化方式上有所不同。为了显示关系,Somers' (D) 可以用 Kendall's (\tau_a) 定义:
[D(Y|X) = \frac{\tau_a(X, Y)}{\tau_a(X, X)}]
假设第一个排名 (X) 有 (r) 个不同的排名,第二个排名 (Y) 有 (s) 个不同的排名。这两个由 (n) 个排名组成的列表也可以看作是一个 (r \times s) 的列联表,其中元素 (i, j) 是排名 (X) 中排名 (i) 和排名 (Y) 中排名 (j) 的对数。因此,somersd 还允许将输入数据提供为单个的二维列联表,而不是两个分开的一维排名。
注意,Somers' (D) 的定义是非对称的:一般来说,(D(Y|X) \neq D(X|Y))。somersd(x, y) 计算的是 Somers' (D(Y|X)):将“行”变量 (X) 视为独立变量,“列”变量 (Y) 视为依赖变量。要计算 Somers' (D(X|Y)),请交换输入列表或转置输入表。
参数:
xarray_like
1D 排名数组,被视为(行)独立变量。或者,一个二维列联表。
yarray_like, optional
如果 x 是一个一维排名数组,y 是相同长度的一维排名数组,被视为(列)依赖变量。如果 x 是二维的,则忽略 y。
alternative{‘two-sided’, ‘less’, ‘greater’}, optional
定义备择假设。默认为 'two-sided'。可用的选项包括:* 'two-sided': 排名相关不为零 * 'less': 排名相关为负(小于零) * 'greater': 排名相关为正(大于零)
返回:
resSomersDResult
一个 SomersDResult 对象,具有以下字段:
statisticfloat
Somers' (D) 统计量。
pvaluefloat
假设检验的 p 值,其零假设是没有关联,即 (D=0)。更多信息请参阅注释。
table2D array
由排名 x 和 y 形成的列联表(或者如果 x 是二维数组,则为提供的列联表)
参见
kendalltau
计算 Kendall's tau,另一种相关度量。
weightedtau
计算 Kendall’s tau 的加权版本。
spearmanr
计算 Spearman 等级相关系数。
pearsonr
计算 Pearson 相关系数。
注意事项
此函数遵循 [2] 和 [3] 的列联表方法。p-值是基于在零假设 (D=0) 下的检验统计分布的渐近逼近计算的。
理论上,基于 Kendall’s (tau) 和 Somers’ (D) 的假设检验应该是相同的。然而,kendalltau 返回的 p-值基于 (X) 和 (Y) 之间独立性的零假设(即从中抽取 (X) 和 (Y) 对的总体包含所有可能对的等数量),这比此处使用的 (D=0) 的零假设更为具体。如果需要独立性的零假设,则可以使用 kendalltau 返回的 p-值和 somersd 返回的统计量,反之亦然。更多信息,请参阅 [2]。
按照 SAS 和 R 使用的约定格式化列联表:第一个提供的排名(x)是“行”变量,第二个提供的排名(y)是“列”变量。这与 Somers 的原始论文的约定相反 [1]。
参考文献
[1]
Robert H. Somers,《用于序数变量的新的非对称关联度量》,《美国社会学评论》,第 27 卷,第 6 期,799–811 页,1962 年。
[2] (1,2)
Morton B. Brown 和 Jacqueline K. Benedetti,《在二维列联表中检验相关性的抽样行为》,《美国统计协会期刊》第 72 卷,第 358 期,309–315 页,1977 年。
[3]
SAS Institute, Inc.,《频数程序(书摘)》,《SAS/STAT 9.2 用户指南,第二版》,SAS Publishing,2009 年。
[4]
Laerd 统计,《使用 SPSS 统计的 Somers’ d》,《SPSS 统计教程和统计指南》,statistics.laerd.com/spss-tutorials/somers-d-using-spss-statistics.php,访问日期为 2020 年 7 月 31 日。
例子
我们为[4]中的示例计算 Somers' D,其中一位酒店连锁店主想要确定酒店房间清洁度与客户满意度之间的关联。自变量酒店房间清洁度在一个有序尺度上进行排名:“低于平均(1)”,“平均(2)”或“高于平均(3)”。因变量客户满意度在第二个尺度上进行排名:“非常不满意(1)”,“中度不满意(2)”,“既不满意也不满意(3)”,“中度满意(4)”,或“非常满意(5)”。共有 189 位顾客参与了调查,结果转化为一个以酒店房间清洁度为“行”变量和客户满意度为“列”变量的列联表。
| 27 | 25 | 14 | 7 | 0 | |
| 7 | 14 | 18 | 35 | 12 | |
| 1 | 3 | 2 | 7 | 17 |
例如,27 位顾客将其房间的清洁度排名为“低于平均(1)”,相应的满意度为“非常不满意(1)”。我们按以下方式进行分析。
>>> from scipy.stats import somersd
>>> table = [[27, 25, 14, 7, 0], [7, 14, 18, 35, 12], [1, 3, 2, 7, 17]]
>>> res = somersd(table)
>>> res.statistic
0.6032766111513396
>>> res.pvalue
1.0007091191074533e-27
Somers' D 统计量的值约为 0.6,表明样本中房间清洁度与客户满意度之间存在正相关关系。 p-value 非常小,表明在零假设下观察到该统计量极端值的概率非常小(我们的样本来自 189 位顾客,整体人群的统计量为零假设)。这支持备择假设,即人群的真实 Somers' D 值不为零。
scipy.stats.siegelslopes
scipy.stats.siegelslopes(y, x=None, method='hierarchical')
对于点集合(x, y),计算Siegel 估计量。
siegelslopes 实现了使用重复中位数进行鲁棒线性回归的方法(参见[1]),以拟合点集(x, y)的直线。该方法对异常值具有 50%的渐近破坏点。
参数:
y 数组型
因变量。
x 数组型或 None,可选
自变量。如果为 None,则使用arange(len(y))代替。
方法{‘层次化’, ‘分离’}
如果是‘层次化’,使用估计的斜率slope估计截距(默认选项)。如果是‘分离’,独立估计截距。详见注释。
返回:
result SiegelslopesResult 实例
返回值是一个具有以下属性的对象:
斜率浮点数
回归线斜率的估计。
截距浮点数
回归线截距的估计。
另请参阅
theilslopes
一种类似的技术,但没有重复中位数
注释
对于n = len(y),将m_j计算为从点(x[j], y[j])到所有其他n-1点的斜率的中位数。然后slope是所有斜率m_j的中位数。可以通过参数method选择两种估计截距的方法。层次化方法使用估计的斜率slope,计算intercept作为y - slope*x的中位数。另一种方法独立估计截距如下:对于每个点(x[j], y[j]),计算通过其余点的所有n-1条线的截距i_j的中位数。intercept是i_j的中位数。
该实现计算大小为n的向量的中位数n次,对于大向量可能较慢。有更高效的算法(参见[2]),此处未实现。
为了与 SciPy 旧版本兼容,返回值行为类似于长度为 2 的namedtuple,包含字段slope和intercept,因此可以继续写:
slope, intercept = siegelslopes(y, x)
参考文献
[1] (1,2)
A. Siegel,“使用重复中位数的鲁棒回归”,Biometrika,Vol. 69,pp. 242-244,1982 年。
[2]
A. Stein 和 M. Werman,“寻找重复中位数回归线”,第三届 ACM-SIAM 离散算法年会论文集,pp. 409-413,1992 年。
示例
>>> import numpy as np
>>> from scipy import stats
>>> import matplotlib.pyplot as plt
>>> x = np.linspace(-5, 5, num=150)
>>> y = x + np.random.normal(size=x.size)
>>> y[11:15] += 10 # add outliers
>>> y[-5:] -= 7
计算斜率和截距。为了比较,还可以使用linregress计算最小二乘拟合:
>>> res = stats.siegelslopes(y, x)
>>> lsq_res = stats.linregress(x, y)
绘制结果。Siegel 回归线以红色显示。绿色线显示最小二乘拟合以供比较。
>>> fig = plt.figure()
>>> ax = fig.add_subplot(111)
>>> ax.plot(x, y, 'b.')
>>> ax.plot(x, res[1] + res[0] * x, 'r-')
>>> ax.plot(x, lsq_res[1] + lsq_res[0] * x, 'g-')
>>> plt.show()
scipy.stats.theilslopes
scipy.stats.theilslopes(y, x=None, alpha=0.95, method='separate')
计算一组点(x, y)的 Theil-Sen 估计器。
theilslopes 实现了一种鲁棒线性回归的方法。它计算斜率作为所有配对值之间斜率的中位数。
参数:
yarray_like
因变量。
xarray_like 或 None,可选
自变量。如果为 None,则使用arange(len(y))。
alphafloat,可选
置信度在 0 到 1 之间,默认为 95% 置信度。请注意,alpha 对称地围绕 0.5,即 0.1 和 0.9 都被解释为“查找 90% 置信区间”。
方法{‘joint’, ‘separate’},可选
用于计算截距估计的方法。支持以下方法,
- ‘joint’: 使用 np.median(y - slope * x) 作为截距。
- ‘separate’: 使用 np.median(y) - slope * np.median(x)
- 作为截距。
默认值为‘separate’。
版本 1.8.0 中的新功能。
返回:
resultTheilslopesResult 实例
返回值是一个具有以下属性的对象:
斜率 float
Theil 斜率。
截距 float
Theil 线的截距。
低斜率 float
斜率置信区间的下限。
高斜率 float
斜率置信区间的上限。
参见
siegelslopes
使用重复中位数的类似技术
注意事项
theilslopes 的实现遵循 [1]。在 [1] 中未定义截距,在这里定义为 median(y) - slope*median(x),这在 [3] 中给出。文献中也有其他截距的定义,例如在 [4] 中的 median(y - slope*x)。确定如何计算截距可以通过参数 method 来确定。由于文献中未涉及,因此没有给出截距的置信区间。
为了与 SciPy 的旧版本兼容,返回值表现得像一个长度为 4 的 namedtuple,具有字段 slope、intercept、low_slope 和 high_slope,因此可以继续写:
slope, intercept, low_slope, high_slope = theilslopes(y, x)
参考文献
[1] (1,2,3)
P.K. Sen, “基于 Kendall's tau 的回归系数估计”, J. Am. Stat. Assoc., Vol. 63, pp. 1379-1389, 1968.
[2]
H. Theil, “一种秩不变的线性和多项式回归分析方法 I, II 和 III”, Nederl. Akad. Wetensch., Proc. 53:, pp. 386-392, pp. 521-525, pp. 1397-1412, 1950.
[3]
W.L. Conover, “实用非参数统计”, 第 2 版, John Wiley and Sons, 纽约, pp. 493.
[4]
zh.wikipedia.org/wiki/Theil%E2%80%93Sen%E5%9B%9E%E5%BD%92
示例
>>> import numpy as np
>>> from scipy import stats
>>> import matplotlib.pyplot as plt
>>> x = np.linspace(-5, 5, num=150)
>>> y = x + np.random.normal(size=x.size)
>>> y[11:15] += 10 # add outliers
>>> y[-5:] -= 7
计算斜率、截距和 90%置信区间。为了比较,还使用 linregress 计算最小二乘拟合:
>>> res = stats.theilslopes(y, x, 0.90, method='separate')
>>> lsq_res = stats.linregress(x, y)
绘制结果。Theil-Sen 回归线显示为红色,虚线红线表示斜率的置信区间(请注意,虚线红线不是回归的置信区间,因为截距的置信区间未包括在内)。绿色线显示最小二乘拟合以便比较。
>>> fig = plt.figure()
>>> ax = fig.add_subplot(111)
>>> ax.plot(x, y, 'b.')
>>> ax.plot(x, res[1] + res[0] * x, 'r-')
>>> ax.plot(x, res[1] + res[2] * x, 'r--')
>>> ax.plot(x, res[1] + res[3] * x, 'r--')
>>> ax.plot(x, lsq_res[1] + lsq_res[0] * x, 'g-')
>>> plt.show()
scipy.stats.page_trend_test
scipy.stats.page_trend_test(data, ranked=False, predicted_ranks=None, method='auto')
执行 Page 的检验,衡量处理之间观察结果的趋势。
Page’s Test(也称为 Page 的(L)检验)在以下情况下很有用:
-
至少有(n \geq 3)个处理,
-
(m \geq 2)个受试者观察每种处理,并且
-
假设观察结果具有特定的顺序。
具体来说,该检验考虑的是零假设,即
[m_1 = m_2 = m_3 \cdots = m_n,]
其中(m_j)是在处理(j)下观察量的平均值,对立假设是
[m_1 \leq m_2 \leq m_3 \leq \cdots \leq m_n,]
其中至少有一处不等式是严格的。
正如[4]所指出的,Page 的(L)检验在趋势差异的替代假设下比 Friedman 检验具有更强的统计功效,因为 Friedman 的检验只考虑观察值的平均值差异而不考虑它们的顺序。而 Spearman 的(\rho)则考虑两个变量(例如燕子的飞行速度与它所携带的椰子的重量)的排名观察之间的相关性,Page 的(L)则关注观察(例如燕子的飞行速度)在几种不同处理(例如携带不同重量的五个椰子)中的趋势,即使在多个受试者(例如一个欧洲燕子和一个非洲燕子)重复观察的情况下也是如此。
参数:
data类似数组
一个(m \times n)数组;第(i)行第(j)列的元素是与受试者(i)和处理(j)对应的观察结果。默认情况下,假设列按预测均值递增的顺序排列。
ranked布尔值,可选
默认情况下,数据被假定为观察值而不是排名;将使用scipy.stats.rankdata沿axis=1对其进行排名。如果数据以排名形式提供,请传递参数True。
predicted_ranks类似数组,可选
列均值的预测排名。如果未指定,默认假设列按预测均值递增的顺序排列,因此默认的predicted_ranks是([1, 2, \dots, n-1, n])。
method{‘auto’, ‘asymptotic’, ‘exact’},可选
选择用于计算p-值的方法。以下选项可用。
-
‘auto’:在合理时间内选择‘exact’和‘asymptotic’之间以获得合理精度的结果(默认)
-
‘asymptotic’:将标准化的检验统计量与正态分布进行比较
-
‘exact’:通过比较所有可能的排名排列(在零假设下,每个排列等可能)来计算精确的p-值
返回:
resPage 趋势检验结果
一个包含属性的对象:
statisticfloat
Page’s (L) 测试统计量。
pvaluefloat
相关 p-值
方法{‘渐近’, ‘精确’}
用于计算 p-值的方法
另见
rankdata, friedmanchisquare, spearmanr
注释
如 [1] 所述,“这里的 (n) ‘处理’ 也可以表示 (n) 个对象、事件、表演、人员或试验,按排名排序。” 同样,(m) ‘主体’ 也可以等同于能力分组、某种控制变量的分组、进行排名的评委或某种随机复制。
计算 (L) 统计量的过程,改编自 [1],如下:
-
“预先用严谨的逻辑确定关于实验结果预测排序的适当假设。如果没有关于任何处理排序的合理依据,那么 (L) 检验不适用。”
-
“与其他实验一样,确定在何种置信水平下你将拒绝零假设,即实验结果与单调假设不一致。”
-
“将实验材料分类为具有 (n) 列(处理、排名对象、条件)和 (m) 行(主体、复制组、控制变量水平)的二向表。”
-
“记录实验观察时,对每行进行排名”,例如
ranks = scipy.stats.rankdata(data, axis=1)。 -
“对每一列中的排名求和”,例如
colsums = np.sum(ranks, axis=0)。 -
“将每个列的排名总和乘以该列的预测排名”,例如
products = predicted_ranks * colsums。 -
“将所有这类乘积求和”,例如
L = products.sum()。
[1] 进一步建议使用标准化统计量
[\chi_L² = \frac{\left[12L-3mn(n+1)²\right]²}{mn²(n²-1)(n+1)}]
“近似服从自由度为 1 的卡方分布。普通使用 (\chi²) 表相当于进行双侧一致性检验。如果需要进行单侧检验,几乎总是如此,则应将卡方表中的概率 减半。”
然而,这种标准化统计量不能区分观察值是与预测排名良好相关还是与预测排名反相关。因此,我们遵循 [2] 并计算标准化统计量
[\Lambda = \frac{L - E_0}{\sqrt{V_0}},]
其中 (E_0 = \frac{1}{4} mn(n+1)²) 和 (V_0 = \frac{1}{144} mn²(n+1)(n²-1)),“这在零假设下渐近地服从正态分布”。
p-值method='exact'是通过将L的观察值与所有*(n!)^m可能的排名排列生成的L*值进行比较而生成的。计算是使用[5]的递归方法执行的。
p-值未针对出现并列的情况进行调整。当存在并列时,报告的 'exact' p-值可能比真实的p-值稍大(即更保守)。然而,'asymptotic' p-值往往比 'exact' p-值小(即不那么保守)。
参考文献
[1] (1,2,3,4)
Ellis Batten Page,“多重处理的有序假设:线性等级的显著性检验”,美国统计协会杂志 58(301),第 216-230 页,1963 年。
[2] (1,2)
Markus Neuhauser,非参数统计检验:计算方法,CRC Press,第 150-152 页,2012 年。
[3] (1,2)
Statext LLC,“Page's L Trend Test - Easy Statistics”,Statext - 统计学习,www.statext.com/practice/PageTrendTest03.php,访问于 2020 年 7 月 12 日。
[4]
“Page's Trend Test”,维基百科,WikimediaFoundation,en.wikipedia.org/wiki/Page%27s_trend_test,访问于 2020 年 7 月 12 日。
[5]
Robert E. Odeh,“两因素布局中 Page's L 统计量的精确分布”,统计学-模拟与计算,6(1),第 49-61 页,1977 年。
示例
我们使用来自[3]的例子:询问 10 名学生对三种教学方法 - 教程、讲座和研讨会 - 进行 1-5 分的评分,其中 1 分是最低分,5 分是最高分。我们已决定以 99%的置信水平拒绝零假设,支持我们的备择假设:研讨会将获得最高评分,而教程将获得最低评分。最初,数据已经列出,每行表示一个学生对三种方法的评分,顺序如下:教程、讲座、研讨会。
>>> table = [[3, 4, 3],
... [2, 2, 4],
... [3, 3, 5],
... [1, 3, 2],
... [2, 3, 2],
... [2, 4, 5],
... [1, 2, 4],
... [3, 4, 4],
... [2, 4, 5],
... [1, 3, 4]]
因为教程被假设为评分最低,教程排名对应的列应该排在第一位;研讨会被假设为评分最高,所以其列应该排在最后。由于这些列已按照预测均值递增的顺序排列,我们可以直接将表传递给page_trend_test。
>>> from scipy.stats import page_trend_test
>>> res = page_trend_test(table)
>>> res
PageTrendTestResult(statistic=133.5, pvalue=0.0018191161948127822,
method='exact')
这个p-值表明,在零假设下,L统计量达到如此极端值的概率为 0.1819%。因为 0.1819%小于 1%,我们有证据拒绝零假设,支持我们的备择假设,在 99%的置信水平下。
(L) 统计量的值为 133.5. 为了手动验证这一点,我们对数据进行排名,使高分对应高排名,并通过平均排名解决并列情况:
>>> from scipy.stats import rankdata
>>> ranks = rankdata(table, axis=1)
>>> ranks
array([[1.5, 3\. , 1.5],
[1.5, 1.5, 3\. ],
[1.5, 1.5, 3\. ],
[1\. , 3\. , 2\. ],
[1.5, 3\. , 1.5],
[1\. , 2\. , 3\. ],
[1\. , 2\. , 3\. ],
[1\. , 2.5, 2.5],
[1\. , 2\. , 3\. ],
[1\. , 2\. , 3\. ]])
我们在每列内添加排名,将总和乘以预测排名,然后求和。
>>> import numpy as np
>>> m, n = ranks.shape
>>> predicted_ranks = np.arange(1, n+1)
>>> L = (predicted_ranks * np.sum(ranks, axis=0)).sum()
>>> res.statistic == L
True
如在 [3] 中所述,p 值的渐近近似是正态分布的生存函数,其在标准化检验统计量处的值:
>>> from scipy.stats import norm
>>> E0 = (m*n*(n+1)**2)/4
>>> V0 = (m*n**2*(n+1)*(n**2-1))/144
>>> Lambda = (L-E0)/np.sqrt(V0)
>>> p = norm.sf(Lambda)
>>> p
0.0012693433690751756
这与上文由 page_trend_test 报告的 p 值不完全匹配。对于 (m \leq 12) 和 (n \leq 8),渐近分布并不准确,也不保守,因此 page_trend_test 根据表格的维度和 Page 原文中的建议选择了 method='exact'。若要覆盖 page_trend_test 的选择,请提供 method 参数。
>>> res = page_trend_test(table, method="asymptotic")
>>> res
PageTrendTestResult(statistic=133.5, pvalue=0.0012693433690751756,
method='asymptotic')
如果数据已经排名,我们可以传入 ranks 而不是 table 来节省计算时间。
>>> res = page_trend_test(ranks, # ranks of data
... ranked=True, # data is already ranked
... )
>>> res
PageTrendTestResult(statistic=133.5, pvalue=0.0018191161948127822,
method='exact')
假设原始数据的制表顺序与预测均值的顺序不同,比如讲座、研讨会、教程。
>>> table = np.asarray(table)[:, [1, 2, 0]]
由于该表格的排列与假定的顺序不一致,我们可以重新排列表格或提供 predicted_ranks。请记住,预计讲座将排在中间位置,研讨会最高,教程最低,我们传递:
>>> res = page_trend_test(table, # data as originally tabulated
... predicted_ranks=[2, 3, 1], # our predicted order
... )
>>> res
PageTrendTestResult(statistic=133.5, pvalue=0.0018191161948127822,
method='exact')
scipy.stats.multiscale_graphcorr
scipy.stats.multiscale_graphcorr(x, y, compute_distance=<function _euclidean_dist>, reps=1000, workers=1, is_twosamp=False, random_state=None)
计算多尺度图相关(MGC)检验统计量。
具体而言,对于每个点,MGC 找到一个属性的k个最近邻(例如云密度),和另一个属性的l个最近邻(例如草湿度)[1]。这对*(k, l)*被称为“尺度”。然而,事先不知道哪些尺度会最具信息性。因此,MGC 计算所有距离对,然后有效地计算所有尺度的距离相关性。局部相关性显示哪些尺度相对于关系是最具信息性的。因此,成功发现和解释不同数据模态之间关系的关键是自适应确定哪些尺度最具信息性,以及最具信息性尺度的几何含义。这不仅提供了是否模态相关的估计,还揭示了如何进行该决定的见解。在高维数据中尤为重要,因为简单的可视化无法揭示关系给肉眼。特别是,这一实现的表征已经从[2]中得出,并在内部进行了基准测试。
参数:
x, y ndarray
如果x和y的形状为(n, p)和(n, q),其中n是样本数,p和q是维度数,则将运行 MGC 独立性检验。另外,如果x和y的形状为(n, n),并且它们是距离或相似性矩阵,则compute_distance必须发送到None。如果x和y的形状为(n, p)和(m, p),则将运行不配对双样本 MGC 检验。
compute_distance可调用对象,可选
计算每个数据矩阵中样本之间的距离或相似性的函数。如果x和y已经是距离矩阵,则设置为None。默认使用欧氏距离度量。如果调用自定义函数,请先创建距离矩阵或创建形如compute_distance(x)的函数,其中x是计算成对距离的数据矩阵。
reps整数,可选
使用排列测试估计零假设时的复制次数。默认为1000。
workers整数或类似映射的可调用对象,可选
如果 workers 是一个整数,那么将人群细分为 workers 部分,并并行评估(使用 multiprocessing.Pool <multiprocessing>)。提供 -1 来使用所有可用于进程的核心。或者提供一个类似映射的可调用对象,例如 multiprocessing.Pool.map 用于并行评估 p 值。此评估作为 workers(func, iterable) 进行。要求 func 可以被 pickle。默认为 1。
is_twosampbool, optional
如果 True,将运行双样本检验。如果 x 和 y 的形状为 (n, p) 和 (m, p),则此选项将被覆盖并设置为 True。如果 x 和 y 都具有形状 (n, p),并且希望运行双样本检验,则设置为 True。默认为 False。请注意,如果输入为距离矩阵,则不会运行此操作。
random_state{None, int, numpy.random.Generator,
numpy.random.RandomState}, optional
如果 seed 为 None(或 np.random),则使用 numpy.random.RandomState 单例。如果 seed 是一个整数,则使用一个新的 RandomState 实例,其种子为 seed。如果 seed 已经是 Generator 或 RandomState 实例,则使用该实例。
返回:
resMGCResult
包含属性的对象:
statisticfloat
样本 MGC 测试统计量位于 [-1, 1]。
pvaluefloat
通过置换获得的 p 值。
mgc_dictdict
包含额外有用结果:
mgc_mapndarray
关系的潜在几何的二维表示。
opt_scale(int, int)
估计的最优尺度为 (x, y) 对。
null_distlist
来自置换矩阵的空分布。
另请参见
pearsonr
Pearson 相关系数和用于测试非相关性的 p 值。
kendalltau
计算 Kendall's tau。
spearmanr
计算 Spearman 秩相关系数。
注释
MGC 过程及其在神经科学数据上的应用的描述可在 [1] 中找到。它通过以下步骤执行:
-
计算并修改为零均值列的两个距离矩阵 (D^X) 和 (D^Y)。这导致两个 (n \times n) 距离矩阵 (A) 和 (B)(中心化和无偏修改) [3]。
-
对于所有的值 (k) 和 (l),从 (1, ..., n),
-
对于每个属性,计算 (k) 近邻图和 (l) 近邻图。这里,(G_k (i, j)) 表示 (A) 的第 (i) 行的 (k) 个最小值,(H_l (i, j)) 表示 (B) 的第 (i) 行的 (l) 个最小值
-
让 (\circ) 表示逐元素矩阵乘积,然后使用以下统计量对局部相关性进行求和和归一化:
-
[c^{kl} = \frac{\sum_{ij} A G_k B H_l} {\sqrt{\sum_{ij} A² G_k \times \sum_{ij} B² H_l}}]
- MGC 测试统计量是 ({ c^{kl} }) 的平滑最优局部相关性。将平滑操作表示为 (R(\cdot))(本质上将所有孤立的大相关性设置为 0,将连接的大相关性保持不变),见[3]。MGC 是,
[MGC_n (x, y) = \max_{(k, l)} R \left(c^{kl} \left( x_n, y_n \right) \right)]
由于归一化,测试统计量返回一个值在 ((-1, 1)) 之间。
返回的 p 值是使用置换检验计算的。这个过程首先通过随机置换 (y) 来估计零分布,然后计算在零分布下观察到的测试统计量至少与观察到的测试统计量一样极端的概率。
MGC 需要至少 5 个样本才能获得可靠的结果。它还可以处理高维数据集。此外,通过操纵输入数据矩阵,双样本检验问题可以简化为独立性检验问题[4]。给定大小为 (p \times n) 和 (p \times m) 的样本数据 (U) 和 (V),可以如下创建数据矩阵 (X) 和 (Y):
[X = [U | V] \in \mathcal{R}^{p \times (n + m)} Y = [0_{1 \times n} | 1_{1 \times m}] \in \mathcal{R}^{(n + m)}]
然后,MGC 统计量可以像平常一样计算。这种方法可以扩展到类似的测试,比如距离相关性[4]。
1.4.0 版本中的新功能。
参考文献
[1] (1,2)
Vogelstein, J. T., Bridgeford, E. W., Wang, Q., Priebe, C. E., Maggioni, M., & Shen, C. (2019). 发现和解读不同数据模态之间的关系。《ELife》。
[2]
Panda, S., Palaniappan, S., Xiong, J., Swaminathan, A., Ramachandran, S., Bridgeford, E. W., … Vogelstein, J. T. (2019). mgcpy:一个全面的高维独立性检验 Python 包。arXiv:1907.02088
[3] (1,2)
Shen, C., Priebe, C.E., & Vogelstein, J. T. (2019). 从距离相关性到多尺度图相关性。《美国统计协会杂志》。
[4] (1,2)
Shen, C. & Vogelstein, J. T. (2018). 距离和核方法在假设检验中的精确等价性。arXiv:1806.05514
示例
>>> import numpy as np
>>> from scipy.stats import multiscale_graphcorr
>>> x = np.arange(100)
>>> y = x
>>> res = multiscale_graphcorr(x, y)
>>> res.statistic, res.pvalue
(1.0, 0.001)
要运行一个不配对的双样本检验,
>>> x = np.arange(100)
>>> y = np.arange(79)
>>> res = multiscale_graphcorr(x, y)
>>> res.statistic, res.pvalue
(0.033258146255703246, 0.023)
或者,如果输入的形状相同,
>>> x = np.arange(100)
>>> y = x
>>> res = multiscale_graphcorr(x, y, is_twosamp=True)
>>> res.statistic, res.pvalue
(-0.008021809890200488, 1.0)
scipy.stats.chi2_contingency
scipy.stats.chi2_contingency(observed, correction=True, lambda_=None)
列联表中变量独立性的卡方检验。
此函数计算列联表中观察频率独立性的卡方统计量和 p 值 [1] observed 的假设检验。基于独立性假设下的边际和计算期望频率;参见 scipy.stats.contingency.expected_freq。自由度的数量使用 numpy 函数和属性表达:
dof = observed.size - sum(observed.shape) + observed.ndim - 1
参数:
observed:array_like
列联表。表中包含每个类别中的观察频率(即出现次数)。在二维情况下,该表通常描述为“R x C 表”。
correction:bool,可选。
如果为 True,并且自由度为 1,则对连续性应用 Yates 校正。校正的效果是将每个观察值调整 0.5 向相应的期望值。
lambda_:float 或 str,可选。
默认情况下,此测试中计算的统计量是皮尔逊卡方统计量 [2]。lambda_ 允许使用来自 Cressie-Read 功率差异族的统计量 [3]。有关详细信息,请参阅 scipy.stats.power_divergence。
返回:
res:Chi2ContingencyResult
一个包含属性的对象:
统计量:float
测试统计量。
p 值:float
测试的 p 值。
dof:int
自由度。
expected_freq:与observed具有相同的形状。
基于表的边际和的期望频率。
另请参见
scipy.stats.contingency.expected_freq
scipy.stats.fisher_exact
scipy.stats.chisquare
scipy.stats.power_divergence
scipy.stats.barnard_exact
scipy.stats.boschloo_exact
注意事项
这种计算的有效性的常被引用的一个准则是,只有在每个单元格中的观察频率和期望频率至少为 5 时,才应使用该测试。
这是对人口不同类别独立性的检验。当观察维度为二或更多时,该检验才有意义。将该检验应用于一维表将导致期望等于观察且卡方统计量等于 0。
由于计算中存在缺失值,此函数不处理掩码数组。
像 scipy.stats.chisquare 一样,此函数计算卡方统计量;此函数提供的便利性在于从给定的列联表中确定预期频率和自由度。如果这些已知,并且不需要 Yates 修正,可以使用 scipy.stats.chisquare。也就是说,如果调用:
res = chi2_contingency(obs, correction=False)
则以下为真:
(res.statistic, res.pvalue) == stats.chisquare(obs.ravel(),
f_exp=ex.ravel(),
ddof=obs.size - 1 - dof)
lambda_ 参数是在 scipy 的版本 0.13.0 中添加的。
参考文献
[1]
“列联表”,zh.wikipedia.org/wiki/%E5%88%97%E8%81%94%E8%A1%A8
[2]
“皮尔逊卡方检验”,zh.wikipedia.org/wiki/%E7%9A%AE%E5%B0%94%E9%80%8A%E5%8D%A1%E6%96%B9%E6%A3%80%E9%AA%8C
[3]
Cressie, N. 和 Read, T. R. C.,“Multinomial Goodness-of-Fit Tests”,J. Royal Stat. Soc. Series B,Vol. 46, No. 3(1984),pp. 440-464。
[4]
Berger, Jeffrey S. 等人。 “Aspirin for the Primary Prevention of Cardiovascular Events in Women and Men: A Sex-Specific Meta-analysis of Randomized Controlled Trials.” JAMA, 295(3):306-313, DOI:10.1001/jama.295.3.306, 2006。
例子
在[4]中,研究了阿司匹林在预防女性和男性心血管事件中的应用。研究显著结论为:
…阿司匹林疗法通过降低女性缺血性中风的风险,从而减少心血管事件的复合风险 [...]
文章列出了各种心血管事件的研究。我们将重点放在女性的缺血性中风上。
下表总结了参与者连续多年定期服用阿司匹林或安慰剂的实验结果。记录了缺血性中风的案例:
Aspirin Control/Placebo
Ischemic stroke 176 230
No stroke 21035 21018
有证据表明阿司匹林减少了缺血性中风的风险吗?我们首先提出一个零假设 (H_0):
阿司匹林的效果等同于安慰剂。
让我们通过卡方检验来评估这一假设的合理性。
>>> import numpy as np
>>> from scipy.stats import chi2_contingency
>>> table = np.array([[176, 230], [21035, 21018]])
>>> res = chi2_contingency(table)
>>> res.statistic
6.892569132546561
>>> res.pvalue
0.008655478161175739
使用 5%的显著水平,我们将拒绝零假设,支持备择假设:“阿司匹林的效果不等同于安慰剂的效果”。因为scipy.stats.contingency.chi2_contingency执行的是双侧检验,备择假设并不指示效果的方向。我们可以使用stats.contingency.odds_ratio来支持结论,即阿司匹林减少缺血性中风的风险。
下面是进一步的示例,展示如何测试更大的列联表。
一个二路示例(2 x 3):
>>> obs = np.array([[10, 10, 20], [20, 20, 20]])
>>> res = chi2_contingency(obs)
>>> res.statistic
2.7777777777777777
>>> res.pvalue
0.24935220877729619
>>> res.dof
2
>>> res.expected_freq
array([[ 12., 12., 16.],
[ 18., 18., 24.]])
使用对数似然比(即“G 检验”)而不是皮尔逊卡方统计量来进行测试。
>>> res = chi2_contingency(obs, lambda_="log-likelihood")
>>> res.statistic
2.7688587616781319
>>> res.pvalue
0.25046668010954165
一个四路示例(2 x 2 x 2 x 2):
>>> obs = np.array(
... [[[[12, 17],
... [11, 16]],
... [[11, 12],
... [15, 16]]],
... [[[23, 15],
... [30, 22]],
... [[14, 17],
... [15, 16]]]])
>>> res = chi2_contingency(obs)
>>> res.statistic
8.7584514426741897
>>> res.pvalue
0.64417725029295503
scipy.stats.fisher_exact
scipy.stats.fisher_exact(table, alternative='two-sided')
在 2x2 列联表上执行 Fisher 精确检验。
零假设是观察到的表的边际必须等于这些总体的边际条件下,真实几率比是一的真实几率比,并且观察是从这些总体中抽取的。返回的统计量是几率比的无条件最大似然估计,p 值是在零假设下获得至少与实际观察到的表格一样极端的概率。与 Fisher 精确检验相关的统计量和双侧 p 值定义还有其他可能的选择,请参阅注释获取更多信息。
参数:
table整数数组
2x2 列联表。元素必须是非负整数。
alternative{‘two-sided’, ‘less’, ‘greater’},可选
定义备择假设。以下选项可用(默认为‘two-sided’):
-
‘two-sided’: 底层总体的几率比不是一
-
‘less’: 底层总体的几率比一小
-
‘greater’: 底层总体的几率比一大
详细信息请参阅注释。
返回:
resSignificanceResult
包含属性的对象:
统计量 float
这是先前的几率比,而不是后验估计。
p 值 float
在零假设下,获得至少与实际观察到的表格一样极端的概率。
另请参见
chi2_contingency
列联表中变量独立性的卡方检验。当表中的数字较大时,可以用作fisher_exact的替代方法。
contingency.odds_ratio
计算 2x2 列联表的几率比(样本或条件极大似然估计)。
barnard_exact
Barnard 精确检验,对于 2x2 列联表来说比 Fisher 精确检验更为强大的替代方法。
boschloo_exact
Boschloo 精确检验,对于 2x2 列联表来说比 Fisher 精确检验更为强大的替代方法。
注释
零假设和 p 值
零假设是观察下层群体的真实比率为一,且这些观察是从这些群体中随机抽样的条件下成立的:结果表的边际必须与观察表的边际相等。等价地,零假设是输入表来自超几何分布,其参数为 (如 hypergeom 中所用) M = a + b + c + d, n = a + b 和 N = a + c,其中输入表为 [[a, b], [c, d]]。这个分布的支持区间为 max(0, N + n - M) <= x <= min(N, n),或者用输入表中的值来说是 min(0, a - d) <= x <= a + min(b, c)。x 可以解释为一个 2x2 表的左上元素,因此分布中的表格形式为:
[ x n - x ]
[N - x M - (n + N) + x]
例如,如果:
table = [6 2]
[1 4]
那么支持区间为 2 <= x <= 7,并且分布中的表格为:
[2 6] [3 5] [4 4] [5 3] [6 2] [7 1]
[5 0] [4 1] [3 2] [2 3] [1 4] [0 5]
每个表格的概率由超几何分布 hypergeom.pmf(x, M, n, N) 给出。例如,这些分别是(精确到三个有效数字):
x 2 3 4 5 6 7
p 0.0163 0.163 0.408 0.326 0.0816 0.00466
可以用以下方式计算:
>>> import numpy as np
>>> from scipy.stats import hypergeom
>>> table = np.array([[6, 2], [1, 4]])
>>> M = table.sum()
>>> n = table[0].sum()
>>> N = table[:, 0].sum()
>>> start, end = hypergeom.support(M, n, N)
>>> hypergeom.pmf(np.arange(start, end+1), M, n, N)
array([0.01631702, 0.16317016, 0.40792541, 0.32634033, 0.08158508,
0.004662 ])
双侧 p 值是,在零假设下,一个随机表的概率等于或小于输入表的概率。对于我们的示例,输入表的概率(其中 x = 6)为 0.0816。概率不超过这个值的 x 值为 2、6 和 7,因此双侧 p 值为 0.0163 + 0.0816 + 0.00466 ~= 0.10256:
>>> from scipy.stats import fisher_exact
>>> res = fisher_exact(table, alternative='two-sided')
>>> res.pvalue
0.10256410256410257
对于 alternative='greater',单侧 p 值是随机表具有 x >= a 的概率,例如在我们的示例中是 x >= 6,或 0.0816 + 0.00466 ~= 0.08626:
>>> res = fisher_exact(table, alternative='greater')
>>> res.pvalue
0.08624708624708627
这相当于在 x = 5 处计算分布的生存函数(从输入表中减去 x,因为我们想要在总和中包括 x = 6 的概率):
>>> hypergeom.sf(5, M, n, N)
0.08624708624708627
对于 alternative='less',单侧 p 值是随机表具有 x <= a 的概率(例如我们的示例中 x <= 6),或 0.0163 + 0.163 + 0.408 + 0.326 + 0.0816 ~= 0.9949:
>>> res = fisher_exact(table, alternative='less')
>>> res.pvalue
0.9953379953379957
这相当于在 x = 6 处计算分布的累积分布函数:
>>> hypergeom.cdf(6, M, n, N)
0.9953379953379957
比率
计算得到的比率与 R 函数 fisher.test 计算的值不同。此实现返回“样本”或“无条件”最大似然估计,而 R 中的 fisher.test 使用条件最大似然估计。要计算比率的条件最大似然估计,请使用 scipy.stats.contingency.odds_ratio.
参考文献
[1]
费舍尔,罗纳德·A,“实验设计:一位女士品茶的数学。” ISBN 978-0-486-41151-4, 1935.
[2]
“费舍尔精确检验”,zh.wikipedia.org/wiki/费舍尔精确检验
[3]
Emma V. Low 等人,“确定乙酰唑胺预防急性高山病的最低有效剂量:系统评价和荟萃分析”,BMJ,345,DOI:10.1136/bmj.e6779,2012 年。
示例
在 3 中,对乙酰唑胺预防急性高山病的有效剂量进行了研究。研究显著结论如下:
每日服用 250 mg、500 mg 和 750 mg 乙酰唑胺都能有效预防急性高山病。有可用证据表明,乙酰唑胺 250 mg 是这一适应症的最低有效剂量。
以下表格总结了实验结果,一些参与者每日服用 250 mg 乙酰唑胺,而其他参与者服用安慰剂。记录了急性高山病的发病情况:
Acetazolamide Control/Placebo
Acute mountain sickness 7 17
No 15 5
有证据表明乙酰唑胺 250 mg 能减少急性高山病的风险吗?我们首先制定一个零假设 (H_0):
使用乙酰唑胺治疗和使用安慰剂的急性高山病发病几率相同。
让我们用费舍尔检验评估这一假设的可信度。
>>> from scipy.stats import fisher_exact
>>> res = fisher_exact([[7, 17], [15, 5]], alternative='less')
>>> res.statistic
0.13725490196078433
>>> res.pvalue
0.0028841933752349743
使用 5%的显著水平,我们会拒绝零假设,支持备择假设:“与安慰剂相比,使用乙酰唑胺治疗的急性高山病发病几率较低。”
注意
因为费舍尔精确检验的零分布假设是在假定行和列的总和都是固定的情况下形成的,所以在行总和不固定的实验中应用时,其结果是保守的。
在这种情况下,列的总和是固定的;每组有 22 名受试者。但是急性高山病的发病例数却不是(也不能在进行实验前被固定)。这是一个结果。
博斯洛检验不依赖于行总和固定的假设,因此在这种情况下提供了更强大的检验。
>>> from scipy.stats import boschloo_exact
>>> res = boschloo_exact([[7, 17], [15, 5]], alternative='less')
>>> res.statistic
0.0028841933752349743
>>> res.pvalue
0.0015141406667567101
我们验证 p 值小于fisher_exact。
scipy.stats.barnard_exact
scipy.stats.barnard_exact(table, alternative='two-sided', pooled=True, n=32)
对一个 2x2 列联表执行 Barnard 精确检验。
参数:
table 整数的 array_like
一个 2x2 列联表。元素应为非负整数。
alternative {'two-sided', 'less', 'greater'},可选
定义零假设和备择假设。默认为“双侧”。请参阅下面注释部分中的解释。
pooled 布尔值,可选
是否计算具有混合方差(如学生 t 检验中)或非混合方差(如韦尔奇 t 检验中)的分数统计。默认为True。
n 整数,可选
用于构建采样方法的采样点数。请注意,由于使用scipy.stats.qmc.Sobol选择样本点,此参数将自动转换为下一个更高的 2 次幂。默认值为 32。必须为正。在大多数情况下,32 个点足以达到良好的精度。更多的点会带来性能成本。
返回:
ber BarnardExactResult
一个结果对象,具有以下属性。
统计值 浮点数
与用户选择的pooled相对应的具有混合或非混合方差的 Wald 统计量。
p 值浮点数
P 值,即在假设原假设为真的情况下,获得至少与实际观察到的分布一样极端的概率。
另请参见
chi2_contingency
列联表中变量独立性的卡方检验。
fisher_exact
一个 2x2 列联表的 Fisher 精确检验。
boschloo_exact
Boschloo 的 2x2 列联表的精确检验,这是比 Fisher 精确检验更强大的替代方法。
注释
Barnard 检验是用于分析列联表的精确检验。它检验两个分类变量的关联,并且对于 2x2 列联表而言,比 Fisher 精确检验更具有力量。
让我们定义 (X_0) 为一个 2x2 矩阵,表示观察样本,其中每列存储二项实验,如下例所示。我们还定义 (p_1, p_2) 为 (x_{11}) 和 (x_{12}) 的理论二项概率。当使用 Barnard 精确检验时,我们可以断言三种不同的零假设:
-
(H_0 : p_1 \geq p_2) 对 (H_1 : p_1 < p_2),其中 alternative = “less”
-
(H_0 : p_1 \leq p_2) 对 (H_1 : p_1 > p_2),其中 alternative = “greater”
-
(H_0 : p_1 = p_2) 对 (H_1 : p_1 \neq p_2),其中 alternative = “two-sided”(默认值)
为了计算 Barnard's 精确检验,我们使用带有汇总或非汇总方差的 Wald 统计量 [3]。在默认假设下,即两个方差相等(pooled = True),统计量计算如下:
[T(X) = \frac{ \hat{p}_1 - \hat{p}_2 }{ \sqrt{ \hat{p}(1 - \hat{p}) (\frac{1}{c_1} + \frac{1}{c_2}) } }]
其中(\hat{p}_1, \hat{p}_2)和(\hat{p})分别是(p_1, p_2)和(p)的估计量,后者是联合概率,假设(p_1 = p_2)。
如果这个假设无效(pooled = False),则统计量为:
[T(X) = \frac{ \hat{p}_1 - \hat{p}_2 }{ \sqrt{ \frac{\hat{p}_1 (1 - \hat{p}_1)}{c_1} + \frac{\hat{p}_2 (1 - \hat{p}_2)}{c_2} } }]
然后计算 p 值:
[\sum \binom{c_1}{x_{11}} \binom{c_2}{x_{12}} \pi^{x_{11} + x_{12}} (1 - \pi)^{t - x_{11} - x_{12}}]
在此和所有 2x2 列联表(X)的和上,如下:当alternative* = "less"时,* T(X) \leq T(X_0) ,当alternative = "greater"时,* T(X) \geq T(X_0) ,或者 * T(X) \geq |T(X_0)| * 当alternative* = "two-sided"。上面,(c_1, c_2)是第 1 和 2 列的和,(t)是总和(4 个样本元素的和)。
返回的 p 值是在烦扰参数(\pi)上取的最大 p 值,其中(0 \leq \pi \leq 1)。
此函数的复杂度为(O(n c_1 c_2)),其中n是样本点的数量。
参考文献
[1]
Barnard, G. A. “2x2 表的显著性检验”。 Biometrika。 34.1/2 (1947): 123-138. DOI:dpgkg3
[2] (1,2)
Mehta, Cyrus R., 和 Pralay Senchaudhuri. “比较两个二项分布的条件与非条件精确检验”。 Cytel Software Corporation 675 (2003): 1-5.
[3]
“Wald 检验”。 维基百科。 en.wikipedia.org/wiki/Wald_test
例子
在[2]中展示了 Barnard's 检验的一个示例。
考虑疫苗有效性研究的以下示例(Chan, 1998)。在一个 30 名受试者的随机临床试验中,15 名接种了重组 DNA 流感疫苗,另外 15 名接种了安慰剂。安慰剂组中的 15 名受试者中有 12 名最终感染了流感,而对于疫苗组,只有 15 名受试者中的 7 名(47%)感染了流感。数据表现为一个 2 x 2 表格:
Vaccine Placebo Yes 7 12 No 8 3
在进行统计假设检验时,通常使用阈值概率或显著水平来决定是否拒绝零假设(H_0)。假设我们选择了常见的显著性水平 5%。
我们的备择假设是,疫苗将降低感染该病毒的概率;即,接种疫苗后感染病毒的概率(p_1)将小于未接种疫苗后感染病毒的概率(p_2)。因此,我们使用barnard_exact选项alternative="less"调用:
>>> import scipy.stats as stats
>>> res = stats.barnard_exact([[7, 12], [8, 3]], alternative="less")
>>> res.statistic
-1.894...
>>> res.pvalue
0.03407...
在零假设下,即疫苗不会降低感染几率的情况下,获得至少与观察数据一样极端的测试结果的概率约为 3.4%。由于这个 p 值小于我们选择的显著性水平,我们有证据来拒绝 (H_0),支持备择假设。
假设我们使用了费舍尔精确检验:
>>> _, pvalue = stats.fisher_exact([[7, 12], [8, 3]], alternative="less")
>>> pvalue
0.0640...
在相同的显著性阈值 5%下,我们无法拒绝零假设,支持备择假设。正如在[2]中所述,巴纳德检验比费舍尔精确检验更具统计功效,因为巴纳德检验不依赖于任何边际条件。费舍尔检验应仅在两组边际都固定的情况下使用。
scipy.stats.boschloo_exact
scipy.stats.boschloo_exact(table, alternative='two-sided', n=32)
在二维列联表上执行 Boschloo 的精确检验。
参数:
table整数数组
一个二维列联表。元素应为非负整数。
alternative{‘two-sided’, ‘less’, ‘greater’},可选
定义了零假设和备择假设。默认为 'two-sided'。请参见下面的注释部分中的解释。
nint,可选
用于构建抽样方法中使用的抽样点数量。请注意,由于使用 scipy.stats.qmc.Sobol 选择抽样点,因此此参数将自动转换为更高的 2 的幂次方。默认为 32。必须为正数。在大多数情况下,32 个点足以达到良好的精度。更多点会带来性能成本。
返回:
berBoschlooExactResult
一个结果对象,具有以下属性。
统计量:浮点数
Boschloo 检验中使用的统计量;即来自 Fisher 精确检验的 P 值。
P 值:浮点数
P 值,即在假设原假设成立的情况下,观察到至少与实际观察到的分布一样极端的分布的概率。
参见
chi2_contingency
二维列联表中变量独立性的卡方检验。
fisher_exact
Fisher 精确检验在二维列联表上的应用。
barnard_exact
Barnard 的精确检验,它是二维列联表中比 Fisher 精确检验更强大的替代方法。
注释
Boschloo 的检验是用于分析列联表的精确检验。它检验两个分类变量之间的关联,并且是二维列联表中比 Fisher 精确检验更具统一更强大的替代方法。
Boschloo 的精确检验使用 Fisher 精确检验的 P 值作为统计量,而 Boschloo 的 P 值是在零假设下观察到这种统计量的极端值的概率。
让我们定义 (X_0) 为一个表示观察样本的二维矩阵,其中每列存储二项式实验,如下例所示。让我们还定义 (p_1, p_2) 为 (x_{11}) 和 (x_{12}) 的理论二项式概率。在使用 Boschloo 精确检验时,我们可以提出三种不同的备择假设:
-
(H_0 : p_1=p_2) 对 (H_1 : p_1 < p_2),alternative = “less”
-
(H_0 : p_1=p_2) 对 (H_1 : p_1 > p_2),alternative = “greater”
-
(H_0 : p_1=p_2) 对 (H_1 : p_1 \neq p_2),alternative = “two-sided”(默认)
当空值分布不对称时,计算双边 p 值的多种约定。在这里,我们应用这样一种约定,即双边检验的 p 值是单边检验 p 值的两倍(截断为 1.0)。请注意,fisher_exact遵循不同的约定,因此对于给定的table,boschloo_exact报告的统计量可能与fisher_exact报告的 p 值不同,当alternative='two-sided'时。
新版本 1.7.0 中的新增内容。
参考文献
[1]
R.D. Boschloo,“在检验两个概率相等时提升 2 x 2 表的条件显著水平”,Statistica Neerlandica,24(1),1970 年
[2]
“Boschloo's test”,维基百科,en.wikipedia.org/wiki/Boschloo%27s_test
[3]
Lise M. Saari 等人,“员工态度和工作满意度”,人力资源管理,43(4),395-407,2004 年,DOI:10.1002/hrm.20032。
示例
在下面的例子中,我们考虑了文章“员工态度和工作满意度”[3],该文章报告了对 63 名科学家和 117 名大学教授进行的调查结果。在 63 名科学家中,有 31 名表示他们对工作非常满意,而在 117 名大学教授中,有 74 名表示他们对工作非常满意。这是否是大学教授比科学家更满意他们的工作的重要证据?下表总结了上述数据:
college professors scientists
Very Satisfied 74 31
Dissatisfied 43 32
在进行统计假设检验时,我们通常会选择一个阈值概率或显著性水平,用来决定是否拒绝零假设(H_0)。假设我们选择常见的 5%显著性水平。
我们的备择假设是大学教授对他们的工作更满意,而不是科学家。因此,我们期望(p_1)非常满意的大学教授的比例要大于(p_2),即非常满意的科学家的比例。因此,我们调用boschloo_exact并选择alternative="greater"选项:
>>> import scipy.stats as stats
>>> res = stats.boschloo_exact([[74, 31], [43, 32]], alternative="greater")
>>> res.statistic
0.0483...
>>> res.pvalue
0.0355...
在零假设下,即科学家比大学教授在工作中更快乐,获得至少与观察数据一样极端测试结果的概率约为 3.55%。由于此 p 值小于我们选择的显著性水平,我们有证据拒绝(H_0),支持备择假设。
scipy.stats.ttest_ind_from_stats
scipy.stats.ttest_ind_from_stats(mean1, std1, nobs1, mean2, std2, nobs2, equal_var=True, alternative='two-sided')
两个独立样本的均值 t 检验,从描述统计学数据。
这是检验两个独立样本具有相同平均(期望)值的零假设的检验。
参数:
mean1array_like
样本 1 的均值。
std1array_like
样本 1 的修正样本标准差(即ddof=1)。
nobs1array_like
样本 1 的观察次数。
mean2array_like
样本 2 的均值。
std2array_like
样本 2 的修正样本标准差(即ddof=1)。
nobs2array_like
样本 2 的观察次数。
equal_varbool, optional
如果为 True(默认),执行假设总体方差相等的标准独立两样本检验[1]。如果为 False,执行不假设总体方差相等的 Welch's t 检验[2]。
alternative{‘two-sided’, ‘less’, ‘greater’}, optional
定义备择假设。以下选项可用(默认为‘two-sided’):
-
‘two-sided’: 分布的均值不相等。
-
‘less’: 第一个分布的平均值小于第二个分布的平均值。
-
‘greater’: 第一个分布的平均值大于第二个分布的平均值。
1.6.0 版中的新功能。
返回:
statisticfloat or array
计算得到的 t 统计量。
pvaluefloat or array
双侧 p 值。
另请参阅
scipy.stats.ttest_ind
注释
统计量计算为(mean1 - mean2)/se,其中se为标准误差。因此,当mean1大于mean2时,统计量为正;当mean1小于mean2时,统计量为负。
此方法不会检查std1或std2的任何元素是否为负数。如果在调用此方法时std1或std2的任何元素为负数,则此方法将返回与分别传递numpy.abs(std1)和numpy.abs(std2)相同的结果;不会抛出异常或警告。
参考文献
[1]
[zh.wikipedia.org/wiki/T 檢定#獨立樣本 t 檢定](zh.wikipedia.org/wiki/T 檢定#獨立樣本 t 檢定)
[2]
[zh.wikipedia.org/wiki/Welch's_t 检验](zh.wikipedia.org/wiki/Welch'… 检验)
示例
假设我们有两个样本的汇总数据如下(其中样本方差为修正的样本方差):
Sample Sample
Size Mean Variance
Sample 1 13 15.0 87.5
Sample 2 11 12.0 39.0
对这些数据应用 t 检验(假设总体方差相等):
>>> import numpy as np
>>> from scipy.stats import ttest_ind_from_stats
>>> ttest_ind_from_stats(mean1=15.0, std1=np.sqrt(87.5), nobs1=13,
... mean2=12.0, std2=np.sqrt(39.0), nobs2=11)
Ttest_indResult(statistic=0.9051358093310269, pvalue=0.3751996797581487)
对比起来,这是摘要统计数据来自的数据。利用这些数据,我们可以使用scipy.stats.ttest_ind计算相同的结果:
>>> a = np.array([1, 3, 4, 6, 11, 13, 15, 19, 22, 24, 25, 26, 26])
>>> b = np.array([2, 4, 6, 9, 11, 13, 14, 15, 18, 19, 21])
>>> from scipy.stats import ttest_ind
>>> ttest_ind(a, b)
Ttest_indResult(statistic=0.905135809331027, pvalue=0.3751996797581486)
假设我们有二进制数据,并希望应用 t 检验来比较两个独立组中 1 的比例:
Number of Sample Sample
Size ones Mean Variance
Sample 1 150 30 0.2 0.161073
Sample 2 200 45 0.225 0.175251
样本均值 (\hat{p}) 是样本中 1 的比例,而二进制观察的方差由 (\hat{p}(1-\hat{p})) 估算。
>>> ttest_ind_from_stats(mean1=0.2, std1=np.sqrt(0.161073), nobs1=150,
... mean2=0.225, std2=np.sqrt(0.175251), nobs2=200)
Ttest_indResult(statistic=-0.5627187905196761, pvalue=0.5739887114209541)
对比起来,我们可以使用 0 和 1 的数组以及scipy.stat.ttest_ind计算 t 统计量和 p 值,就像上面一样。
>>> group1 = np.array([1]*30 + [0]*(150-30))
>>> group2 = np.array([1]*45 + [0]*(200-45))
>>> ttest_ind(group1, group2)
Ttest_indResult(statistic=-0.5627179589855622, pvalue=0.573989277115258)
scipy.stats.poisson_means_test
scipy.stats.poisson_means_test(k1, n1, k2, n2, *, diff=0, alternative='two-sided')
执行泊松均值检验,又称“E-测试”。
This is a test of the null hypothesis that the difference between means of two Poisson distributions is diff. The samples are provided as the number of events k1 and k2 observed within measurement intervals (e.g. of time, space, number of observations) of sizes n1 and n2.
Parameters:
k1int
Number of events observed from distribution 1.
n1: float
Size of sample from distribution 1.
k2int
Number of events observed from distribution 2.
n2float
Size of sample from distribution 2.
difffloat, default=0
The hypothesized difference in means between the distributions underlying the samples.
alternative{‘two-sided’, ‘less’, ‘greater’}, optional
Defines the alternative hypothesis. The following options are available (default is ‘two-sided’):
- ‘two-sided’: the difference between distribution means is not equal to diff
- ‘less’: the difference between distribution means is less than diff
- ‘greater’: the difference between distribution means is greater than diff
Returns:
statisticfloat
测试统计量(见[1] 方程式 3.3)。
pvaluefloat
The probability of achieving such an extreme value of the test statistic under the null hypothesis.
Notes
Let:
[X_1 \sim \mbox{Poisson}(\mathtt{n1}\lambda_1)]
be a random variable independent of
[X_2 \sim \mbox{Poisson}(\mathtt{n2}\lambda_2)]
and let k1 and k2 be the observed values of (X_1) and (X_2), respectively. Then poisson_means_test uses the number of observed events k1 and k2 from samples of size n1 and n2, respectively, to test the null hypothesis that
[H_0: \lambda_1 - \lambda_2 = \mathtt{diff}]
A benefit of the E-test is that it has good power for small sample sizes, which can reduce sampling costs [1]. It has been evaluated and determined to be more powerful than the comparable C-test, sometimes referred to as the Poisson exact test.
References
[1] (1,2)
Krishnamoorthy, K., & Thomson, J. (2004). A more powerful test for comparing two Poisson means. Journal of Statistical Planning and Inference, 119(1), 23-35.
[2]
Przyborowski, J., & Wilenski, H. (1940). Homogeneity of results in testing samples from Poisson series: With an application to testing clover seed for dodder. Biometrika, 31(3/4), 313-323.
Examples
假设一个园艺师希望测试从种子公司购买的苜蓿种子袋中的病草(杂草)种子数量。先前已经确定苜蓿中病草种子的数量服从泊松分布。
从袋子中取出 100 克样本,并在运送给园丁之前进行分析。样本经分析后发现不含有爬根藤种子;也就是说,k1为 0。然而,园丁到货后又从袋中取出 100 克样本。这次,在样本中发现了三颗爬根藤种子;也就是说,k2为 3。园丁想知道这种差异是否显著且不是由于偶然因素引起的。零假设是两个样本之间的差异仅仅是由于偶然因素引起的,即 (\lambda_1 - \lambda_2 = \mathtt{diff}),其中 (\mathtt{diff} = 0)。备择假设是差异不是由偶然因素引起的,即 (\lambda_1 - \lambda_2 \ne 0)。园丁选择了 5%的显著水平,以拒绝零假设,支持备择假设[2]。
>>> import scipy.stats as stats
>>> res = stats.poisson_means_test(0, 100, 3, 100)
>>> res.statistic, res.pvalue
(-1.7320508075688772, 0.08837900929018157)
P 值为 0.088,表明在零假设下观察到测试统计量的值的几率接近 9%。这超过了 5%,因此园丁不拒绝零假设,因为在这个水平上不能认为差异是显著的。
scipy.stats.ttest_ind
原文链接:
docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.stats.ttest_ind.html#scipy.stats.ttest_ind
scipy.stats.ttest_ind(a, b, axis=0, equal_var=True, nan_policy='propagate', permutations=None, random_state=None, alternative='two-sided', trim=0, *, keepdims=False)
计算两个独立样本的均值的 T 检验。
这是一个测试,用于检验两个独立样本的均值(期望值)是否相同的空假设。该测试默认假设总体具有相同的方差。
参数:
a, b:数组类型
数组必须具有相同的形状,除了与 axis 对应的维度(默认为第一维)。
axis:整数或 None,默认为 0
如果是整数,则计算输入的轴(例如行)上的统计量。输入的每个轴切片的统计量将出现在输出的相应元素中。如果为 None,则在计算统计量之前会对输入进行展平。
equal_var:布尔值,可选
如果为 True(默认),执行假设两个独立样本具有相等总体方差的标准独立 2 样本测试 [1]。如果为 False,则执行威尔奇 t 检验,该检验不假设相等的总体方差 [2]。
自版本 0.11.0 新增。
nan_policy:{‘propagate’, ‘omit’, ‘raise’}
定义如何处理输入中的 NaN。
-
propagate:如果在计算统计量的轴切片(例如行)中存在 NaN,则输出的相应条目将为 NaN。 -
omit:在执行计算时将 NaN 剔除。如果在计算统计量的轴切片中剩余的数据不足,则输出的相应条目将为 NaN。 -
raise:如果存在 NaN,则会引发ValueError。
permutations:非负整数、np.inf 或 None(默认),可选
如果为 0 或 None(默认),则使用 t 分布计算 p 值。否则,permutations 是用于使用排列检验估计 p 值的随机排列次数。如果 permutations 等于或超过池化数据的不同分区数,则执行精确测试(即每个不同分区仅使用一次)。有关详细信息,请参阅注释。
自版本 1.7.0 新增。
random_state:{None, 整数, numpy.random.Generator
如果 seed 为 None(或 np.random),则使用 numpy.random.RandomState 单例。如果 seed 是整数,则使用新的 RandomState 实例,并使用 seed 进行种子化。如果 seed 已经是 Generator 或 RandomState 实例,则使用该实例。
用于生成排列的伪随机数生成器状态(仅在 permutations 不为 None 时使用)。
1.7.0 版本中的新功能。
alternative{‘two-sided’,‘less’,‘greater’},可选
定义了备择假设。以下选项可用(默认为‘双侧’):
-
‘two-sided’:样本底层分布的均值不相等。
-
‘less’:第一个样本潜在分布的均值小于第二个样本潜在分布的均值。
-
‘greater’:第一个样本潜在分布的平均值大于第二个样本潜在分布的平均值。
1.6.0 版本中的新功能。
trimfloat,可选
如果非零,执行修剪(Yuen’s)t 检验。定义从输入样本的每端修剪的元素比例。如果为 0(默认),则不会从任何一侧修剪元素。每个尾部修剪元素的数量是修剪次数乘以元素数量的地板值。有效范围为 0,.5)。
1.7 版本中的新功能。
keepdimsbool,默认值:False
如果设置为 True,则减少的轴作为维度大小为一的结果保留在结果中。使用此选项,结果将正确广播到输入数组。
返回:
result[TtestResult
具有以下属性的对象:
statisticfloat 或 ndarray
t 统计量。
pvaluefloat 或 ndarray
与给定备择假设相关联的 p 值。
dffloat 或 ndarray
用于计算 t 统计量的自由度数。对于排列 t 检验,此值始终为 NaN。
1.11.0 版本中的新功能。
该对象还具有以下方法:
confidence_interval(confidence_level=0.95)
为给定置信水平计算两总体均值差异的置信区间。置信区间以namedtuple返回,具有low和high字段。执行排列 t 检验时,不计算置信区间,low和high字段包含 NaN。
1.11.0 版本中的新功能。
注释
假设我们观察到两个独立样本,例如花瓣长度,并且我们正在考虑这两个样本是否来自同一总体(例如同一种花或具有相似花瓣特征的两种物种)或两个不同的总体。
t 检验量化两样本算术均值之间的差异。p 值量化在假设空白假设为真的情况下观察到的或更极端值的概率,即样本来自具有相同总体均值的总体。大于所选阈值的 p 值(例如 5%或 1%)表示我们的观察不太可能是偶然发生的。因此,我们不拒绝等同总体均值的空白假设。如果 p 值小于我们的阈值,则我们有反对等同总体均值的空白假设的证据。
默认情况下,通过将观察数据的 t 统计量与理论 t 分布进行比较来确定 p 值。当 1 < permutations < binom(n, k) 时,其中
-
k是 a 中的观察次数, -
n是 a 和 b 中的总观察数, -
binom(n, k)是二项式系数(n选k),
数据被汇总(连接起来),随机分配到 a 组或 b 组,并计算 t 统计量。这个过程重复进行(permutation 次),生成零假设下 t 统计量的分布,将观察数据的 t 统计量与此分布进行比较,以确定 p 值。具体来说,报告的 p 值是在 [3] 第 4.4 节中定义的 “实现显著性水平”(ASL)。请注意,还有其他使用随机置换测试估计 p 值的方法;有关其他选项,请参见更一般的 permutation_test。
当 permutations >= binom(n, k) 时,进行精确检验:数据按每种不同方式精确分组一次。
置换检验可能计算成本高,并且不一定比分析检验更准确,但它不对基础分布的形状做出强烈假设。
常见的修剪方法被称为修剪 t 检验。有时被称为尤恩 t 检验,这是 Welch t 检验的扩展,区别在于在方差计算中使用修剪平均数以及在统计量计算中使用修剪样本大小。如果基础分布呈长尾分布或受离群值污染,建议使用修剪方法 [4]。
统计量计算为 (np.mean(a) - np.mean(b))/se,其中 se 是标准误差。因此,当 a 的样本均值大于 b 的样本均值时,统计量为正;当 a 的样本均值小于 b 的样本均值时,统计量为负。
从 SciPy 1.9 开始,将不建议使用的 np.matrix 输入转换为 np.ndarray 后再执行计算。在这种情况下,输出将是一个适当形状的标量或 np.ndarray,而不是 2D 的 np.matrix。类似地,虽然忽略掩码数组的掩码元素,但输出将是一个标量或 np.ndarray,而不是具有 mask=False 的掩码数组。
参考文献
[1]
en.wikipedia.org/wiki/T-test#Independent_two-sample_t-test
[2]
en.wikipedia.org/wiki/Welch%27s_t-test
[3]
- Efron 和 T. Hastie. Computer Age Statistical Inference. (2016).
[4]
Yuen, Karen K. “不等总体方差的两样本修剪 t。”Biometrika,vol. 61,no. 1,1974 年,pp. 165-170。JSTOR,www.jstor.org/stable/2334… 年 3 月 30 日。
[5]
Yuen, Karen K.和 W.J. Dixon. “两样本修剪 t 的近似行为和性能。”Biometrika,vol. 60,no. 2,1973 年,pp. 369-374。JSTOR,www.jstor.org/stable/2334… 年 3 月 30 日。
示例
>>> import numpy as np
>>> from scipy import stats
>>> rng = np.random.default_rng()
具有相同均值的样本的检验:
>>> rvs1 = stats.norm.rvs(loc=5, scale=10, size=500, random_state=rng)
>>> rvs2 = stats.norm.rvs(loc=5, scale=10, size=500, random_state=rng)
>>> stats.ttest_ind(rvs1, rvs2)
Ttest_indResult(statistic=-0.4390847099199348, pvalue=0.6606952038870015)
>>> stats.ttest_ind(rvs1, rvs2, equal_var=False)
Ttest_indResult(statistic=-0.4390847099199348, pvalue=0.6606952553131064)
ttest_ind 低估了不等方差情况下的 p 值:
>>> rvs3 = stats.norm.rvs(loc=5, scale=20, size=500, random_state=rng)
>>> stats.ttest_ind(rvs1, rvs3)
Ttest_indResult(statistic=-1.6370984482905417, pvalue=0.1019251574705033)
>>> stats.ttest_ind(rvs1, rvs3, equal_var=False)
Ttest_indResult(statistic=-1.637098448290542, pvalue=0.10202110497954867)
当n1 != n2时,等方差 t 统计量不再等于不等方差 t 统计量:
>>> rvs4 = stats.norm.rvs(loc=5, scale=20, size=100, random_state=rng)
>>> stats.ttest_ind(rvs1, rvs4)
Ttest_indResult(statistic=-1.9481646859513422, pvalue=0.05186270935842703)
>>> stats.ttest_ind(rvs1, rvs4, equal_var=False)
Ttest_indResult(statistic=-1.3146566100751664, pvalue=0.1913495266513811)
具有不同均值、方差和 n 的 t 检验:
>>> rvs5 = stats.norm.rvs(loc=8, scale=20, size=100, random_state=rng)
>>> stats.ttest_ind(rvs1, rvs5)
Ttest_indResult(statistic=-2.8415950600298774, pvalue=0.0046418707568707885)
>>> stats.ttest_ind(rvs1, rvs5, equal_var=False)
Ttest_indResult(statistic=-1.8686598649188084, pvalue=0.06434714193919686)
在执行排列测试时,更多的排列通常会产生更准确的结果。使用np.random.Generator来确保可重复性:
>>> stats.ttest_ind(rvs1, rvs5, permutations=10000,
... random_state=rng)
Ttest_indResult(statistic=-2.8415950600298774, pvalue=0.0052994700529947)
取这两个样本,其中一个有一个极端的尾部。
>>> a = (56, 128.6, 12, 123.8, 64.34, 78, 763.3)
>>> b = (1.1, 2.9, 4.2)
使用trim关键字执行修剪(Yuen)t 检验。例如,使用 20%修剪,trim=.2,测试将从样本a的每个尾部减少一个元素(np.floor(trim*len(a)))。它对样本b没有影响,因为np.floor(trim*len(b))为 0。
>>> stats.ttest_ind(a, b, trim=.2)
Ttest_indResult(statistic=3.4463884028073513,
pvalue=0.01369338726499547)