Реализация теста Колмогорова Смирнова в Python Scipy

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? Если это так, как я могу проверить, что мои данные распределены по гауссовски, но с другими mu и sigma?)

Answer 1

Ваши данные были сгенерированы с 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

Answer 2

Обновление ответа unutbu:

Для распределений, которые зависят только от местоположения и масштаба, но не имеют параметра формы, распределения нескольких статистических данных критерия соответствия не зависят от значений местоположения и масштаба. Распределение является нестандартным, однако его можно табулировать и использовать с любым местоположением и масштабом базового распределения.

Тест Колмогорова-Смирнова для нормального распределения с расчетным местоположением и масштабом также называют тестом Лиллифорса .

Теперь доступно в 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)

Большинство исследований Монте-Карло показывают, что критерий Андерсона-Дарлинга является более мощным, чем критерий Колмогорова-Смирнова. Он доступен в scipy.stats с критическими значениями и в statsmodels с приблизительными значениями p:

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

Ни один из тестов не отвергает нулевую гипотезу о том, что выборка распределена нормально. В то время как kstest в вопросе отвергает нулевую гипотезу о том, что выборка стандартная нормальная распределенная.

Answer 3

Вы также можете рассмотреть возможность использования теста Шапиро-Уилка, который «проверяет нулевую гипотезу о том, что данные получены из нормального распределения». Это также реализовано в 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, если вы предпочитаете - это не имеет значения), у нас есть веская причина отклонить нулевую гипотезу о том, что эти данные были получены из нормального распределения. </p>

Answer 4

В качестве дополнения к ответу @unutbu вы также можете указать параметры распределения для тестового распределения в kstest.Предположим, что у нас было несколько выборок из переменной (и назвали их datax), и мы хотели проверить, не могли ли эти выборки быть получены из логнормального, униформного или нормального.Обратите внимание, что для статистики scipy способ ввода входных параметров для каждого распределения немного варьируется.Теперь, благодаря "args" (кортежу или последовательности) в kstest, можно предоставить аргументы для дистрибутива scipy.stats, с которым вы хотите проверить.

:) Я также добавил возможность использовать тест с двумя выборками, на случай, если вы захотите сделать это любым из способов:

import numpy as np
from math import sqrt
from scipy.stats import kstest, ks_2samp, lognorm
import scipy.stats

def KSSeveralDists(data,dists_and_args,samplesFromDists=100,twosampleKS=True):
    returnable={}
    for dist in dists_and_args:
        try:
            if twosampleKS:
                try:
                    loc=dists_and_args[dist][0]
                    scale=dists_and_args[dist][1]
                    expression='scipy.stats.'+dist+'.rvs(loc=loc,scale=scale,size=samplesFromDists)'
                    sampledDist=eval(expression)
                except:
                    sc=dists_and_args[dist][0]
                    loc=dists_and_args[dist][1]
                    scale=dists_and_args[dist][2]
                    expression='scipy.stats.'+dist+'.rvs(sc,loc=loc,scale=scale,size=samplesFromDists)'
                    sampledDist=eval(expression)
                D,p=ks_2samp(data,sampledDist)
            else:
                D,p=kstest(data,dist,N=samplesFromDists,args=dists_and_args[dist])
        except:
            continue
        returnable[dist]={'KS':D,'p-value':p}
    return returnable

a=lambda m,std: m-std*sqrt(12.)/2.
b=lambda m,std: m+std*sqrt(12.)/2.
sz=2000

sc=0.5 #shape 
datax=lognorm.rvs(sc,loc=0.,scale=1.,size=sz)
normalargs=(datax.mean(),datax.std())

#suppose these are the parameters you wanted to pass for each distribution
dists_and_args={'norm':normalargs,
               'uniform':(a(*normalargs),b(*normalargs)),
               'lognorm':[0.5,0.,1.]
              }
print "two sample KS:"
print KSSeveralDists(datax,dists_and_args,samplesFromDists=sz,twosampleKS=True)
print "one sample KS:"
print KSSeveralDists(datax,dists_and_args,samplesFromDists=sz,twosampleKS=False)

, который дает в качестве вывода что-то вроде:

два образца KS: {'lognorm': {'KS': 0.023499999999999965, 'p-значение': 0.63384188886455217}, 'норма': {'KS': 0.10600000000000004, 'p-значение': 2.918766666723155e-10},'iform ': {' KS ': 0.15300000000000002,' p-value ': 6.443660021191129e-21}}

один образец KS: {' lognorm ': {' KS ': 0.01763415915126032,' p-значение ': 0,56275820961065193},' норма ': {' KS ': 0,10792612430093562,' p-значение ': 0,0},' униформа ': {' KS ': 0,14910036159697559,' p-значение ': 0,0}}

Примечание: Для равномерного распределения scipy.stats a и b принимаются как a = loc и b = loc + scale (см. документация ).