1. Наряду с проверкой различия параметров положения двух генеральных совокупностей часто требуется проверить различие двух параметров масштаба.
Пусть– случайная выборка объёма
с плотностью распределения
и
– случайная выборка объёма
с плотностью распределения
, где
– симметричная функция плотности распределения случайной величины,
.
Отметим, что в оговоренных условиях, имеем
откуда, в частности, следует
2. Интересующая нас гипотеза имеет вид
где. Далее рассматривается, основанная на рангах, процедура для проверки указанной гипотезы – тест Флигнера-Киллина (Fligner-Killeen test), а также строятся точечная оценка и доверительный интервал для параметра
.
3. Так как тест для проверки гипотезыдолжен не зависеть от параметров положения, то следует центрировать имеющиеся наблюдения. В параметрическом F-тесте (
var.test {stats}
), основанном на отношении выборочных дисперсий, для этой цели используются выборочные средние; вместо них мы применим, равные им в силу симметрии, медианы:
Исключив из дальнейшего рассмотрения возможные нулевые значения, найдём логарифмы абсолютных величин центрированных наблюдений. В результате получим две выборкии
, разница в параметрах положения которых (ключевой момент!) равна
.
Таким образом, задача о различии параметров масштаба сведена к задаче о различии параметров положения. Применим для её решения анализ, основанный на score функциях, разбор которого начат в предшествующей статье «Статистически устойчивый анализ данных: тест Манна-Уитни-Уилкоксона и Score-функции».
4. Так как логарифмическая функция монотонно возрастает, то ранги исходных наблюдений и ранги их логарифмов совпадают. Это позволяет сформировать следующую статистику
где ранги рассчитываются по объединенной выборке. Данная статистика может быть использована для проверки исследуемой нулевой гипотезыо величине параметра
[см. пункты 6, 7 прошлой статьи].
Напомним, что в условиях нулевой гипотезы распределение статистикиасимптотически нормально, с нулевым математическим ожиданием
и дисперсией
Это позволяет использовать стандартизованную статистику
чтобы отвергнуть гипотезуна уровне значимости
, если
, где
– квантиль уровня
стандартного нормального распределения.
Решение вопроса о том, какую score функцию использовать при расчете статистикизависит от вида функции
. В распространенном случае, когда
имеет нормальное распределение, оптимальной является квадратичная нормальная score функция (squared-normal score) вида
где– cdf стандартного нормального распределения, FK – аббревиатура, соответствующая работе Fligner и Killeen (1976) [1], где данная функция была предложена.
5. Для демонстрации изложенной выше методики рассмотрим данные конкретного эксперимента (Doksum and Sievers (1976) [2]) о влиянии озона на прирост в весе у крыс. Контрольная группа крысв течение
дней находилась в обычной среде, в то время как экспериментальная группа крыс
в течение 7 дней находилась в среде с повышенной концентрацией озона. В конце недели изменение в весе каждой крысы были взяты в качестве отклика.
![«Ящик-с-усами» исходных данных. «Ящик-с-усами» исходных данных.](https://habrastorage.org/getpro/habr/upload_files/137/616/10a/13761610ab5b4c91c2723c533b7a8ae5.png)
Следующая программа на зыке R содержит алгоритм расчета статистикии соответствующего значения p-value для проверки двухсторонней гипотезы о различии параметров масштаба двух выборок.
> x <-
+ c(
+ 41.0, 38.4, 24.4, 25.9, 21.9,
+ 18.3, 13.1, 27.3, 28.5, -16.9,
+ 26.0, 17.4, 21.8, 15.4, 27.4,
+ 19.2, 22.4, 17.7, 26.0, 29.4,
+ 21.4, 26.6, 22.7)
> y <-
+ c(
+ 10.1, 6.1, 20.4, 7.3, 14.3,
+ 15.5, -9.9, 6.8, 28.2, 17.9,
+ -9.0, - 12.9, 14.0, 6.6, 12.1,
+ 15.7, 39.9, -15.9, 54.6, -14.7,
+ 44.1, -9.0)
>
> # Зададим squared-normal score функцию
> fkscores = new(
+ "scores",
+ phi = function(u) {
+ (qnorm((u + 1) / 2)) ^ 2 - 1
+ },
+ Dphi = function(u) {
+ qnorm((u + 1) / 2) / dnorm(qnorm((u + 1) / 2))
+ }
+ )
> # Проверим гипотезу о равенстве 0 параметра Delta
> zed = c(abs(x - median(x)), abs(y - median(y)))
> n1 = length(x)
> n2 = length(y)
> n = n1 + n2
>
> # В расчётe используется стандартизированная score функция
> rz = rank(zed)/(n + 1)
> afk = Rfit::getScores(fkscores, rz)
>
> # Расчет статистики Sfk, zfk, p-value
> Sfk = sum(afk[(n1+1):n])
>
> varfk = ((n1 * n2)/(n*(n-1))) * sum(afk^2)
> zfk = Sfk/sqrt(varfk)
> see = "two.sided"
> if (see == "two.sided") {
+ pvalue = 2 * (1 - pnorm(abs(zfk)))
+ }
> if (see == "greater") {
+ pvalue = 1 - pnorm(zfk)
+ }
> if (see == "less") {
+ pvalue = pnorm(zfk)
+ }
>
> # Вывод результатов
> res <- list(statistic = zfk, p.value = pvalue)
> with(res, cat("statistic = ", statistic, ", p-value = ", p.value, "\n"))
statistic = 2.095976 , p-value = 0.03608434
Стандартизованная статистика теста Флигнера-Киллина имеет значениеc p-value
, то есть можно сказать, что наличие озона влияет на вариабельность в приросте веса у крыс.
6. Для того чтобы получить точечную оценку и доверительный интервал для параметра сформулируем нашу задачу в виде лог-модели. Пусть
– вектор-индикатор, с первыми
элементами равными нулю и последними
элементами равными единице. Тогда регрессионная задача [см. пункт 8 прошлой статьи] имеет вид
Применяя алгоритм построения ранговой регрессии (rfit {Rfit}
) с использованием квадратичной нормальной score функцийможно найти оценку
параметра
. Откуда оценка для параметра
получается в виде
.
Далее, используя найденную в ходе построения регрессии величину стандартной ошибки, можно построить, например, приблизительный
доверительный интервал для
, основанный на квантиле уровня
t-распределения с
степенями свободы
где– объем выборки
,
. Рассчитанные границы также преобразуются в соответствующие границы доверительного интервала для
:
. Отметим, что найденные таким образом доверительные границы всегда положительны.
7. Следующая программа на языке R содержит алгоритм расчета точечной оценки и доверительного интервала для параметра.
> # Подготовим данные
> xs = abs(x - median(x))
> xs = xs[xs != 0]
> xstarl = log(xs)
> ys = abs(y - median(y))
> ys = ys[ys != 0]
> ystarl = log(ys)
> n1s = length(xs)
> n2s = length(ys)
> ns = n1s + n2s
>
> # Построим уравнение регрессии
> cvec = c(rep(0, n1s), rep(1, n2s))
> zed = c(xstarl, ystarl)
> fitz = rfit(zed ~ cvec, scores = fkscores)
> (sumf = summary(fitz))
Call:
rfit.default(formula = zed ~ cvec, scores = fkscores)
Coefficients:
Estimate Std. Error t.value p.value
(Intercept) 1.42288 0.38905 3.6573 0.0007044 ***
cvec 0.86586 0.42783 2.0238 0.0493802 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Multiple R-squared (Robust): 0.07815811
Reduction in Dispersion Test: 3.56096 p-value: 0.06608
> delta = coef(sumf)[2, 1]
>
> # Оценка отношения параметров масштаба
> (eta = exp(delta))
[1] 2.377049
>
> # Построим доверительный интервал
> conf.level = 0.95
> se = coef(sumf)[2, 2]
> tc = abs(qt(((1 - conf.level)/2), ns - 2))
> lb = delta - tc * se
> ub = delta + tc * se
> conf.int = c(exp(lb), exp(ub))
> cat(100 * conf.level, " percent confidence interval:\n", conf.int)
95 percent confidence interval:
1.002462 5.636488
Таким образом, ранговая оценка параметраравна
c
доверительным интервалом
. Это также подтверждает вывод о том, что нулевую гипотезу
о равенстве параметров масштаба следует отвергнуть.
8. В заключение отметим, что рассматриваемая в статье ранговая процедура на основе предложенной Флигнером и Киллингом квадратичной нормальной score функции позволяет успешно проверять гипотезы о различии параметров масштаба, в двух важных случаях:
когда плотность распределение функции ошибки
симметрична, то есть процедура обобщает чувствительный к отклонениям от нормальности F-критерий;
когда распределение генеральной совокупности относится к классу засорённых нормальных распределений (skewed contaminated normal distributions). В этом случае случайная величина получается по закону
, где
и
имеют
и
распределения соответственно,
– имеет распределение Бернулли с вероятностью успеха
и
,
и
независимы в совокупности.
Ссылки на литературу:
Fligner, M. A. and Killeen, T. J. (1976), Distribution-free two-sample test for scale, Journal of the American Statistical Association, 71, 210-213.
Doksum, K. A. and Sievers, G. L. (1976), Plotting with confidence: Graphical comparisons of two populations, Biometrika, 63, 421-434.