Mittearvuliste tunnuste vahelised seosed

Tänases praktikumis

  • hii-ruut-test
  • seose tugevuse hindamine
  • Fisheri täpne test
  • Monte-Carlo simulatsioon

Hii-ruut-test

Hii-ruut-test on statistiline test, millega saab hinnata kahe kategoriaalse tunnuse sõltuvust. Test töötab risttabeliga, milles on lihtsagedused (mitte suhtelised ega normaliseeritud sagedused).

Tegelikult on hii-ruut-testid terve testide pere, mille arvutamiseks kasutatakse sama valemit, ent mis erinevad selle poolest, milline peab olema valim, mis on nullhüpotees ja kuidas testi tulemust tõlgendada.

Franke, Todd Michael, Timothy Ho & Christina A. Christie 2012. The Chi-Square Test: Often Used and More Often Misinterpreted. American Journal of Evaluation 33(3), 448-458.

Kooskõlamäära testis (goodness of fit) vaatleme mingi kategoriaalse tunnuse jaotumist ühes valimis ning testime, kas see valim võiks tulla populatsioonist, mille sama tunnuse jaotumise kohta on meil juba teadmine olemas. See testi liik on niisiis analoogne nt ühe valimiga t-testile, kus võrdleme valimi mingit arvulist näitajat (nt laste pikkused) populatsiooni omaga, otsustamaks, kui tõenäoline on, et valim on samast populatsioonist. Arvulise näitaja asemel on siin aga mingi kategoriaalse tunnuse väärtuste sagedused.

Sõltumatuse testis (independence) uurime kahe tunnuse vahelist seost samas valimis, homogeensuse testis (homogeneity) aga uurime, kas ühe tunnuse kategooriate proportsioonid on kahes sõltumatus valimis erinevad. See eristus võib aga kohati olla keeruline, kuna alati pole selge, mis hetkest peaks rääkima eraldi populatsioonidest/gruppidest. Testi tegemine käib mõlemal juhul samamoodi, ent erinevused tulenevad valimi moodustamisest ja testi tõlgendamisest.

Selles praktikumis keskendume sõltumatuse testile, kus uurime kahe tunnuse omavahelise seose olemasolu: kas ühe tunnuse väärtuste teadmine aitab ennustada teise tunnuse väärtust (chi-square test of independence)?

Vaatleme näitena andmeid Ann Siimani artiklist “Ainsuse sisseütleva vormi valiku seos morfosüntaktiliste ja semantiliste tunnustega - materjali ning meetodi sobivus korpusanalüüsiks” (2015), kus ta kasutas hii-ruut-testi, et selgitada välja, millest sõltub sisseütleva käände kasutuses see, kas kasutatakse lühikest sisseütleva käände vormi ehk aditiivi (majja) või pikka sisseütleva käände vormi ehk illatiivi (majasse).

Uurime siin esmalt, kas vormi valik sõltub sellest, kas tegemist on pärisnime või üldnimega (nt Tartusse vs. ajajärku). Teeme seega sõltumatuse testi, kuna iga valimisse sattunud sõna kohta on kodeeritud kaks seda sõna iseloomustavat tunnust (sisseütleva vorm ja nimelisus).

Siiman 2015: 220 

# Koostame artiklis esitatud andmete põhjal
# esmalt risttabeli.
tab <- matrix(c(43, 77, 377, 343), # koosesinemise sagedused
              nrow = 2, # mitmesse ritta sagedused jaotada?
              dimnames = list(c("Lühike", "Pikk"), # tabeli ridade nimed
                              c("Pärisnimi", "Üldnimi"))) # tabeli tulpade nimed

tab
##        Pärisnimi Üldnimi
## Lühike        43     377
## Pikk          77     343

Kui meil oleks vormikasutuse kohta põhjalikum andmestik, saaks risttabelit saaks mõistagi teha ka andmestiku tunnustest, nt table(andmed$tunnus1, andmed$tunnus2).

library(tidyverse)

# teeme esmalt table-objektist data.frame-objekti
df_pikk <- tab %>% # võtame table-tüüpi objekti
  as.data.frame() %>% # teisendame selle data.frame-tüüpi objektiks
  rownames_to_column(var = "Vorm") %>% # paneme pakutud reanimed tulpa "Vorm"
  pivot_longer(cols = c("Pärisnimi", "Üldnimi"), # teeme laiast tabelist "pika", nii et võtame pärisnime ja üldnime tulbad,
               names_to = "Liik", # paneme need tulbanimed ("Pärisnimi", "Üldnimi") pika tabeli tulpa "Liik" 
               values_to = "Sagedus") %>% # tulpades olevad sagedused pika tabeli tulpa "Sagedus"
  slice(rep(1:n(), times = Sagedus)) %>% # korda igat rida nii mitu korda kui on väärtus sama rea tulbas "Sagedus"
  select(-Sagedus) # eemalda tulp "Sagedus"
df_pikk

table(df_pikk$Vorm, df_pikk$Liik)
## # A tibble: 840 × 2
##    Vorm   Liik     
##    <chr>  <chr>    
##  1 Lühike Pärisnimi
##  2 Lühike Pärisnimi
##  3 Lühike Pärisnimi
##  4 Lühike Pärisnimi
##  5 Lühike Pärisnimi
##  6 Lühike Pärisnimi
##  7 Lühike Pärisnimi
##  8 Lühike Pärisnimi
##  9 Lühike Pärisnimi
## 10 Lühike Pärisnimi
## # ℹ 830 more rows
##         
##          Pärisnimi Üldnimi
##   Lühike        43     377
##   Pikk          77     343

Visualiseerime sagedusi.

par(mfrow = c(1, 2)) # näita jooniseid ühes reas kahes tulbas (= kõrvuti)
# esimene joonis:
barplot(tab, # sisendiks tavaline risttabel
        col = c("gold1", "grey50"), 
        legend.text = rownames(tab), # legendi väärtused võta tabeli reanimedest
        beside = TRUE, # tulbad kõrvuti
        args.legend = list(x = "topleft")) # legend ülemisse vasakusse nurka

# teine joonis:
barplot(prop.table(tab, margin = 2), # sisendiks suhteliste sagedustega risttabel
        col = c("gold1", "grey50"))

par(mfrow = c(1, 1)) # vii sätted uuesti algsele kujule

Kui meil on kasutada ainult sagedustabel (st meil pole algandmeid, kus ridades on vaatlused ja tulpades on tunnused, nt vorm ja nimelisus), siis peame ggploti kasutamiseks sagedustabeli esmalt data.frame tüüpi tabeliks teisendama.

df <- tab %>% # võtame table-tüüpi objekti
  as.data.frame() %>% # teisendame selle data.frame-tüüpi objektiks
  rownames_to_column(var = "Vorm") %>% # paneme pakutud reanimed tulpa "Vorm"
  pivot_longer(cols = c("Pärisnimi", "Üldnimi"), # teeme laiast tabelist "pika", nii et võtame pärisnime ja üldnime tulbad,
               names_to = "Liik", # paneme need tulbanimed ("Pärisnimi", "Üldnimi") pika tabeli tulpa "Liik" 
               values_to = "Sagedus")
# esimene joonis:
ggplot(df) +
  geom_col(aes(x = Liik, y = Sagedus, fill = Vorm), position = "dodge") +
  scale_fill_manual(values = c("gold1", "grey50")) +
  theme(legend.position = "top") -> p1

# teine joonis:
ggplot(df) +
  geom_col(aes(x = Liik, y = Sagedus, fill = Vorm), position = "fill", show.legend = F) +
  scale_fill_manual(values = c("gold1", "grey50")) +
  labs(y = "Suhteline sagedus") -> p2

# kuva jooniseid kõrvuti
# install.packages("patchwork")
library(patchwork)
p1 + p2

Hii-ruut-test põhineb tegelike ja teoreetiliste sageduste võrdlemisel.

Tegelikud/vaadeldud (Observed) sagedused on tegelikud sagedused igas tabeli lahtris (nt 43, 77 jne).
Teoreetilised/oodatud (Expected) sagedused leitakse valemiga reasumma*veerusumma/kogusumma.

Teoreetilised sagedused väljendavad kahe tunnuse kategooriate koosesinemiste sagedusi nullhüpoteesi korral ehk juhul, kui kahe vaadeldava tunnuse erinevad väärtused esineksid koos täiesti juhuslikult ning nende kasutuses ei esineks mingeid mustreid ega seoseid. Sellisel juhul peaks saama tabeli lahtrite väärtused tuletada puhtalt ridade ja veergude kogusagedustest: kui üht kategooriat (nt üldnimesid) esineb andmestikus rohkem, esineb teda ühtviisi rohkem nii lühikeste kui ka pikkade sisseütleva vormide korral.

Mida suurem on tegelike sageduste erinevus teoreetilistest sagedustest, seda suurem on ka hii-ruut-statistiku väärtus.

(tab2 <- addmargins(tab)) # lisame tabelile rea- ja veerusummad
##        Pärisnimi Üldnimi Sum
## Lühike        43     377 420
## Pikk          77     343 420
## Sum          120     720 840
# Näiteks lühikese sisseütlevaga pärisnimede (kuhu? - Tartu) teoreetiline sagedus oleks
# (reasumma*veerusumma/kogusumma)
tab2["Lühike","Sum"]*tab2["Sum","Pärisnimi"]/tab2["Sum","Sum"]
## [1] 60

Teoreetilisi sagedusi ei pea aga ise käsitsi arvutama, need saab kätte ka hii-ruut-testi väljundist.

chisq.test(tab)$expected # teoreetilised sagedused
##        Pärisnimi Üldnimi
## Lühike        60     360
## Pikk          60     360



Hii-ruut-testi eeldused:

  • Valim on üldkogumist juhuslikult valitud.
  • Vaatlused on sõltumatud.
  • Iga vaatlus peab sobima risttabelis täpselt ainult ühte lahtrisse (st iga uuritava tunnuse ainult ühte kategooriasse).

Uuritava andmestiku puhul ei ole tegelikult vaatluste sõltumatuse nõue täidetud, kuna korpuses ei ole igalt kõnelejalt ainult üks mingi sisseütlevas käändes sõna kasutuskord. See tähendab, et kõneleja/kirjutaja, kes on rohkem sisseütlevas käändes sõnu tarvitanud, mõjutab enam ka testi tulemusi. Keelekorpustega, kus pahatihti on ka keeruline tuvastada, kes täpselt kõneleja/kirjutaja on, vaatluste sõltumatuse nõuet sageli eiratakse, ent sellisel juhul võib alati küsida, kas ikka kirjeldame keelt üldiselt või mõne domineeriva keelekasutaja idiolekti.

Sageli viidatakse hii-ruut-testi eeldusena ka nn Cochrani reeglile, mille põhjal ei tohiks hii-ruut-testi tegemisel ükski tegelik sagedus olla null ega ükski teoreetiline sagedus olla väiksem kui 5. Nendele tingimustele on esitatud ka mööndusi, eriti keeleteaduses (nt Stefanowitsch 2020: 177 või Wallis 2021: 117), kuna reeglit peetakse liiga konservatiivseks. Klassikaliselt soovitatakse siiski sellistel juhtudel kasutada mõnd muud sarnast testi, nt Fisheri testi või G-testi.



Sõnastame hüpoteesid.

H0: Sõna nimelisus ja tema sisseütleva käände vorm on teineteisest sõltumatud tunnused ja ühe tunnuse väärtuse teadmine ei aitaks ennustada teise tunnuse väärtust (tabeli lahtrites olevad väärtused on tuletatavad tabeli marginaalsagedustest).
H1: Sõna nimelisus ja sisseütleva käände vorm on teineteisest sõltuvad tunnused ja ühe tunnuse väärtuse teadmine aitab ennustada teise tunnuse väärtust (tabeli lahtrites olevad väärtused ei ole tuletatavad tabeli marginaalsagedustest).

Teeme nüüd päriselt hii-ruut-testi.

chisq.test(tab, correct = FALSE)
## 
##  Pearson's Chi-squared test
## 
## data:  tab
## X-squared = 11.239, df = 1, p-value = 0.000801

Hii-ruut-statistik leitakse valemiga \(x^2 = \sum \frac{(o-e)^2}{e}\), kus o (observed) on tegelik sagedus ja e (expected) on teoreetiline sagedus. See tähendab, et iga risttabeli lahtri sagedusest lahutatakse selle teoreetiline sagedus (o-e), saadud tehte tulemused võetakse ruutu ((o-e)^2), ruutu võetud erinevused jagatakse läbi lahtrite teoreetiliste sagedustega ((o-e)^2/e) ning kõik selle tulemusel lahtritesse saadud väärtused liidetakse kokku.

# salvestame testi tulemuse objekti
t <- chisq.test(tab, correct = FALSE)

# leiame hii-ruut-statistiku käsitsi
sum(((t$observed - t$expected)^2)/t$expected)
## [1] 11.23889

Mida suurem erinevus tegelike ja oodatud sageduste vahel on, seda suurem on statistiku väärtus, ehkki konkreetne väärtus sõltub oluliselt marginaalsagedustest.

Kuna p-väärtus on siin väga väike (0.000801), on väga väike tõenäosus (täpsemalt 0,0801%), et valimis täheldatud erinevused oodatud ja tegelike sageduste vahel on juhuslikud ja esineksid ka siis, kui tegelikult kehtiks nullhüpotees ja seost kahe tunnuse vahel ei oleks. Seega on alust nullhüptees hüljata ja võtta vastu sisukas hüpotees, mille kohaselt kaks tunnust on omavahel sõltuvuses.

Kasutasime praegu testi tegemisel argumenti correct = FALSE. Vaikimisi tehakse R-is 2x2 tabelite puhul test nn Yatesi pidevusparandusega (Yates’ continuity correction). Pidevusparanduse käigus lahutatakse iga lahtri sageduse ja selle oodatud sageduse vahe absoluutväärtusest veel 0,5. Nii väheneb hii-ruut-statistiku väärtus ja suureneb p-väärtus. Tegemist on niisiis pisut konservatiivsema taktikaga, mida võiks rakendada juhul, kui 2x2 tabelis on mõni teoreetiline sagedus < 5.

chisq.test(tab)
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  tab
## X-squared = 10.587, df = 1, p-value = 0.001139

Pidevusparanduse käigus muutuvad ainult statistik ja p-väärtus, oodatud sagedusi ja jääke pidevusparanduse kasutamine kuidagi ei mõjuta.

Pane tähele, et kahe tunnuse sõltuvust testival hii-ruut-testil on vabadusastmete arv enam mitte üks vähem kui vaatluste arv, vaid (ridade arv - 1) * (tulpade arv - 1). Vabadusastmed iseloomustavad siin nende risttabeli lahtrite arvu, mille väärtusi saaks vabalt muuta, ilma et ridade ja tulpade marginaalsagedused või üldine sageduste summa muutuks. 2x2 risttabelis on selliseid lahtreid ainult 1, sest kohe, kui oleme ühte lahtrisse mingi sageduse kirjutanud, on teiste lahtrite sagedused marginaalsageduste põhjal juba ette kindlaks määratud.

Testi p-väärtus näitab niisiis, kui tõenäoline on selle vabadusastmete arvu ja määratud olulisusnivoo korral (nt df = 1, alfa = 0.05), et saaksime sellise hii-ruut-statistiku ka siis, kui kahe tunnuse vahel poleks tegelikult mingit seost. Hii-ruut-jaotuse kriitilisi väärtusi, millest alates saame lugeda seose oluliseks, saab vaadata vastavatest tabelitest, nt siit.

Tuleb aga silmas pidada, et mida suurem on vabadusastmete arv ehk mida rohkem kategooriaid kahel tunnusel kokku on, seda tõenäolisemalt hii-ruut-test sealt mõne olulise seose tuvastab. Samuti kipuvad andmed sellistes tabelites kategooriate vahel hajuma ning tegelike sageduste hulka satub palju nulle. Seepärast ei soovitata teha hii-ruut-testi suurtel tabelitel, kus on palju ridu ja/või tulpasid.

Harjutus 1

Laadi kursusel kasutatud küsimustiku andmestik ning vaata esmalt, kas leiad statistiliselt olulise seose lemmikjoogi ja kvantitatiivsete meetodite kogemuse vahel ning programmeerimisoskuse ja kvantitatiivsete meetodite kogemuse vahel. Kui kasutada on terve andmetabel, kus ridades on vaatlused ja tulpades tunnused, siis ei pea hii-ruut-testi sisendiks andma tingimata risttabelit, vaid piisab ka lihtsalt tunnuste vektoritest.

chisq.test(table(andmed$tunnus1, andmed$tunnus2))
chisq.test(andmed$tunnus1, andmed$tunnus2)


Seejärel loe sisse 2023. aasta valimiste andmestik, jäta välja üksikkandidaatide andmed ning kontrolli, kas on statistiliselt oluline seos soo ja valimisnimekirja vahel.

Hii-ruut-testi jäägid

Hii-ruut-test ise ütleb lihtsalt, et kahe kategoriaalse tunnuse vahel on statistiliselt oluline seos, aga ei ütle midagi selle kohta, milline see seos on. Seepärast oleks kasulik vaadata Pearsoni jääkide abil, millised erinevused tegelike ja oodatud sageduste vahel on suuremad. Pearsoni jääke arvutatakse valemiga \(\frac{(o-e)}{\sqrt{e}}\).

# Näiteks lühikese sisseütlevaga pärisnimede jääk on
(t$observed["Lühike", "Pärisnimi"]-t$expected["Lühike", "Pärisnimi"])/sqrt(t$expected["Lühike", "Pärisnimi"])
## [1] -2.194691
# Jäägid saab kätte ka testi väljundist
t$residuals
##        Pärisnimi    Üldnimi
## Lühike -2.194691  0.8959787
## Pikk    2.194691 -0.8959787

Kui jääk on positiivne, on tegelik sagedus suurem kui oodatud ja konkreetne tunnuste kombinatsioon esineb sagedamini kui marginaalsageduste põhjal võiks oletada. Kui jääk on negatiivne, on tegelik sagedus väiksem kui oodatud. Näiteks pärisnimed esinevad oodatust harvemini lühikeses sisseütlevas käändes (nt läksin Tartu).

Teeme jääkide graafiku.

barplot(t$residuals, 
        beside = TRUE, 
        legend = TRUE,
        ylab = "Pearsoni jäägid")

# ggplotiga
jaak_pikk <- t$residuals %>%
  as.data.frame() %>%
  rownames_to_column(var = "Vorm") %>%
  pivot_longer(cols = c("Pärisnimi", "Üldnimi"), names_to = "Liik", values_to = "Jääk")

ggplot(jaak_pikk) +
  geom_col(aes(x = Liik, fill = Vorm, y = Jääk), position = "dodge")

Mida suurem on Pearsoni jäägi absoluutväärtus, seda enam on tõenäoliselt konkreetne lahter statistiliselt olulise hii-ruut-statistiku väärtusesse panustanud.

Selleks aga, et näha, millised jäägid esindavad olulist erinevust oodatud sagedustest, vaatame standardiseeritud jääke. Standardiseeritud jäägid R-i hii-ruut-testi kontekstis (ka Pearsoni jäägid on põhimõtteliselt standardiseeritud) väljendavad erinevust tegeliku ja oodatud sageduse vahel, mis on jagatud mitte ruutjuurega oodatud sagedustest (nagu tavaliste Pearsoni jääkide puhul), vaid ruutjuurega avaldisest, kus korrutatakse oodatud sagedused, 1 - oodatud sageduste osakaal ridade summast ja 1 - oodatud sageduste osakaal tulpade summast. Need on tuntud ka kui adjusted residuals.

Näiteks lühikeses sisseütlevas käändes pärisnime standardiseeritud jääk oleks

(teg <- tab["Lühike", "Pärisnimi"]) # tegelik sagedus
## [1] 43
(ood <- sum(tab["Lühike",]) * sum(tab[,"Pärisnimi"]) / sum(tab)) # oodatud sagedus
## [1] 60
(std <- ood * (1 - ood/sum(tab["Lühike",])) * (1 - ood/sum(tab[,"Pärisnimi"]))) # korrutis
## [1] 25.71429
(teg-ood)/sqrt(std)
## [1] -3.352445
t$stdres
##        Pärisnimi   Üldnimi
## Lühike -3.352445  3.352445
## Pikk    3.352445 -3.352445

Kui olulisuse nivoo on 0.05, siis on statistiliselt olulised need jäägid, mis on suuremad kui 1.96 või väiksemad kui -1.96 (tihti ümardatakse ka 2 ja -2). Vastava kriitilise väärtuse võib leida nt funktsiooniga qnorm(1-olulisuse_nivoo/2).

barplot(t$stdres, 
        beside = TRUE, 
        legend = TRUE,
        ylab = "Kohandatud Pearsoni jäägid")
abline(h = 1.96, col = "red", lty = "dashed")
abline(h = -1.96, col = "red", lty = "dashed")

# ggplotiga
t$stdres %>%
  as.data.frame() %>%
  rownames_to_column(var = "Vorm") %>%
  pivot_longer(., cols = c("Pärisnimi", "Üldnimi"), names_to = "Liik", values_to = "Std_jääk") %>%
  ggplot() +
  geom_col(aes(x = Liik, fill = Vorm, y = Std_jääk), position = "dodge") +
  geom_hline(aes(yintercept = 1.96), col = "red", linetype = "dashed") +
  geom_hline(aes(yintercept = -1.96), col = "red", linetype = "dashed")

Seose tugevuse hindamine

Kategoriaalsete tunnuste vahelise seose tugevust saab mõõta näiteks Craméri V või φ-kordaja abil.

# install.packages("vcd")
library(vcd)
## Loading required package: grid
assocstats(tab)
##                     X^2 df   P(> X^2)
## Likelihood Ratio 11.373  1 0.00074531
## Pearson          11.239  1 0.00080101
## 
## Phi-Coefficient   : 0.116 
## Contingency Coeff.: 0.115 
## Cramer's V        : 0.116

Sümmeetrilise 2x2 risttabeli puhul on Craméri V (ruutjuur hii-ruut-statistiku ja vaatluste arv * vabadusastmete arv jagatisest) ja φ-kordaja (ruutjuur hii-ruut-statistiku ja vaatluste arvu jagatisest) samad ja varieeruvad vahemikus nullist üheni.

  • Nõrk seos: 0.1-0.3
  • Keskmine seos: 0.3-0.5
  • Tugev seos: > 0.5

Seose tugevust seostatakse vahel ka vabadusastmete arvuga (mida rohkem on vabadusastmeid, seda madalam on tugeva seose nivoo).

https://peterstatistics.com/CrashCourse/3-TwoVarUnpair/NomNom/NomNom3c.html

NB! Seose tugevuse hindamine ei käi samadel alustel nagu korrelatsioonikordajate puhul. Samuti ei näita Craméri V ega φ-kordaja seose suunda.

Suuremate (mitteruut)tabelite puhul võivad hakata kordajate väärtused erinema. Craméri V jääb alati vahemikku 0-1, φ-kordaja väärtus võib minna ka üle 1.

Tavaliselt kasutatakse 2x2 tabelist suuremate tabelite puhul Craméri V-d, sest see on võrreldav ka eri suurustega tabelite korral.

Kuidas hii-ruut-testi tulemust raporteerida

Hii-ruut-testi tegemisel on tavaks esitada hii-ruut-statistiku väärtus koos vabadusastmete arvu ja vaatluste arvuga ning p-väärtus. Võib lisada ka info mõju suuruse kohta.

## 
##  Pearson's Chi-squared test
## 
## data:  tab
## X-squared = 11.239, df = 1, p-value = 0.000801

Hii-ruut-testi tulemus on statistiliselt oluline: X2(1, N = 840) = 11.239, p < 0.001. Craméri V on 0.12, mille põhjal on tegemist siiski nõrga seosega.

Testi võimsus

Ka hii-ruut-testi jaoks võib arvutada välja testi võimsuse, mis ütleb, kui suure tõenäosusega leiame sellise valimi suuruse, vabadusastmete arvu, mõju suuruse ja olulisuse nivoo korral olulise seose juhul, kui see tõesti olemas on.

library(pwr)
w <- assocstats(tab)$cramer # mõju suurus (siin Craméri V) 
N <- sum(tab) # vaatluste arv
df <- t$parameter # hii-ruut-testi vabadusastmete arv (t <- chisq.test(tab))

pwr.chisq.test(w = w, N = N, df = df, sig.level = 0.05)
## 
##      Chi squared power calculation 
## 
##               w = 0.1156703
##               N = 840
##              df = 1
##       sig.level = 0.05
##           power = 0.9181117
## 
## NOTE: N is the number of observations

Siin näiteks on tõenäosus, et me ei tee teist tüüpi viga (ütleme, et seost ei ole, kui tegelikult on), ~92%. See on päris hea.

Harjutus 2

Leia, millised soo ja valimisnimekirja kombinatsioonid panustavad kõige enam hii-ruut-testi tulemusse (standardiseeritud jäägid), kui tugev on seos soo ja valimisnimekirja vahel ning milline on testi võimsus.

Kontrolli seejärel, kas äkki tuli küsimustiku andmete põhjal seos lemmikjoogi ja kvantitatiivsete meetodite kogemuse vahel ebaoluline seetõttu, et testi võimsus polnud piisav.

Fisheri täpne test ja G-test

Mida teha siis, kui andmestikus on teoreetilised sagedused, mis on < 5 või lausa 0? Tihti räägitakse, et sellistel juhtudel hii-ruut-test ei sobi.

Alternatiivid on näiteks Fisheri täpne test või G-test.

Vaatleme nüüd käändevormi sõltuvust sellest, millise sõnaliigiga on tegemist.

Siiman 2015: 215 

tab2 <- matrix(c(4, 4, 2, 6, 389, 391, 25, 19), 
               nrow = 2, ncol = 4, 
               dimnames = list(c("Lühike", "Pikk"), 
                               c("Arvsõna", "Asesõna", "Nimisõna", "Omadussõna")))

tab2
##        Arvsõna Asesõna Nimisõna Omadussõna
## Lühike       4       2      389         25
## Pikk         4       6      391         19
(t2 <- chisq.test(tab2))
## Warning in chisq.test(tab2): Chi-squared approximation may be incorrect
## 
##  Pearson's Chi-squared test
## 
## data:  tab2
## X-squared = 2.8233, df = 3, p-value = 0.4197
t2$expected
##        Arvsõna Asesõna Nimisõna Omadussõna
## Lühike       4       4      390         22
## Pikk         4       4      390         22

Neli oodatud sagedust on < 5.

fisher.test(tab2)
## 
##  Fisher's Exact Test for Count Data
## 
## data:  tab2
## p-value = 0.4329
## alternative hypothesis: two.sided

Fisheri test arvutab täpset tõenäosust kohata valimis esinevat tegelike sageduste jaotust. Tõenäosuse arvutamiseks loob ta ridade ja tulpade marginaalsageduste põhjal kõik võimalikud viisid, kuidas sagedused võiksid olla lahtrites jaotunud. See aga teeb Fisheri testi arvutuslikult väga kompleksseks ja ressursiahneks, mistõttu ei tööta test väga suurte andmestikega.
Argumendiga workspace võib anda algoritmile juurde lisamälu, ent see võib teha aeglaseks kõik muud protsessid teie arvutis.

valimised <- read.delim("data/kandidaadid_2023.csv")
fisher.test(valimised$nimekiri, valimised$sugu, workspace = 100000000000)

G-test on hii-ruut-testile väga sarnane, ent selle valem statistiku arvutamiseks on pisut erinev: \(G = 2\sum_i O_i log(\frac{O_i}{E_i})\).

DescTools::GTest(tab2)
## 
##  Log likelihood ratio (G-test) test of independence without correction
## 
## data:  tab2
## G = 2.9189, X-squared df = 3, p-value = 0.4043

Vaata ka teisi kategoriaalsete andmete teste.

Monte-Carlo simulatsioon

Veel üks variant seoste hindamiseks tabelites, kus on ka väikeseid sagedusi, mille puhul tavaline hii-ruut-test muutub pirtsakaks, on kasutada p-väärtuse leidmiseks nn Monte-Carlo simulatsiooni. Selle käigus seotakse testi p-väärtus lahti eeldusest, et hii-ruut-statistik on pärit pidevast hii-ruut-jaotusest, kus kahe tunnuse vaheline seos muutub statistiliselt oluliseks kindla vabadusastmete arvu ja sellele vastava statistiku väärtuse korral.

chisq.test(tab2, simulate.p.value = TRUE)
## 
##  Pearson's Chi-squared test with simulated p-value (based on 2000
##  replicates)
## 
## data:  tab2
## X-squared = 2.8233, df = NA, p-value = 0.4288

Monte-Carlo simulatsioonis tehakse risttabeli oodatud sageduste põhjal mitu (nt 2000) juhuslikku “tegelike sagedustega” tabelit, varieerides iga kord pisut tõenäosust, millega mingi tunnuste kombinatsioon esineda võiks. “Tegelike sagedustega” juhuslikud tabelid on aga ikkagi küllalt lähedal oodatud sagedustega tabelile (ehk tabelile, mille põhjal kahe tunnuse vahel ei tohiks olla mingit seost).

Näiteks kui oodatud sagedustega tabel näeb välja selline:

chisq.test(tab2)$expected
## Warning in chisq.test(tab2): Chi-squared approximation may be incorrect
##        Arvsõna Asesõna Nimisõna Omadussõna
## Lühike       4       4      390         22
## Pikk         4       4      390         22

siis kolm juhuslikult genereeritud tabelit võiksiv välja näha umbes sellised:

(ts1 <- matrix(c(3,6,1,4,368,410,19,29), nrow= 2, 
               dimnames = list(c("Lühike", "Pikk"), 
                               c("Arvsõna", "Asesõna", "Nimisõna", "Omadussõna"))))
##        Arvsõna Asesõna Nimisõna Omadussõna
## Lühike       3       1      368         19
## Pikk         6       4      410         29
(ts2 <- matrix(c(8,5,3,7,383,381,25,28), nrow= 2, 
               dimnames = list(c("Lühike", "Pikk"),
                               c("Arvsõna", "Asesõna", "Nimisõna", "Omadussõna"))))
##        Arvsõna Asesõna Nimisõna Omadussõna
## Lühike       8       3      383         25
## Pikk         5       7      381         28
(ts3 <- matrix(c(3,2,3,1,394,395,21,21), nrow= 2, 
               dimnames = list(c("Lühike", "Pikk"), 
                               c("Arvsõna", "Asesõna", "Nimisõna", "Omadussõna"))))
##        Arvsõna Asesõna Nimisõna Omadussõna
## Lühike       3       3      394         21
## Pikk         2       1      395         21
# jne

Juhuslikult genereeritud tabelite põhjal leitakse hii-ruut-statistikud (kui on 2000 kordust, siis 2000 statistiku väärtust), kasutades esialgse tabeli oodatud sagedusi.

hiid <- c(
  sum(((ts1-t2$expected)^2)/t2$expected),
  sum(((ts2-t2$expected)^2)/t2$expected),
  sum(((ts3-t2$expected)^2)/t2$expected)
)
hiid
## [1] 8.403030 9.128788 3.946037

Lõpuks arvutatakse nende kordade osakaal simuleeritud tabelitest, mille puhul leitud hii-ruut-statistik oli sama suur või suurem kui algse (ilma pidevusparanduseta) testi statistik 2.82331. See osakaal ongi testi p-väärtus.

chisq.test(tab2, correct = FALSE)
## Warning in chisq.test(tab2, correct = FALSE): Chi-squared approximation may be
## incorrect
## 
##  Pearson's Chi-squared test
## 
## data:  tab2
## X-squared = 2.8233, df = 3, p-value = 0.4197
set.seed(1)
chisq.test(tab2, simulate.p.value = TRUE, B = 3000)
## 
##  Pearson's Chi-squared test with simulated p-value (based on 3000
##  replicates)
## 
## data:  tab2
## X-squared = 2.8233, df = NA, p-value = 0.4459

Järgmisel korral

Järgmisel korral läheme kahe tunnuse vaheliste seoste hindamise juurest regressioonimudelite juurde, mis võimaldavad ühe tunnuse muutumist kirjeldada enamate tunnuste kaudu. NB! Järgmine praktikum alles ülejärgmisel nädalal (01.04)!

Funktsioonide sõnastik

  • chisq.test() - tee hii-ruut-test, sisendiks sagedustega risttabel või andmetabeli kahe tulba nimed
  • DescTools::Gtest() - tee G-test, sisendiks sagedustega risttabel või andmetabeli kahe tulba nimed
  • fisher.test() - tee Fisheri täpne test, sisendiks sagedustega risttabel või andmetabeli kahe tulba nimed
  • pwr::pwr.chisq.test() - leia hii-ruut-testi võimsus
  • tidyverse::rownames_to_column() - muuda tabeli või data.frame-objekti reanimed eraldi tunnuseks/tulbaks andmestikus
  • vcd::assocstats() - leia kategoriaalsete tunnuste vahelise seose tugevuse mõõdikud