Testid jaotuse kohta,
usaldusintervallid (keskväärtusele ja mediaanile)
Pidev tunnus


Normaaljaotuse testimine

Shapiro-Wilk'i testi abil saab kontrollida, kas uuritava tunnuse jaotuseks on normaaljaotus (H0: tegemist on normaaljaotusega). Erinevalt teistest alternatiivsetest testidest ei pea me eelnevalt teadma, millise normaaljaotusega (millise keskväärtuse ja dispersiooniga) on tegemist:

shapiro.test(pikkus)

Tulemus:

        Shapiro-Wilk normality test

data:  pikkus
W = 0.9729, p-value = 1.034e-09

Näeme, et tudengite pikkus pole normaaljaotusega: p-väärtus on 1.034e-09 = 0,000000001034, mis on märksa väiksem kui harilikult kasutatav olulisuse nivoo 0,05. Seega ei saa nullhüpotees õige olla: uuritav tunnus (tudengite kaalud) pole normaaljaotusega.

Vaatame, kas naistudengite (sugu=1) pikkuste jaotus võiks olla normaaljaotus:

shapiro.test(pikkus[sugu==1])

Tulemus:

        Shapiro-Wilk normality test

data:  pikkus[sugu == 1]
W = 0.9954, p-value = 0.1337
Näeme, et naistudengite pikkuste jaotuseks võib olla normaaljaotus: Shapiro-Wilk'i testi p-väärtus on suhteliselt suur, suurem kui 0,05.


Kolmogorov-Smirnovi test
Kas uuritava tunnuse jaotus on just selline .... jaotus?

Kas uuritava tunnuse jaotus võiks olla justnimelt see etteantud pidev jaotus (Eksponentjaotus parameetriga 3 või normaaljaotus keskväärtusega 2 ja standardhälbega 3 või ...)? Sellisele küsimusele võime vastust otsida Kolmogorov-Smirnovi testi abil.

Kolmogorov-Smirnovi test võrdleb teoreetilist jaotusfunktsiooni empiirilise jaotusfunktsiooniga ja leiab nende kahe jaotusfunktsiooni vahelise kauguse. Kontrollib seejärel, kas nähtud erinevus võib olla tingitud juhususest (juhtumisi saime halva valimi?).

Näide 1
Soovime kontrollida, kas naistudengite (sugu=1) pikkuste jaotuseks võiks olla normaaljaotus keskväärtusega 170 ja standardhälbega 6, ehk kas naistudengite jaoks kehtib nullhüpotees H0: pikkus~N(170; 6)?

ks.test(pikkus[sugu==1], pnorm, mean=170, sd=6)

Selle käsu tulemus:

One-sample Kolmogorov-Smirnov test

data: pikkus[sugu == 1]
D = 0.1716, p-value = 1.62e-13
alternative hypothesis: two-sided

Warning message:
In ks.test(pikkus[sugu == 1], pnorm, mean = 170, sd = 6) :
ties should not be present for the Kolmogorov-Smirnov test

Väljundi selgitused:

Näeme, et naiste pikkused ei ole meie poolt väljapakutud normaaljaotusega: p-väärtus on äärmiselt pisike (1.62e-13 = 0,000000000000162), seega nullhüpotees ei saa kehtida (meie poolt väljapakutud jaotus ei sobi - naiste pikkused ei ole sellise normaaljaotusega, mille keskväärtus on 170 ja standardhälve 6).

Tuleb meeles pidada, et Kolmogorov-Smirnovi test kontrollib samal ajal nii seda, kas uuritava tunnuse jaotuseks on justnimelt normaaljaotus, kas uuritava tunnuse keskväärtus on 170 ja kas uuritava tunnuse standardhälve on 6. Kui ükski neist eeldustest ei pea paika, võib Kolmogorov-Smirnovi test otsustada alternatiivse hüpoteesi kasuks. Seega antud tulemuse korral ei oska me ka kindlalt öelda, miks nullhüpoteesi väide kummutati: kas sellepärast, et väljapakutud jaotuste pere - normaaljaotus - on vale, või on vale meie poolt pakutud keskväärtus või standardhälve.

Hoiatus (Warning message: ... ties should not be present for the Kolmogorov-Smirnov test) juhib tähelepanu kahele asjale korraga. Esmalt ei tohiks pideva jaotuse korral põhimõtteliselt esineda täpselt ühte ja sama väärtust kaks korda. Antud juhul tekib probleem sellest, et pikkused ei ole mõõdetud täpselt - parimal juhul on mõõtmistäpsuseks 1mm. Seega võiks öelda, et pikkused on ümmardatud täismillimeetriteks ja mitte kirja pandud täpselt. Normaaljaotus teoreetilise jaotusena ei luba aga ümmardamisvigasid.

Antud hoiatus juhib tähelepanu ka sellele, et kui valim on väike ja kui kokkulangevaid väärtuseid peaks olema palju, siis võib Kolmogorov-Smirnovi testi p-väärtus osutuda ka ebatäpselt arvutatuks, raporteeritud p-väärtus võib olla vale. Antud juhul ei tohiks p-väärtuse arvutusega siiski midagi väga valesti olla, sest vaatluste arv on suur ja kokkulangevaid väärtuseid pole siiski hirmus palju - seega võime loota, et p-väärtuse arvutamisel tehtav viga on väga väike, praktikas arvatavasti ignoreeritav.

Näide 2
Kontrollime kas tunnuse pvalues jaotuseks võiks olla ühtlane jaotus:

pvalues=c(0.001, 0.01, 0.02, 0.03, 0.2, 0.5, 0.6)

ks.test(pvalues, punif)

Tulemuseks saame:

One-sample Kolmogorov-Smirnov test

data: pvalues
D = 0.5414, p-value = 0.01897
alternative hypothesis: two-sided

Näeme, et vektorisse pvalues kogutud vaatluste jaotuseks ei saa olla ühtlane jaotus U(0,1), sest Kolmogorov-Smirnovi testi p-väärtus (0,01897) on väiksem kui tavapärane olulisuse nivoo 0,05.

R tunneb paljude jaotuste jaotusfunktsioone. Vajaliku jaotusfunktsiooni võid leida näiteks tippides R-is käsu ?Distributions. Sealt võid edasi uurida ka seda, milliseid jaotusparameetreid üks või teine jaotusfunktsioon lubab määrata (ja millised on nende jaotusparameetrite vaikimisi väärtused).

Näide 3
Keerukamal juhul võime kontrollitava jaotuse jaotusfunktsiooni ka ise kirja panna. Näiteks tahame kontrollida, kas uuritava tunnuse jaotus võiks olla järgmise jaotusfunktsiooniga:

   F(x)=0, kui x<0
   F(x)=x2/2, kui 0≤ x<1
   F(x)=1-(2-x)2/2, kui 1≤x<2
   F(x)=1, kui x≥2.

Sellise jaotuse testimiseks peame jaotusfunktsiooni väärtuseid arvutava funktsiooni esmalt ise kirja panema:

F=function(x){
    tulemus=rep(0, length(x))
    tulemus[0<=x & x<1]=x[0<=x & x<1]^2/2
    tulemus[1<=x & x<2]=1-(2-x[1<=x & x<2])^2/2
    tulemus[x>=2]=1
    return(tulemus)
}

Peale jaotusfunktsiooni väärtuseid arvutava funktsiooni tekitamist võime juba kontrollima hakata, kas meie vaatlused on justnimelt meie poolt kirjapandud jaotusega:

# Tekitame algandmed - antud juhul laseme arvutil vaatlused genereerida: vaatlused=runif(1000)+runif(1000)

# Kontrollime Kolmogorov-Smirnovi testi abil jaotuse sobivust:
ks.test(vaatlused, F)

Tulemus:

        One-sample Kolmogorov-Smirnov test

data:  tunnus
D = 0.0075, p-value = 0.622
alternative hypothesis: two-sided

Näeme, et antud juhul peame jääma nullhüpoteesi juurde - meie vaatlused võivad tõepoolest olla pärit ülaltpool kirjapandud jaotusest.


Hii-ruut test
Jaotuse sobivuse testimine - alternatiiv Kolmogorov-Smirnovi testile.

Kolmogorov-Smirnovi testi võimsus osutub vahel väga väikeseks - mõnikord ei suuda antud test aru saada isegi üsnagi suurte valimite korral, et uuritava tunnuse tegelikuks jaotuseks pole nullhüpoteesi poolt väidetav jaotus. Praktikas ettetulevais olukordades on sageli suurema võimsusega ehk terasemaks testiks hii-ruut test.

Hii-ruut testi puhul tuleb ette anda mingite sündmuste toimumiste sagedused - enamasti leitakse need käsu table abil - ja nende sündmuste toimumise teoreetilised tõenäosused (kui tõenäoliselt peaks üks või teine sündmus aset leidma nullhüpoteesi H0 kehtides, vaata ka näidet siit. Antud juhul vaatame sündmustena seda, kas juhuslik suurus sattus mingisse konkreetsesse väärtuste vahemikku. Vahemikud valime selliselt, et kui nullhüpotees kehtiks, siis peaks igasse vahemikku sattuma ootuspäraselt sama arv juhuslikke suuruseid.

Näide 1
Alljärgnev programm kontrollib, kas 15-ne statistilise testi poolt raporteeritud p-väärtuste jaotuseks võiks olla ühtlane jaotus.

Kuna paljude testide korral (näiteks t-test) on nullhüpoteesi kehtides p-väärtuste jaotuseks ühtlane jaotus, siis kontrollib alltoodud hii-ruut test sisuliselt järgmiseid hüpoteese:
H0: kõigi teostatud testide korral kehtis H0
H1: leidus teste, mille korral kehtis H1.

# Algandmed
pvaartused=c(0.06, 0.09, 0.12, 0.14, 0.2, 0.23, 0.234, 0.25,
    0.3, 0.31, 0.45, 0.55, 0.68, 0.78, 0.9)

# Mitmeks vahemikuks jagame
mitmeks=3

# valime lõikepiirid nii, et H0 kehtides peaksid
# tekkima võrdtõenäolised vahemikud - igasse vahemikku
# sattumise tõenäosus peaks olema sama.
# Lõikepiirid valime praegu ühtlase jaotuse kvantiile kasutades.
piirid=qunif(seq(0,1,length=mitmeks+1))

# Teeme hii-ruut testi.
# Igasse vaadeldavasse vahemikku sattumise tõenäosus peaks H0 korral
# olema samasuur: 1/mitmeks.
chisq.test(table(cut(pvaartused, piirid)), p=rep(1/mitmeks, mitmeks))

Ülaltoodu testi tulemuseks saame:


        Chi-squared test for given probabilities

data:  table(cut(pvaartused, piirid))
X-squared = 7.6, df = 2, p-value = 0.02237

Hii-ruut testi tulemusena näeme, et vektorisse pvaartused kirjutatud arvud ei ole pärit ühtlasest jaotusest (p-väärtus 0,02237 on väiksem kui 0,05. Seega peame nullhüpoteesi kummutama.

Järeldus: Antud 15 testi seas peab olema selliseid teste, mille puhul kehtis alternatiivne hüpotees. Kuigi me antud juhul ei saa täie kindlusega näidata milliste testide puhul just kehtis alternatiivne hüpotees, siiski aga saame öelda - nende testide seas leidus selliseid, kus õige pidi olema alternatiivne hüpotees (või mille puhul testi eeldused olid rikutud).

Näide 2
Kontrollime hii-ruut testi abil, kas meestudengite pikkus võiks olla normaaljaotusega (sellise hüpoteesi puhul oleks ehk targem kasutada Shapiro-Wilk'i testi, aga ka hii-ruut test on ehk vastuvõetav alternatiiv...)

Esimene katse - ära kasuta - tulemuseks saadud p-väärtus on vale!!!

# Võtame vaatluse alla meestudengite (sugu=2) pikkused.
# Eemaldame puuduvad väärtused käsuga na.omit.
tunnus=na.omit(pikkus[sugu==2])

# Mitmeks vahemikuks jagame
mitmeks=20

# valime lõikepiirid nii, et H0 kehtides peaksid
# tekkima võrdtõenäolised vahemikud - igasse vahemikku
# sattumise tõenäosus peaks olema sama.
piirid=qnorm(seq(0,1,length=mitmeks+1), mean=mean(tunnus), sd=sd(tunnus))

# Esimene test - esialgu vale tulemus!!!
chisq.test(table(cut(tunnus, piirid)), p=rep(1/mitmeks, mitmeks))

Tulemuseks saame:

        Chi-squared test for given probabilities

data:  table(cut(tunnus, piirid))
X-squared = 26.0541, df = 19, p-value = 0.1287

Miks tulemus ei sobi? Antud tulemus oleks õige, kui oleksime püstitanud enne andmete kogumist hüpoteesi: meestudengite pikkuste jaotuseks on normaaljaotus keskväärtusega 181,9878 ja standardhälbega 6,820965. Paraku tulime sellise keskväärtuse ja standardhälbe peale alles peale olemasolevate vaatlustega tutvumist - tegemist on valimi põhjal hinnatud suurustega. Seega oleme hinnanud (täiendavalt) kaks parameetrit olemasolevate vaatluste põhjal ja hii-ruut testi vabadusastmete arv peaks seega olema mitte 19 vaid kahe võrra väiksem, 19-2=17.

Õige p-väärtuse saamiseks peaksime seega muutma hii-ruut testi vabadusastmete arvu. Selleks salvestame leitud hii-ruut statistiku väärtuse ja algse vabadusastmete arvu ning vähendame vabadusastmete arvu kahe võrra (kaks vaatluste pealt hinnatud parameetrit). Seejärel küsime, kui tõenäoline on näha sedavõrd suurt (või veel suuremat) hii-ruut statistiku väärtust siis, kui kehtiks nullhüpotees (leiame p-väärtuse) kasutades juba õiget vabadusastmete arvu:

# Salvestame esialgse hii-ruut testi tulemused
abi=chisq.test(table(cut(tunnus, piirid)), p=rep(1/mitmeks, mitmeks))

# Leiame korrigeeritud vabadusastmete arvu kasutades p-väärtuse:
1-pchisq(abi$statistic, df=abi$parameter-2)

Tulemuseks saame:

 X-squared 
0.07348134 
Seega on hii-ruut testi korrigeeritud p-väärtus 0,07348... ehk ikkagi peame jääma nullhüpoteesi juurde: on võimalik, et meestudengite pikkuse jaotuseks on normaaljaotus.


Usaldusintervall keskväärtusele, hüpoteesid keskväärtuse kohta
T-test hüpoteeside kontrollimiseks uuritava tunnuse keskväärtuse kohta

Vahel soovime kirjeldada, kui täpselt me ikkagi teame uuritava tunnuse keskväärtust. Valimi keskmine on ju kõigest hinnang (ja enamasti vähemalt veidikene vale hinnang) uuritava populatsiooni keskväärtusele. Võimalikku vahemikku, kuhu võiks jääda uuritava tunnuse keskväärtus (keskmine üle kogu uuritava populatsiooni) kirjeldab usaldusintervall.

Näide. Leiame 95%-usaldusintervalli (95%-Confidence Interval, 95%-CI) tudengite pikkuste keskväärtusele:

t.test(pikkus)

Käsu tulemuseks saame:

        One Sample t-test

data:  pikkus
t = 523.2705, df = 659, p-value < 2.2e-16
alternative hypothesis: true mean is not equal to 0
95 percent confidence interval:
 170.4746 171.7588
sample estimates:
mean of x 
 171.1167 
Näeme, et 95%-usaldusintervall keskväärtusele on 170,474..171,759. Usaldusintervalli on viisakas ümmardada laiemaks (seega usaldusintervalli alumise piirina läks kirja 170,474, mitte aga 170,475 nagu tavaliste ümmardamisreeglite puhul).

Kui soovid usaldusintervalli salvestada (näiteks vektorisse UI) hilisemaks kasutamiseks (näiteks joonisele kandmiseks), siis seda saab teha järgmise käsuga:

UI=t.test(pikkus)$conf.int

Hiljem saab salvestatud väärtuseid kasutada üsna tavapäraselt:

> UI
[1] 170.4746 171.7588
attr(,"conf.level")
[1] 0.95
> UI[1]
[1] 170.4746
> UI[2]
[1] 171.7588

Kui vaja oleks leida 90%-usaldusintervalli, siis sobiv käsk selleks on järgmine:

t.test(pikkus, conf.level=0.9)

T-test
Vaikimisi kontrollib t.test käsu poolt tehtav t-test hüpoteesi, kas uuritava tunnuse keskväärtus võiks olla 0 (antud juhul siis H0: tudengite pikkuste keskväärtus on 0). Hüpotees on üsna jabur, õnneks lükkab ka t-test selle jabura hüpoteesi ümber. Kui aga sooviksime kontrollida mingit muud ja mõistlikumat hüpoteesi, siis peame nullhüpoteesi väite käsule ette andma lisaparameetri mu= abil. Näiteks kui soovime kontrollida, kas tudengite pikkuste keskväärtus võiks olla 170, saame seda kontrollida nii:

t.test(pikkus, mu=170)

        One Sample t-test

data:  pikkus
t = 3.4147, df = 659, p-value = 0.0006776
alternative hypothesis: true mean is not equal to 170
95 percent confidence interval:
 170.4746 171.7588
sample estimates:
mean of x 
 171.1167 

Tulemustest näeme, et p-väärtus tuli 0,0006776 ehk kõvasti väiksem kui 0,05 (väiksem tavapärasest olulisuse nivoost) - seega peame nullhüpoteesi kummutama (tudengite pikkuste keskväärtus ei saa olla 170).

T-testi tegemisel tasub piiluda ka vabadusastmete arvu (df=659). Siit näeme, et t-testi tegemiseks (ja usaldusintervalli leidmiseks) on kasutatud 659+1=660 tudengi andmeid. Kui algandmetega on midagi valesti, siis vale vabadusastmete arv võib olla üheks väga heaks võimaluseks tuvastada algandmetes esinevat (või nende väljanoppimise käigus tekkivat) viga.

Eeldustest
Leitud usaldusintervall ja t-testi p-väärtus on korrektselt arvutatud kui on täidetud järgmised nõuded:

  • Tegemist on juhusliku esindava valimiga, vaatlused on teineteisest sõltumatud (kes järgmisena valimisse sattub ei sõltu sellest, kes juba valimis on);
     
  • Uuritav tunnus on normaaljaotusega

    või

    uuritav tunnus pole küll normaaljaotusega, kuid valim on suur ja uuritava tunnuse keskväärtus eksisteerib (näiteks teame, et uuritava tunnuse võimalikud väärtused peavad asuma mingis kindlas vahemikus - näiteks inimeste pikkused jäävad kindlasti vahemikku 0cm..500cm).

    Kui uuritava tunnuse jaotus pole normaaljaotusega, kuid valim suhteliselt suur, siis tulemused on usaldusväärsed siis, kui kasutame tavapäraseid olulisuse nivoosid (0,1...0,01) ja vaatame klassikalisi usaldusintervalle (95%-usaldusintervalli; 90%-usaldusintervall vms). Kui soovime kasutada ekstreemseid usaldusnivoosid (näiteks seepärast, et teeme sadu tuhandeid teste ja soovime kasutada Bonferroni korrektsiooni), siis ei pruugi isegi väga suurtest valimitest piisata ja ka väikseimaki erinevused normaaljaotusest võivad tekitada probleeme.


Usaldusintervall mediaanile, hüpoteesid mediaani kohta
Märgitest hüpoteeside kontrollimiseks mediaani kohta ja samal testil põhinev usaldusintervall.

Keskväärtus pakub enam huvi sümmeetriliste jaotuste puhul. Ebasümmetriliste jaotuste korral - näiteks siis, kui esineb üksikuid väga suuri väärtuseid - ei iseloomusta keskväärtus kuigivõrd hästi uuritava tunnuse vaatluste paiknemist. Mediaan võib sellistel juhtudel osutuda hoopis informatsioonirikkamaks statistikuks. Aga kui täpselt me ikkagi teame populatsiooni mediaani? Kui soovime leida usaldusintervalli populatsiooni mediaanile, siis võime kasutada lisamoodulis BSDA paikneva funktsiooni SIGN.test abi.

Kui kasutatavale R-le pole veel lisamoodulit BSDA installeeritud, siis tee seda kas käsuga

install.packages("BSDA")

Mainitud käsku tuleb anda vaid kord - edaspidi on vajalik lisamoodul ja temas sisalduvad täiendavad funktsioonid sul kasutamiseks kenasti olemas.

Edasi tuleb lisamoodul BSDA ka kasutusse võtta (seda tuleb teha kord peale R-programmi käivitamist. Kui sulged R-i ja taaskäivitad ta, siis pead selle käsu uuesti andma):

library(BSDA)

Nüüd on R'i tekkinud lisaks terve hulk kasulikke funktsioone, nende seas ka funktsioon SIGN.test. Viimatimainitud funktsiooni kasutades leiamegi usaldusintervalli meestudengite kaalule:

SIGN.test(kaal[sugu==2], conf.level = 0.95)

Tulemus:

        One-sample Sign-Test

data:  kaal[sugu == 2]
s = 149, p-value < 2.2e-16
alternative hypothesis: true median is not equal to 0
95 percent confidence interval:
 73.04881 77.95119
sample estimates:
median of x 
         75 

                  Conf.Level  L.E.pt  U.E.pt
Lower Achieved CI     0.9289 74.0000 77.0000
Interpolated CI       0.9500 73.0488 77.9512
Upper Achieved CI     0.9511 73.0000 78.0000
Kuidas saadud tulemust interpreteerida, milleks sedavõrd palju usaldusintervalle? Usaldusintervalle trükitakse välja palju seepärast, et täpselt 95%-usaldusintervalli leida ei saa. Küll aga saab leida 92,89%-usaldusintervalli mediaanile ja 95,11%-usaldusintervalli mediaanile. Ohutum on raporteerida neist rangemat usaldusintervalli, seega peaksime 95%-usaldusintervallina pigem raporteerima vahemikku 73..78 (meestudengite kaalu mediaan asub mainitud vahemikus 95%-se kindlusega, täpsemalt öeldes isegi 95,11%-se kindlusega).

Testid mediaani kohta
Kui soovime testida, kas uuritava tunnuse mediaan on näiteks 70, siis seda saame teha märgitesti abil. Testi läbiviimiseks R-is on kaks võimalust. Esiteks lisamooduli BSDA funktsiooni SIGN.test abil:

SIGN.test(kaal[sugu==2], md=70)

Tulemus:

        One-sample Sign-Test

data:  kaal[sugu == 2]
s = 97, p-value = 5.815e-06
alternative hypothesis: true median is not equal to 70
95 percent confidence interval:
 73.04881 77.95119
sample estimates:
median of x 
         75 

                  Conf.Level  L.E.pt  U.E.pt
Lower Achieved CI     0.9289 74.0000 77.0000
Interpolated CI       0.9500 73.0488 77.9512
Upper Achieved CI     0.9511 73.0000 78.0000

Näeme, et märgitesti p-väärtus on 5,815e-06 = 0,000005815 ehk märgatavalt väiksem tavapärasest olulisuse nivoost 0,05-st, seega võime väita et uuritava tunnuse (meestudengite pikkuse) mediaan ei saa olla 70.

Variant 2
Teine võimalus märgitesti teha on kasutada R'i enda funktsiooni binom.test. Sisuliselt kontrollib ju märgitest seda, kas mediaanist suurema väärtuse saamise tõenäosus on 0.5. Hüpoteese tõenäosuse kohta saab aga kontrollida binom.test-käsu abil, vaata näiteid binom.test-käsu kohta ka siit. Kontrollimaks, kas tunnuse x mediaan uuritavas populatsioonis on 70, tuleb anda järgmine käsk:

binom.test(table(x-70>0))

Hoiatus! Binom.test-käsku kasuta vaid siis, kui uuritaval tunnusel ei ole korduvaid väärtuseid!!!


Märgitesti eeldustest:
Märgitest eeldab esindava valimi olemasolu; vaatlused peaksid olema sõltumatud.
Uuritava tunnuse jaotuseks võib olla suvaline pidev jaotus. Korduvate väärtuste (kui sama uuritava tunnuse väärtust esineb andmestikus palju kordi) on taunitav - vaata kas mõistad, mis on märgitesti nullhüpoteesi väide sellisel juhul.