Оригинал: Contingency Tables in R
Как создать таблицу
Во-первых, давайте получим некоторые данные. Пакет MASS содержит данные о 93 автомобилях, поступивших в продажу в США в 1993 году. Они хранятся в объекте Cars93 и включают 27 признаков для каждого автомобиля, некоторые из которых являются категориальными. Итак, давайте загрузим пакет MASS и посмотрим на тип транспортных средств, включенных в cars93:
library(MASS)
Cars93$Type
## [1] Small Midsize Compact Midsize Midsize Midsize Large Large
## [9] Midsize Large Midsize Compact Compact Sporty Midsize Van
## [17] Van Large Sporty Large Compact Large Small Small
## [25] Compact Van Midsize Sporty Small Large Small Small
## [33] Compact Sporty Sporty Van Midsize Large Small Sporty
## [41] Sporty Small Compact Small Small Sporty Midsize Midsize
## [49] Midsize Midsize Midsize Large Small Small Compact Van
## [57] Sporty Compact Midsize Sporty Midsize Small Midsize Small
## [65] Compact Van Midsize Compact Midsize Van Large Sporty
## [73] Small Compact Sporty Midsize Large Compact Small Small
## [81] Small Compact Small Small Sporty Midsize Van Small
## [89] Van Compact Sporty Compact Midsize
## Levels: Compact Large Midsize Small Sporty Van
У нас там 6 типов автомобилей. Функция table сообщает, сколько у нас их каждого типа:
table(Cars93$Type)
##
## Compact Large Midsize Small Sporty Van
## 16 11 22 21 14 9
prop.table преобразует их в доли:
prop.table(table(Cars93$Type))
prop.table(table(Cars93$Type))
##
## Compact Large Midsize Small Sporty Van
## 0.17204301 0.11827957 0.23655914 0.22580645 0.15053763 0.09677419
То же самое с происхождением автомобилей:
table(Cars93$Origin)
##
## USA non-USA
## 48 45
prop.table(table(Cars93$Origin))
##
## USA non-USA
## 0.516129 0.483871
Как создать таблицу сопряженности
Отлично, мы увидели, что наш набор данных содержит одинаковое количество автомобилей из США и не из США, и что наиболее распространенными типами являются Midsize и Small. Однако, может быть, США и неамериканцы различаются по типам?
Давайте посмотрим на типы автомобилей с точки зрения их происхождения. Мы можем снова использовать table, но теперь с двумя аргументами. Первый станет переменной строки, а второй станет переменной столбца:
table(Cars93$Type, Cars93$Origin)
##
## USA non-USA
## Compact 7 9
## Large 11 0
## Midsize 10 12
## Small 7 14
## Sporty 8 6
## Van 5 4
Теперь мы увидели то, что все знают и так: американцы любят большие машины! В таблице выше показано совместное распределение двух категориальных переменных (Type и Origin). Такие таблицы называются таблицами сопряженности.
Как получить маржинальную форму таблицы сопряженности
Функции rowSums и colSums действительно говорят сами за себя
(tab1<-table(Cars93$Type, Cars93$Origin))
##
## USA non-USA
## Compact 7 9
## Large 11 0
## Midsize 10 12
## Small 7 14
## Sporty 8 6
## Van 5 4
rowSums(tab1)
## Compact Large Midsize Small Sporty Van
## 16 11 22 21 14 9
colSums(tab1)
## USA non-USA
## 48 45
Как получить проценты из таблицы сопряженности
prop.table, вложенная в table, дает частоты:
prop.table(table(Cars93$Type, Cars93$Origin))
##
## USA non-USA
## Compact 0.07526882 0.09677419
## Large 0.11827957 0.00000000
## Midsize 0.10752688 0.12903226
## Small 0.07526882 0.15053763
## Sporty 0.08602151 0.06451613
## Van 0.05376344 0.04301075
Преобразование в проценты - это просто умножение на 100:
prop.table(table(Cars93$Type, Cars93$Origin))*100
##
## USA non-USA
## Compact 7.526882 9.677419
## Large 11.827957 0.000000
## Midsize 10.752688 12.903226
## Small 7.526882 15.053763
## Sporty 8.602151 6.451613
## Van 5.376344 4.301075
Обратите внимание, что это совместное распределение вероятностей, из которого мы видим, например, что около 7,5% автомобилей являются небольшими и имеют американское происхождение.
Чаще всего нас интересует распределение одной переменной внутри групп, созданных другой. Здесь интересным представляется распределение типов автомобилей среди американских и (отдельно) неамериканских автомобилей. Чтобы получить это, мы используем аргумент margin для функции prop.table. Он сообщает, где в строках (margin = 1) или в столбцах (margin = 2) переменная группировки:
prop.table(table(Cars93$Type, Cars93$Origin), margin=2)*100
##
## USA non-USA
## Compact 14.583333 20.000000
## Large 22.916667 0.000000
## Midsize 20.833333 26.666667
## Small 14.583333 31.111111
## Sporty 16.666667 13.333333
## Van 10.416667 8.888889
Теперь мы можем легко увидеть, что маленькие автомобили встречаются за пределами США в два раза чаще, чем в США.
Также обратите внимание, что проценты составляют до 100 в столбцах, тогда как в объединенной таблице распределения (таблица без аргумента margin) 100 представляла собой сумму всей таблицы.
(tab2<-prop.table(table(Cars93$Type, Cars93$Origin), margin=2)*100)
##
## USA non-USA
## Compact 14.583333 20.000000
## Large 22.916667 0.000000
## Midsize 20.833333 26.666667
## Small 14.583333 31.111111
## Sporty 16.666667 13.333333
## Van 10.416667 8.888889
colSums(tab2)
## USA non-USA
## 100 100
Критерий хи-квадрат
Наиболее распространенный вопрос, который возникает в таблицах сопряженности, заключается в том, являются ли переменные строк и столбцов независимыми. Самый простой способ ответить на него - запустить тест хи-квадрат. Он подробно описан в этом руководстве. Давайте проверим, являются ли Type и Origin независимыми:
chisq.test(Cars93$Type, Cars93$Origin)
## Warning in chisq.test(Cars93$Type, Cars93$Origin): Chi-squared
## approximation may be incorrect
##
## Pearson's Chi-squared test
##
## data: Cars93$Type and Cars93$Origin
## X-squared = 14.08, df = 5, p-value = 0.01511
По-видимому, это не так, но мы также получили предупреждение хи-квадрат. Это связано с тем, что статистика хи-квадрат следует распределению хи-квадрат только приблизительно. Чем больше наблюдений мы имеем, тем лучше приближение. Функция chisq.test выдает вышеупомянутое предупреждение всякий раз, когда одно из ожидаемых значений меньше 5 (значение «ожидаемое количество» см. в руководстве, указанном выше).
Точный тест Фишера
Точный критерий Фишера является альтернативой критерию хи-квадрат, используемой в основном, когда аппроксимация хи-квадрат не является удовлетворительной. Давайте запустим его:
fisher.test(Cars93$Type, Cars93$Origin)
##
## Fisher's Exact Test for Count Data
##
## data: Cars93$Type and Cars93$Origin
## p-value = 0.007248
## alternative hypothesis: two.sided
Результаты очень похожи на результаты из хи-квадрат, но это не всегда так. Одним из основных недостатков точного теста Фишера является то, что для больших таблиц (или больших выборок) он становится неэффективным в вычислительном отношении.
G-тест
Другой альтернативой является так называемый G-тест. Его статистика также приблизительно распределена по хи-квадрат, но для небольших выборок это приближение лучше, чем при использовании теста хи-квадрат. Для G-теста мы можем использовать функцию GTest из пакета DescTools. Результаты снова очень похожи на два предыдущих теста: Type и Origin не являются независимыми.
library(DescTools)
GTest(Cars93$Type, Cars93$Origin)
##
## Log likelihood ratio (G-test) test of independence without
## correction
##
## data: Cars93$Type and Cars93$Origin
## G = 18.362, X-squared df = 5, p-value = 0.002526
Поправка Йейтса
В таблицах сопряженности 2x2 критерий хи-квадрат может быть улучшен с помощью поправки непрерывности Йейтса. Она просто вычитает 0,5 из каждого члена | Наблюдаемый - Ожидаемый | в статистике хи-квадрат. Снова обратитесь к этому руководству, если вы не понимаете, о чем идет речь. Более того, R применяет поправку Йейтса автоматически, когда это необходимо. Давайте посмотрим на наличие версий с механической коробкой передач для автомобилей США и неамериканских:
(tab3<-table(Cars93$Man.trans.avail, Cars93$Origin))
##
## USA non-USA
## No 26 6
## Yes 22 39
chisq.test(tab3)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tab3
## X-squared = 15.397, df = 1, p-value = 8.712e-05
Поправка Йейтса была применена автоматически, и R сказал нам об этом.
Таблица размерности 3 (или больше)
Давайте посмотрим на Man.trans.avail, Origin и Type одновременно. table разбивает набор данных в соответствии с переменной, представленной как третья:
table(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)
## , , = Compact
##
##
## USA non-USA
## No 2 0
## Yes 5 9
##
## , , = Large
##
##
## USA non-USA
## No 11 0
## Yes 0 0
##
## , , = Midsize
##
##
## USA non-USA
## No 9 4
## Yes 1 8
##
## , , = Small
##
##
## USA non-USA
## No 0 0
## Yes 7 14
##
## , , = Sporty
##
##
## USA non-USA
## No 0 0
## Yes 8 6
##
## , , = Van
##
##
## USA non-USA
## No 4 2
## Yes 1 2
ftable дает более компактный вид:
ftable(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)
## Compact Large Midsize Small Sporty Van
##
## No USA 2 11 9 0 0 4
## non-USA 0 0 4 0 0 2
## Yes USA 5 0 1 7 8 1
## non-USA 9 0 8 14 6 2
Тест Кокрана-Мантеля-Хензеля
Из приведенных выше результатов должно быть ясно, что зависимость между Origin и Man.trans.avail отличается от Type. Например, для небольших и спортивных автомобилей ассоциация отсутствует: каждый американский, а также любой неамериканский маленький или спортивный автомобиль выпускается с ручной коробкой передач. С другой стороны, большинство автомобилей среднего размера в США не имеют такой версии, в то время как большинство неамериканских автомобилей этого типа имеют. Чтобы учесть возможные различия, можно использовать критерий Кохрана-Мантеля-Хензеля:
mantelhaen.test(Cars93$Man.trans.avail, Cars93$Origin, Cars93$Type)
##
## Mantel-Haenszel chi-squared test with continuity correction
##
## data: Cars93$Man.trans.avail and Cars93$Origin and Cars93$Type
## Mantel-Haenszel X-squared = 8.0153, df = 1, p-value = 0.004638
## alternative hypothesis: true common odds ratio is not equal to 1
## 95 percent confidence interval:
## 2.226531 76.891307
## sample estimates:
## common odds ratio
## 13.08438
Третий аргумент, переданный mantelhaen.test, определяет страты. Сравните приведенные выше результаты с результатами без стратификации (пример поправки Йейтса выше). Ассоциация все еще присутствует, но доказательства для этого гораздо слабее.
Меры ассоциации
Как только мы обнаружили некоторые связи между переменными, пришло время измерить их силу. Существует множество различных типов возможных мер. Многие из них описаны здесь. Сейчас я остановлюсь на двух наиболее широко используемых.
V-коэффициент Крамера
V основан на статистике хи-квадрат:
где:
N - общий итог таблицы сопряженности (сумма всех ее ячеек),
C - количество столбцов
R - количество строк.
V [0; 1]. Чем больше V, тем сильнее связь между переменными. V = 0 можно интерпретировать как независимость (поскольку V = 0 тогда и только тогда, когда χ2 = 0). Основным недостатком V является отсутствие точной интерпретации. Является ли V = 0,6 сильной, умеренной или слабой ассоциацией?
Функция CramerV от DescTools может рассчитать это для нас:
CramerV(Cars93$Type, Cars93$Origin)
## [1] 0.3890967
Еще раз: это сильная, умеренная или слабая ассоциация?
Лямбда Гудмена — Краскэла
Лямбда Гудмена — Краскэла - пример меры, основанной на пропорциональном уменьшении вариаций. Эти виды мер направлены на имитацию R2 - коэффициента детерминации линейной регрессии:
они принимают значения из [0, 1]
они выражают долю вариаций, объясняемых независимой переменной.
Давайте использовать столбец-переменную как независимую. Тогда формула лямбда Гудмена — Краскэла:
где:
Lj - сумма немодальных частот в j-м столбце,
L - сумма немодальных частот в столбце «total»
Это лучше всего объяснить на примере. Давайте посмотрим на Type и Origin, последняя будет независимой:
table(Cars93$Type, Cars93$Origin)
##
## USA non-USA
## Compact 7 9
## Large 11 0
## Midsize 10 12
## Small 7 14
## Sporty 8 6
## Van 5 4
Мода столбца US - 11, поэтому L1 = 7 + 10 + 7 + 8 + 5 = 37. Мода столбца non-US - 14, поэтому L2 = 9 + 0 + 12 + 6 + 4 = 31.
Столбец «Total»
rowSums(table(Cars93$Type, Cars93$Origin))
## Compact Large Midsize Small Sporty Van
## 16 11 22 21 14 9
Мода равна 22, поэтому L = 16 + 11 + 21 + 14 + 9 = 71 и
λ=(71−(37+31))/71=0.042
Вместо этого мы можем использовать функцию Lambda из пакета DescTools:
Lambda(Cars93$Type, Cars93$Origin, direction='row')
## [1] 0.04225352
Параметр direction определяет, где (в строке или в столбце) находится зависимая переменная.
Origin объясняет только около 4% изменений в Type. Обратите внимание на формулу, что лямбда определяет вариацию как дихотомию между принадлежностью и не принадлежностью к самой большой группе.
Стоит отметить, что лямбда равна нулю, когда модальная категория каждого столбца одинакова. Посмотрите на эту таблицу:
(lambdaTab<-cbind(c(0,100), c(49,51)))
## [,1] [,2]
## [1,] 0 49
## [2,] 100 51
chisq.test(lambdaTab)
##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: lambdaTab
## X-squared = 62.279, df = 1, p-value = 2.981e-15
Модальная категория одинакова для каждого столбца (это вторая строка). Таким образом, лямбда должна быть нулевой, несмотря на значительную и видимую связь:
Lambda(lambdaTab, direction='row')
## [1] 0
Комментариев нет:
Отправить комментарий