Python SciPy stats.quantile_test用法及代碼示例

本文簡要介紹 python 語言中 scipy.stats.quantile_test 的用法。

用法:

scipy.stats.quantile_test(x, *, q=0, p=0.5, alternative='two-sided')#

執行分位數測試並計算分位數的置信區間。

此函數測試原假設，即 q 是與樣本 x 的總體概率 p 相關的分位數值。例如，使用默認參數，它測試 x 背後的總體中位數是否為零。該函數返回一個對象，其中包括檢驗統計量、p 值以及計算分位數周圍置信區間的方法。

參數 ：：

x： array_like

一維樣本。

q： 浮點數，默認值：0

分位數的假設值。

p： 浮點數，默認值：0.5

與分位數相關的概率；即小於 q 的人口比例為 p。必須嚴格介於 0 和 1 之間。

alternative： {‘雙麵’，‘less’, ‘greater’}，可選

定義備擇假設。可以使用以下選項(默認為“雙麵”)：

• “雙麵”：與概率 p 相關的分位數不是 q。

• ‘less’：與概率 p 小於 q 相關的分位數。

• ‘greater’：與概率 p 大於 q 相關的分位數。

返回 ：：

result： QuantileTestResult

具有以下屬性的對象：

統計 浮點數

可用於分位數檢驗的兩個檢驗統計量之一。第一個檢驗統計量，T1, 是樣本的比例x小於或等於假設的分位數q。第二個檢驗統計量，T2, 是樣本的比例x嚴格小於假設的分位數q.

當 alternative = 'greater' 、 T1 用於計算 p 值且 statistic 設置為 T1 時。

當 alternative = 'less' 、 T2 用於計算 p 值且 statistic 設置為 T2 時。

當 alternative = 'two-sided' 時，同時考慮 T1 和 T2，並使用導致最小 p 值的那個。

statistic_type int

任何一個1或者2取決於哪一個T1或者T2用於計算 p 值。

p值 浮點數

與給定替代方案相關的 p 值。

該對象還具有以下方法：

confidence_interval(confidence_level=0.95)

計算與概率相關的總體分位數周圍的置信區間p。置信區間以namedtuple帶字段低的和高的。值為楠當沒有足夠的觀測值來計算所需置信度的置信區間時。

注意：

該檢驗及其計算置信區間的方法是非參數的。當且僅當觀察值獨立同分布時它們才有效。

測試的實施遵循 Conover [1]。考慮兩個檢驗統計量。

T1： 中的觀測值數量x小於或等於q.

T1 = (x <= q).sum()

T2： 中的觀測值數量x嚴格小於q.

T2 = (x < q).sum()

有必要使用兩個檢驗統計量來處理 x 是從離散分布或混合分布生成的可能性。

檢驗的原假設是：

H0: The \(p^{\mathrm{th}}\) population quantile is q.

每個測試統計量的零分布是 \(\mathrm{binom}\left(n, p\right)\) 。當 alternative='less' 時，備擇假設為：

H1: The \(p^{\mathrm{th}}\) population quantile is less than q.

p 值是二項式隨機變量的概率

\[Y \sim \mathrm{binom}\left(n, p\right)\]

大於或等於觀測值T2。

當 alternative='greater' 時，備擇假設為：

H1: The \(p^{\mathrm{th}}\) population quantile is greater than q

p 值是二項式隨機變量 Y 小於或等於觀測值 T1 的概率。

當 alternative='two-sided' 時，備擇假設為

H1: q is not the \(p^{\mathrm{th}}\) population quantile.

對於 'less' 和 'greater' 情況，p 值是較小 p 值的兩倍。對於相同的數據，這兩個 p 值都可能超過 0.5，因此該值被剪裁到區間 \([0, 1]\) 中。

置信區間的方法歸功於 Thompson [2]，後來被證明適用於任何 i.i.d 集合。樣品[3]。該計算基於以下觀察：分位數 \(q\) 大於任何觀察值 \(x_m (1\leq m \leq N)\) 的概率可以計算為

\[\mathbb{P}(x_m \leq q) = 1 - \sum_{k=0}^{m-1} \binom{N}{k} q^k(1-q)^{N-k}\]

默認情況下，置信區間是針對 95% 的置信水平計算的。對 95% 置信區間的常見解釋是，如果 i.i.d.從同一總體中重複抽取樣本，每次都會形成置信區間，置信區間將包含大約 95% 的試驗中指定分位數的真實值。

QuantileNPCI R 包 [4] 中提供了類似的函數。基礎是相同的，但它通過在樣本值之間進行插值來計算置信區間邊界，而該函數僅使用樣本值作為邊界。因此，quantile_test.confidence_interval 返回更保守的間隔(即更大)。

confintr 包 [5] 中包含相同的分位數置信區間計算。

不能保證兩側置信區間是最優的；即，可能存在更緊密的區間，其中可能包含概率大於置信水平的感興趣分位數。在不對樣本進行進一步假設(例如，基礎分布的性質)的情況下，單側間隔是最佳的緊的。

參考：

[1]

23. 10. 科諾弗。實用非參數統計，第三版。 1999 年。

[2]

W. R. Thompson，“關於未知分布形式總體的中位數和其他期望分布的置信範圍”，《數理統計年鑒》，卷。 7、沒有。 3，第 122-128 頁，1936 年，訪問時間：2019 年 9 月 18 日。[在線]。可用：https://www.jstor.org/stable/2957563。

[3]

H. A. David 和 H. N. Nagaraja，“訂單統計中的非參數推理中的訂單統計”，John Wiley & Sons, Ltd，2005 年，第 159-170 頁。可用的：https://onlinelibrary.wiley.com/doi/10.1002/0471722162.ch7.

[4]

N. Hutson、A. Hutson、L. Yan，“QuantileNPCI：分位數的非參數置信區間”，R 包，https://cran.r-project.org/package=QuantileNPCI

[5]

M. Mayer，“confintr：置信區間”，R 包，https://cran.r-project.org/package=confintr

例子：

假設我們希望檢驗總體中位數等於 0.5 的原假設。我們選擇99%的置信度；也就是說，如果 p 值小於 0.01，我們將拒絕原假設並支持替代假設。

當測試中位數為 0.5 的標準均勻分布的隨機變量時，我們期望數據在大多數情況下與原假設一致。

>>> import numpy as np
>>> from scipy import stats
>>> rng = np.random.default_rng()
>>> rvs = stats.uniform.rvs(size=100, random_state=rng)
>>> stats.quantile_test(rvs, q=0.5, p=0.5)
QuantileTestResult(statistic=45, statistic_type=1, pvalue=0.36820161732669576)

正如預期的那樣，p 值不低於我們的閾值 0.01，因此我們不能拒絕原假設。

當測試中位數為 0 的標準正態分布數據時，我們預計零假設會被拒絕。

>>> rvs = stats.norm.rvs(size=100, random_state=rng)
>>> stats.quantile_test(rvs, q=0.5, p=0.5)
QuantileTestResult(statistic=67, statistic_type=2, pvalue=0.0008737198369123724)

事實上，p 值低於我們的閾值 0.01，因此我們拒絕零假設，轉而支持默認的 “two-sided” 替代方案：總體中位數不等於 0.5。

然而，假設我們要針對總體中位數大於 0.5 的片麵替代方案檢驗原假設。由於標準正態分布的中位數小於 0.5，因此我們預計原假設不會被拒絕。

>>> stats.quantile_test(rvs, q=0.5, p=0.5, alternative='greater')
QuantileTestResult(statistic=67, statistic_type=1, pvalue=0.9997956114162866)

毫不奇怪，當 p 值大於我們的閾值時，我們不會拒絕零假設而支持所選的替代方案。

分位數檢驗可用於任何分位數，而不僅僅是中位數。例如，我們可以測試樣本分布的第三個四分位數是否大於 0.6。

>>> rvs = stats.uniform.rvs(size=100, random_state=rng)
>>> stats.quantile_test(rvs, q=0.6, p=0.75, alternative='greater')
QuantileTestResult(statistic=64, statistic_type=1, pvalue=0.00940696592998271)

p 值低於閾值。我們拒絕原假設，轉而采用另一種假設：樣本分布的第三個四分位數大於 0.6。

quantile_test 還可以計算任何分位數的置信區間。

>>> rvs = stats.norm.rvs(size=100, random_state=rng)
>>> res = stats.quantile_test(rvs, q=0.6, p=0.75)
>>> ci = res.confidence_interval(confidence_level=0.95)
>>> ci
ConfidenceInterval(low=0.284491604437432, high=0.8912531024914844)

當測試單邊替代方案時，置信區間包含所有觀察結果，如果作為 q 傳遞，測試的 p 值將大於 0.05，因此不會拒絕原假設。例如：

>>> rvs.sort()
>>> q, p, alpha = 0.6, 0.75, 0.95
>>> res = stats.quantile_test(rvs, q=q, p=p, alternative='less')
>>> ci = res.confidence_interval(confidence_level=alpha)
>>> for x in rvs[rvs <= ci.high]:
...     res = stats.quantile_test(rvs, q=x, p=p, alternative='less')
...     assert res.pvalue > 1-alpha
>>> for x in rvs[rvs > ci.high]:
...     res = stats.quantile_test(rvs, q=x, p=p, alternative='less')
...     assert res.pvalue < 1-alpha

此外，如果為隨機樣本重複生成 95% 置信區間，則置信區間將在大約 95% 的重複中包含真實分位數值。

>>> dist = stats.rayleigh() # our "unknown" distribution
>>> p = 0.2
>>> true_stat = dist.ppf(p) # the true value of the statistic
>>> n_trials = 1000
>>> quantile_ci_contains_true_stat = 0
>>> for i in range(n_trials):
...     data = dist.rvs(size=100, random_state=rng)
...     res = stats.quantile_test(data, p=p)
...     ci = res.confidence_interval(0.95)
...     if ci[0] < true_stat < ci[1]:
...         quantile_ci_contains_true_stat += 1
>>> quantile_ci_contains_true_stat >= 950
True

隻要樣本是獨立同分布的，這適用於任何分布和任何分位數。