Изменения

Перейти к: навигация, поиск
м
Сравнительная таблица реализации критериев в пакетах
На практике мы встречаемся с двумя вариантами задач по проверке принадлежности распределения нормальному закону: для одномерного и многомерного распределения.
 
{{Pkg-req-notice}}
== Одномерное нормальное распределение ==
В качестве Нулевой гипотезой (<math>H_0</math> ) для всех нижеприведённых критериев является предположение, что «случайная величина <math>X</math> распределена нормально». Для демонстрации работы функций, реализующих различные критерии проверки принадлежности распределения нормальному закону сгенерируем вектор случайных чисел, имеющих стандартное нормальное распределение: {{r-code|code=<nowiki>> x <- rnorm(n = 1000)</nowiki>}} === Статистические критерии === В R реализовано множество критериев проверки соответствия распределения нормальному закону. ==== Сравнительная таблица реализации критериев в пакетах ==== {| class="wide wikitable sortable" style="text-align: center"! Критерии !! {{r-package|stats|core=true}} !! {{r-package|nortest}} !! {{r-package|moments}} !! {{r-package|fBasics}} !! {{r-package|tseries}} !! {{r-package|lawstat}}|-| style="text-align: left" | Критерий Шапиро - Уилка || <code>shapiro.test</code> || - || - || <code>shapiroTest</code> || - || -|-| style="text-align: left" | Критерий Колмогорова - Смирнова || <code>ks.test</code><ref>Для оценки нормальности вызов выглядит следующим образом: <code>ks.test(x, y = "pnorm")</code>.</ref> || - || - || <code>ksnormTest</code> || - || -|-| style="text-align: left" | Критерий Андерсона - Дарлинга || - || <code>ad.test</code> || - || <code>adTest</code> || - || -|-| style="text-align: left" | Критерий Крамера - фон Мизеса || - || <code>cvm.test</code> || - || <code>cvmTest</code> || - || -|-| style="text-align: left" | Критерий Лиллиефорса || - || <code>lillie.test</code> || - || <code>lillieTest</code> || - || -|-| style="text-align: left" | Критерий <math>\chi^2</math> Пирсона || - || <code>pearson.test</code> || - || <code>pchiTest</code> || - || -|-| style="text-align: left" | Критерий Шапиро - Франчия || - || <code>sf.test</code> || - || <code>sfTest</code> || - || -|-| style="text-align: left" | Критерий Д'Агостино || - || - || <code>agostino.test</code> || <code>dagoTest</code> || - || -|-| style="text-align: left" | Критерий Бонетта – Сайера || - || - || <code>bonett.test</code> || - || - || -|-| style="text-align: left" | Критерий Жарка - Бера || - || - || <code>jarque.test</code> || <code>jarqueberaTest</code> || <code>jarque.bera.test</code> || <code>rjb.test</code>|} Пакет {{r-package|fBasics}} содержит также функцию <code>normalTest()</code>, которая является «обёрктой» для ряда функций из того же пакета. Необходимый критерий можно задать с помощью аргумента <code>method</code>. Доступны следующие критерии: * <code>sw</code> - критерий Шапиро - Уилка* <code>jb</code> - критерий Жарка-Бера* <code>ks</code> - критерий Колмогорова - Смирнова* <code>da</code> - критерий Д'Агостино* <code>ad</code> - критерий Андерсона - Дарлинга. Пример вызова данной функции: {{r-code|code=<nowiki>> normalTest(x, method = "sw") Title: Shapiro - Wilk Normality Test Test Results: STATISTIC: W: 0.9831 P VALUE: 0.2301  Description: Fri Feb 14 19:59:59 2014 by user:</nowiki>}} Пакет <code>lawstat</code> содержит также функцию <code>sj.test()</code>, которая является реализацией рабастного критерия нормальности, созданного на основа критерия Шапиро - Уилка. Пакет <code>TeachingDemos</code> содержит функцию <code>SnowsPenultimateNormalityTest()</code>, реализующую неописанный в литературе критерий. Данная функция возвращает только уровень статистической значимости, свидетельствующий об отклонения распределения от нормального закона. ==== Маленькие хитрости ==== ===== Применение функций к нескольким переменным ===== С помощью <code>apply</code>-функций можно последовательно применить функцию к вектору, списку или массиву. Прежде чем всего нам необходимо сформировать таблицу данных. С помощью функции <code>replicate()</code> сгенерируем 10 переменных, имеющих стандартное нормальное распределение, которые объединяются в класс <code>data.frame</code>.
Для демонстрации работы функций{{r-code|code=<nowiki>> DF <- data.frame(replicate(n = 10, реализующий различные критерий проверки принадлежности распределения нормальному закону сгенерируем вектор случайных чисел, имеющих стандартное нормальное распределение:rnorm(n = 100)))</nowiki>}}
<syntaxhighlight lang="rsplus">> x <- rnorm(n = 100)</syntaxhighlight>Структура сгенерированной таблицы выглядит следующим образом:
{{r-code|code=== Пакет <codenowiki>stats> str(DF)'data.frame': 100 obs. of 10 variables: $ X1 : num 1.051 1.08 -0.477 -1.396 3.423 ... $ X2 : num -0.602 2.29 -0.758 -1.615 -0.364 ... $ X3 : num 0.0559 -1.0117 0.5242 0.4105 -0.3191 ... $ X4 : num -0.0965 0.2006 0.29 0.7702 -0.0182 ... $ X5 : num -0.7074 -1.6111 0.3478 0.2504 0.0609 ... $ X6 : num -1.432 0.535 -0.932 0.581 -1.606 ... $ X7 : num -1.42407 -0.31827 -2.04648 -0.19856 0.00301 ... $ X8 : num 0.511 0.192 0.467 -1.308 2.496 ... $ X9 : num -0.8508 0.4481 -0.2828 -0.5464 0.0605 ... $ X10: num 1.421 0.408 1.254 -0.956 -1.91 ...</codenowiki> ===}}
==== Для решения поставленной задачи можно воспользоваться функцией <code>sapply()</code>. Но прежде, нам необходимо немного отформатировать формат вывода результатов нашей функции: нам нужно извлечь значения критерия и его уровень значимости, т.к. результат функции <code>shapiro.test()</code> - критерий содержит также информацию, которая не подлежит включению в итоговую таблицу, например, информация об используемом методе (критерии) и уточнение характера альтернативной гипотезы. Вывод результатов тест Шапиро - Уилка ====выглядит следующим образом:
<syntaxhighlight lang{{r-code|code="rsplus"><nowiki>> shapiro.test(x)
Shapiro-Wilk normality test
data: x
W = 0.98319903, p-value = 0.23016882</syntaxhighlightnowiki>}}
Структура результата применения функции <code>shapiro.test()</code> представлена ниже:
{{r-code|code==== <codenowiki>>ksstr(shapiro.test</code> (x))List of 4 $ statistic: Named num 0.99 ..- критерий Колмогорова attr(*, "names")= chr "W" $ p.value : num 0.688 $ method : chr "Shapiro- Смирнова ===Wilk normality test" $ data.name: chr "x" - attr(*, "class")=chr "htest"</nowiki>}}
<syntaxhighlight lang="rsplus">> ks.test(xКак видим, y = "pnorm")помимо значений критерия и уровня значимости здесь содержится информация о применяемом методе. Мы можем отфильтровать вывод следующим образом:
One{{r-sample Kolmogorovcode|code=<nowiki>> normTest <-Smirnov function (x) {+ res <- shapiro.test(x)+ return(c(res$statistic, p.value = res$p.value))+ }</nowiki>}}
dataРезультат теперь будет выглядеть следующим образом: xD = 0.0534, p-value = 0.938alternative hypothesis: two-sided</syntaxhighlight>
{{r-code|code=== Пакет <codenowiki>nortest> normTest(x) W p.value 0.9903 0.6882</codenowiki> ===}}
Перед использованием функций из данного пакета, его необходимо предварительно установить и загрузить:Теперь можно использовать данную функцию при обработке столбцов нашей таблицы.
<syntaxhighlight lang{{r-code|code="rsplus"><nowiki> install.packages(pkgs = "nortest")> libraryt(package = "nortest"sapply(DF, normTest)) W p.valueX1 0.9831 0.2301X2 0.9936 0.9213X3 0.9800 0.1333X4 0.9829 0.2219X5 0.9874 0.4625X6 0.9862 0.3874X7 0.9839 0.2617X8 0.9833 0.2360X9 0.9915 0.7834X10 0.9808 0.1531</syntaxhighlightnowiki>}}
==== Того же результата можно добиться и с помощью функции <code>ad.testlapply()</code> - критерий Андерсона - Дарлинга ====<ref>По результатам сравнения производительности, данный вариант оказался чуть быстрее предыдущего.</ref>:
<syntaxhighlight lang{{r-code|code="rsplus"><nowiki> ad> do.testcall(xrbind, lapply(DF, normTest)) W p.valueX1 0.9831 0.2301X2 0.9936 0.9213X3 0.9800 0.1333X4 0.9829 0.2219X5 0.9874 0.4625X6 0.9862 0.3874X7 0.9839 0.2617X8 0.9833 0.2360X9 0.9915 0.7834X10 0.9808 0.1531</nowiki>}}
Anderson-Darling normality test===== Применение функций к нескольким по группам =====
dataДобавим к нашей таблице группы испытуемых: xA = 0.4181, p-value = 0.323</syntaxhighlight>
{{r-code|code==== <codenowiki>cvm.test> DF$GRP </code> - критерий Крамера - фон Мизеса ==factor(sample(LETTERS[1:3], size =100, replace =TRUE))</nowiki>}}
<syntaxhighlight lang="rsplus">> cvm.test(x)Состав групп получился следующим:
Cramer{{r-von Mises normality testcode|code=<nowiki>> table(DF$GRP) A B C 38 25 37 </nowiki>}}
data: xW = 0.048, pРассчитаем значения критерия Шапиро -value = 0.5359</syntaxhighlight>Уилка для первого столбца для каждоый группы испытуемых:
{{r-code|code==== <codenowiki>lillie> do.call(rbind, tapply(DF$X1, DF$GRP, normTest)) W p.valueA 0.9522 0.13281B 0.9607 0.28697C 0.9410 0.test07256</codenowiki> - критерий Лиллиефорса ====}}
<syntaxhighlight lang="rsplus">> lillie.test(x)== Графические методы ===
Lilliefors (KolmogorovМногие исследователи также используют графические методы для определения степени отклонения распределения от нормального закона. В R реализована возможность построения Q-Smirnov) normality testQ графиков, гистограмм и кривых распределения плотностей вероятности.
data: xD = 0.0521, p-value = 0.7258</syntaxhighlight>== Гистограмма ====
Гистограмма представляет собой графическое изображение зависимости частоты попадания элементов выборки от соответствующего интервала группировки. Построить гистограмму в R можно с помощью следующей команды:
{{r-code|code==== <codenowiki>pearson.test</code> - критерий <math>\chi^2hist(x)</mathnowiki> Пирсона ====}}
<syntaxhighlight lang="rsplus">> pearson[[Файл:Graphics-hist.test(x)svg|400px|центр]]
Pearson chi-square normality testНа гистограмме изображены абсолютные частоты. Также можно построить гистограмму, отражающую плотности вероятностей:
data: xP = 17.78, p{{r-value code|code= 0.05879<nowiki>> hist(x, freq = FALSE)</syntaxhighlightnowiki>}}
==== <code>sf.test</code> [[Файл:Graphics- критерий Шапиро hist- Франчия ====probs.svg|400px|центр]]
<syntaxhighlight lang="rsplus">> sf.test(x)=== График плотностей вероятности ====
Shapiro-Francia normality test===== Пакет <code>stats</code> =====
data: xW = 0.9781, p{{r-value code|code= 0.08574<nowiki>> plot(density(x))</syntaxhighlightnowiki>}}
[[Файл:Stats-density.svg|400px|центр]]
===== Пакет <code>momentscar</code> =====
Перед использованием функций из данного пакета, его необходимо предварительно установить и загрузить:{{r-code|code=<nowiki>> densityPlot(x)</nowiki>}}
<syntaxhighlight lang="rsplus">> install[[Файл:Car-densityPlot.packages(pkgs = "moments")> library(package = "moments")</syntaxhighlight>svg|400px|центр]]
==== <code>agostino.test</code> - критерий Д'Агостино Гистограммы с наложением графика плотностей вероятнотси ====
<syntaxhighlight lang="rsplus"==== Пакет <code>stats</code> agostino.test(x)=====
D'Agostino skewness test{{r-code|code=<nowiki>> hist(x, freq = FALSE)> lines(density(x))</nowiki>}}
data[[Файл: xskew = 0.1087, z = 0.3106, pStats-value = 0hist-density.7561alternative hypothesis: data have a skewness</syntaxhighlight>svg|400px|центр]]
==== <code>bonett.test</code> - критерий Бонетта – Сайера ====Теперь наложим на наш график кривую плотностей вероятности для нормального распределения:
{{r-code|code=<syntaxhighlight langnowiki>> xfit <- seq(min(x), max(x), length ="rsplus"100) # Координаты по оси X>yfit <- dnorm(xfit, mean = mean(x), sd = sd(x)) # Вычисление координат по оси Y> bonett.testhist(x, freq = FALSE)> lines(density(x), col = "red") # Накладываем кривую плотностей вероятности> lines(xfit, yfit, col = "blue") # Накладываем «нормальную» кривую</nowiki>}}
Bonett[[Файл:Stats-Seier test for Geary kurtosisdensity-compare.svg|400px|центр]]
data: xtau = 0.8158, z = -0.6097, p-value = 0.5421alternative hypothesis: kurtosis is not equal to sqrt(2/pi)== Пакет <code>gamlss</syntaxhighlightcode>=====
==== Более простой способ сравнение графиков плотностей вероятности представлен в функции <code>jarque.testhistDist</code> - критерий Жарка-Бера ====из пакета <code>gamlss</code>:
<syntaxhighlight lang{{r-code|code="rsplus"><nowiki> jarque.test> histDist(x, family = "NO", density = TRUE)
Jarque-Bera Normality TestFamily: c("NO", "Normal") Fitting method: "nlminb"
dataCall: xJB gamlssML(y = 1.356y, p-value family = 0.5076alternative hypothesis: greater</syntaxhighlight>"NO", formula = x)
=== Пакет <code>fBasics</code> ===Mu Coefficients:[1] -0.0462Sigma Coefficients:[1] 0.023
В данном пакете не предлагается никакой оригинальной реализации критериев - код в основном заимствован из пакетов <code>stats</code>, <code>nortest</code>, <code>moments</code> Degrees of Freedom for the fit: 2 Residual Deg. Данный пакет предлагает альтернативный вывод результатов в виде объектов S4-класса <code>fHTESTof Freedom 998 Global Deviance: 2884 AIC: 2888 SBC: 2898</codenowiki>, в том время как все предыдущие функции использовал S3-класс <code>htest</code>.}}
Перед использованием функций из данного пакета, его необходимо предварительно установить и загрузить[[Файл:Gamlss-histdist.svg|400px|центр]]
С помощью аргумента <syntaxhighlight lang="rsplus"code>family</code> install.packages(pkgs = "fBasics")можно задать семейство распределений для подгонки и сравнения<ref>Более подробную информацию о доступных семействах распределений можно получить с помощью команды <code> libraryhelp(package = "fBasicsgamlss.family")</syntaxhighlightcode>.</ref>.
Функция <code>normalTest()</code> является «обёрктой» для остальных функций из пакета <code>fBasics</code>. Задать необходимый критерий моно с помощью аргумента <code>method</code>. Доступны следующие критерии:==== Q-Q график ====
* <code>sw</code> Q- критерий Шапиро Q график (Q - Уилка* <code>jb</code> - критерий Жарка-Бера* <code>ks</code> - критерий Колмогорова - Смирнова* <code>da</code> - критерий Д'Агостино* <code>ad</code> - критерий Андерсона - Дарлингаквантиль) — это график, на котором квантили из двух распределений расположены относительно друг друга. Чем ближе точки на графике к диагональной прямой, тем ближе распределение исследуемой переменной к нормальному закону.
==== <code>shapiroTest</code> - критерий Шапиро - Уилка ====Построение квантильных графиков в R реализовано в нескольких пакетах.
<syntaxhighlight lang="rsplus"==== Пакет <code>stats</code> shapiroTest(x)=====
TitleПостроение Q–Q plot с помощью пакета <code>stats</code> выглядит следующим образом: Shapiro - Wilk Normality Test
Test Results:{{r-code|code= STATISTIC:<nowiki>> qqnorm(x) W: 0.9831> qqline(x)</nowiki> P VALUE: 0.2301 }}
Description[[Файл: Fri Feb 14 17:27:33 2014 by user: </syntaxhighlight>Stats-qqnorm.svg|400px|центр]]
==== = Пакет <code>ksnormTestQTLRel</code> - критерий Колмогорова - Смирнова =====
Построение Q–Q plot с помощью пакета <syntaxhighlight lang="rsplus"code>QTLRel</code> ksnormTest(x)выглядит следующим образом:
Title: One{{r-sample Kolmogorov-Smirnov testcode|code=<nowiki>> qqPlot(x, x = "norm")</nowiki>}}
Test Results[[Файл: STATISTIC: D: 0.0534 P VALUE: Alternative TwoQtlrel-Sided: 0.938 Alternative Less: 0.7049 Alternative Greater: 0qqplot.5654 svg|400px|центр]]
Description===== Пакет <code>car</code> ===== Альтернативный вариант реализован в функции <code>qqPlot()</code> из пакета <code>car</code>: Fri Feb 14 17{{r-code|code=<nowiki>> qqPlot(x, distribution = "norm")</nowiki>}} [[Файл:25:27 2014 by userCar-qqPlot.svg|400px|центр]] ===== Пакет <code>e1071</code> ===== Построение Q-Q plot можно осуществить с помощью функции <code>probplot</code> из пакета <code>e1071</code>:  {{r-code|code=<nowiki>> probplot(x, qdist = qnorm)</syntaxhighlightnowiki>}} [[Файл:E1071-probplot.svg|400px|центр]]
== Многомерное нормальное распределение ==
 
Перед началом обзора функций, реализующий критерии проверки многомерной нормальности, сгенерируем массив данных. Сделать это можно при помощью следующих функций
* <code>mvrnorm</code> из пакета <code>MASS</code>
* <code>rmvnorm</code> из пакета <code>mvtnorm</code>
* <code>rmnorm</code> из пакета <code>mnormt</code>
 
Вот пример кода, генерирующего массив данных, имеющих многомерное нормальное распределение:
 
{{r-code|code=
> means <- c(0, 0, 0, 0) # средние для переменных
> sigmas <- diag(length(means)) # ковариационная матрица
> mx <- rmvnorm(100, mean = means, sigma = sigmas)
}}
 
Пакет {{r-package|mvnormtest}} реализует модификацию критерия Шапиро - Уилка для многомерных данных - функция <code>mshapiro.test()</code><ref>В качестве аргумента необходимо передать транспонированную матрицу: <code>mshapiro.test(t(mx))</code>.</ref>.
 
Пакет {{r-package|ICS}} предлагает реализацию критериев эксцесса и асимметрии для многомерных данных: <code>mvnorm.kur.test()</code>, <code>mvnorm.skew.test()</code>.
 
Пакет {{r-package|energy}} реализует E-статистики для сравнения распределений. Критерия для проверки гипотезы о соответствия распределения многомерной переменной многомерному нормальному распределению предлагается функция <code>mvnorm.etest()</code><ref>Для вычисления уровня значимости критерия используется метод бутстрепа (bootstrap). Число итераций для бутстрепа можно задать с помощью аргумента <code>R</code>.</ref>.
 
== Ссылки ==
 
* Juergen Gross and bug fixes by Uwe Ligges (2012). nortest: Tests for Normality. R package version 1.0-2.
*: http://CRAN.R-project.org/package=nortest
* Lukasz Komsta and Frederick Novomestky (2012). moments: Moments, cumulants, skewness, kurtosis and related tests. R package version 0.13.
*: http://CRAN.R-project.org/package=moments
* Diethelm Wuertz, Rmetrics core team members, uses code builtin from the following R contributed packages: gmm from Pierre Chauss, gld from Robert King, gss from Chong Gu, nortest from Juergen Gross, HyperbolicDist from David Scott, sandwich from Thomas Lumley, Achim Zeileis, fortran/C code from Kersti Aas and akima from Albrecht Gebhardt (2013). fBasics: Rmetrics - Markets and Basic Statistics. R package version 3010.86.
*: http://CRAN.R-project.org/package=fBasics
* Adrian Trapletti and Kurt Hornik (2013). tseries: Time Series Analysis and Computational Finance. R package version 0.10-32.
*: http://CRAN.R-project.org/package=tseries
* Joseph L. Gastwirth; Yulia R. Gel <ygl@math.uwaterloo.ca>; W. L. Wallace Hui <wlwhui@uwaterloo.ca>; Vyacheslav Lyubchich <vlyubchich@uwaterloo.ca>; Weiwen Miao <miao@macalester.edu>; Kimihiro Noguchi <kinoguchi@ucdavis.edu> (2013). lawstat: An R package for biostatistics, public policy, and law. R package version 2.4.1.
*: http://CRAN.R-project.org/package=lawstat
* John Fox and Sanford Weisberg (2013). car: Companion to Applied Regression. R package version 2.0-19/r346.
*: http://R-Forge.R-project.org/projects/car/
* Mikis Stasinopoulos, Bob Rigby with contributions from Calliope Akantziliotou and Vlasios Voudouris (2014). gamlss: Generalised Additive Models for Location Scale and Shape. R package version 4.2-7.
*: http://CRAN.R-project.org/package=gamlss
* Riyan Cheng (2013). QTLRel: Tools for mapping of quantitative traits of genetically related individuals and calculating identity coefficients from a pedigree. R package version 0.2-14.
*: http://CRAN.R-project.org/package=QTLRel
* David Meyer, Evgenia Dimitriadou, Kurt Hornik, Andreas Weingessel and Friedrich Leisch (2014). e1071: Misc Functions of the Department of Statistics (e1071), TU Wien. R package version 1.6-2.
*: http://CRAN.R-project.org/package=e1071
 
== Примечания ==
<references />
[[Категория:R]]
[[Категория:Статистический анализ]]
[[Категория:Статистические критерии]]
[[Категория:Проверка статистических гипотез]]

Навигация