在Python scipy中实现Kolmogorov Smirnov检验

Question

我有一组N个数字的数据集，想要测试它们是否符合正态分布。 我知道scipy.stats有一个kstest函数， 但是没有找到关于如何使用它和如何理解结果的例子。 这里有没有人熟悉这个，可以给我一些建议？

根据文档，使用kstest会返回两个数字，分别是KS检验统计量D和p值。 如果p值大于显著性水平（比如5%），那么我们就不能拒绝数据来自于某个特定分布的假设。

当我从一个正态分布中抽取10000个样本并进行高斯性测试时：

import numpy as np
from scipy.stats import kstest

mu,sigma = 0.07, 0.89
kstest(np.random.normal(mu,sigma,10000),'norm')

我得到了以下输出：

(0.04957880905196102, 8.9249710700788814e-22)

p值小于5%，这意味着我们可以拒绝数据是正态分布的假设。 但是这些样本确实是从正态分布中抽取的啊！

有没有人能理解并解释一下这里的矛盾？

（测试正态性是否假设均值mu = 0和标准差sigma = 1？如果是这样，我该如何测试我的数据是否是高斯分布，但均值和标准差不同呢？）

Answer 1

你可能还想考虑使用Shapiro-Wilk检验，这个检验是用来“测试数据是否来自正态分布”的。这个检验也可以在scipy库中找到：

http://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.shapiro.html

你需要把你的数据直接传递给这个函数。

import scipy

W, p = scipy.stats.shapiro(dataset)
print("Shapiro-Wilk test statistic, W:", W, "\n", "p-value:", p)

这个函数会返回类似这样的结果：

 Shapiro-Wilk test statistic, W: 0.7761164903640747 
 p-value: 6.317247641091492e-37

如果p值小于0.01（或者0.05，随你喜欢 - 这没关系），那么我们有充分的理由拒绝“这些数据来自正态分布”的假设。

Answer 2

关于unutbu的回答的更新：

对于那些只依赖于位置和规模，但没有形状参数的分布来说，几个适合度检验统计量的分布与位置和规模的值是独立的。虽然这种分布不是标准的，但可以列出表格，并且可以用于任何基础分布的具体位置和规模。

对于正态分布，使用估计的位置信息和规模的Kolmogorov-Smirnov检验也被称为Lilliefors检验。

现在在statsmodels中可以找到这个检验，并且提供了相关决策范围的近似p值。

>>> import numpy as np
>>> mu,sigma = 0.07, 0.89
>>> x = np.random.normal(mu, sigma, 10000)
>>> import statsmodels.api as sm
>>> sm.stats.lilliefors(x)
(0.0055267411213540951, 0.66190841161592895)

大多数蒙特卡洛研究表明，Anderson-Darling检验比Kolmogorov-Smirnov检验更有效。这个检验在scipy.stats中可以找到，并且有临界值，在statsmodels中也有近似的p值：

>>> sm.stats.normal_ad(x)
(0.23016468240712129, 0.80657628536145665)

这两个检验都没有拒绝样本是正态分布的零假设。而问题中的kstest则拒绝了样本是标准正态分布的零假设。

Answer 3

你的数据是用均值（mu）为0.07和标准差（sigma）为0.89生成的。你正在将这些数据与一个均值为0、标准差为1的正态分布进行比较。

原假设（H0）是指你的数据所代表的分布和标准正态分布是一样的，标准正态分布的均值是0，标准差是1。

小的p值表示，像D这样大的测试统计量出现的概率是p值。

换句话说，（p值大约是8.9e-22）这意味着原假设（H0）很可能不成立。

这很合理，因为均值和标准差并不匹配。

你可以将你的结果与以下内容进行比较：

In [22]: import numpy as np
In [23]: import scipy.stats as stats
In [24]: stats.kstest(np.random.normal(0,1,10000),'norm')
Out[24]: (0.007038739782416259, 0.70477679457831155)

为了测试你的数据是否符合高斯分布，你可以对其进行平移和缩放，使其变成均值为0、标准差为1的正态分布：

data=np.random.normal(mu,sigma,10000)
normed_data=(data-mu)/sigma
print(stats.kstest(normed_data,'norm'))
# (0.0085805670733036798, 0.45316245879609179)

---

警告：（感谢用户user333700（也就是scipy的开发者Josef Perktold））如果你不知道mu和sigma，那么估算这些参数会使得p值无效：

import numpy as np
import scipy.stats as stats

mu = 0.3
sigma = 5

num_tests = 10**5
num_rejects = 0
alpha = 0.05
for i in xrange(num_tests):
    data = np.random.normal(mu, sigma, 10000)
    # normed_data = (data - mu) / sigma    # this is okay
    # 4915/100000 = 0.05 rejects at rejection level 0.05 (as expected)
    normed_data = (data - data.mean()) / data.std()    # this is NOT okay
    # 20/100000 = 0.00 rejects at rejection level 0.05 (not expected)
    D, pval = stats.kstest(normed_data, 'norm')
    if pval < alpha:
        num_rejects += 1
ratio = float(num_rejects) / num_tests
print('{}/{} = {:.2f} rejects at rejection level {}'.format(
    num_rejects, num_tests, ratio, alpha))     

打印输出

20/100000 = 0.00 rejects at rejection level 0.05 (not expected)

这表明，如果样本使用样本的均值和标准差进行归一化，stats.kstest可能不会拒绝预期的原假设数量。

normed_data = (data - data.mean()) / data.std()    # this is NOT okay