19  Тесты хи-квадрат

Тесты хи-квадрат – это целый набор статистических методов для исследования связи между переменными. Значения делятся по группам

\[\chi^2 = \sum^n_{i = 1} \frac{{(O_i - E_i)^2}}{{E_i}}\] \(E_i\) – ожидаемая частота, \(O_i\) – наблюдаемая частота.

Соответственно, чем больше расхождения между ожидаемым распределением по группам и реальным, тем больше статистика \(\chi^2\) и тем меньше p-value.

19.0.1 Тест хи-квадрат Пирсона на независимость

Тест хи-квадрат Пирсона на независимость – это статистический метод, используемый для определения того, есть ли статистическая связь между категориальными переменными.

Например, с помощью этого теста мы можем посмотреть, различается ли распределение по полу супергероев у двух основных издателей супергероики – Marvel Comics и DC Comics.

library(tidyverse)

Warning: package 'dplyr' was built under R version 4.2.3

Warning: package 'stringr' was built under R version 4.2.3

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.4     ✔ readr     2.1.4
✔ forcats   1.0.0     ✔ stringr   1.5.1
✔ ggplot2   3.4.4     ✔ tibble    3.2.1
✔ lubridate 1.9.3     ✔ tidyr     1.3.0
✔ purrr     1.0.2     
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

heroes <- read_csv("https://raw.githubusercontent.com/Pozdniakov/tidy_stats/master/data/heroes_information.csv",
                   na = c("-", "-99", "NA", " "))

New names:
• `` -> `...1`

Warning: One or more parsing issues, call `problems()` on your data frame for details,
e.g.:
  dat <- vroom(...)
  problems(dat)

Rows: 734 Columns: 11
── Column specification ────────────────────────────────────────────────────────
Delimiter: ","
chr (8): name, Gender, Eye color, Race, Hair color, Publisher, Skin color, A...
dbl (3): ...1, Height, Weight

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

Для начала нам нужно создать таблицу сопряженности (contingency table) для двух переменных с помощью встроенной функции table(). Эта функция похожа на count() из пакета {dplyr}.

gender_publisher <- heroes %>%
  drop_na(Gender) %>%
  filter(Publisher %in% c("Marvel Comics", "DC Comics")) %>%
  select(Gender, Publisher)

gender_publisher

# A tibble: 577 × 2
   Gender Publisher    
   <chr>  <chr>        
 1 Male   Marvel Comics
 2 Male   DC Comics    
 3 Male   Marvel Comics
 4 Male   Marvel Comics
 5 Male   Marvel Comics
 6 Male   DC Comics    
 7 Female Marvel Comics
 8 Male   Marvel Comics
 9 Male   Marvel Comics
10 Male   Marvel Comics
# ℹ 567 more rows

table(gender_publisher)

        Publisher
Gender   DC Comics Marvel Comics
  Female        61           111
  Male         153           252

С помощью дженерик-функции summary() примененной на построенной таблице сопряженности можно посчитать тест на независимость хи-квадрат.

summary(table(gender_publisher))

Number of cases in table: 577 
Number of factors: 2 
Test for independence of all factors:
    Chisq = 0.27673, df = 1, p-value = 0.5988

То же самое можно сделать и с помощью специальной функции chisq.test():

chisq.test(table(gender_publisher))

    Pearson's Chi-squared test with Yates' continuity correction

data:  table(gender_publisher)
X-squared = 0.18649, df = 1, p-value = 0.6659

Как можно заметить, в этом случае p-value несколько отличается. Дело в том, что p-value для хи-квадрат рассчитывается по непрерывному хи-квадрат распределению, тогда как хи-квадрат статистики распределены таким образом только в случае достаточно большой выборке. Чтобы скорретировать на малый размер групп используется поправка на непрерывность Йейтса.

chisq.test(table(gender_publisher), correct = FALSE)

    Pearson's Chi-squared test

data:  table(gender_publisher)
X-squared = 0.27673, df = 1, p-value = 0.5988

Другой подход для более точного подсчета p-value – это использование точного теста Фишера:

fisher.test(table(gender_publisher))

    Fisher's Exact Test for Count Data

data:  table(gender_publisher)
p-value = 0.6381
alternative hypothesis: true odds ratio is not equal to 1
95 percent confidence interval:
 0.6120819 1.3321610
sample estimates:
odds ratio 
 0.9053112