Porovnáváme obdržené počty v diskrétních kategoriích s očekávanými hodnotami.
Jsou odlišnosti od H0 jen dílem náhody?
H0: ano
H1: není to náhoda
Příklad: Na otázku „váš nejoblíbenější předmět“ odpovědělo \(n = 70\) žáků následovně:
Na hladině významnosti \(\alpha = 0.05\) rozhodněte mezi H0: předměty jsou stejně oblíbené, H1: mezi předměty jsou výrazné rozdíly.
Vzhledem k H0, která předpokládá rovnoměrné zastoupení, jsou očekávané hodnoty počtů všech skupin stejné. Můžeme tedy použít defaultní nastavení funkce chisq.test.
chisq.test(c(20, 21, 29))##
## Chi-squared test for given probabilities
##
## data: c(20, 21, 29)
## X-squared = 2.0857, df = 2, p-value = 0.3524P-hodnota není menší než 0.05, takže nezamítáme H0 na hladině významnosti \(\alpha = 0.05\). Pokud máme pocit, že tělocvik má navrch, musíme sehnat více dat a uvidíme, jak test dopadne. chisq.test(c(77, 37), p = c(1147/12, 1145/12), rescale.p = TRUE)
Očekávané počty můžeme specifikovat buď formou očekávaných pravděpodobností (parametr p), a nebo stejným parametrem zadat vektor očekávaných četností pro H0 a přidat parametr rescale.p = TRUE, který řekne funkci, aby si uvnitř tyto hodnoty přepočetla na pravděpodobnosti s celkovým součtem 1.
Příklad: V období duben až říjen proběhlo v městě 78 svateb, zatímco od listopadu do března jen 36. H0: svatby jsou napříč těmito obdobími rovnoměrně rozděleny, H1: nejsou rovnoměrně rozděleny.
Rok má 12 měsíců a pro očekávané četnosti je tedy dobré zohlednit nestejnoměrnou velikost sledovaných období. Duben až říjen je 7 měsíců, listopad až březen je 5 měsíců.
svatby <- c(78, 36)
celkem <- sum(svatby)
chisq.test(svatby, p = c(celkem*7/12, celkem*5/12), rescale.p = TRUE)##
## Chi-squared test for given probabilities
##
## data: svatby
## X-squared = 4.7729, df = 1, p-value = 0.02891P-hodnota je nižší než 0.05, proto přijímáme hypotézu H1 na hladině významnosti \(\alpha = 0.05\).
Výše jsme pracovaly s vektory souhrnných hodnot. Pokud máme ale data.frame (získaný např. načtením tabulky ze souboru) s odpověďmi jednotlivých dotázaných, můžeme souhrn počtů podle požadovaných kategorií vyrobit funkcí table().
data(survey, package="MASS") # interní příklad v R, dotazník 237 studentů statistiky na University of Adelaide
head(survey)## Sex Wr.Hnd NW.Hnd W.Hnd Fold Pulse Clap Exer Smoke Height M.I
## 1 Female 18.5 18.0 Right R on L 92 Left Some Never 173.00 Metric
## 2 Male 19.5 20.5 Left R on L 104 Left None Regul 177.80 Imperial
## 3 Male 18.0 13.3 Right L on R 87 Neither None Occas NA <NA>
## 4 Male 18.8 18.9 Right R on L NA Neither None Never 160.00 Metric
## 5 Male 20.0 20.0 Right Neither 35 Right Some Never 165.00 Metric
## 6 Female 18.0 17.7 Right L on R 64 Right Some Never 172.72 Imperial
## Age
## 1 18.250
## 2 17.583
## 3 16.917
## 4 20.333
## 5 23.667
## 6 21.000Zaměřme se na otázku kouření (sloupec Smoke).
kuraci <- table(survey$Smoke)
kuraci##
## Heavy Never Occas Regul
## 11 189 19 17Řešme otázku, že H0: mezi studenty statistiky je 85 % nekuřáků a zbylé kategorie jsou rovnoměrně zastoupeny, každá 5 %. H1: rozdělení je jiné.
chisq.test(kuraci, p=c(.05, .85, .05, .05)) # pravděpodobnosti musejí odpovídat pořadí hodnot v tabulce##
## Chi-squared test for given probabilities
##
## data: kuraci
## X-squared = 7.4098, df = 3, p-value = 0.05992Na hladině významnosti \(\alpha = 0.05\) nezamítáme H0.
Očekávané četnosti můžeme zjistit přidáním expected a pomocí residuals vypíšeme hodnoty \(\chi\) bez umocnění, takže můžeme sledovat, kde jsou největší odchylky od očekávaných hodnot podle H0.
vysledek <- chisq.test(kuraci, p=c(.05, .85, .05, .05))
vysledek$expected## Heavy Never Occas Regul
## 11.8 200.6 11.8 11.8vysledek$residuals##
## Heavy Never Occas Regul
## -0.2328890 -0.8190163 2.0960010 1.5137785Použití testu \({\chi ^2}\) je jen aproximací problému. Funguje dobře, pokud jsou očekávané četnosti větší nebo rovny 5 v 80 % kategorií a žádná nemá očekávaný počet 0. Jinak je dobré kategorie s malým zastoupením sloučit do společných skupin. Pro malé četnosti je také možné použít Fisher’s exact test – funkce fisher.test().
Kontingenční tabulka – tabulka „se dvěma vstupy“. Je souvislost mezi oběma jevy?
Příklad: Dotazníku se zúčastnilo \(n = 400\) studentů. Mimo jiné odpovídali na otázky 1. „Byl(a) jste loni ubytován(a) na kolejích?“ (byl(a)/nebyl(a)) a 2. „Jaký byl váš loňský průměrný prospěch?“ (A: < 1.6, B: 1.6 – 2.1, C: > 2.1).
H0: ubytování na koleji s prospěchem nesouvisí, H1: souvisí.
Na \({\chi ^2}\) testu je dobré, že není nutné, aby byly kategorie rovnoměrně zastoupeny (např. kolej ano/ne). Stačí, když je v každé kategorii dostatečné zastoupení.
Jednou z možností, jak počty do R zadat, je vytvořit si dílčí vektory po řádcích, a ty pak spojit to jednoho data.frame. Pojmenování sloupců je sice hezké, ale pro test nemá význam.
kolejAno = c(39, 107, 93)
kolejNe = c(41, 73, 47)
tabulka <- as.data.frame(rbind(kolejAno, kolejNe)) # spojování vektorů jako řádky, z jmen vektorů se stanou jména řádků
tabulka## V1 V2 V3
## kolejAno 39 107 93
## kolejNe 41 73 47names(tabulka) <- c("A", "B", "C") # pojmenování sloupců - není nutné, ale je to hezké :-)
tabulka## A B C
## kolejAno 39 107 93
## kolejNe 41 73 47chisq.test(tabulka)##
## Pearson's Chi-squared test
##
## data: tabulka
## X-squared = 6.6286, df = 2, p-value = 0.03636fisher.test(tabulka) # přesnější, pokud jde použít (nenapíše error :-) )##
## Fisher's Exact Test for Count Data
##
## data: tabulka
## p-value = 0.03761
## alternative hypothesis: two.sidedNa hladině významnosti \(\alpha = 0.05\) přijímáme H1.
Často se setkáváme s úlohami následujícího typy (úloha i čísla jsou naprosto smyšlená, sledujme však princip). Výzkumný pracovník chce srovnat, zda o sobě autoři diplomových prací v textu mluví v jednotném či množném čísle (“provedli jsme měření…” vs. “změřila jsem…”). K dispozici má tabulku, ve které jsou uvedeny souhrnné počty pozorovaných případů napříč mnoha diplomovými pracemi, přičemž hodnoty jsou rozděleny do skupin dle pohlaví. Hypotézou je, že mezi pohlavími existuje rozdíl ve využívání množného čísla, protože z poměrů se zdá, že muži používají množné číslo častěji.
## Jednotné číslo Množné číslo
## Ženy 191 263
## Muži 138 211Výzkumný pracovník však nezná, kolik mužů a žen (tedy kolik diplomových prací) bylo pro analýzu použito a jediné, co ví, je, že byly započítány všechny výskyty v rámci celých textů. Je tedy možné, že jeden autor diplomové práce bude do tabulky přispívat 15 počty (např. 13x množné číslo a 2x jednotné), další třeba 23 počty (vše v množném čísle) a jiný 19 počty (18x jednotné číslo a 1x množné). Tato konkrétní čísla však k dispozici výzkumník nemá.
Hned na začátku musíme konstatovat, že takováto data jsou k ničemu a žádná statistická metoda nám k analýze nepomůže. Jisté je, že jsou zde replikace (více položek od jednoho člověka) a není to tak, že co +1 v tabulce, to jeden náhodný nezávislý výběr z populace.
Hrůzu situace si naplno představíme analogickým příkladem, že známe výsledek předvolebního průzkumu mezi dvěma kandidáty na prezidenta, ovšem data sbíralo dohromady mnoho brigádníků a běžně se stávalo, že jeden člověk byl anonymně dotázán v průběhu času vícekrát, a protože se jednalo o delší časové období, mohl i postupně svůj názor měnit. Takže se klidně stalo, že 3x odpověděl jasnou volbu prvního kandidáta, ale napočtvrté už odpovídal preferenci druhého kandidáta. Další respondent byl podobně dotázán dokonce sedmkrát, zatímco jiný respondent jen jednou. Nyní je snad zřejmé, že takováto data neznamenají 12 nezávislých zástupců populace, ale jsou to pořád jen 3 lidé, přičemž někteří jsou zastoupeni vícekrát a překlápějí tak svými více hlasy nesprávně váhu na svou stranu.
Co s tím? U dat je potřeba znát identifikátor (klidně v anonymizované formě) respondenta (v našem případě autora diplomové práce), aby bylo zřejmé, od koho pocházejí které položky. A pak je nutné data v rámci každé osoby normalizovat počtem položek. Pokud je u jednoho poměr 3:7, do tabulky je nutné zanést +0.3 a +0.7 (celkem = 1 člověk). Pokud je u jiného poměr 2:0, do tabulky se zanese +1 a +0 (celkem = 1 člověk). Zatím máme celkovou velikost výběrového souboru 2 a i součet všech počtů v tabulce je 2. A vůbec nevadí, že v buňkách nejsou celá čísla.
Na takto normalizovaná data můžeme použít Fisher’s exact test (viz sekce výše) a nedopustíme se problému pseudoreplikace a máme správně vyřešeno to, že každý respondent byl v datech zastoupen různým počtem položek.
Mysleme na způsob analýz již při sběru dat! Pokud jsou data zaznamenaná špatně, na výběr vhodné metody může být už pozdě.
Podobný příklad je vysvětlen zde (doporučuji první odpověď s nejvíce hlasy): https://stats.stackexchange.com/questions/123609/exact-two-sample-proportions-binomial-test-in-r-and-some-strange-p-values
Využijme interní data.frame mtcars a vytáhneme souhrnné počty testovaných automobilů podle počtu válců (cyl) a převodovky (am: 0 = automatická / 1 = manuální), pro více info o datech viz ?mtcars.
head(mtcars)## mpg cyl disp hp drat wt qsec vs am gear carb
## Mazda RX4 21.0 6 160 110 3.90 2.620 16.46 0 1 4 4
## Mazda RX4 Wag 21.0 6 160 110 3.90 2.875 17.02 0 1 4 4
## Datsun 710 22.8 4 108 93 3.85 2.320 18.61 1 1 4 1
## Hornet 4 Drive 21.4 6 258 110 3.08 3.215 19.44 1 0 3 1
## Hornet Sportabout 18.7 8 360 175 3.15 3.440 17.02 0 0 3 2
## Valiant 18.1 6 225 105 2.76 3.460 20.22 1 0 3 1tab <- table(mtcars$cyl, mtcars$am)
tab##
## 0 1
## 4 3 8
## 6 4 3
## 8 12 2Nová kontingenční tabulka podle vybraných faktorů z již existující větší kontingenční tabulky (s více faktory)
data(UCBAdmissions) # data o přijatých a nepřijatých uchazečích podle pohlaví a oddělení
DF <- as.data.frame(UCBAdmissions)
head(DF)## Admit Gender Dept Freq
## 1 Admitted Male A 512
## 2 Rejected Male A 313
## 3 Admitted Female A 89
## 4 Rejected Female A 19
## 5 Admitted Male B 353
## 6 Rejected Male B 207tab <- xtabs(Freq ~ Gender + Admit, data = DF) # kontingenčí tabulka hodnot Freq dle Gender a Admit (Department nás tedy v tuto chvíli nezajímá)
tab## Admit
## Gender Admitted Rejected
## Male 1198 1493
## Female 557 1278(Převzato z http://www.r-tutor.com/elementary-statistics/goodness-fit/chi-squared-test-independence)
Vraťme se k příkladu se studenty statistiky na University of Adelaide. Zajímá nás, zda je vztah mezi kuřáctvím a frekvencí sportování.
data(survey, package="MASS") # interní příklad v R, dotazník 237 studentů statistiky na University of
tbl <- table(survey$Smoke, survey$Exer) # vytvoření kontingenční tabulky závislosti 2 věcí
tbl##
## Freq None Some
## Heavy 7 1 3
## Never 87 18 84
## Occas 12 3 4
## Regul 9 1 7H0: kouření nemá souvislost s frekvencí sportování, H1: souvislost je.
chisq.test(tbl)## Warning in chisq.test(tbl): Chi-squared approximation may be incorrect##
## Pearson's Chi-squared test
##
## data: tbl
## X-squared = 5.4885, df = 6, p-value = 0.4828P-hodnota sice ukazuje na nezamítnutí H0, ale zároveň dostáváme varování o malém zastoupení kategorií.
Zkusme tedy některé sloučit kategorie None a Some do jedné společné.
ctbl <- cbind(tbl[,"Freq"], tbl[,"None"] + tbl[,"Some"])
ctbl## [,1] [,2]
## Heavy 7 4
## Never 87 102
## Occas 12 7
## Regul 9 8chisq.test(ctbl)##
## Pearson's Chi-squared test
##
## data: ctbl
## X-squared = 3.2328, df = 3, p-value = 0.3571P-hodnota je sice nižší, ale stále nezamítáme H0 na hladině významnosti \(\alpha = 0.5\). A co je hlavní, nyní máme mnohem lepší pocit, že se nedopouštíme statistického zločinu, protože funkce již žádné varování nehlásí.
Situace: máme „metodu“, zkoumáme úspěch před a po.
H0: metoda nemá vliv, neboli úspěch->neúspěch == neúspěch->úspěch, H1: metoda má vliv (ať kladný či záporný).
Poznámka: McNemarův test je pouze aproximací pro tuto úlohu a zvláště pro nízké počty může být nepřesný. Zde ho tedy uvádíme spíše z historických důvodů. Vhodnější je pro tyto úlohy použít tzv. oboustranný znaménkový test popsaný v kapitole 16. Testy středních hodnot.
mat <- as.table(rbind(c(38, 14),
c( 22, 36) ))
colnames(mat) <- rownames(mat) <- c("úspěch", "neúspěch")
names(dimnames(mat)) <- c("Před", "Po")
mat## Po
## Před úspěch neúspěch
## úspěch 38 14
## neúspěch 22 36mcnemar.test(mat, correct = TRUE)##
## McNemar's Chi-squared test with continuity correction
##
## data: mat
## McNemar's chi-squared = 1.3611, df = 1, p-value = 0.2433Na hladině významnosti \(\alpha = 0.05\) nezamítáme H0.
Poznámka: pro testy s více situacemi (více než 2) se používá Cochran’s Q Test.