Изменения

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

> x <- rnorm(n = 100)

=== Статистические критерии ===

Пример вызова данной функции:

> normalTest(x, method = "sw")

Description:

Fri Feb 14 19:59:59 2014 by user:

Помимо функции <code>normalTest()</code> данный пакет включает в себя следующие функции:

С помощью <code>apply</code>-функций можно последовательно применить функцию к вектору, списку или массиву. Прежде чем всего нам необходимо сформировать таблицу данных. С помощью функции <code>replicate()</code> сгенерируем 10 переменных, имеющих стандартное нормальное распределение, которые объединяются в класс <code>data.frame</code>.

> DF <- data.frame(replicate(n = 10, rnorm(n = 100)))

Структура сгенерированной таблицы выглядит следующим образом:

> str(DF)

'data.frame': 100 obs. of 10 variables:

$ X9 : num -0.8508 0.4481 -0.2828 -0.5464 0.0605 ...

$ X10: num 1.421 0.408 1.254 -0.956 -1.91 ...

Для решения поставленной задачи можно воспользоваться функцией <code>sapply()</code>. Но прежде, нам необходимо немного отформатировать формат вывода результатов нашей функции: нам нужно извлечь значения критерия и его уровень значимости, т.к. результат функции <code>shapiro.test()</code> содержит также информацию, которая не подлежит включению в таблицу.

> shapiro.test(x)

data: x

W = 0.9903, p-value = 0.6882

Структура результата применения функции <code>shapiro.test()</code> представлена ниже:

> str(shapiro.test(x))

List of 4

$ data.name: chr "x"

- attr(*, "class")= chr "htest"

Как видим, помимо значений критерия и уровня значимости здесь содержится информация о применяемом методе. Мы можем отфильтровать вывод следующим образом:

> normTest <- function (x) {

+ res <- shapiro.test(x)

+ return(c(res$statistic, p.value = res$p.value))

+ }

Результат теперь будет выглядеть следующим образом:

> normTest(x)

W p.value

0.9903 0.6882

Теперь можно использовать данную функцию при обработке столбцов нашей таблицы.

> t(sapply(DF, normTest))

W p.value

X9 0.9915 0.7834

X10 0.9808 0.1531

Того же результата можно добиться и с помощью функции <code>lapply()</code><ref>По результатам сравнения производительности, данный вариант оказался чуть быстрее предыдущего.</ref>:

> do.call(rbind, lapply(DF, normTest))

W p.value

X9 0.9915 0.7834

X10 0.9808 0.1531

===== Применение функций к нескольким по группам =====

Добавим к нашей таблице группы испытуемых:

> DF$GRP <- factor(sample(LETTERS[1:3], size = 100, replace = TRUE))

Состав групп получился следующим:

> table(DF$GRP)

A B C

38 25 37

Рассчитаем значения критерия Шапиро - Уилка для первого столбца для каждоый группы испытуемых:

> do.call(rbind, tapply(DF$X1, DF$GRP, normTest))

W p.value

B 0.9607 0.28697

C 0.9410 0.07256

=== Графические методы ===

Построение Q–Q plot с помощью пакета <code>stats</code> выглядит следующим образом:

> qqnorm(x)

> qqline(x)

[[Файл:Stats-qqnorm.svg|400px|центр]]

Построение Q–Q plot с помощью пакета <code>QTLRel</code> выглядит следующим образом:

> qqPlot(x, x = "norm")

[[Файл:Qtlrel-qqplot.svg|400px|центр]]

Альтернативный вариант реализован в функции <code>qqPlot()</code> из пакета <code>car</code>:

> qqPlot(x, distribution = "norm")

[[Файл:Car-qqplot.svg|400px|центр]]

Построение P-P plot можно осуществить с помощью функции <code>probplot</code> из пакета <code>e1071</code>:

> probplot(x, qdist = qnorm)

[[Файл:E1071-probplot.svg|400px|центр]]

Ещё один интересный способ графического анализа представлен функцией <code>histDist</code> из пакета <code>gamlss</code>:

> histDist(x, family = "NO", density = TRUE)

AIC: 307.414

SBC: 312.624

[[Файл:gamlss-histdist.svg|400px|центр]]