分析单个样本

首先，我们创建一些随机变量。我们设置一个种子，以便在每次运行中都能获得相同的结果。例如，我们从 Student 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

因为正态性被强烈拒绝，我们可以检查正态性检验是否适用于其他情况

>>> 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 分布和正态分布的随机变量。