分析一个样本

首先，我们创建一些随机变量。我们设置一个种子，以便在每次运行中获得相同的结果以供查看。例如，我们从学生 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，这意味着，例如，在 alpha 误差为 10% 的情况下，我们不能拒绝样本均值等于零的假设，即标准 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 分布和正态分布的随机变量，在小样本中。