分析單一樣本

首先，我們建立一些隨機變數。我們設定一個種子，以便每次執行都能得到相同的結果來查看。作為範例，我們從學生 t 分佈中取樣。

>>> import numpy as np
>>> import scipy.stats as stats
>>> x = stats.t.rvs(10, size=1000)

在這裡，我們將 t 分佈所需的形狀參數（在統計學中對應於自由度）設定為 10。使用 size=1000 表示我們的樣本包含 1000 個獨立抽取的（偽）隨機數。由於我們沒有指定關鍵字參數 loc 和 scale，因此它們被設定為預設值零和一。

描述性統計

x 是一個 numpy 陣列，我們可以直接存取所有陣列方法，例如：

>>> print(x.min())   # equivalent to np.min(x)
-3.78975572422  # random
>>> print(x.max())   # equivalent to np.max(x)
5.26327732981  # random
>>> print(x.mean())  # equivalent to np.mean(x)
0.0140610663985  # random
>>> print(x.var())   # equivalent to np.var(x))
1.28899386208  # random

樣本屬性與其理論對應物相比如何？

>>> m, v, s, k = stats.t.stats(10, moments='mvsk')
>>> n, (smin, smax), sm, sv, ss, sk = stats.describe(x)

>>> sstr = '%-14s mean = %6.4f, variance = %6.4f, skew = %6.4f, kurtosis = %6.4f'
>>> print(sstr % ('distribution:', m, v, s ,k))
distribution:  mean = 0.0000, variance = 1.2500, skew = 0.0000, kurtosis = 1.0000  # random
>>> print(sstr % ('sample:', sm, sv, ss, sk))
sample:        mean = 0.0141, variance = 1.2903, skew = 0.2165, kurtosis = 1.0556  # random

注意：stats.describe 使用變異數的無偏估計量，而 np.var 是有偏估計量。

對於我們的樣本，樣本統計量與其理論對應物之間只有些微差異。

T 檢定和 KS 檢定

我們可以使用 t 檢定來檢驗樣本的平均值是否與理論期望值存在顯著的統計差異。

>>> print('t-statistic = %6.3f pvalue = %6.4f' %  stats.ttest_1samp(x, m))
t-statistic =  0.391 pvalue = 0.6955  # random

p 值為 0.7，這表示在例如 10% 的 alpha 誤差下，我們無法拒絕樣本平均值等於零（標準 t 分佈的期望值）的假設。

作為練習，我們也可以直接計算 t 檢定，而無需使用提供的函數，這應該會給我們相同的答案，而且確實如此

>>> tt = (sm-m)/np.sqrt(sv/float(n))  # t-statistic for mean
>>> pval = stats.t.sf(np.abs(tt), n-1)*2  # two-sided pvalue = Prob(abs(t)>tt)
>>> print('t-statistic = %6.3f pvalue = %6.4f' % (tt, pval))
t-statistic =  0.391 pvalue = 0.6955  # random

Kolmogorov-Smirnov 檢定可用於檢驗樣本是否來自標準 t 分佈的假設

>>> print('KS-statistic D = %6.3f pvalue = %6.4f' % stats.kstest(x, 't', (10,)))
KS-statistic D =  0.016 pvalue = 0.9571  # random

同樣地，p 值夠高，我們無法拒絕隨機樣本確實是根據 t 分佈分佈的假設。在實際應用中，我們不知道底層分佈是什麼。如果我們對樣本執行針對標準常態分佈的 Kolmogorov-Smirnov 檢定，那麼我們也無法拒絕樣本是由常態分佈產生的假設，因為在這個範例中，p 值幾乎為 40%。

>>> print('KS-statistic D = %6.3f pvalue = %6.4f' % stats.kstest(x, 'norm'))
KS-statistic D =  0.028 pvalue = 0.3918  # random

然而，標準常態分佈的變異數為 1，而我們的樣本的變異數為 1.29。如果我們將樣本標準化並針對常態分佈進行檢定，那麼 p 值再次夠大，我們無法拒絕樣本來自常態分佈的假設。

>>> d, pval = stats.kstest((x-x.mean())/x.std(), 'norm')
>>> print('KS-statistic D = %6.3f pvalue = %6.4f' % (d, pval))
KS-statistic D =  0.032 pvalue = 0.2397  # random

注意：Kolmogorov-Smirnov 檢定假設我們針對具有給定參數的分佈進行檢定，因為在最後一種情況下，我們估計了平均值和變異數，因此此假設被違反，並且檢定統計量的分佈（p 值基於此分佈）不正確。

分佈的尾部

最後，我們可以檢查分佈的上尾部。我們可以使用百分點函數 ppf（它是 cdf 函數的反函數）來獲得臨界值，或者更直接地，我們可以使用生存函數的反函數

>>> crit01, crit05, crit10 = stats.t.ppf([1-0.01, 1-0.05, 1-0.10], 10)
>>> print('critical values from ppf at 1%%, 5%% and 10%% %8.4f %8.4f %8.4f' % (crit01, crit05, crit10))
critical values from ppf at 1%, 5% and 10%   2.7638   1.8125   1.3722
>>> print('critical values from isf at 1%%, 5%% and 10%% %8.4f %8.4f %8.4f' % tuple(stats.t.isf([0.01,0.05,0.10],10)))
critical values from isf at 1%, 5% and 10%   2.7638   1.8125   1.3722

>>> freq01 = np.sum(x>crit01) / float(n) * 100
>>> freq05 = np.sum(x>crit05) / float(n) * 100
>>> freq10 = np.sum(x>crit10) / float(n) * 100
>>> print('sample %%-frequency at 1%%, 5%% and 10%% tail %8.4f %8.4f %8.4f' % (freq01, freq05, freq10))
sample %-frequency at 1%, 5% and 10% tail   1.4000   5.8000  10.5000  # random

在這三種情況下，我們的樣本在上尾部中的權重都比底層分佈更大。我們可以簡單地檢查一個更大的樣本，看看是否能更接近。在這種情況下，經驗頻率非常接近理論機率，但如果我們重複多次，波動仍然相當大。

>>> freq05l = np.sum(stats.t.rvs(10, size=10000) > crit05) / 10000.0 * 100
>>> print('larger sample %%-frequency at 5%% tail %8.4f' % freq05l)
larger sample %-frequency at 5% tail   4.8000  # random

我們也可以將其與常態分佈的尾部進行比較，常態分佈在尾部中的權重較小

>>> print('tail prob. of normal at 1%%, 5%% and 10%% %8.4f %8.4f %8.4f' %
...       tuple(stats.norm.sf([crit01, crit05, crit10])*100))
tail prob. of normal at 1%, 5% and 10%   0.2857   3.4957   8.5003

卡方檢定可用於檢驗對於有限數量的組距，觀察到的頻率是否與假設分佈的機率顯著不同。

>>> quantiles = [0.0, 0.01, 0.05, 0.1, 1-0.10, 1-0.05, 1-0.01, 1.0]
>>> crit = stats.t.ppf(quantiles, 10)
>>> crit
array([       -inf, -2.76376946, -1.81246112, -1.37218364,  1.37218364,
        1.81246112,  2.76376946,         inf])
>>> n_sample = x.size
>>> freqcount = np.histogram(x, bins=crit)[0]
>>> tprob = np.diff(quantiles)
>>> nprob = np.diff(stats.norm.cdf(crit))
>>> tch, tpval = stats.chisquare(freqcount, tprob*n_sample)
>>> nch, npval = stats.chisquare(freqcount, nprob*n_sample)
>>> print('chisquare for t:      chi2 = %6.2f pvalue = %6.4f' % (tch, tpval))
chisquare for t:      chi2 =  2.30 pvalue = 0.8901  # random
>>> print('chisquare for normal: chi2 = %6.2f pvalue = %6.4f' % (nch, npval))
chisquare for normal: chi2 = 64.60 pvalue = 0.0000  # random

我們看到標準常態分佈被明確拒絕，而標準 t 分佈則無法被拒絕。由於我們樣本的變異數與這兩個標準分佈都不同，我們可以再次重新進行檢定，並將尺度和位置的估計值納入考量。

分佈的 fit 方法可用於估計分佈的參數，並且使用估計分佈的機率重複檢定。

>>> tdof, tloc, tscale = stats.t.fit(x)
>>> nloc, nscale = stats.norm.fit(x)
>>> tprob = np.diff(stats.t.cdf(crit, tdof, loc=tloc, scale=tscale))
>>> nprob = np.diff(stats.norm.cdf(crit, loc=nloc, scale=nscale))
>>> tch, tpval = stats.chisquare(freqcount, tprob*n_sample)
>>> nch, npval = stats.chisquare(freqcount, nprob*n_sample)
>>> print('chisquare for t:      chi2 = %6.2f pvalue = %6.4f' % (tch, tpval))
chisquare for t:      chi2 =  1.58 pvalue = 0.9542  # random
>>> print('chisquare for normal: chi2 = %6.2f pvalue = %6.4f' % (nch, npval))
chisquare for normal: chi2 = 11.08 pvalue = 0.0858  # random

考量到估計的參數，我們仍然可以拒絕樣本來自常態分佈的假設（在 5% 水準下），但同樣地，在 p 值為 0.95 的情況下，我們無法拒絕 t 分佈。

常態分佈的特殊檢定

由於常態分佈是統計學中最常見的分佈，因此有幾個額外的函數可用於檢驗樣本是否可能來自常態分佈。

首先，我們可以檢驗樣本的偏度和峰度是否與常態分佈的偏度和峰度顯著不同

>>> print('normal skewtest teststat = %6.3f pvalue = %6.4f' % stats.skewtest(x))
normal skewtest teststat =  2.785 pvalue = 0.0054  # random
>>> print('normal kurtosistest teststat = %6.3f pvalue = %6.4f' % stats.kurtosistest(x))
normal kurtosistest teststat =  4.757 pvalue = 0.0000  # random

這兩個檢定在常態性檢定中結合在一起

>>> print('normaltest teststat = %6.3f pvalue = %6.4f' % stats.normaltest(x))
normaltest teststat = 30.379 pvalue = 0.0000  # random

在這三個檢定中，p 值都非常低，我們可以拒絕樣本具有常態分佈的偏度和峰度的假設。

由於樣本的偏度和峰度基於中心動差，因此如果我們檢定標準化樣本，我們會得到完全相同的結果

>>> print('normaltest teststat = %6.3f pvalue = %6.4f' %
...       stats.normaltest((x-x.mean())/x.std()))
normaltest teststat = 30.379 pvalue = 0.0000  # random

由於常態性被如此強烈地拒絕，我們可以檢查 normaltest 在其他情況下是否給出合理的結果

>>> print('normaltest teststat = %6.3f pvalue = %6.4f' %
...       stats.normaltest(stats.t.rvs(10, size=100)))
normaltest teststat =  4.698 pvalue = 0.0955  # random
>>> print('normaltest teststat = %6.3f pvalue = %6.4f' %
...              stats.normaltest(stats.norm.rvs(size=1000)))
normaltest teststat =  0.613 pvalue = 0.7361  # random

當檢驗 t 分佈觀測值的小樣本和常態分佈觀測值的大樣本的常態性時，在這兩種情況下我們都無法拒絕樣本來自常態分佈的虛無假設。在第一種情況下，這是因為檢定力不足以區分小樣本中的 t 分佈和常態分佈隨機變數。