scipy.stats.

page_trend_test#

scipy.stats.page_trend_test(data, ranked=False, predicted_ranks=None, method='auto')[source]#

执行 Page 检验,这是一种衡量处理间观察趋势的方法。

Page 检验(也称为 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布尔值,可选

默认情况下,假定 data 是观察结果而不是等级;它将使用 scipy.stats.rankdata 沿 axis=1 进行排序。如果以等级的形式提供 data,则传递参数 True

predicted_ranks类数组,可选

列平均值的预测等级。如果未指定,则假定这些列按预测平均值递增的顺序排列,因此默认的 predicted_ranks\([1, 2, \dots, n-1, n]\)

method{‘auto’、‘asymptotic’、‘exact’},可选

选择用于计算 *p*- 值的方法。以下选项可用。

  • ‘auto’:在 ‘exact’ 和 ‘asymptotic’ 之间进行选择,以在合理的时间内获得相当准确的结果(默认)

  • ‘asymptotic’:将标准化检验统计量与正态分布进行比较

  • ‘exact’:通过将观察到的 \(L\) 统计量与等级的所有可能排列实现的统计量进行比较来计算精确的 *p*- 值(在每个排列同样可能的零假设下)

返回:
resPageTrendTestResult

包含属性的对象

statistic浮点数

Page 的 \(L\) 检验统计量。

pvalue浮点数

关联的 *p*- 值

method{‘asymptotic’、‘exact’}

用于计算 *p*- 值的方法

注释

正如 [1] 中指出的那样,“\(n\) 个‘处理’可以同样代表 \(n\) 个被排序的对象或事件或表现或人或试验。” 类似地,\(m\) 个‘对象’可以同样代表 \(m\) 个“按能力或其他控制变量分组的组,或进行排序的法官,或某种随机复制。”

计算 \(L\) 统计量的过程,改编自 [1],如下:

  1. “用严谨的逻辑预先确定有关实验结果的预测排序的适当假设。如果不知道对任何处理进行排序的合理依据,则 \(L\) 检验是不合适的。”

  2. “与其他实验一样,确定您将在哪个置信度上拒绝实验结果与单调假设没有一致性的零假设。”

  3. “将实验材料放入一个包含 \(n\) 列(处理、排序的对象、条件)和 \(m\) 行(对象、复制组、控制变量的水平)的二维表中。”

  4. “当记录实验观察结果时,对每行中的它们进行排序”,例如 ranks = scipy.stats.rankdata(data, axis=1)

  5. “将每列中的等级相加”,例如 colsums = np.sum(ranks, axis=0)

  6. “将每个等级和乘以同一列的预测等级”,例如 products = predicted_ranks * colsums

  7. “将所有这些乘积相加”,例如 L = products.sum()

[1] 继续建议使用标准化统计量

\[\chi_L^2 = \frac{\left[12L-3mn(n+1)^2\right]^2}{mn^2(n^2-1)(n+1)}\]

“它近似地服从自由度为 1 的卡方分布。通常使用 \(\chi^2\) 表相当于对一致性进行双侧检验。如果需要单侧检验,*几乎总是这样*,则卡方表中发现的概率应*减半*。”

但是,此标准化统计量无法区分观察到的值与预测等级高度相关还是与预测等级*反*相关。相反,我们遵循 [2] 并计算标准化统计量

\[\Lambda = \frac{L - E_0}{\sqrt{V_0}},\]

其中 \(E_0 = \frac{1}{4} mn(n+1)^2\)\(V_0 = \frac{1}{144} mn^2(n+1)(n^2-1)\),“它在零假设下是渐近正态的”。

代码 method='exact' 的 *p*- 值是通过将观察到的 \(L\) 值与为所有 \((n!)^m\) 个可能的等级排列生成的 \(L\) 值进行比较生成的。计算是使用 [5] 的递归方法执行的。

*p*- 值未针对领带的可能性进行调整。当存在领带时,报告的 'exact' *p*- 值可能比真实的 *p*- 值 [2] 稍微大一些(即,更保守)。但是,'asymptotic'` *p*- 值往往小于(即,不太保守)'exact' *p*- 值。

参考

[1] (1,2,3,4)

Ellis Batten Page,“Ordered hypotheses for multiple treatments: a significant test for linear ranks”,*Journal of the American Statistical Association* 58(301), p. 216–230, 1963。

[2] (1,2)

Markus Neuhauser,《非参数统计检验:一种计算方法》,CRC Press,p. 150–152, 2012。

[3] (1,2)

Statext LLC,“Page’s L Trend Test - Easy Statistics”,*Statext - Statistics Study*,https://www.statext.com/practice/PageTrendTest03.php,于 2020 年 7 月 12 日访问。

[4]

“Page’s Trend Test”,*维基百科*,WikimediaFoundation,https://en.wikipedia.org/wiki/Page%27s_trend_test,于 2020 年 7 月 12 日访问。

[5]

Robert E. Odeh,“The exact distribution of Page’s L-statistic in the two-way layout”,*Communications in Statistics - Simulation and Computation*,6(1), p. 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 原始论文 [1] 中的建议使用 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')