10 Dvouvýběrové parametrické testy V předchozí kapitole jsme se zaměřili na situace, kdy jsme porovnávali vybraný parametr jednoho náhodného výběru (ať už parametr p nebo a2, za předpokladu, že náhodný výběr pochází z normálního rozdělení, parametr p za předpokladu, že náhodný výběr pochází z dvourozměrného normálního rozdělení, nebo parametr p, za předpokladu, že data pochází z alternativního rozdělení) s konkrétní hodnotou získanou například z literatury. V této kapitole se zaměříme na situaci, kdy vzájemně porovnáváme dva různé navzájem nezávislé náhodné výběry, popisující stejný znak (výšku, hmotnost novorozence, délku femuru, největší šířku mozkovny, apod). K tomu nám, za splnění určitých předpokladů, poslouží dvouvýběrové parametrické testy. Vzhledem k tomu, že náhodný výběr je reprezentantem vybrané populace, umožňují nám dvouvýběrové parametrické testy porovnat prostřednictvím náhodných výběrů a procesem testování hypotéz navzájem dvě populace. Konkrétním příkladem na využití dvouvýběrových testů je porovnávání znaků pohlavního dimorfismus (viz příklad ??. Příklad 10.1. Příklad aplikace dvouvýběrového parametrického testu Předpokládejme, že chceme vytvořit studii porovnávající výšku mužů a žen. V rámci studie tedy chceme porovnávat dvě populace, a sice populaci mužů a populaci žen. Sledovaným znakem (vlastností) je výška postavy. Pro první populaci vytvoříme reprezentativní vzorek změřením výšky postavy ri\ mužů. Získáme tak první náhodný výběr Xn,..., Xini. Pro druhou populaci vytvoříme reprezentativní vzorek změřením výšky postavy žen. Získáme tak druhý náhodný výběr X21,..., X'2n2- Poznamenejme, že rozsahy náhodných výběrů ri\ a se mohou lišit, tj. můžeme mít jiný počet naměřených výšek mužů a jiný počet naměřených výšek žen. Naším cílem by potom mohlo být porovnání střední hodnoty výšky mužů se střední hodnotou výšky žen pomocí procesu testování hypotéz, Dalším příkladem využití dvouvýběrových testů je například porovnání téhož znaku (šířka nadočnicového oblouku, délka levé klíční kosti, interorbitální šířka, apod.) u dvou populací (například malajské a čínské populace, německé a bantuské populace, či řecké a indické populace, apod.), nebo u dvou různých sexuálních orientací (heterosexuální populace a homosexuální populace). Metody uvedené v této kapitole jsou převážně rozšířením teorie uvedené v kapitole ?? na situaci dvou náhodných výběrů. Konkrétně si zde představíme test o podílu rozptylů o\ja\ dvou náhodných výběrů pocházejících z normálních rozdělení N(pi,a'f) a N(p2,c'Í), test o rozdílu středních hodnot p\ — P2 dvou náhodných výběrů pocházejících z normálních rozdělení N(pi,a'() a 7V(/i2,o"2)> kde rozptyly a\ a a\ jsou neznámé, ale shodné, tj. a'i = a'Í = o"2, a test o rozdílu středních hodnot p\ — P2 dvou náhodných výběrů pocházejících z normálních rozdělení N(pi,a'f) a N(p2, 2, > 2. Na hladině významnosti a testujeme jednu z následujících tří hypotéz oproti příslušné alternativní hypotéze. Hqi : crjVcrf = °o oproti Hn : <72/a'i 7^ tjg (pravostranná alt.) H03 '■ a\la\ ^ oproti 7íi3 : o\ja\ < tjg (levostranná alt.) Test nazýváme dvouvýběrový F-test o podílu rozptylů o\la\. Testovací statistika má tvar Fw = f|, (10.1) kde S'f je výběrový rozptyl prvního náhodného výběru a S% je výběrový rozptyl druhého náhodného výběru. Za platnosti nulové hypotézy pochází statistika Fyy z Fisherova-Snedecorova F-rozdělení o ri\ — 1 a ni — 1 stupních volnosti, tj. p _ ^1 v° p ^If — "^2 ~ ^ni-l,n2-l-^2 Kritický obor podle zvolené alternativní hypotézy má tvar Hn : a2/a2 ^ a2 W = (0; Fni_i,„2_i(a/2)> U (Fni_i,„2_i(l - a/2); 00) H12 : af/a'i > W = (Fni_iiTl2_i(l - a); 00) H13 : < ^0 W = (0; Fni_i,„2_i(a)> kde F„1_ii„2_i(a/2), Fni_iiTl2_i(l - a/2), Fni_iiII2_i(a), Fni_iiTl2_i(l - a) jsou kvantily F rozdělení o nx - 1 a ri2 — 1 stupních volnosti, jejichž hodnoty získáme pomocí softwaru W a implementované funkce qf(). Interval spolehlivosti má podle zvolené alternativní hypotézy jeden z následujících tvarů 4/4 . 4/4 Hn-.af/al^ai (d,h) - , V^m-i^-lU- - a/2) i'„1_ii„2_i(a/2) #12 : aí/A > (d, 00) = f--^4-T 5 00 v-Frii-l,7ij-l(l — a) H13 : < 4 (0, Ä) = (o; V'2 , , V -rn1-l,n2-llQíJ Poznámka: Protože parametry tr2 i a\ jsou z definice větší než 0, je i jejich podíl o\ja\ vždy větší než 0. Proto pravostranný interval spolehlivosti omezíme zdola hodnotou 0, namísto mínus nekonečnem. p-hodnota má v závislosti na zvolené alternativní hypotéze jeden z následujících tvarů Hn : a\la\ ^ a'l p-hodnota = 2min{Pr(rV < fw), Pr(Fw > fw)} H12 ■ v'llvl > ao p-hodnota = Yr(Fw > fw) = 1 - Pr(Fw < fw) H13 : O1/02 < °o p-hodnota = Pr(Fw < fw) kde Fyy je náhodná veličina, fw je realizace testovací statistiky Fyy (viz vzorec 10.1), tedy konkrétní číslo, a Pr(i c'i/c'i ao> ^e a'o = 1 (oboustranná alternativa) 2. Volba hladiny významnosti • Hladinu významnosti volíme v souladu se zadáním jako a = 0.05. 3. Testování kritickým oborem • Testovací statistika S'f _ 10.142622 _ 102.8727 ň Fw = -4 =-- =-= 3.60482 = 3.6048 W Sš 5.3420552 28.53755 • Kritický obor W=(0; F„1_li„2_1(a/2))U(FII1_liII2_1(l-a/2); 00) = (0; F10,22(0.025)) U (F10,22(0.975); 00) = (0 ; 0.2950131) U (2.699813 ; 00) 6 9 alpha <- 0.05 10 sl <- sd(bigo.Wms) 11 s2 <- sd(bigo.Wfs) 12 Fw <- sl"2 / s2"2 # 3.60482 13 qf(alpha/2, nl -1, n2 - 1) # 0.2950131 14 qf(l - alpha/2, nl -1, n2 - 1) # 2.699813 • Závěr testování Protože realizace testovací statistiky fw = 3.6048 náleží do kritického oboru, tj. fw e W, Hq zamítáme na hladině významnosti a = 0.05. 4. Testování intervalem spolehlivosti • Interval spolehlivosti (d, h) 4/4 4/4 Fni_i,„2_i(l-Q!/2) ' Fni_i,„2_i(Q!/2) 10.142622/5.3420552 10.142622/5.3420552 ^106,214(0.975) ' F106,2i4(0.025) 3.60482 3.60482 \ 2.699813 ' 0.2950131 / = (1.335211; 12.21919) = (1.3352; 12.2192) • Závěr testování Protože ctq = 1 nenáleží do Waldova 95% empirického oboustranného intervalu spolehlivosti, tj. tjg = 1 IS, Hq zamítáme na hladině významnosti a = 0.05. 5. Testování p-hodnotou • p-hodnota p-hodnota = 2min{Pr(Fw < fw), Pr(Fw > fw)} = 2min{Pr(FvK < 3.60482), 1 - Pr(rV < 3.60482)} = 2 min{0.9941918, 0.005808181} = 2 x 0.005808181 = 0.01161636 = 0.01162 15 p.val <- 2*min(pf(Fw, nl - 1, n2 - 1), 1 - pf(Fw, nl - 1, n2 - 1)) # 0.01161636 • Závěr testování Protože p-hdonota = 0.01162 je menší než a = 0.05, Hq zamítáme na hladině významnosti a = 0.05. 6. Interpretace výsledků: Na základě všech tří způsobů testování zamítáme hypotézu o shodě rozptylů a\ a o\. Mezi rozptylem šířky dolní čelisti mužů orientovaných jinak než heterosexuálně a rozptylem šířky dolní čelisti žen orientovaných jinak než heterosexuálě existuje statisticky významný rozdíl. Poznámka: K otestování nulové hypotézy o podílu rozptylů můžeme využít funkci var.test(). Vstupními parametry budou nejprve dva vektory reprezentující náhodné výběry, tj. bigo.Wms a bigo.Wfs, dále hodnota hladiny významnosti a zadaná prostřednictvím koeficientu spolehlivosti 1 — a nastavením hodnoty argumentu conf.level = 0.95 a nakonec typ zvolené alternativní hypotézy (oboustranná), zadaný pomocí argumentu alternativě == 'two.sideď. 7 16 var.test(bigo.Wms, bigo.Wfs, conf.level = 0.95, alternative = 'two.sided') 17 F test to compare two variances 18 19 data: bigo.Wins and bigo.Wfs 20 F = 3.6048, mil df = 10, denom df = 22, p-value = 0.01162 21 alternative hypothesis: true ratio of variances is not equal to 1 22 95 percent confidence interval: 23 1.335211 12.219186 24 sample estimates: 25 ratio of variances 26 3.60482 27 Součástí výstupu je hodnota testovací statistiky F = 3.6048, počty stupňů volnosti Fisherova rozdělení num df = 10 a denom df = 22, hranice intervalu spolehlivosti 1.335211 a 12.219186 a p-hodnota p-value = 0.01162. Jediné, co musíme stanovit zvlášť, jsou dolní a horní hranice kritického oboru. ★ 8 Příklad 10.3. Test o podílu rozptylů Mějme datový soubor ll-two-samples-means-skull.txt a proměnnou (skuli.H), popisující basion-bregmatickou výšku lebky v mm (viz sekce ??). Předpokládejme, že náhodná veličina X, popisující výšku lebky u žen, pochází z normálního rozdělení N(fii,a'f), a že náhodná veličina Y, popisující výšku lebky u mužů, pochází z normálního rozdělení N(p2, &%)• Na hladině významnosti a = 0.01 otestujte nulovou hypotézu, že rozptyl výšky lebky u mužů je větší nebo roven výšce u žen. Řešení příkladu 10.3 Pomocí příkazu read.delim() načteme datový soubor a příkazem na.omit() odtraníme ze souboru NA hodnoty. Pomocí operátoru [] vybereme z datové tabulky údaje o basion-bregmatické výšce lebky (skuli.H) u mužů (sex == 'm'), resp. u žen (sex == T). Dále zjistíme rozsahy obou náhodných výběrů. 28 data <- read.delim('00-Data//11-two-samples-means-skuli.txt') 29 data <- na.omit(data) 30 #head(data, n = 3) 31 skuli.Hm <- data[data$sex == 'm', 'skuli.H'] 32 skuli.Hf <- data[data$sex == 'f, 'skuli.H'] 33 34 nl <- length(skuli.Hm) # 215 35 n2 <- length(skull.Hf) # 107 Datový soubor obsahuje celkem 215 údajů o basion-bregmatické výšce lebky u mužů a 107 údajů o basion-bregmatické výšce lebky u žen. Řešení příkladu vede na test o podílu rozptylů. Před testováním hypotézy ze zadání musíme ověřit předpoklad normality pro oba náhodné výběry. Jelikož není uvedeno jinak, zvolíme pro test normality hladinu významnostia = 0.05. Nejprve testujeme nulovou hypotézu Hqi: Náhodný výběr výšek lebky mužů pochází z normálního rozdělení, oproti alternativní hypotéze Hn: Náhodný výběr výšek lebky mužů nepochází z normálního rozdělení. Protože rozsah náhodného výběru n = 215 je větší než 30, zvolíme na test normality Lillieforsův test. Normalitu ověříme graficky pomocí kvantilového di-gramu a histogramu superponovaného křivkou normálního rozdělení (viz obrázek ??). Datový soubor rozdělíme na základě Sturgerova pravidla do devíti ekvidistatních intervalů s šířkou 3 mm prostřednictvím stanovených hranic 119, 121,..., 146. 120.5 126.5 132.5 138.5 144.5 -3 -2 -1 0 1 2 3 teoreticky kvantu vyska lebky (v mm) Lillieforsův test: p-hodnota =0.1263 Obrázek 8: Histogram a kvantilový diagram výšky lebky mužů Protože p-hodnota = 0.1263 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr výšek lebky mužů tedy pochází z normálního rozdělení. Analogicky testujeme na hladině významnosti a = 0.05 nulovou hypotézu H02: Náhodný výběr výšek lebky žen pochází z normálního rozdělení, oproti alternativní hypotéze Hyi- Náhodný výběr výšek lebky žen nepochází z 9 normálního rozdělení Protože rozsah náhodného výběru n = 107 je větší než 30, zvolíme na test normality opět Lil-lieforsův test. V rámci histogramu rozdělíme soubor do osmi ekvidistatních intervalů s šířkou 3 mm prostřednictvím stanovených hranic 114, 117,..., 138 (viz obrázek 9). Obrázek 9: Histogram a kvantilový diagram výšky lebky žen Protože p-hodnota = 0.0989 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr výšek lebky žen tedy pochází z normálního rozdělení. Jelikož oba náhodné výběry pochází z normálních rozdělení, můžeme k testování hypotézy ze zadání použít parametrický test o podílu rozptylů. Zde je vhodné upozornit, že v zadání máme uvedené přesné znění nulové hypotézy. Zbývá tedy vhodně zvolit hypotézu alternativní. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Rozptyl výšky lebky mužů je větší nebo roven rozptylu výšky lebky žen. Hi : Rozptyl výšky lebky mužů je menší než rozptyl výšky lebky žen. • matematická formulace nulové a alternativní hypotézy H0:a'i>a'i o\jal > crg, kde crg = 1 H\ : a\ < a2 —> c'i/c'i < °o> kde a'o = 1 (levostranná alternativa) 2. Volba hladiny významnosti • Hladinu významnosti volíme ze zadání jako a = 0.01. 3. Testování kritickým oborem • Testovací statistika S_l 4-835494^ = = = W Si 4.6532562 21.65279 36 alpha <- 0.01 37 sl <- sd(skull.Hm) 38 s2 <- sd(skull.Hf) 39 Fw <- sl"2 / s2"2 # 1.079861 10 • Kritický obor W=(0;Fni_hn2_1(a)) = (0; F2i4,i07(0.01)> = (0; 0.683192) = (0; 0.6832) 40 qf(alpha, nl - 1, n2 - 1) # 0.683192 • Závěr testování Protože realizace testovací statistiky fyy = 1.0799 nenáleží do kritického oboru, tj. fyy £ W, Hq nezamítáme na hladině významnosti a = 0.01. 4. Testování intervalem spolehlivosti • Interval spolehlivosti (0,h)=(0; 4/4 ) 4.8354942/4.6532562\ ^214,106(0.01) ) _ ŕ 1.079861 \ ~ V ' 0.683192 / = (0; 1.580611) = (0; 1.5806) 41 HH <- (sl"2 / s2"2) / qf(alpha, nl - 1, n2 - 1) # 1.580611 • Závěr testování Protože ctq = 1 náleží do Waldova 99% empirického pravostranného intervalu spolehlivosti, tj. 0g = 1 € IS, Hq nezamítáme na hladině významnosti a = 0.01. 5. Testování p-hodnotou • p-hodnota p-hodnota = Pr(Fw < fw) = Pr(Fw < 0.92060453) = 0.3165 42 p.val <- pf(Fw, nl - 1, n2 - 1) # 0.3165 • Závěr testování Protože p-hdonota = 0.3165 je větší než a = 0.01, Hq nezamítáme na hladině významnosti a = 0.01. 6. Interpretace výsledků: Rozptyl výšky lebky mužů není statisticky významně vyšší než rozptyl výšky lebky žen. Poznámka: K otestování nulové hypotézy o podílu rozptylů můžeme využít funkci var.test(), kde hodnotu argumentu conf.level nastavíme na hodnotu 0.99 a typ alternativní hypotézy zvolíme pomocí argumentu alternativě = 'less' jako levostranný. 11 43 var.test(skull.Hm, skull.Hf, conf.level = 0.99, alternative = 'less') 44 F test to compare two varianc es 45 46 data: skull.Hm and skull.Hf 47 F = 1.0799, mil df = 214, denom df = 106, p-value = 0.6685 48 alternative hypothesis: true ratio of variances is less than 1 49 99 percent confidence interval: 50 0.000000 1.580611 51 sample estimates: 52 ratio of variances 53 1.079861 54 Součástí výstupu je hodnota testovací statistiky F = 1.0799, počty stupňů volnosti Fisherova rozdělení num df = 214 a denom df = 106, hranice 99% Waldova empirického pravostranného intervalu spolehlivosti 0 a 1.5806 a p-hodnota p-value = 0.6685. Jediné, co musíme stanovit zvlášť, je horní hranice kritického oboru. 12 Příklad 10.4. Test o podílu rozptylů Mějme datový soubor 18-more-samples-variances-clavicle.txt a proměnnou cla.L popisující největší délku klíční kosti z pravé strany v mm (viz sekce ??). Předpokládejme, že náhodná veličina X, popisující největší délku klíční kosti z pravé strany u anglické populace, pochází z normálního rozdělení N(fii,a'f), a že náhodná veličina Y, popisující největší délku klíční kosti z pravé strany u indické populace z Amritsaru, pochází z normálního rozdělení N(p2, Na hladině významnosti a = 0.10 zjistěte, zda je rozptyl největší délky klíční kosti z pravé strany u anglické populace statisticky významně větší než u indické populace z Amritsaru. Řešení příkladu 10.4 Pomocí příkazu read.delim() načteme datový soubor a příkazem na.omit() odtraníme ze souboru NA hodnoty. Pomocí operátoru [] vybereme z datové tabulky údaje o největší délce klíční kosti z pravé strany (cla.L) u jedinců anglické populace (pop == 'eng'), resp. indické populace z Amritsaru (pop == 'indl'). Nakonec zjistíme rozsahy obou náhodných výběrů. 55 data <- read.delim('00-Data//18-more-samples-variances-clavicle.txt') 56 data <- na.omit(data) 57 #head(data, n = 3) 58 cla.Le <- data[data$pop == 'eng', 'cla.L'] 59 cla.Li <- data[data$pop == 'indl', 'cla.L'] 60 61 nl <- length(cla.Le) # 70 62 n2 <- length(cla.Li) # 120 Datový soubor obsahuje celkem 70 údajů o největší délce klíční kosti z pravé strany u jedinců anglické populace a 120 údajů o největší délce klíční kosti z pravé strany u jedinců indické populace z Amritsaru. Řešení příkladu vede na test o podílu rozptylů. Před testováním nulové hypotézy ze zadání musíme ověřit normalitu obou náhodných výběrů. Hladinu významnosti a pro test normality stanovíme štandartne, t. a = 0.05. Nejprve testujeme hypotézu Hq\: Náhodný výběr největších délek klíčních kostí z pravé strany u jedinců anglické populace pochází z normálního rozdělení, oproti alternativní hypotéze Hn: Náhodný výběr největších délek klíčních kostí z pravé strany u jedinců anglické populace nepochází z normálního rozdělení. Protože rozsah náhodného výběru je větší než 30, ověříme předpoklad normality Lillieforsovým testem. Grafické ověření provedeme na základě kvantilového diagramu a histogramu superponovaného křivkou normálního rozdělení (viz obrázek ??). Datový soubor rozdělíme do devíti ekvi-distatních intervalů s šířkou 3 mm stanovením hranic 126,133,..., 175. 129.5 143.5 157.5 171.5 -2 -1 0 1 2 teoreticky kvantil délka pravé klicni kosti (v mm) Lillieforsuv test: p-hodnota =0.1958 Obrázek 10: Histogram a kvantilový diagram délky pravé klíční kosti u mužů anglické populace Jelikož p-hodnota = 0.1958 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný největších délek klíčních kostí z pravé strany u jedinců anglické populace pochází z normálního rozdělení. 13 Analogicky testujeme na hladine významnosti a = 0.05 nulovou hypotézu H02: Náhodný výběr největších délek klíčních kostí z pravé strany u jedinců indické populace z Amritsaru pochází z normálního rozdělení, oproti alternativní hypotéze H\2'- Náhodný výběr největších délek klíčních kostí z pravé strany u jedinců indické populace z Amritsaru nepochází z normálního rozdělení. Protože rozsah náhodného výběru n = 120 je větší než 30, ověříme předpoklad normality opět Lillieforsovým testem. V rámci histogramu rozdělíme data do osmi ekvidistatních intervalů s šířkou 8 mm prostřednictvím stanovených hranic 122,130,..., 170 (viz obrázek E10-Test-ss-cla.LAI). 126 134 142 150 158 166 -2 -1 0 1 2 teoreticky kvantil délka pravé klicni kosti (v mm) Lillieforsuv test: p-hodnota =0.0956 Obrázek 11: Histogram a kvantilový diagram délky pravé klíční kosti u mužů indické populace Protože p-hodnota = 0.0956 je větší než 0.05, nulovou hypotézu o normalitě dat nezamítáme na hladině významnosti a = 0.05. Náhodný výběr největších délek klíčních kostí z pravé strany u jedinců indické populace z Amritsaru pochází z normálního rozdělení. Jelikož oba náhodné výběry pochází z normálních rozdělení, můžeme k testování hypotézy ze zadání použít parametrický test o podílu rozptylů. Zde je vhodné upozornit, že v zadání příkladu se nepíše nic o znění nulové hypotézy. Ze zadání víme, že se snažíme prokázat, že rozptyl u anglické populace je větší než rozptyl u indické populace z Armitsaru. Jak bylo zmíněno v úvodní části kapitoly ??, tvrzení, jehož platnost se snažíme dokázat, je vždy součástí alternativní hypotézy. Ze zadání tedy známe tvar alternativní hypotézy, zatímto znění nulové hypotézy vhodně doplníme. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Rozptyl největší délky klíční kosti z pravé strany u anglické populace je menší nebo roven rozptylu u indické populace z Amritsaru. Hi : Rozptyl největší délky klíční kosti z pravé strany u anglické populace je větší než rozptyl u indické populace z Amritsaru. • matematická formulace nulové a alternativní hypotézy H0:af a2 —> o"?/0"! > °oj kde tjg = 1 (pravostranná alternativa) 2. Volba hladiny významnosti • Hladinu významnosti volíme jako a = 0.10 (viz zadání příkladu). 3. Testování kritickým oborem • Testovací statistika S'{ _ 10.510072 _ 110.4616 Fw = -^ = ———^ =_______= 1.448243 = 1.4482 '-'2 8.7334322 76.27283 14 63 alpha <- 0.10 64 sl <- sd(cla.Le) 65 s2 <- sd(cla.Li) 66 Fw <- sl"2 / s2"2 # 1.448242 • Kritický obor W = (Fni_iiTl2_i(l - a); 00) = (r69,H9(0.90); oo> = (1.307459oo) = (1.3075; 00) 67 qf(l - alpha, nl - 1, n2 - 1) # 1.307459 • Závěr testování Protože realizace testovací statistiky fw = 1.4482 náleží do kritického oboru, tj. fw e W, Hq zamítáme na hladině významnosti a = 0.10. 4. Testování intervalem spolehlivosti • Interval spolehlivosti (eř, 00) 4/4 Fni-i,„2-i(a); 00 10.510072/8.7334322 ^69,119(0.90); 00 1.448243 1.307459 ' (1.107678; 00) = (1.1077; 00) 68 (DH <- (sl"2 / s2"2) / qf(1 - alpha, nl - 1, n2 - 1)) # 1.107677 [1] 1.107677 69 • Závěr testování Protože ctq = 1 nenáleží do Waldova 90% empirického levostranného intervalu spolehlivosti, tj. 0g = 1 € IS, Hq zamítáme na hladině významnosti a = 0.10. 5. Testování p-hodnotou • p-hodnota p-hodnota = Pr(rV > fw) = 1 - Pr(Fw < fw) = 1 - Pr(rV < 1-307459) = 0.03845745 = 0.0385 • Závěr testování Protože p-hdonota = 0.0385 je menší než a = 0.10, Hq zamítáme na hladině významnosti a = 0.1. 6. Interpretace výsledků: Za základě všech tří typů testování zamítáme nulovou hypotézu na hladině významnosti a = 0.10. Rozptyl největší délky klíční kosti z pravé strany u anglické populace je statisticky významně větší než rozptyl u indické populace z Amritsaru. 15 70 p.val <- 1 - pf(Fw, nl - 1, n2 - 1) # 0.03845745 Poznámka: K otestování nulové hypotézy o podílu rozptylů můžeme využít funkci var.test(), kde hodnotu argumentu conf.level nastavíme na hodnotu 0.90 a typ alternativní hypotézy zvolíme pomocí argumentu alternativě = 'greater' jako pravostranný. 71 var.test(cla.Le, cla.Li, conf.level = 0.90, alternativě = 'greater') F test to compare two variances data: cla.Le and cla.Li F = 1.4482, num df = 69, denom df = 119, p-value = 0.03846 alternative hypothesis: true ratio of variances is greater than 1 90 percent confidence interval: 1.107677 Inf sample estimates: ratio of variances 1.448242 Součástí výstupu je hodnota testovací statistiky F = 1.4482, počty stupňů volnosti Fisherova rozdělení num df = 69 a denom df = 119, hranice 90% Waldova empirického levostranného intervalu spolehlivosti 1.1077 a Inf a p-hodnota p-value = 0.0385. Jediné, co musíme stanovit zvlášť, je dolní hranice kritického oboru. 16 Příklad 10.5. Test o podílu rozptylů Mějme datový soubor 19-more-samples-correlations-skull.txt a proměnnou intorb.B popisující interorbitální šířku v mm (viz sekce ??). Na hladině významnosti a = 0.01 zjistěte, zda je rozptyl interorbitální šířky mužů malajské populace statisticky významně větší než rozptyl interorbitální šířky mužů čínské populace. Řešení příkladu 10.5 Pomocí příkazu read.delim() načteme datový soubor a odstraníme ze souboru chybějící hodnoty. Pomocí operátoru [] vybereme z datové tabulky údaje o interorbitální šířce (intorb.B) u mužů malajské populace (pop == 'mal'), resp. u mužů čínské populace (pop == 'cin'). Nakonec zjistíme rozsahy obou náhodných výběrů. 83 data <- read.delim('00-Data//19-more-samples-correlations-skuli.txt') 84 data <- na.omit(data) 85 intorb.BM <- data[data$pop == 'mal', 'intorb.B'] 86 intorb.BC <- data[data$pop == 'cin', 'intorb.B'] 87 88 nl <- length(intorb.BM) # 73 89 n2 <- length(intorb.BC) # 19 Datový soubor obsahuje celkem údaje o interorbitální šířce 73 mužů malajské populace a 19 mužů čínské populace. Řešení příkladu vede na test o podílu rozptylů. Před testováním hypotézy ze zadání musíme ověřit normalitu obou náhodných výběrů. Hladinu významnosti zvolíme pro tento účel a = 0.05. Protože rozsah náhodného výběru interorbitálních šířek mužů malajské populace je větší než 30, použijeme na ověření předpokladu normality Lillieforsův test. Rozdělení datového souboru vizualizujeme pomocí kvantilového diagramu a histogramu (viz obrázek 12). Datový soubor rozdělíme do sedmi ekvidistatních intervalů s šířkou 2 mm stanovením hranic 16,18,..., 30 mm. 17 19 21 23 25 27 29 -2 -1 0 1 2 teoreticky kvantil interorbitalni sirka (v mm) Lillieforsův test: p-hodnota = 9e-04 Obrázek 12: Histogram a kvantilový diagram interorbitální šířky mužů malajské populace Jelikož p-hodnota = 9 x 10-4 je menší než 0.05, nulovou hypotézu zamítáme na hladině významnosti a = 0.05. Takéž histogram ukazuje porušení předpoklaud normality zejména vyšikmením naměřených hodnot doleva oproti předpokládanému rozdělení. V kvantilovém diagramu vidíme odklon hodnot od referenční přímky v pravém horním rohu. Náhodný výběr interorbitálních šířek mužů malajské populace nepochází z normálního rozdělení. Dále otestujeme hypotézu o normalitě náhodného výběru interorbitálních šířek mužů čínské populace. Protože rozsah náhodného výběru n = 19 je menší než 30, ověříme předpoklad normality Shapiro-Wilkovým testem. V rámci histogramu rozdělíme data do pěti ekvidistatních intervalů s šířkou 2 mm prostřednictvím stanovených hranic 18, 20,..., 28 (viz obrázek 13). Protože p-hodnota = 0.3475 je větší než 0.05, nulovou hypotézu o normalitě dat nezamítáme na hladině významnosti a = 0.05. Z obrázku 13 vidíme, že histogram věrně kopíruje tvar křivky normálního rozdělení. Kvantilový diagram ukazuje těsnou příchylnost bodů k referenční přímce vyjma krajného levého bodu. Odchýlení tohoto bodů však 17 19 21 23 25 27 interorbitalni sirka (v mm) -2-10 1 2 teoreticky kvantil Shapiro-Wilkuv test: p-hodnota = 0.3475 Obrázek 13: Histogram a kvantilový diagram interorbitalni šířky mužů čínské populace podle Shapiro-Wilkova testu nenarušuje fatálně předpoklad normality náhodného výběru. Náhodný výběr interor-bitálních šířek u mužů čínské populace pochází z normálnho rozdělení. Jelikož náhodný výběr interorbitálních šířek mužů malajské populace nesplňuje předpoklad normality, nemůžeme provést parametrický test o podílu rozptylů. 18 10.2 Klasický dvouvýběrový ŕ—test o rozdílu středních hodnot \i\ — /i2 Nechť .Xii,..., X\ni je náhodný výběr z N(p,i, a'f), a X21, ■ ■ ■, X'2n2 Je na něm nezávislý náhodný výběr z rozdělení N(fi2, o"2), přičemž ri\ > 2, n2 > 2 a a\ a o\ jsou neznámé ale shodné rozptyly, tj. a'f = a\. Nechť /uq je konstanta. Na hladině významnosti a testujeme jednu z následujících tří hypotéz oproti příslušné alternativní hypotéze. #01 : Mi ~~ M2 = Mo oproti Třn : Mi — M2 7^ Mo (oboustranná alt.) #02 : Mi ~~ M2 < Mo oproti 7íi2 : Mi — M2 > Mo (pravostranná alt.) 7ío3 : Mi — M2 > Mo oproti .H13 : Mi — M2 < Mo (levostranná alt.) Test nazýváme klasický dvouvýběrový í-test o rozdílu středních hodnot pi — M2- Testovací statistika má tvar TW = ^-M^-^. (10.2) S** — + — kde Mi je výběrový průměr a ri\ je rozsah prvního náhodného výběru, M2 je výběrový průměr a n2 je rozsah druhého náhodného výběru, dále S* je aritmetický průměr výběrových rozptylů S'f a 5| (viz kapitola ??), kde 5f je výběrový rozptyl prvního náhodného výběru a S% je výběrový rozptyl druhého náhodného výběru a konečně Mo je konstanta z nulové hypotézy. Za platnosti nulové hypotézy pochází statistika Tw ze Studentova rozdělení o ni + n2 — 2 stupních volnosti, tj. (Mi - M2) - mo Ho . jvk — -; ~ ~ *n1+n2-2- Kritický obor podle zvolené alternativní hypotézy má tvar Hn : Mi - M2 7^ Mo W = (-00; í„1+„2_2(a/2)) U (í„1+„2_2(l - a/2); 00) #12 : Mi - M2 > Mo = (íru+n2-2(l - a); 00) #13 : Mi - M2 < Mo = (-00; í„1+„2_2(a)) kde íij1+IJ2_2(qí/2), íij1+IJ2_2(1 —a/2), íij1+IJ2_2(qí) a íIJ1+IJ2_2(l —a) jsou kvantity Studentova rozdělení o ri\ +n2 — 2 stupních volnosti, jejichž hodnoty získáme pomocí 'Q- a implementované funkce qt(). Interval spolehlivosti má podle zvolené alternativní hypotézy jeden z následujících tvarů #11 : Mi -M2 7^ Mo (d,h) = [mx - m2 - s*y^ + ^í„1+„2_2(l - a/2); mx - m2 - + ^í„1+„2_2(a/2)) #12 : Mi - M2 > Mo (eř,oo) = (mi - m2 - s*^^ + ^í„1+„2_2(l - a); 00) #13 : Mi - M2 < Mo (-00, h) = (-00; mi - m2 - s* y^T^í„i+rl2-2(a)) p-hodnota má v závislosti na zvolené alternativní hypotéze jeden z následujících tvarů Tin : fi jí fi0 p-hodnota = 2min{Pr(7V < tw), Pr(7V > i^)} #12 : M > Mo p-hodnota = Pr(7V > = 1 — Pr(?V < tw) H13 : m < Mo p-hodnota = Pr(2V < tw) kde Tvk je náhodná veličina, tw je realizace testovací statistiky Tw (viz vzorec 10.2), tedy konkrétní číslo, a Pr(7V < tw) je distribuční funkce Studentova rozdělení o ri\ + n2 — 2 stupních volnosti, jejíž hodnotu získáme pomocí a implementované funkce pt(). 19 Příklad 10.6. Test o rozdílu středních hodnot Mějme datový soubor 01-one-sample-mean-skull-mf.txt a proměnnou skuli.B popisující největší šířku mozkovny v mm (viz sekce ??). Předpokládejme, že náhodná veličina X, popisující největší šířku mozkovny u mužů, pochází z normálního rozdělení N(fii,a'f), a že náhodná veličina Y, popisující největší šířku mozkovny u žen, pochází z normálního rozdělení N(p2, c'Í)- Na hladině významnosti a = 0.05 otestujte hypotézu o shodě středních hodnot největší šířky mozkovny u mužů a u žen. Řešení příkladu 10.6 Příkazem read.delim() načteme datový soubor a příkazem na.omit() z něj odtraníme NA hodnoty. Pomocí operátoru [] vybereme z datové tabulky údaje o největší šířce mozkovny (skuli.B) u mužů (sex == 'm'), resp. u žen (sex == 'f) a zjistíme rozsahy obou náhodných výběrů. 90 data <- read.delim('00-Data//01-one-sample-mean-skuli-mf.txt') 91 data <- na.omit(data) 92 #head(data, n = 3) 93 skuli.Bm <- data[data$sex == 'm', 'skuli.B'] 94 skuli.Bf <- data[data$sex == 'f, 'skuli.B'] 95 96 nl <- length(skuli.Bm) # 216 97 n2 <- length(skull.Bf) # 109 Datový soubor obsahuje celkem 216 údajů o největší šířce mozkovny u mužů a 107 údajů o největší šířce mozkovny u žen. Řešení příkladu vede na test o rozdílu středních hodnot. Před testováním nulové hypotézy ze zadání musíme ověřit splnění předpokladu normalitu obou náhodných výběrů, abychom zjistili, zda použít parametrický nebo neparametrický test. V případě splnění předpokladu normality musíme provést test o shodě rozptylů, abychom určili, zda máme použít klasický dvouvýběrový parametrický í-test nebo Welchův dvouvýběrový í-test. Jelikož není stanoveno jinak, zvolíme pro testy normality i případný test o podílu rozptylů hladinu významnosti a = 0.05. Nejprve testujeme hypotézu H^i'- Náhodný výběr největších šířek mozkovny u mužů pochází z normálního rozdělení. oproti alternativní hypotéze Hn: Náhodný výběr největší šířce mozkovny u mužů nepochází z normálního rozdělení. Jelikož náhodný výběr sestává z 215 pozorování, což je více než 30, provedeme test normality Lillieforsovým testem. Normalitu ověříme též graficky vykreslením kvantilového diagramu (viz sekce ??) a pomocí histogramu superponovaného křivkou normálního rozdělení (viz sekce ??). Datový soubor rozdělíme na základě Sturgerova pravidla (viz sekce ??) do devíti ekvidistatních intervalů s šířkou 3 mm prostřednictvím stanovených hranic 123,126,..., 150 (viz obrázek 14). 124.5 130.5 136.5 142.5 148.5 -3 -2 -1 0 1 2 3 teoreticky kvantu nejvetsi sirka mozkovny (v mm) Lillieforsuv test: p-hodnota =0.0766 Obrázek 14: Histogram a kvantilový diagram největší šířky mozkovny mužů Protože p-hodnota = 0.0766 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. 20 Náhodný výběr nej větších šířek mozkovny mužů pochází z normálního rozdělení. Vykreslený histogram i kvantilový diagram výsledek testu normality podporují. Histogram kopíruje tvar křivky normálního rozdělení a ve vykresleném kvantilovém diagramu se body pohybují okolo referenční čáry. Analogicky testujeme hypotézu H02: Náhodný výběr největších šířek mozkovky žen pochází z normálního rozdělení. oproti alternativní hypotéze H12: Náhodný výběr největších šířek mozkovny žen nepochází z normálního rozdělení. Jelikož rozsah náhodného výběru je dostatečně vysoký, použijeme na otestování hypotézy o normalitě Lillieforsův test. V rámci histogramu rozdělíme soubor do osmi ekvidistatních intervalů s šířkou 4mm prostřednictvím stanovených hranic 116,120 ..., 148 (viz obrázek 15). 118 126 134 142 -2 -1 0 1 2 teoreticky kvantil nejvetsi sirka mozkovny (v mm) Lillieforsův test: p-hodnota =0.0638 Obrázek 15: Histogram a kvantilový diagram největší šířky mozkovny žen Protože p-hodnota = 0.0638 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr největší šířky mozkovny žen také pochází z normálního rozdělení. Vykreslený histogram i kvantilový diagram podporují výsledek testování. Jelikož oba náhodné výběry splňují předpoklad normality, můžeme k testování hypotézy ze zadání použít parametrický test. Zda je vhodné zvolit klasický dvouvýběrový í-test nebo Welchův dvouvýběrový í-test rozhodneme na základě výsledku testu o shodě rozptylů obou náhodných výběrů. Na hladině významnosti a = 0.05 testujeme hypotézu Hq : uf/a'^ = 1 (rozptyly a\ a o-\ jsou shodné) oproti alternativní hypotéze H\ : o-f/a'^ 7^ 1 (rozptyly a\ a a\ nejsou shodné). Nulovou hypotézu otestujeme pomocí p-hodnoty získané pomocí funkce var.test(). Kompletní postup testování by byl analogický postupu uvedenému v příkladu 10.2. var.test(skuli.Bm , skuli.Bf, alternativě = 'two.sided')$p.val [1] 0.761025 Jelikož p-hodnota = 0.761025 je větší než 0.05, nulovou hypotézu o shodě rozptylů nezamítáme na hladině významnosti a = 0.05. Na základě výsledků testování tedy předpokládáme, že oba rozptyly jsou shodné. Protože oba výběry pochází z normálních rozdělení, jejichž rozptyly jsou shodné, použijeme na otestování hypotézy ze zadání klasický dvouvýběrový í-test o rozdílu středních hodnot pi — P2. Řešení si nyní uvedeme v posloupnosti sedmi kroků. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Střední hodnota největší šířky mozkovny mužů a žen jsou shodné. H\ : Střední hodnota největší šířky mozkovny mužů a žen nejsou shodné. 21 • matematická formulace nulové a alternativní hypotézy Ho ■ Mi = M2 Mi - M2 = Mo, kde /i0 = 0 Hi : Hij^ H2 Mi - M2 7^ Mo, kde /i0 = 0 (oboustranná alternativa) Volba hladiny významnosti • Hladinu významnosti volíme ze zadání a = 0.05. Testování kritickým oborem • Aritmetický průměr výběrových rozptylů S* (ni _ i)S'f + (n2 - 1)S22 ni + n2 (216 - 1) x 4.8246422 + (109 - 1) x 4.6959912 216 + 109-2 215 x 23.27717+ 108 x 22.05233 323 7386.243 = ^22.86763 = 4.7820109 = 4.7820 ÓZlÓ 100 sl <- sd(skull.Bm) 101 s2 <- sd(skull.Bf) 102 sh <- sqrt(((nl - 1) * sl"2 + (n2 - 1) * s2"2) / (nl + n2 - 2)) # 4.7820109 • Testovací statistika (Mi - M2) - Mo J-w —-, V 711 712 _ (137.1852 - 134.1468) -0 4.7820109^ + ^ _ 3.038396 ~ 4.7820109 x 0.1174902 3.038396 5.40794 = 5.4079 0.5618394 103 alpha <- 0.05 104 muO <- 0 105 ml <- mean(skull.Bm) 106 m2 <- mean(skull.Bf) 107 tw <- ((ml - m2) - muO) / (sh * sqrt(l / nl + 1 / n2)) # 5.40794 Kritický obor W = (-00; í„1+„2_2(a/2)) U (í„1+„2_2(l - a/2); oo) = (-oo; í2i6+io9-2(0.05/2)) U (í2i6+i09-2(l " 0.05/2); oo) = (-oo; í323(0.025)) U (í323(0.975); oo) = (-oo; -1.967336) U (1.967336; oo) 22 108 qt(alpha/2, nl + n2 - 2) # -1.967336 109 qt(l-alpha/2, nl + n2 - 2) # 1.967336 • Závěr testování Protože realizace testovací statistiky tyy = 5.4079 náleží do kritického oboru, tj. tyy e W, Hq zamítáme na hladině významnosti a = 0.05. 4. Testování intervalem spolehlivosti • Interval spolehlivosti (d,h) = [m1-m2-stj\ + — í„1+„2_2(l - a/2); m1 - m2 - s*W \ + — í„1+„2_2(a/2) 1 " nl n2 \ nl n2 ((137.1852 - 134.1468) - 4.7820109a/-|- + -^323(0.975); (137.1852 - 134.1468) - 4.7820109\/ V 215 109 V (3.038396 - 4.7820109 x 0.1174902 x 1.967336; 3.038396 - 4.7820109 x 0.1174902 x (-1.967336)) (1.933069; 4.143723) = (1.9331; 4.1437) 110 dh <- (ml - m2) - sh * sqrt(l / nl + 1 / n2) * qt(1 - alpha / 2, nl + n2 - 2) # 1.93307 111 hh <- (ml - m2) - sh * sqrt(l / nl + 1 / n2) * qt(alpha / 2, nl + n2 - 2) # 4.143723 • Závěr testování Protože jj,q = 0 nenáleží do Waldova 95% empirického oboustranného intervalu spolehlivosti, tj. jj,q = 0 ^ IS, Hq zamítáme na hladině významnosti a = 0.05. 5. Testování p-hodnotou • p-hodnota p-hodnota = 2min{Pr(Tw < tw), Pr(Tw > tw)} = 2min{Pr(TvK < 5.40794), 1 - Pľ(Tw < 5.40794)} = 2 min{0.9999999, 6.212613 x 10"8} = 2 x (6.212613 x 10"8) = 1.242523 x 10"7 = 1.2425 x 10"7 112 p.val <- 2*min(pt(tw, nl + n2 - 2) , 1 - pt (tw, nl + n2 - 2)) # 1. 242523e-07 • Závěr testování Protože p-hodnota = 1.2425 x 10-7 je menší než a = 0.05, Hq zamítáme na hladině významnosti a = 0.05. 6. Grafická vizualizace výsledků testování Vhodný grafem porovnávajícím střední hodnoty obou náhodných výběrů je krabicový diagram (obrázek 16). 7. Interpretace výsledků: Na základě všech tří způsobů testování zamítáme hypotézu o shodě středních hodnot největší šířky mozkovny mužů a žen. Mezi největší šířkou mozkovny mužů a žen existuje statisticky významný rozdíl. 1 1 215 + 109 23 s s o > o M N O S 155 - 150 - 145 - 140 - 135 130 - 125 120 - prumer median muzi zeny pohlaví Obrázek 16: Krabicový diagram největší šířky mozkovny mužů a žen 113 t.test(skull.Bm, skull.Bf, conf.level = 0.95, alternative = 'two.sided', var.equal = T) Poznámka: K otestování nulové hypotézy o rozdílu středních hodnot můžeme využít funkci t.test(). Vstupními parametry budou nejprve dva vektory reprezentující náhodné výběry, tj. skuli.Bm a skuli.Bf, dále hodnota hladiny významnosti a zadaná prostřednictvím koeficientu spolehlivosti 1 — a nastavením hodnoty argumentu conf.level = 0.95, typ zvolené alternativní hypotézy (oboustranná) zadaný pomocí argumentu alternativě = 'two.sideď a nakonec argument var.equal == T volbu klasického í-testu založeného na předpokladu shody rozptylů a'f a M2 Mi - M2 > Mo, kde Mo = 0 (pravostranná alternativa) Volba hladiny významnosti • Hladinu významnosti volíme ze zadání a = 0.05. Testování kritickým oborem • Aritmetický průměr výběrových rozptylů S* (ni - 1)S? + (n2 - 1)S22 ni + n2 (18 - 1) x 4.5632812 + (69 - 1) x 4.9524592 18 + 69-2 17 x 20.82353 + 68 x 24.52685 85 ^^ = 721^8619 85 4.877109 = 4.8771 135 sl <- sd(upface.Hc) 136 s2 <- sd (upf ace . Hm) 137 sh <- sqrt(((nl - 1) * sl"2 + (n2 - 1) * s2"2) / (nl + n2 - 2)) # 4.877109 • Testovací statistika (Mi - M2) - mo J-w —-/ (72 - 70.13043) - 4.877109,/jig + -1 69 1.86957 4.877109 x 0.2646664 1.86957 1.290807 1.44837 = 1.4484 138 alpha <- 0.05 139 muO <- 0 140 ml <- mean(upface.Hc) 141 m2 <- mean ( upf ace . Hm) 142 tw <- ((ml - m2) - muO) / (sh * sqrt(l / nl + 1 / n2)) # 1.448369 • Kritický obor W = (í„1+„2_2(l - a); 00) = (ŕ18+69_2(l - 0.05); 00) = (í85(0.95); 00) = (1.662978; 00) 27 143 qt(l - alpha, nl + n2 - 2) # 1.662978 • Závěr testování Protože realizace testovací statistiky tyy = 1.44837 nenáleží do kritického oboru, tj. tyy £ W, Hq nezamítáme na hladině významnosti a = 0.05. 4. Testování intervalem spolehlivosti • Interval spolehlivosti (eř,oo) = \ mi-m2- stJ \ + — í„i+„2_2(l - a); oo nl ri2 (72 - 70.13043) - 4.877109^ + ^í85(0.95); oo (1.869565 - 4.877109 x 0.2646664 x 1.662978 ; oo) (-0.277018; oo) = (-0.2770; oo) 144 DH <- (ml - m2) - sh * sqrt(l / nl + 1 / n2) * qt(1 - alpha, nl + n2 - 2) # -0.2770188 • Závěr testování Protože jj,0 = 0 náleží do Waldova 95% empirického levostranného intervalu spolehlivosti, tj. p0 = 0 G IS, Hq nezamítáme na hladině významnosti a = 0.05. Testování p-hodnotou • p-hodnota p-hodnota = Pr(Tw >tw) = l- V?{TW < tw) = 1 - Pr(Tw < 1.44837) = 0.9999 145 p.val <- 1 - pt(tw, nl + n2 -2) # 0.07559631 • Závěr testování Protože p-hodnota = 0.07556 je větší než a = 0.05, Hq nezamítáme na hladině významnosti a = 0.05. 6. Grafická vizualizace výsledků testování Vhodný grafem porovnávajícím střední hodnoty obou náhodných výběrů je krabicový diagram (viz obrázek 19). 7. Interpretace výsledků: Na základě všech tří způsobů testování nezamítáme hypotézu o tom, že střední hodnota výšky horní části tváře čínské populace je menší nebo rovna střední hodnotě výšky horní části tváře malajské populace. Testování tedy neprokázalo, že by výška horní části tváře čínské populace byla statisticky významně vyšší než u malajské populace. Poznámka: K otestování nulové hypotézy o rozdílu středních hodnot můžeme využít funkci t.test() s argumenty conf.level = 0.95 pro volbu hladiny významnosti a = 0.05, alternativě = 'greater' pro volbu pravostranné alternativní hypotézy a var.equal = T pro volbu klasického í-testu. 28 • prumer - median O O o 1-1- cinska malajska populace Obrázek 19: Krabicový diagram výšky horní části tváře mužů čínské a malajské populace 146 t.test(upface.He, upface.Hm, conf.level = 0.95, alternative = 'greater', var.equal = T) o > 90 H 85 80 75 70 65 H 60 55 H Two Sample t-test data: upface.He and upface.Hm t = 1.4484, df = 85, p-value = 0.0756 alternative hypothesis: true difference in means is greater than 0 95 percent confidence interval: -0.2770188 Inf sample estimates: mean of x mean of y 72.00000 70.13043 29 Příklad 10.8. Test o rozdílu středních hodnot Mějme datový soubor 19-more-samples-correlations-skull.txt a proměnnou nose.B popisující šířku nosu v mm (viz sekce ??). Předpokládejme, že náhodná veličina X, popisující šířku nosu nemecké populace, pochází z normálního rozdělení N(/j,\,&'(), a že náhodná veličina Y, popisující šířku nosu bantuské populace, pochází z normálního rozdělení N(p2, )• Na hladině významnosti a = 0.01 zjistěte, zda šířka nosu německé populace je statisticky významně menší než šířka nosu bantuské populace. Řešení příkladu 10.8 Příkazem read.delim() načteme datový soubor a příkazem na.omit() odstraníme NA hodnoty. Pomocí operátoru [] vybereme z datové tabulky údaje o šířce nosu (nose.B) německé (pop == 'nem'), resp. bantuské (pop == 'ban') populace. 158 data <- read.delim('00-Data//19-more - samples-correlations-skull.txt') 159 data <- na.omit(data) 160 #head(data, n = 3) 161 nose.Bn <- data[data$pop == 'nem', ' nose.B'] 162 nose.Bb <- data [data$pop == 'ban', ' nose.B'] 163 164 nl <- length(nose.Bn) # 20 165 n2 <- length(nose.Bb) # 14 Datový soubor obsahuje naměřené hodnoty šířek nosu 20 jedinců německé populace a 14 jedinců bantuské populace. Nejprve otestujeme předpoklad normality obou náhodných výběrů. Na hladině významnosti a testujeme hypotézu Hqi: Náhodný výběr šířky nosu německé populace pochází z normálního rozdělení, oproti alternativní hypotéze Hn: Náhodný výběr šířky nosu německé populace nepochází z normálního rozdělení. K otestování předpokladu normality použijeme Shapiro-Wilkův test, kvantilový diagram a histogram. Datový soubor rozdělíme do pěti ekvidistatních intervalů s šířkou 3 mm prostřednictvím stanovených hranic 17, 20,..., 32 (viz obrázek 20). 18.5 21.5 24.5 27.5 30.5 -2 -1 0 1 2 teoreticky kvantil sirka nosu (v mm) Shapiro-Wilkuv test: p-hodnota =0.5899 Obrázek 20: Histogram a kvantilový diagram šířky nosu mužů německé populace Protože p-hodnota = 0.5899 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr šířky nosu německé populace pochází z normálního rozdělení. Vykreslený histogram dostatečným způsobem kopíruje tvar křivky normálního rozdělení. Splnění předpokladu normality je tentokrát viditelné na kvan-tilovém diagramu, kde se všechny body pohybují okolo referenční přímky. Dále testujeme hypotézu H02: Náhodný výběr šířky nosu bantuské populace pochází z normálního rozdělení, oproti alternativní hypotéze H12: Náhodný výběr šířky nosu bantuské populace nepochází z normálního rozdělení. Předpoklad normality otestujeme Shapiro-Wilkovým testem, kvantilovým diagramem a histogramem. Datový soubor rozdělíme do pěti ekvidistatních intervalů s šířkou 2 mm prostřednictvím stanovených hranic 21, 23 ..., 31 (viz obrázek 21). 30 Obrázek 21: Histogram a kvantilový diagram šířky nosu mužů bantuské populace Protože p-hodnota = 0.1511 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr šířky nosu bantuské populace pochází z normálního rozdělení. Při pohledu na histogram bychom o normálním rozdělení naměřených hodnot mohli pochybovat. Nesmíme však zapomínat, že pro takto nízký počet pozorování může i nepříliš příznivě vypadající histogram reflektovat normálně rozdělená data. Kvantilový diagram naproti tomu normlání charakter dat popisuje uspokojivě. Zbývá ověřit potenciální shodu rozptylů obou náhodných výběrů. Na hladině významnosti a = 0.05 testujeme hypotézu Hq : o\ja\ = 1 oproti alternativní hypotéze H\ : o\ja\ ^ 1. Nulovou hypotézu testujeme pomocí funkce var.testQ na základě p-hodnoty (viz příklad 10.2). 166 var.test(nose.Bn, nose.Bb, alternativě = 'two.sided')$p.val [1] 0.2877778 167 Jelikož p-hodnota = 0.2878 je větší než 0.05, nulovou hypotézu o shodě rozptylů nezamítáme na hladině významnosti a = 0.05. Na základě výsledků testování tedy předpokládáme, že rozptyly obou výběrů jsou shodné. Protože oba náhodné výběry pochází z normálních rozdělení, o jejichž rozptylech předpokládáme, že jsou shodné, použijeme na otestování hypotézy ze zadání klasický dvouvýběrový í-test o rozdílu středních hodnot pi — P2-V zadání není uvedeno znění nulové hypotézy. Zadaným úkolem je prověřit, zda je šířka nosu německé populace statisticky významně menší, než šířka nosu malajské populace. Toto tvrzení je zněním alternativní hypotézy. Nulová hypotéza tvoří potom doplněk tohoto tvrzení. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Střední hodnota šířky nosu německé populace je větší nebo rovna střední hodnotě šířky nosu bantuské populace. H\ : Střední hodnota šířky nosu německé populace je menší než střední hodnota šířky nosu bantuské populace. • matematická formulace nulové a alternativní hypotézy Ho ■ Pi > P2 Pi - P2 > Po, kde p0 = 0 Hi : pi < P2 pi - P2 < Po, kde p0 = 0 (levostranná alternativa) 2. Volba hladiny významnosti • Hladina významnosti a = 0.01 (viz zadání příkladu). 31 3. Testování kritickým oborem • Aritmetický průměr výběrových rozptylů (ni - l)S'f + (n2 - 1)S22 ni + n2 - 2 (20 - 1) x 2.4192212 + (14 - 1) x 3.1492192 20 + 14-2 19 x 5.85263 + 13 x 9.91758 32 240.1285 32 2.739346 = 2.7393 ^7.504016 168 si <- sd(nose.Bn) 169 s2 <- sd(nose.Bb) 170 sh <- sqrt(((nl - 1) * sl"2 + (n2 - 1) * s2"2) / (nl + n2 - 2)) # 2.739346 • Testovací statistika (Mi - m2) - /íq J-w —-, V 711 712 _ (24.2 - 26.92857) 2.739346^ + ^ -2.72857 2.739346 x 0.348466 -2.72857 0.9545689 -2.85843 = -2.8584 171 alpha <- 0.01 172 muO <- 0 173 ml <- mean(nose.Bn) 174 m2 <- mean(nose.Bb) 175 tw <- ((ml - m2) - muO) / (sh * sqrt(l / nl + 1 / n2)) # -2.858432 • Kritický obor W = (-oo; i„1+„2_2(a)) = (-oo; Í20+i4-2(0.01)> = (-oo; i32(0.01)> = (-oo; -2.448678) 176 qt(alpha, nl + n2 - 2) # -2.448678 • Závěr testování Protože realizace testovací statistiky tyy = —2.8584 náleží do kritického oboru, tj. tyy e W, Hq zamítáme na hladině významnosti a = 0.01. 32 4. Testování intervalem spolehlivosti • Interval spolehlivosti (d, oo) = ( -oo ; m1 - m2 - st\j -^ + ^-í„i+„2_2(l - a); oo -oo ; (24.2 - 26.92857) - 2.739346^/ ^ + -^í32(0.01) (-oo; -2.72857 - 2.739346 x 0.348466 x (-2.448678)) (-oo; -0.39113) = (-oo; -0.3911) 177 HH <- (ml - m2) - sh * sqrt(l / nl + 1 / n2) * qt(alpha, nl + n2 - 2) # -0.3911394 • Závěr testování Protože jj,q = 0 nenáleží do Waldova 99% empirického pravostranného intervalu spolehlivosti, tj. po = 0 ^ IS, Hq zamítáme na hladině významnosti a = 0.01. 5. Testování p-hodnotou • p-hodnota p-hodnota = ~Pr(Tyy < tyy) = Pr(Tw < -2.858432) = 0.003714569 = 0.0037 178 p.val <- pt(tw, nl + n2 -2) # 0. 003714569 • Závěr testování Protože p-hodnota = 0.0037 je menší než a = 0.01, Hq zamítáme na hladině významnosti a = 0.01. 6. Grafická vizualizace výsledků testování Výsledek testování vizualizujeme pomocí krabicového diagramu (viz obrázek 22). 7. Interpretace výsledků: Na základě všech tří způsobů testování zamítáme hypotézu o tom, že střední hodnota šířky nosu německé populace je větší nebo rovna střední hodnotě šířky nosu bantuské populace. Testování tedy prokázalo, že šířka nosu německé populace je statisticky významně nižší než šířka nosu bantuské populace. Poznámka: K otestování nulové hypotézy o rozdílu středních hodnot můžeme využít funkci t.test() s argumenty conf.level = 0.99 pro volbu hladiny významnosti a = 0.01, alternativě = 'less' pro volbu levostranné alternativní hypotézy a var.equal = T pro volbu klasického í-testu. 179 t.test(nose.Bn, nose.Bb, conf.level = 0.99, alternativě = 'less', var.equal = T) Two Sample t-test data: nose.Bn and nose.Bb t = -2.8584, df = 32, p-value = 0.003715 alternative hypothesis: true difference in means is less than 0 99 percent confidence interval: -Inf -0.3911394 sample estimates: mean of x mean of y 24.20000 26.92857 180 181 182 183 184 185 186 187 188 189 190 33 35 30 - 25 - 20 - prumer median -1-1- německa bantuska populace Obrázek 22: Krabicový diagram šířky nosu mužů německé a bantuské popul 34 Příklad 10.9. Test o rozdílu středních hodnot Mějme datový soubor 10-two-samples-means-birth.txt a proměnnou weight.W popisující porodní hmotnost novorozence v g (viz sekce ??). Na hladině významnosti a = 0.05 zjistěte, zda se porodní porodní hmotnost novorozence bez staršího sourozence a porodní hmotnost novorozence s jedním starším sourozencem liší. Řešení příkladu 10.9 Příkazem read.delim() načteme datový soubor a odstraníme chybějící hodnoty. Pomocí operátoru [] vybereme z datové tabulky údaje o porodní hmotnosti (weight.W) novorozenců bez staršího sourozence (o.sib.N == 0), resp. novorozenců s jedním starším sourozencem (o.sib.N == 1) populace. 191 data <- read.delim('00-Data//10-two-samples-means-birth.txt') 192 data <- na.omit(data) 193 #head (data, n = 3) 194 weight.WO <- data [data$o . sib . N == 0, 'birth.W] 195 weight.Wl <- data [data$o . sib . N == 1, 'birth.W] 196 197 nl <- length(weight.W0) # 297 198 n2 <- length(weight.Wl) # 276 Datový soubor obsahuje naměřené porodní hmotnosti 297 novorozenců bez staršího sourozence a 276 novorozenců s jedním starším sourozencem. Nejprve otestujeme předpoklad normality obou náhodných výběrů (a = 0.05). Protože rozsah náhodného výběru porodních hmotností novorozenců bez staršího sourozence je větší než 30, použijeme na ověření předpokladu normality náhodného výběru Lillieforsův test v kombinaci s kvantilovým diagramem a histogramem (viz obrázek 23). Datový soubor rozdělíme do devíti ekvidistatních intervalů s šířkou 453 mm prostřednictvím stanovených hranic 886,1339,..., 4064 g. Obrázek 23: Histogram a kvantilový diagram porodní hmotnosti novorozenců bez staršího sourozence Protože p-hodnota = 0.0014 je menší než 0.05, nulovou hypotézu zamítáme na hladině významnosti a = 0.05. Vykreslený histogram je oproti křivce hustoty normálního rozdělení posunutý doprava. Taktéž kvantilový diagram ukazuje odlehlost vykreslených bodů od referenční přímky. Náhodný výběr porodní hmotnost novorozenců bez staršího sourozence nepochází z normálního rozdělení. Dále testujeme hypotézu o normalitě náhodného výběru porodních hmotností novorozenců s jedním starším sourozencem. Předpoklad normality otestujeme opět Lillieforsovým testem, kvantilovým diagramem a histogramem. Datový soubor rozdělíme do devíti ekvidistatních intervalů s šířkou 443 mm prostřednictvím stanovených hranic 986,1429,..., 4973 g (viz obrázek 24). Error in quantile(y, probs , names = FALSE, type = qtype, na.rm = TRUE): object Jweigth.WlJ not f ound 199 35 Obrázek 24: Histogram a kvantilový diagram porodní hmotnosti novorozenců s jedním starším sourozencem Protože p-hodnota = 1 x 10-6 je menší než 0.05, nulovou hypotézu zamítáme na hladině významnosti a = 0.05. Při pohledu na histogram vidíme jednak posunutí naměřených hodnot oproti křivce normálního rozdělení směrem doprava a větší špičatost. Kvantilový digram ukazuje odlehlost hodnot od referenční křivky. Náhodný výběr porodních hmotností novorozenců s jedním starším sourozencem nepochází z normálního rozdělení. Jelikož ani jeden z náhodných výběrů nepochází z normálního rozdělení, nemůžeme na hypotézu ze zadání použít parametrický test o rozdílu středních hodnot pí — Hypotézu ze zadání bychom testovali některou z metod uvedených v sekci ??. 36 10.3 Welchův dvouvýběrový ŕ-test o rozdílu středních hodnot \i\ — /i2 Nechť Xn,..., X\ni je náhodný výběr z N(p,i, a'f), a X21, ■ ■ ■, X'ín% je na něm nezávislý náhodný výběr z rozdělení -/V(p2, o"2), přičemž ri\ > 2, n2 > 2 a a\ a a\ jsou neznámé různé rozptyly, tj. a\ 7^ a\. Nechť po je konstanta. Na hladině významnosti a testujeme jednu z následujících tří hypotéz oproti příslušné alternativní hypotéze. Hqi : pi — p2 = Mo oproti Tin : Mi — M2 7^ Mo (oboustranná alt.) #02 : Mi — M2 < Mo oproti 7ři2 : Mi — M2 > Mo (pravostranná alt.) 7řo3 : Mi — M2 > Mo oproti .H13 : Mi — M2 < Mo (levostranná alt.) Test nazýváme Welchův dvouvýběrový í-test o rozdílu středních hodnot pi — p2. Testovací statistika má tvar TW = ^-M^-^. (10.3) kde Mi je výběrový průměr a ri\ je rozsah prvního náhodného výběru, M2 je výběrový průměr a n2 je rozsah druhého náhodného výběru, S* je aritmetický průměr výběrových rozptylů S'f a S% (viz kapitola ??) a /uq je konstanta z nulové hypotézy. Za platnosti nulové hypotézy pochází statistika Tyy ze Studentova rozdělení o df stupních volnosti, tj- (Mi - M2) - p0 Ho Jw =-— ~ *#> Jl 1 ^2. ni 712 kde (ěL + s|\2 \ ni \ 1 J V 2 / 711 —1 712 — 1 je Welchova aproximace stupňů volnosti Studentova rozdělení navržená tak, aby zohledňovala různost rozptylů af a a2. Kritický obor podle zvolené alternativní hypotézy má tvar 2 11 H H13 M ^ Mo W = (-00 ; t#(a/2)) U (íd/(l - a/2); 00) M > Mo = (íd/(l - a); 00) M < Mo W = (-00; tdf (a)) kde tdf(a/2), í Mo (rf, 00) = \ nii - m2 - \J^ + ^tdf(l - a); 00 #13 : Mi - M2 < Mo (-00, /i) = ^-00 ; mi - m2 - \[^+^Uf p-hodnota má v závislosti na zvolené alternativní hypotéze jeden z následujících tvarů 11 H- H13 M 7^ Mo p-hodnota = 2min{Pr(7V < , Pr(?V > tw)} M > Mo p-hodnota = Pr(7V > tw) = 1 — Pr(?V < p < po p-hodnota = Pr(2V < íw) kde Tvk je náhodná veličina, tyy je realizace testovací statistiky Tyy (viz vzorec 10.3), tedy konkrétní číslo, a Pr(7V < tw) je distribuční funkce Studentova rozdělení o df stupních volnosti, jejíž hodnotu získáme pomocí a implementované funkce pt(). 37 Příklad 10.10. Test o rozdílu středních hodnot Mějme datový soubor 16-anova-head.txt a proměnnou bigo.W popisující šířku dolní čelisti v mm (viz sekce ??). Předpokládejme, že náhodná veličina X, popisující šířku dolní čelisti u žen orientovaných heterosexuálně, pochází z normálního rozdělení N({ii,a'(), a že náhodná veličina Y, popisující šířku dolní čelisti žen orientovaných jinak než výhradně heterosexuálně, pochází z normálního rozdělení N(p2, Na hladině významnosti a = 0.05 zjistěte, zda šířka dolní čelisti žen orientovaných heterosexuálně je statisticky významně menší než šířka dolní čelisti žen orientovaných jinak než výhradně heterosexuálně. Řešení příkladu 10.10 Příkazem read.delim() načteme datový soubor a příkazem na.omit() odstraníme NA hodnoty. Pomocí operátoru [] vybereme z datové tabulky údaje o šířce dolní čelisti (bigo.W) u žen (sex == 'f') orientovaných heterosexuálně (sexor == 'op'), resp. orientovaných jinak než heterosexuálně (sexor == 'sa') . 200 data <- read.delim('00-Data//16-anova-head.txt') 201 data <- na.omit(data) 202 #head (data, n = 3) 203 bigo.Wfo <- data [data$sex == 'f' k data$sexor =■■ 204 bigo.Wfs <- data [data$sex == 'f' k data$sexor =■■ 205 206 nl <- length(bigo.Wfo) # 62 207 n2 <- length(bigo.Wfs) # 23 Datový soubor obsahuje hodnoty šířek dolní čelisti u 62 žen orientovaných heterosexuálně a 23 žen orientovaných jinak než heterosexuálně. Řešení příkladu začneme testováním předpokladu normality obou náhodných výběrů. Na hladině významnosti a testujeme hypotézu Hqi: Náhodný výběr šířky dolní čelisti u žen orientovaných heterosexuálně pochází z normálního rozdělení, oproti alternativní hypotéze Hn: Náhodný výběr šířky dolní čelisti u žen orientovaných heterosexuálně nepochází z normálního rozdělení. K otestování normality náhodného výběru použijeme Lillieforsův test, kvantilový diagram a histogram (viz obrázek 25). Datový soubor rozdělíme na základě Sturgesova pravidla do sedmi třídicích intervalů s ekvidistantní šířkou 3 mm prostřednictvím stanovených hranic 90, 93,..., 111. 91.5 97.5 103.5 109.5 -2 -1 0 1 2 teoreticky kvantil sirka dolni čelisti (v mm) Lillieforsův test: p-hodnota =0.4252 Obrázek 25: Histogram a kvantilový diagram šířky dolní čelisti heterosexuálně orientovaných žen Protože p-hodnota = 0.4252 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr šířky dolní čelisti heterosexuálně orientovaných žen pochází z normálního rozdělení. Histogram, stejně jako kvantilový diagram, předpoklad normality též potvrzují. Dále testujeme hypotézu Hq2'- Náhodný výběr šířky dolní čelisti žen orientovaných jinak než heterosexuálně pochází z normálního rozdělení, oproti alternativní hypotéze H12: Náhodný výběr šířky dolní čelisti žen orientovaných jinak než heterosexuálně nepochází z normálního rozdělení. Předpoklad normality otestujeme Shapiro-Wilkovým testem, 'op', 'bigo.W'] 'sa', 'bigo .W] 38 kvantilovým diagramem a histogramem. Datový soubor rozdělíme do pěti ekvidistatních intervalů s šířkou 4 mm prostřednictvím stanovených hranic 91, 95 ..., 111 (viz obrázek 26). 0.00 H i—i—i—i—i—r 101 105 109 sirka dolni čelisti (v mm) > > a -O > 110 - 105 - 100 - 95 ii r -2-10 1 2 teoreticky kvantil Shapiro-Wilkuv test: p-hodnota =0.1672 Obrázek 26: Histogram a kvantilový diagram šířky dolní čelisti žen orientovaných jinak než heterosexuálně Protože p-hodnota = 0.1672 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr šířky dolní čelisti žen orientovaných jinak než heterosexuálně pochází z normálního rozdělení. Zbývá ověřit předpoklad shody rozptylů obou náhodných výběrů. Na hladině významnosti a = 0.05 testujeme nulovu hypotézu Hq : o"2/o"2 = 1 oproti hypotéze H\ : cr'i/a'2 7^ 1. Nulovou hypotézu testujeme pomocí funkce var.test() na základě p-hodnoty (viz příklad 10.2). 208 var.test(bigo.Wfo, bigo.Wfs, alternativě = 'two.sided')$p.val [1] 0.0336774 209 Jelikož p-hodnota = 0.03368 je menší než 0.05, nulovou hypotézu o shodě rozptylů zamítáme na hladině významnosti a = 0.05. Rozptyly obou náhodných výběrů se statisticky významně liší. Oba náhodné výběry pochází z normálního rozdělení. Na otestování hypotézy ze zadání tedy můžeme použít parametrický test. Protože však není splněn předpoklad shody rozptylů obou náhodných výběrů, použijeme na otestování hypotézy ze zadání namísto klasického dvou výběrového í-tesu Welchův dvou výběrový í-test. Zadaným úkolem je prověřit, zda je šířka dolní čelisti heterosexuálně orientovaných žen statisticky významně menší, než šířka dolní čelisti žen orientovaných jinak než heterosexuálně. Toto tvrzení je zněním alternativní hypotézy a nulová hypotéza tvoří doplněk tohoto tvrzení. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Střední hodnota šířky dolní čelisti heterosexuálně orientovaných žen je větší nebo rovna střední hodnotě šířky dolní čelisti žen orientovaných jinak než heterosexuálně. H\ : Střední hodnota šířky dolní čelisti heterosexuálně orientovaných žen je menší než střední hodnota šířky dolní čelisti žen orientovaných jinak než heterosexuálně. • matematická formulace nulové a alternativní hypotézy Ho ■ Pi > P2 Pi - P2 > Po, kde p0 = 0 Hi : pi < P2 pi - P2 < Po, kde p0 = 0 (levostranná alternativa) 2. Volba hladiny významnosti 39 • Hladina významnosti a = 0.05 (viz zadání príkladu). 3. Testování kritickým oborem • Testovací statistika (Mi - M2) - po w tl\ ri2 (99.5 - 103.08695) - 0 3.7667662 , 5.3420552 62 23 3.58696 ^1.46961 -3.58696 1.212275 -2.958867 = -2.9589 210 alpha <- 0.05 211 muO <- 0 212 ml <- mean(bigo.Wfo) 213 m2 <- mean(bigo.Wfs) 214 tw <- ((ml - m2) - muO) / (sqrt(sl " 2 / nl + s2 " 2 / n2)) # -2.958864 • Welchova aproximace počtu stupňů volnosti df + n\—\ 3.7667662 , 5.3420552 X 2 62 23 / 3.7667662 ^ / 5.3420552 ^ z \ 62 ) I V 23 ) 62-1 23-1 (0.2288473 + 1.240763)2 0.22884732 , 1.2407632 61 22 2.159754 0.07083549 30.48972 215 df <- (sl " 2 / nl + s2 " 2 / n2)"2 / ((sl " 2 / nl) " 2 / (nl - 1) + (s2 " 2 / n2) " 2 / (n2 - 1)) # 30.48972 • Kritický obor W=(-oo;tdf(a)) = (-OO; Í30.48972(0.05)) = (-oo; -2.448678) 216 qt(alpha, df) # -1.696393 40 • Závěr testování Protože realizace testovací statistiky tyy = —2.9589 náleží do kritického oboru, tj. tyy e W, Hq zamítáme na hladině významnosti a = 0.05. 4. Testování intervalem spolehlivosti • Interval spolehlivosti (d, oo) = —oo ; m\ — m2 — \--1--trjf (1 — a); oo ni říz (-oo ; (99.5 - 103.08695) - f-^^- + ^§^30.48972(0.05) (-oo; -3.586957 - 1.212275 x (-1.696393)) (-oo; -1.530462) = (-oo; -1.5305) 217 HH <- (ml - m2) - sqrt(sl"2 / nl + s2"2 / n2) * qt(alpha, df) # -1.530462 • Závěr testování Protože jj,q = 0 nenáleží do Waldova 95% empirického pravostranného intervalu spolehlivosti, tj. po = 0 IS, Hq zamítáme na hladině významnosti a = 0.05. 5. Testování p-hodnotou • p-hodnota p-hodnota = Pr(7V < tw) = Pr(Tw < -2.958864) = 0.002961084 = 0.002961 218 p.val <- pt(tw, df) # 0. 002961084 • Závěr testování Protože p-hodnota = 0.002961 je menší než a = 0.05, H0 zamítáme na hladině významnosti a = 0.01. 6. Grafická vizualizace výsledků testování Výsledek testování vizualizujeme pomocí krabicového diagramu (viz obrázek 27). 7. Interpretace výsledků: Na základě všech tří způsobů testování zamítáme hypotézu o tom, že střední hodnota šířky dolní čelisti heterosexuálně orientovaných žen je větší nebo rovna střední hodnotě šířky dolní čelisti žen orientovaných jinak než heterosexuálně. Šířka dolní čelisti žen orientovaných heterosexuálně je statisticky významně menší než šířka dolní čelisti žen orientovaných jinak než heterosexuálně. Poznámka: K otestování nulové hypotézy o rozdílu středních hodnot můžeme využít funkci t.test(). Vstupními parametry budou nejprve dva vektory reprezentující náhodné výběry, tj. bigo.Wfo a bigo.Wfs, dále hodnota hladiny významnosti a zadaná prostřednictvím koeficientu spolehlivosti 1 — a nastavením hodnoty argumentu conf.level = 0.95, typ zvolené alternativní hypotézy (levostranná) zadaný pomocí argumentu alternativě == 'less' a nakonec argument var.equal == F pro volbu Welchova dvouvýběrového í-testu založeného na předpokladu rozdílu mezi rozptyly a\ a a\. 41 115 - 110 - prumer median 'S 105 - "o •a 100 - 95 - 90 - -1-1- heterosexuální výhradne j ina orientace Obrázek 27: Krabicový diagram šířky dolní čelisti žen orientovaných heterosexuálně a žen orientovaných jinak než heterosexuálně 219 t.test(bigo.Wfo, bigo.Wfs, conf.level = 0.95, alternative = 'less', var.equal = F) Welch Two Sample t-test data: bigo.Wfo and bigo.Wfs t = -2.9589, df = 30.49, p-value = = 0.002961 alternative hypothesis: true diff srence in means is less than 0 95 percent confidence interval: -Inf -1.530462 sample estimates: mean of x mean of y 99.500 103.087 220 221 222 223 224 225 226 227 228 229 230 Součástí výstupu je hodnota testovací statistiky t = -2.9589, Welchova aproximace počtu stupňů volnosti Studentova rozdělení num df = 30.49, hranice intervalu spolehlivosti -Inf a -1.530462 a p-hodnota p-value = 0.002961. Jediné, co musíme stanovit zvlášť, je dolní hranice kritického oboru. 42 Příklad 10.11. Test o rozdílu středních hodnot Mějme datový soubor 13-two-samples-correlations-trunk.txt a proměnnou tru.L popisující délku trupu v mm (viz sekce ??). Předpokládejme, že náhodná veličina X, popisující délku trupu mužů, pochází z normálního rozdělení N({ii,a'(), a že náhodná veličina Y, popisující délku trupu žen, pochází z normálního rozdělení JV(/lí2, c^). Na hladině významnosti a = 0.01 zjistěte, zda existuje statisticky významný rozdíl mezi délkou trupu mužů a žen. Řešení příkladu 10.11 Příkazem read.delim() načteme datový soubor a příkazem na.omit() odstraníme případné NA hodnoty. Pomocí operátoru [] vybereme z datové tabulky údaje o délce trupu (tru.L) mužů (sex == 'm'), resp. žen (sex == 'f) a zjistíme rozsahy obou náhodných výběrů. 231 data <- read.delim('00-Data//13-two-samples-correlations-trunk.txt') 232 data <- na.omit(data) 233 #head (data, n = 3) 234 tru.Lm <- data[data$sex == 'm', 'tru.L'] 235 tru.Lf <- data[data$sex == 'f', 'tru.L'] 236 237 nl <- length(tru.Lm) # 75 238 n2 <- length(tru.Lf) # 100 Datový soubor obsahuje údaje o délkách trupu 75 mužů a 100 žen. Řešení příkladu vede na test o rozdílu středních hodnot. Nejprve otestujeme normalitu obou náhodných výběrů. Hladinu významnosti pro testy normality zvolíme a = 0.05. Nejprve testujeme hypotézu Hqi: Náhodný výběr délky trupu mužů pochází z normálního rozdělení, oproti alternativní hypotéze Hn: Náhodný výběr délky trupu mužů nepochází z normálního rozdělení. Normalitu otestujeme Lillieforsovým testem, kvantilovým diagramem a histogramem (viz obrázek 28). Naměřené hodnoty rozdělíme do sedmi ekvidistatních intervalů s šířkou 22mmm prostřednictvím stanovených hranic 397, 419,..., 154. Obrázek 28: Histogram a kvantilový diagram délky trupu mužů Protože p-hodnota = 0.8055 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr naměřených délek trupu mužů pochází z normálního rozdělení. Stejně tak je normalita zřetelná z histogramu, kde data kopírují tvar křivky hustoty a kvantilového diagramu, kde se vykreslené body zdržují okolo referenční přímky. Analogicky testujeme hypotézu Hqi: Náhodný výběr délky trupu žen pochází z normálního rozdělení, oproti alternativní hypotéze Hi2'. Náhodný výběr délky trupu žen nepochází z normálního rozdělení. Normalitu otestujeme Lillieforsovým testem. V rámci histogramu rozdělíme soubor do osmi ekvidistatních intervalů s šířkou 22 mm prostřednictvím stanovených hranic 228, 250 ..., 495 (viz obrázek 29). 43 330 374 418 462 délka trupu (v mm) -2-10 1 2 teoreticky kvantil Lillieforsuv test: p-hodnota =0.5386 Obrázek 29: Histogram a kvantilový diagram délky trupu žen Protože p-hodnota = 0.5386 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběrdélky trupu žen pochází z normálního rozdělení. Vykreslený histogram i kvantilový diagram podporují výsledek testování. Dále ověříme, zda rozptyly obou výběrů jsou shodné. Na hladině významnosti a = 0.05 testujeme hypotézu Hq : o\ja\ = 1 oproti alternativní hypotéze H\ : o\ja\ ^ 1. Nulovou hypotézu otestujeme pomocí p-hodnoty získané příkazem var.test() (viz příklad 10.2). var.test(tru.Lm, tru.Lf, alternativě = ' two.sided')$p.val # 0.0^256397 [1] 0.04256397 Jelikož p-hodnota = 0.04256 je menší než 0.05, nulovou hypotézu o shodě rozptylů zamítáme na hladině významnosti a = 0.05. Na základě výsledků testování tedy předpokládáme, že oba rozptyly jsou shodné. Protože oba náhodné výběry pochází z normálního rozdělení, použijeme na otestování hypotézy ze zadání parametrický test. Jelikož není splněn předpoklad shody rozptylů obou náhodných výběrů, zvolíme na test hypotézy ze zadání Welchův dvouvýběrový í-test. Zadaným úkolem je zjistit, zda existuje statisticky významný rozdíl mezi délkou trupu mužů a žen. Toto tvrzení je zněním alternativní hypotézy, neboť v nulové hypotéze nikdy nepředpokládáme nerovnost. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Střední hodnota délky trupu mužů a střední hodnota délky trupu žen jsou shodné. H\ : Střední hodnota délky trupu mužů a střední hodnota délky trupu žen nejsou shodné. • matematická formulace nulové a alternativní hypotézy Ho ■ Pi = P2 Pi - P2 = Po, kde p0 = 0 Hi : pi ^ P2 pi - P2 7^ Mo, kde p0 = 0 (oboustranná alternativa) 2. Volba hladiny významnosti • Hladinu významnosti volíme ze zadání a = 0.01. 3. Testování kritickým oborem 44 • Testovací statistika Tw (Mi - M2) - fi0 Ěl + Ěl ni 712 (457.1867 - 423.17)-0 27.155552 , 34.022882 75 "r" 100 34.0167 ^21.40788 34.0167 4.626865 7.35199 = 7.35199 241 alpha <- 0.01 242 muO <■ - 0 243 ml <- mean(tru.Lm) 244 m2 <- mean(tru.Lf) 245 sl <- sd(tru.Lm) 246 s2 <- sd(tru.Lf) 247 248 tw <- ((ml - m2) - muO) / (sqrt(sl ' ~ 2 , / nl + s2 " 2 / n2)) # 7.351989 • Welchova aproximace počtu stupňů volnosti df = $1 i $2_ ni 712 2 M i \ n2 + ni—1 T12 — 1 27.155552 , 34.022882 X 2 75 100 's^O^^SS^.^2 ( 34.022SS2 75 ) i l 100 100-1 2 (9.832319 + 11.57556) 9.8323192 i 11.575562-74 99 458.2973 2.659883 172.2998 249 df <- (sl " 2 / nl + s2 " 2 / n2)"2 / ((sl " 2 / nl) " 2 / (nl - 1) + (s2 " 2 / n2) / (n2 - 1)) # 172.2998 Kritický obor W = (-oo; tdf(a/2)) U (íd/(l - a/2); oo) = (-00 ; Í172.2998(0.005)) U (Í172.2998(0.995) ; oo) = (-oo; -2.604664) U (2.604664; oo) 250 qt(alpha/2, df) # -2.6O4664 251 qt(l-alpha/2, df) # 2.6O4664 45 • Závěr testování Protože realizace testovací statistiky tyy = 7.35199 náleží do kritického oboru, tj. tyy e W, Hq zamítáme na hladině významnosti a = 0.01. 4. Testování intervalem spolehlivosti • Interval spolehlivosti s2 s? / s2 s22 (d, h) = [m1-m2- \\^ + ^tdf^ ~ a/2)' mi ~ m'2 ~ ]l ~[ + ~tdf(a/'2) «457.1867 - 423.17) - ,/« + ^hlMm(9m); (457.1867 - 423.17) - + *f (34.0167 - 4.626865 x 2.604664; 34.0167 - 4.626865 x (-2.604664)) (21.96527; 46.06813) = (21.9653; 46.06813) 252 dh <- (ml - m2) - sqrt(sl"2 / nl + s2"2 / n2) * qt (1 - alpha / 2, nl + n2 - 2) # 21.96578 253 hh <- (ml - m2) - sqrt(sl"2 / nl + s2"2 / n2) * qt(alpha / 2, nl + n2 - 2) # 46.06755 • Závěr testování Protože jj,0 = 0 nenáleží do Waldova 99% empirického oboustranného intervalu spolehlivosti, tj. jj,0 0 ^ IS, Hq zamítáme na hladině významnosti a = 0.01. 5. Testování p-hodnotou • p-hodnota p-hodnota = 2min{Pr(!ZV < tw), Pr(7V > tw)} = 2min{Pr(!ZV < 7.35199), 1 - Pr(!ZV < 7.35199)} = 2 min{0.9999999, 3.774536 x 10"12} = 2 x (3.774536 x 10"12) = 7.549072 x 10"12 = 1.5490 x 10"12 254 p.val <- 2*min(pt(tw, df) , 1 - pt (tw, df)) # 7.549072e-12 • Závěr testování Protože p-hodnota = 1.5490 x 10-12 je menší než a = 0.01, Hq zamítáme na hladině významnosti a = 0.01. 6. Grafická vizualizace výsledků testování Rozdíl mezi středními hodnotami populací mužů a žen vizualizujeme pomocí krabicového diagramu (viz obrázek 30). 7. Interpretace výsledků: Na základě všech tří způsobů testování zamítáme hypotézu o shodě středních hodnot. Mezi délkou trupu mužů a žen existuje statisticky významný rozdíl. Poznámka: K otestování nulové hypotézy o rozdílu středních hodnot můžeme využít funkci t.test() s argumenty conf.level = 0.99 pro volbu hladiny významnosti a = 0.01, alternativě = 'two.sideď pro volbu oboustranné alternativní hypotézy a var.equal = F pro volbu Welchova dvouvýběrového í-testu. 46 s s •a 550 - 500 - 450 - 400 - 350 prumer median • o zeny pohlaví Obrázek 30: Krabicový diagram délky trupu mužů a žen 255 t.test(tru.Lm, tru.Lf, conf.level = 0.95, alternative = 'two.sided', var.equal = F) Welch Two Sample t-test data: tru.Lm and tru.Lf t = 7.352, df = 172.3, p-value = 7.549e-12 alternative hypothesis: true difference in means is not equal to 0 95 percent confidence interval: 24.88403 43.14930 sample estimates: mean of x mean of y 457.1867 423.1700 256 257 258 259 260 261 262 263 264 265 266 Součástí výstupu je hodnota testovací statistiky t = 7.352, Welchova aproximace počtu stupňů volnosti Studentova rozdělení num df = 172.3, hranice intervalu spolehlivosti 24.88403 a 43.14930 a p-hodnota p-value = 7.549e-12. Jediné, co musíme stanovit zvlášť, jsou dolní a horní hranice kritického oboru. 47 Příklad 10.12. Test o rozdílu středních hodnot Mějme datový soubor 19-more-samples-correlations-skull.txt a proměnnou nose.B popisující šířku nosu v mm (viz sekce ??). Předpokládejme, že náhodná veličina X, popisující šířku nosu bantuské populace, pochází z normálního rozdělení N(/j,\, a'f), a že náhodná veličina Y, popisující šířku nosu čínské populace, pochází z normálního rozdělení iV(/Lt2, 02). Na hladině významnosti a = 0.05 testujte nulovou hypotézu o tom, že šířka nosu čínské populace je menší nebo rovna šířce nosu bantuské populace. Řešení příkladu 10.12 Příkazem read.delim() načteme datový soubor a příkazem na.omit() z něj odtraníme NA hodnoty. Pomocí operátoru [] vybereme z datové tabulky naměřené hodnoty šířky nosu (nose.B) čínské populace (pop == 'cin'), resp. bantuské populace (pop == 'ban') a zjistíme rozsahy obou náhodných výběrů. 267 data <- read.delim('00-Data//19-more - samples-correlations-skull.txt') 268 data <- na.omit(data) 269 #head(data, n = 3) 270 nose.Be <- data[data$pop == 'cin', ' nose.B'] 271 nose.Bb <- data [data$pop == 'ban', ' nose.B'] 272 273 nl <- length(nose.Be) # 19 274 n2 <- length(nose.Bb) # 14 Datový soubor obsahuje naměřené hodnoty šířky nosu u 19 jedinců čínské populace a 14 jedinců bantusképopulace. Řešení příkladu vede na test o rozdílu středních hodnot. Nejprve otestujeme normalitu obou náhodných výběrů. Na hladině významnosti a testujeme hypotézu Hqi: Náhodný výběr naměřených šířek nosu čínské populace pochází z normálního rozdělení, oproti alternativní hypotéze Hn: Náhodný výběr naměřených šířek nosu čínské populace nepochází z normálního rozdělení. Vzhledem k rozsahu náhodného výběru použijeme k otestování normality Shapiro-Wilkův test, kvantilový graf a histogram (viz obrázek 31). Datový soubor rozdělíme do pěti ekvidistatních intervalů s šířkou 1 mm prostřednictvím stanovených hranic 23, 24..., 28. 23.5 24.5 25.5 26.5 27.5 -2 -1 0 1 2 teoreticky kvantil sirka nosu (v mm) Shapiro-Wilkuv test: p-hodnota =0.1173 Obrázek 31: Histogram a kvantilový diagram šířky nosu mužů čínské populace Protože p-hodnota = 0.1173 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr naměřených šířek nosu čínské populace pochází z normálního rozdělení. Vykreslený kvantilový graf výsledek testování podporuje, histogram víceméně také. Analogicky testujeme hypotézu H02: Náhodný výběr naměřených šířek nosu bantuské populace pochází z normálního rozdělení, oproti alternativní hypotéze H12: Náhodný výběr naměřených šířek nosu bantuské populace nepochází z normálního rozdělení. K otestování předpokladu normality použijeme Shapiro-Wilkův test, kvantilový diagram a 48 histogram. Datový soubor rozdělíme do pěti ekvidistatních intervalů s šířkou 2 mm prostřednictvím stanovených hranic 21, 23,..., 31 (viz obrázek 32). Obrázek 32: Histogram a kvantilový diagram šířky nosu mužů bantuské populace Protože p-hodnota = 0.1511 je větší než 0.05, nulovou hypotézu nezamítáme na hladině významnosti a = 0.05. Náhodný výběr šířek nosu bantuské populace pochází z normálního rozdělení. Ačkoli pohled na histogram nás o normalitě dat příliš nepřesvědčil, nemůžeme případnou normalitu naměřených hodnot navrhnout kvůli malému rozsahu náhodného výběru. Kvantilový graf normalitu naměřených hodnot podporuje. Proto se kloníme k závěru Shapiro-Wilkova testu a data považujeme za normálně rozdělená. Nyní zbývá rozhodnout o splnění nebo nesplnění předpokladu o shodě rozptylů obou náhodných výběrů. Na hladině významnosti a = 0.05 testujeme hypotézu Hq : o\ja\ = 1 (rozptyly a\ a a\ jsou shodné) oproti alternativní hypotéze H\ : o\ja\ ^ 1 (rozptyly a\ a a\ nejsou shodné). Nulovou hypotézu otestujeme pomocí funkce var.test() na základě p-hodnoty (viz příklad 10.2). 275 var.test(nose.Bc, nose.Bb, alternativě = 'two.sided')$p.val [1] 0.01257845 276 Jelikož p-hodnota = 0.01258 je menší než 0.05, nulovou hypotézu o shodě rozptylů zamítáme na hladině významnosti a = 0.05. Testování tedy prokázalo statisticky významný rozdíl v rozptylech obou náhodných výběrů. Náhodné výběry pochází z normálních rozdělení s vzájemně odlišnými rozptyly, proto na otestování hypotézy ze zadání použijeme Welchův dvouvýběrový í-test o rozdílu středních hodnot pi — p2. V zadání příkladu máme výslovně definované znění hypotézy Hq, zbývá nám tedy dodefinovat znění alternativní hypotézy H\. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Střední hodnota šířky nosu čínské populace je menší nebo rovna střední hodnotě šířky nosu bantuské populace. H\ : Střední hodnota šířky nosu čínské populace je větší než střední hodnota šířky nosu bantuské populace. • matematická formulace nulové a alternativní hypotézy Ho ■ Pi < P2 Pi - P2 < Po, kde p0 = 0 Hi : pi > P2 pi - P2 > Po, kde p0 = 0 (pravostranná alternativa) 2. Volba hladiny významnosti 49 • Hladinu významnosti volíme ze zadání a = 0.05. 3. Testování kritickým oborem • Testovací statistika Tw (Mi - M2) - p0 q2 02 °1 I £2. ni H2 (25.21053 - 26.92857) - 0 í1 ,6525722 1 3.1492192 19 14 -1.71804 ^0.8521351 -1.71804 0.9231117 ~ -1.86114 -1.8611 277 alpha <- 0.05 278 muO <- 0 279 ml <- mean(nose.Bc) 280 m2 <- mean(nose.Bb) 281 sl <- sd (nose . Bc ) 282 s2 <- sd(nose.Bb) 283 284 tw <- ((ml - m2) - muO) / (sqrt(sl ' • 2 / nl + s2 ' " 2 / n2)) # -I.86II45 • Welchova aproximace počtu stupňů volnosti (si.siý \ ni no ) d} = 2 \ 2 /„2\2 + n\— 1 712 — 1 1.6525722 1 3.1492192 X " 19 ~r 14 ( 1.6525722 > 2 í 3.14921a2 \ 2 \ 19 ) , V 14 ) 19-1 14-1 (0.1437365 + 0.7083986)2 0.14373652 i 0.70839862 18 13 0.7261342 0.03974999 = 18.26753 = 18.2675 285 df <- (sl " 2 / nl + s2 " 2 / n2)"2 / ((sl " 2 / nl) " 2 / (nl - 1) + (s2 " 2 / n2) " 2 / (n2 - 1)) # 18.26753 • Kritický obor W = {tdf(l -a); oo) = (íi8.26753(l-0.05);oo) = (1.732689; oo) 50 286 qt(l - alpha, df) # 1.732689 • Závěr testování Protože realizace testovací statistiky tyy = —1.861145 nenáleží do kritického oboru, tj. tyy ^ W, Hq nezamítáme na hladině významnosti a = 0.05. 4. Testování intervalem spolehlivosti • Interval spolehlivosti „2 2 bl , Ä2, (d, oo) = \m1-m2-\l— + — tdf (1 - a); oo (^(25.21053 - 26.92857) - ^lm2^ + ^^í18,26753(0.95); oo (-1.71804 - xO.9231117 x 1.732689; oo) (-3.31751; oo) = (-3.3175; oo) 287 DH <- (ml - m2) - sqrt(sl"2 / nl + s2"2 / n2) * qt(1 - alpha, df) # -3.317511 • Závěr testování Protože jj,q = 0 náleží do Waldova 95% empirického levostranného intervalu spolehlivosti, tj. po = 0 G IS, Hq nezamítáme na hladině významnosti a = 0.05. Testování p-hodnotou • p-hodnota ^hodnota = Pr(7V > tw) = 1 - Pr(Tw < tw) = 1 - Pr(7V < -1.861145) = 0.960552 = 0.9606 288 p.val <- 1 - pt(tw, df) # 0.960552 • Závěr testování Protože p-hodnota = 0.9606 je větší než a = 0.05, Hq nezamítáme na hladině významnosti a = 0.05. 6. Grafická vizualizace výsledků testování Vhodný grafem porovnávajícím střední hodnoty obou náhodných výběrů je krabicový diagram viz obrázek 33). 7. Interpretace výsledků: Na základě všech tří způsobů testování nezamítáme hypotézu o tom, že střední hodnota šířky nosu čínské populace je menší nebo rovna střední hodnotě šířky nosu bantuské populace. Nemáme tedy dostatek indicií k zamítnutí tvrzení, že šířka nosu čínské populace je statisticky významně menší nebo rovná šířce nosu bantuské populace. Poznámka: K otestování nulové hypotézy o rozdílu středních hodnot můžeme využít funkci t.test() s argumenty conf.level = 0.95 pro volbu hladiny významnosti a = 0.05, alternativě = 'greater' pro volbu pravostranné alternativní hypotézy a var.equal = F pro volbu Welchova dvouvýběrového í-testu. 51 32 - 30 - 28 - ö 26 24 - 22 - prumer median -1-1- cinska bantuska populace Obrázek 33: Krabicový diagram šířky nosu mužů čínské a bantuské populace 289 t.test(nose.Bb, nose.Be, conf.level = 0.95, alternative = 'greater', var.equal = F) Welch Two Sample t-test data: nose.Bb and nose.Be t = 1.8611, df = 18.268, p-value = 0.03945 alternative hypothesis: true difference in means is greater than 0 95 percent confidence interval: 0.1185797 Inf sample estimates: mean of x mean of y 26.92857 25.21053 290 291 292 293 294 295 296 297 298 299 300 52 10.4 Test o rozdílu korelačních koeficientů p\ — p2 Nechť (Xn, Yn)T ... (Xlni, Ylni)T je náhodný výběr z 7V2(//1;Si) a (X2i,Y2i)T,... (X2n2,Y2n2)T je na něm nezávislý náhodný výběr z N2([i2,'S2). Nechť po Je konstanta. Na hladině významnosti a testujeme jednu z následujících tří hypotéz oproti příslušné alternativní hypotéze. Hqi : pi — p2 = Po oproti Hn : pi — p2 7^ po (oboustranná alt.) Ho2 ■ pi — p2 < Po oproti Hi2 : pi — p2 > po (pravostranná alt.) Hqz : pi — p2 > po oproti H13 : pi — p2 < po (levostranná alt.) Test nazýváme dvouvýběrovým Z-testem o rozdílu korelačních koeficientů pi — p2. Testovací statistika má tvar Zi - Z2 - £0 Zw = — = , (10.5) 1 n\ — 3 kde Z\ = ^ ln , Z2 = ^ ln jsou Fisherovy Z-transformace výběrových korelačních koeficientů R\ a R2 a £0 = i ln \"^_p°0 je Fisherova Z-transformace konstanty p0 z nulové hypotézy, ii\ je rozsah prvního náhodného výběru a n2 je rozsah druhého náhodného výběru. Testovací statistika Zw pochází ze standardizovaného normálního rozdělení, tj. + n\— 3 712— 3 Kritický obor podle zvolené alternativní hypotézy má tvar Hu : pi- p2^ po W = (-00 ; ua/2) U (u1_a/2 ; 00) Hvl : pi - p2 > po W = ; 00) H13 : pi - p2 < po W = (-00; ua) kde «„/2, Ui-a/2, ua, «i-a jsou kvantily standardizovaného normálního rozdělení, jejichž hodnoty získáme pomocí w a implementované funkce qnorm(). Interval spolehlivosti má podle zvolené alternativní hypotézy jeden z následujících tvarů Hn : pi - p2 ^ po (d, h) = (tanh (zx - z2 - sgu1_a/2) ; tanh {zi - z2 - sgua/2)) HV2 : pi - p2 > po (d, 2) = (tanh (zi - z2 - sg«i_a) ; 2) Hiz ■ pi - p2 < Po (-2, h) = (-2 ; tanh (zi - z2 - sgua)) kde sg = yjn^_3 + a tanh je hyperbolický tangens, jehož hodnotu získáme pomocí Ol a implementované funkce tanh(). Poznámka: Dá se ukázat, že mezi korelačním koeficientem po a Co platí vztah po = tanh(£o), který je inverzí ke vztahu £0 = \ ln x_p° • Všimněme si, že z\ — Z2 — sg«i-a/2 a zi ~ z2 ~ sgua/2 Jsou dolní a horní hranice intervalu spolehlivosti pro Z-transformaci £oj neboť z\ i Z2 jsou Z-transformacemi korelačních koeficientů. Hyperbolický tangens potom funguje jako zpětná transformace, která převede hranice intervalu spolehlivosti pro £0 zpátky na hranice intervalu spolehlivosti pro korelační koeficient po. Poznámka: Protože parametry pi i p2 jsou korelační koeficienty, platí, že pi € ( —1; 1) a p2 G ( — 1,1). Rozdíl pi — p2 nabývá tedy hodnoty z intervalu (—2, 2). Proto levostranný interval spolehlivosti omezíme shora hodnotou 2, namísto nekonečnem, a pravostranný interval spolehlivosti omezíme zdola hodnotou -2, namísto mínus nekonečnem. p-hodnota má v závislosti na zvolené alternativní hypotéze jeden z následujících tvarů ifn : pi - p2 7^ po p-hodnota = 2min{Pr(Zw < zw), Pr(Zw > zw)} H12 ■ pi - p2 > Po p-hodnota = Pr(Zw > zw) = 1 - Pr(Zw < zw) H\z ■ pi - p2 < Po p-hodnota = Pr(Zw < zw) 53 kde Zw je náhodná veličina, zyy je realizace testovací statistiky Zyy (viz vzorec 10.5), tedy konkrétní číslo, a Pr(Zw < zyy) je distribuční funkce standardizovaného normálního rozdělení, jejíž hodnotu získáme pomocí a implementované funkce pnorm(). Příklad 10.13. Test o rozdílu dvou korelačních koeficientů pi — p2 Mějme datový soubor 19-more-samples-correlations-skull.txt, proměnnou nose.H popisující délku nosu a proměnnou nose.B popisující šířku nosu (viz sekce ??). Na hladině významnosti a = 0.05 zjistěte, zda se korelační koeficient délky a šířky nosu mužů bantuské populace významně liší od korelačního koeficientu délky a šířky nosu mužů peruánské populace. Řešení příkladu 10.13 Datový soubor načteme příkazem read.delim() a odstraníme z něj chybějící hodnoty příkazem na.omit(). Operátorem [] vybereme z tabulky nejprve údaje o délce nosu (nose.H), resp. šířce nosu nose.B u mužů bantuské populace pop == 'ban' a následně údaje o délce nosu (nose.H), resp. šířce nosu nose.B u mužů peruánské populace pop == 'per'. Pomocí příkazu length() dále zjistíme rozsahy obou náhodných výběrů a pomocí příkazu range() zjistíme rozsahy naměřených hodnot každé proměnné. Datový soubor obsahuje údaje o výšce a šířce nosu u 14 mužů bantuské 301 data <- read.delim('00-Data//19-more-samples-correlations-skuli.txt') 302 data <- na.omit(data) 303 nose.HB <- data [data$pop == 'ban 304 nose.BB <- data [data$pop == 'ban 305 nose.HP <- data [data$pop == 'per 306 nose.BP <- data[data$pop == 'per 307 308 nl <- length(nose.HB) # 14 309 n2 <- length(nose.HP) # 46 310 range (nose . HB) # 4IS8 311 range(nose.BB) # 22-31 312 range(nose.HP) # 44-51 313 range(nose.BP) # 19-26 populace, přičemž naměřené výšky nosu nabývají hodnot v rozmezí 41-58 ,mm, naměřené šířky nosu nabývají hodnot v rozmezí 22-31 mm. Soubor dále obsahuje údaje o výšce a šířce nosu u 46 mužů peruánské populace, přičemž naměřené výšky nosu nabývají hodnot v rozmezí 44-57 mm a naměřené šířky nosu nabývají hodnot v rozmezí 19-26 mm. Abychom mohli hypotézu ze zadání otestovat pomocí parametrického testu, musí oba datové soubory splňovat předpoklad dvourozměrné normality. Před samotným testováním hypotézy uvedené v zadání je potřeba ověřit, zda oba výběry tento předpoklad splňují. Závěr o dvourozměrné normalitě obou náhodných výběrů stanovíme na základě Mardiova testu v kombinaci s grafickou vizualizací dat pomocí 3D grafu a tečkového diagramu superponovaného 95 % elipsou spolehlivosti, analogicky, jako je uvedeno v sekci ??. Hladinu významnosti zvolíme a = 0.05. Začneme s náhodným výběrem mužů bantuské populace. 314 MVN::mvn(cbind(nose.HB, nose.BB), mvnTest = 'mardia')$multivariateNormality 315 # sikmost: 0.9610029 # spicatost : 0.4619885 Protože p-hodnota testu o nevýznamnosti koeficientu šikmosti, tj. 0.9610, je větší než 0.05, hypotézu o ne-významnosti koeficientu šikmosti nezamítáme na hladině významnosti a = 0.05. Dále protože p-hodnota testu o nevýznamnosti koeficientu špičatosti, tj. 0.3198, je větší než 0.05, nezamítáme hypotézu o nevýznamnosti koeficientu špičatosti. Protože náhodný výběr nevykazuje statisticky významné známky zešikmení ani zešpičatění, nezamítáme na hladině významnsoti a = 0.05 hypotézu o dvourozměrné normalitě náhodného výběru výšek a šířek nosu u mužů bantuské populace. Ke stejnému závěru bychom došli také použitím Henze-Zirklerova testu (p-hodnota = 0.7110 > 0.05) i Roystonova testu (p-hodnota = 0.3198 > 0.05). Nyní se podíváme na grafickou vizualizaci náhodného výběru (viz graf 34). Z 3D grafu i tečkového diagramu je zřejmá dvourozměrná normalita náhodného výběru. 3D graf nám ukazuje kopcovitý tvar náhodného výběru mužů ', 'nose.H'] ', 'nose.B'] ', 'nose.H'] ', 'nose.B'] 54 bantuské populace. Při pohledu na tečkový diagram vidíme, že všechny body leží uvnitř 95 % elipsy spolehlivosti. vyska nosu (v mm) Obrázek 34: 3D graf a tečkový diagram s 95% elipsou spolehlivosti pro výšku nosu a šířku nosu mužů bantuské populace (v mm) Dále se zaměříme na ověření předpoklad normality výšky a šířky nosu u mužů peruánské populace. 316 MVN::mvn(cbind(nose.HP, nose.BP), mvnTest = 'mardia')$multivariateNormality # sikmost: 0.409373 # spicatost: 0.2925019 Jelikož p-hodnota testu o nevýznamnosti koeficientu šikmosti, tj. 0.4094, je větší než 0.05, hypotézu o nevýznamnosti koeficientu šikmosti nezamítáme na hladině významnosti a = 0.05. Dále jelikož p-hodnota testu o nevýznamnosti koeficientu špičatosti, tj. 0.2925, je větší než 0.05, nezamítáme hypotézu o nevýznamnosti koeficientu špičatosti. Protože náhodný výběr nevykazuje statisticky významné zešikmení ani zešpičatění, nezamítáme hypotézu o dvourozměrné normalitě náhodného výběru výšek a šířek nosu u mužů bantuské populace na hladině významnsoti a = 0.05. Stejný závěr byhcom získali také použitím Henze-Zirklerova testu (p-hodnota = 0.1306 > 0.05). N azákladě Roystonova testu bychom hypotézu o dvourozměrné normalitě náhodného výběru hraničně zamítli (p-hodnota = 0.4543 < 0.05). Podíváme se tedy, co nám na dvourozměrnou normalitu náhodného výběru poví grafická vizualizace (viz obrázek 35). 3D graf i tečkový diagram ukazují na dvourozměrné normální rozdělení s jedním odlehlým pozorováním, které se nevešlo do 95 % elipsy spolehlivosti. Vzhledem k tomu, že rozsah náhodného výběru je 46, požadujeme, aby elipsa spolehlivosti pokrývala alepsoň 44 bodů. Zbylé dva body mohou ležet mimo elipsu spolehlivosti. V našem případě leží mimo elipsu spolehlivosti pouze jeden bod, proto se kloníme k závěru Mardiova a Henze-Zirklerova testu o dvourozměrné normalitě náhodného výběru výšek a šířek nosu u mužů peruánské populace. Protože oba výběry pochází z dvourozměnrého normálního rozdělení, můžeme hypotézu ze zadání otestovat parametrickým testem o rozdílu korelačních koeficientů p\ — pí. Testování provedeme v posloupnosti šesti kroků. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Korelační koeficient proměnných výška nosu a šířka nosu u mužů bantuské a korelační koeficient proměnných výška nosu a šířka nosu u mužů peruánské populace jsou shodné. Hi : Korelační koeficient proměnných výška nosu a šířka nosu u mužů bantuské a korelační koeficient proměnných výška nosu a šířka nosu u mužů peruánské populace nesou shodné. • matematická formulace nulové a alternativní hypotézy H0 : pi= pí pi- pí = po, kde p0 = 0 Hi'.pi^pi pi - pí ^ po, kde p0 = 0 (oboustranná alternativa) 55 2. Volba hladiny významnosti • Hladinu významnosti volíme v souladu se zadáním a = 0.05. 3. Testování kritickým oborem • Fisherovy Z-transformace výběrových korelačních koeficientů K výpočtu Fisherových Z-transformací potřebujeme znát hodnoty výběrových korelačních koeficientů. Tyto hodnoty získáme pomocí příkazu cor() s nastaveným argumentem method == 'pearson'. Výběrový korelační koeficient výšky a šířky nosu mužů bantuské populace r\ = 0.6719, výběrový korelační koeficient výšky a šířky nosu peruánské populace vyšel r2 = 0.1371. Nyní již můžeme vypočítat Fisheroy Z-transformace těchto korelačních koeficientů i Fisherovu Z-tranformaci konstanty po = 0 z nulové hypotézy. 7 Xl 1 + Rl 1 1 + 0.6719065 1 = 2lni- 0.6719065 = 2 ln5-°95816 = 0.5 x 1.62842 = 0.81421 = 0.8142 Zo = - ln- = llnl + 0.J370691 = lln 2 1 - 0.1370691 2 = 0.5 x 0.2758746 = 0.1379373 llnl+Po 2 1 - Po 1, 1 + 0 - ln-= 2 1-0 0.5 x 0 = 0 ■lni 56 317 alpha <- 0.05 318 ksiO <- 0 319 320 rl <- cor(nose.HB, nose.BB) # 0.6719065 321 r2 <- cor(nose.HP, nose.BP) # 0.1370691 322 323 zl <- 1 / 2 * log((l + rl) / (1 - rl)) # 0.8142106 324 z2 <- 1 / 2 * log((l + r2) / (1 - r2)) # 0.1379373 • Testovací statistika Zi - Z2 - Co — , V ^7=3 + 7^7=3 _ 0.8142106 - 0.1379373 -0 y 14-3 ~ 46-3 __0.6762733_ ~ ^0.09090909 + 0.0232558Ť _ 0.6762733 ~ ^0.1141649 _ 0.6762733 ~ 0.337883 = 2.001501 325 zw <- (zl - z2 - ksiO) / sqrt(l / (nl - 3) + 1 / (n2 - 3)) # 2.001502 • Kritický obor W = (-00 ; ua/2) U (u1_a/2 ; 00) = (-00 ; «0.05/2) U (mi-o.05/2 5 00) = (-00 ; m0.025) U («0.975 ; 00) = (-00; -1.959964) U (1.959964; 00) 326 qnorm(alpha / 2) # -1.959964 327 qnormd - alpha / 2) # 1.959964 • Závěr testování Protože realizace testovací statistiky zyy = 2.001501 náleží do kritického oboru, tj. zyy e W, Hq zamítáme na hladině významnosti a = 0.05. 4. Testování intervalem spolehlivosti • Interval spolehlivosti Pro výpočet intervalu spolehlivosti si nejdříve vypočítáme hodnotu sg. 1 1 + Tli — 3 «2—3 1 1 + ——- = V0.09090909 + 0.02325581 14 - 3 46 -3 ^0.1141649 = 0.337883 = 0.3379 57 (d, h) = (tanh (zx - z2 - sgu1_a/2) ; tanh {zi - z2 - sgua/2)) = (tanh (0.8142106 - 0.1379373 - 0.337883 x «1-0.05/2) ; tanh (0.8142106 - 0.1379373 - 0.337883 x «0.05/2)) = (tanh (0.6762733 - 0.337883 x 1.959964) ; tanh (0.6762733 - 0.337883 x (-1.959964))) = (tanh (0.01403478) ; tanh (1.338512)) = (0.01403386; 0.8713144) 328 sg <- sqrt(l / (nl - - 3) + 1 / (n2 - 3)) # 0 ,337883 329 dh <- tanh(zl - - z2 - * qnorm(l - alpha / 2)) # 0.01403392 330 hh <- tanh(zl - - z2 - " sg * qnorm(alpha / 2)) # 0.8713143 • Závěr testování Protože po = 0 nenáleží do Waldova 95% empirického oboustranného intervalu spolehlivosti, tj. po = 0 ^ IS, Hq zamítáme na hladině významnosti a = 0.05. 5. Testování p-hodnotou • p-hodnota p-hodnota = 2min{Pr(Zíil < zw),Pr(Zyy > zw)} = 2 min{Pr(Zíil < zw), 1 - Pr(Zw < zw)} = 2min{Pr(Zul < 2.001502), 1 - Pr(Zw < 2.001502)} = 2 min{0.9773308, 0.02266918} = 2 x 0.02266918 = 0.04533837 331 2 * min(pnorm(zw), 1 - pnorm(zw)) # 0.04533837 • Závěr testování Protože p-hodnota = 0.04533837 je menší než a = 0.05, Hq zamítáme na hladině významnosti a = 0.05. 6. Grafická vizualizace výsledků testování Vhodný grafem porovnávajícím míru závislosti v obou náhodných výběrech je tečkový diagram superponovaný lineárními regresními přímkami prokládajícími zobrazené body (viz obrázek 36). Nejprve příkazem plot() vykreslíme body výšky a šířky nosu mužů bantuské populace. Koeficienty lineární regresní přímky získáme pomocí funkce lm(), jejímž jedniným argumentem bude vztah nose.BB nose.HB, tj. vztah vyjadřující závislost mezi proměnnou nose.BB na ose y a proměnnou nose. H B na ose x. Koeficienty regresní přímky, které jsou vloženy v položce coeffcieints, získáme z výstupu funkce Im pomocí odkazu $coef. Dále vytvoříme posloupnost tisíce bodů x v rozsahu hodnot proměnné nose. H B a vypočítáme hodnoty regresní přímky v bodech posloupnosti x, které vložíme do proměnné y. Nakonec vykreslíme lineární regresní přímku v bodech x,y příkazem lines(). Do grafu dále přikreslíme body výšky a šířky nosu peruánské populace příkazem points(), které opět superponujeme lineární regresní přímkou získanou funkcí lm(), jejímž vstupním argumentem bude vztah nose.BP nose. H P budou proměnné nose. HP, tj. vztah vyjadřující závislost mezi proměnnou nose.BP na ose y a proměnnou nose. H P na ose x. Koeficienty regresní přímky opět získáme z výstupu funkce Im pomocí odkazu $coef. Dále vytvoříme posloupnost tisíce bodů x v rozsahu hodnot proměnné nose.HP, vypočítáme hodnoty regresní přímky v bodech posloupnosti x a vložíme je do proměnné y. Vykreslíme lineární regresní přímku v bodech x,y příkazem lines(). Nakonec do grafu doplníme pod osu x popisek obsahující hodnoty výběrových korelačních koeficientů (příkaz mtextQ) a legendu (příkaz legendQ). 58 332 par(mar = c(5, 4, 1, 1), family = 'Times') 333 plot(nose.HB, nose.BB, xlab = '', ylab = 'sirka nosu (v mm)', 334 xlim = c (40, 60), las = 1, 335 pen = 21, col = 'darkgreen', bg = 'olivedrab2') 336 k <- lm(nose.BB " nose.HB)$coef 337 x <- seq(min(nose.HB) , max(nose.HB) , length = 1000) 338 y <- k[l] + x * k [2] 339 lines(x, y , col = 'darkgreen' , lwd = 2) 340 341 points(nose.HP, nose.BP, pch = 21, col = 'tomato4', bg = 'yellow') 342 k <- lm(nose.BP " nose.HP)$coef 343 x <- seq(min(nose.HP) , max(nose.HP) , length = 1000) 344 y <- k[l] + x * k [2] 345 lines(x, y, col = 'tomato4', lwd = 2) 346 347 rl <- round(rl, digit = 4) 348 r2 <- round(r2, digit = 4) 349 mtext('vyska nosu (v mm)', side = 1, line = 2.3) 350 mtext (bquote (paste (rho [1] == . (rl), ', ', rho [2] == .(r2))), side = 1, line = 3.7) 351 legend('topleft' , pch = 21, pt.bg = c ( ' olivedrab2' , 'yellow'), 352 col = c('darkgreen', 'tomato4'), legend = c('bantuska p.', 'peruánská p.'), 353 bty = 'n') 7. Interpretace výsledků: Na základě všech tří způsobů testování zamítáme hypotézu o shodě korelačních koeficientů pro výšku nosu a šířku nosu bantuské a peruánské populace. Mezi korelačním koeficientem výšky a šířky nosu u mužů bantuské populace a u mužů peránské populace existuje statisticky významný rozdíl. Mezi výškou a šikou nosu u mužů bantuské populace existuje význačný stupeň přímé lineární závislosti (pi = 0.6719). Mezi výškou a šířkou nosu u mužů peruánské populace existuje nízký stupeň přímé lineární závislosti (p2 = 0.1379) (viz stupnice míry závislosti pro Pearsonův korelační koeficient, kapitola ??). ★ 59 Obrázek 36: Tečkový diagram závislosti výšky nosu a šířky nosu u mužů bantuské a peruánské populace 60 Příklad 10.14. Test o rozdílu dvou korelačních koeficientů p\ — p2 Mějme datový soubor 16-anova-head.txt, proměnnou head.W popisující šířku hlavy a proměnnou bizyg.W popisující šířku tváře (viz sekce ??). Na hladině významnosti a = 0.05 testujte hypotézu o shodě korelačního koeficientu šířky hlavy a šířky tváře mužů a korelačního koeficientu šířky hlavy a šířky tváře žen. Řešení příkladu 10.14 Datový soubor načteme příkazem read.delim() a odstraníme z něj chybějící hodnoty příkazem na.omit(). Operátorem [] vybereme z tabulky údaje o šířce hlavy (head.W) a šířce tváře (bizyg.W) u mužů (sex == 'm'), resp. žen (sex == T). Dále zjistíme rozsahy obou náhodných výběrů příkazem length() a rozsahy naměřených hodnot každé proměnné příkazem range(). Datový soubor obsahuje naměřené hodnoty šířky hlavy a šířky tváře u 63 mužů, 354 data <- read.delim('00-Data//16-anova-head.txt') 355 data <- na.omit(data) 356 head.WM <- data[data$sex == 'm', 'head.W'] 357 bizyg.WM <- data[data$sex == 'm', 'bizyg.W'] 358 head.WF <- data[data$sex == 'f', 'head.W'] 359 bizyg.WF <- data[data$sex == 'f', 'bizyg.W'] 360 361 nl <- length(head.WM) # 63 362 n2 <- length(head.WF) # 85 363 range(head.WM) # 142-170 364 range(bizyg.WM) # 113-155 365 range(head.WF) # 135-162 366 range(bizyg.WF) # 120-151 přičemž naměřené šířky hlavy nabývají hodnot v rozmezí 142-170 mm, naměřené šířky tváře nabývají hodnot v rozmezí 113-155 mm. Soubor dále obsahuje údaje o šířce hlavy a šířce tváře 85 žen, přičemž naměřené šířky hlavy nabývají hodnot v rozmezí 135-162 mm a naměřené šířky tváře nabývají hodnot v rozmezí 120-151 mm. Aby bylo možné otestovat hypotézu ze zadání pomocí parametrického testu o rozdílu dvou korelačních koeficientů, musí oba datové soubory pocházet z dvourozměrného normálního rozdělení. Před samotným testem hypotézy ze zadání je tedy potřeba ověřit, zda oba výběry tento předpoklad splňují. Závěr o dvourozměrné normalitě obou náhodných výběrů stanovíme na základě Roystonova testu (a = 0.05) v kombinaci s grafickou vizualizací dat 367 MVN::mvn(cbind(head.WM, bizyg.WM), mvnTest = 'royston')$multivariateNormality # 0.003572831 Protože p-hodnota Roystonova testu, tj. 0.003573 je menší než 0.05, hypotézu o dvourozměrné normalitě náhodného výběru šířky lebky a šířky tváře u mužů zamítáme na hladině významnosti a = 0.05. Ke stejnému závěru dojdeme také použitím Henze-Zirklerova testu (p-hodnota = 0.002498 < 0.05) i Mardiova testu (p-hodnota pro koeficient šikmosti = 0.006636 < 0.05, p-hodnota pro koeficient špičatosti = 0.2778 > 0.05). Nyní se podíváme na grafickou vizualizaci náhodného výběru (viz obrázek 37). Loading required package: viridisLite 368 3D graf nám nabízí pohled na dvourozměrné normální rozdělení, tvořené pospolitým kopcem s jednou výrazně odlehlou hodnotou a několika méně odlehlými hodnotami. Aby data pocházela z dvourozměrného normálního rozdělení, je potřeba, aby alespoň 60 bodů leželo uvnitř 95% elipsy spolehlivosti. Zbylé 3 body se mohou realizovat mimo elipsu spolehlivosti. Z grafu vidíme, že mimo elipsu spolehlivosti leží právě 3 body, což ukazuje na dvourozměrnou normalizu náhodného výběru. V tomto případě tedy nastává rozkol mezi grafickou vizualizací a závěry testování dvourozměrné normality. Protože však všechny tři testy ukazují na porušení předpokladu normality, předpokládáme, že datový soubor šířek hlavy a šířek tváře mužů nepochází z normálního rozdělení. Nyní prozkoumáme předpoklad normality šířky hlavy a šířky tváře u žen. Jelikož p-hodnota Roystonova testu, tj. 0.02646746, je menší než 0.05, hypotézu o dvourozměrné normalitě náhodného výběru šířky hlavy a šířky tváře žen zamítáme na hladině významnosti a = 0.05. Ke stejnému závěru 61 sirka hlavy (v mm) Obrázek 37: 3D graf a tečkový diagram s 95% elipsou spolehlivosti pro šířku hlavy a šířku tváře mužů (v mm) 369 MVN::mvn(cbind(head.WF, bizyg.WF), mvnTest = 'royston')$multivariateNormality # 0.02646746 bychom došli také Mardiovým testem (p-hodnota pro koeficient šikmosti = 0.00289277 > 0.05, p-hodnota pro koeficient špičatosti = 0.10385298 > 0.05). Na základě Henze-Zirklerova testu bychom hypotézu o dvourozměrné normalitě náhodného výběru nezamítli (p-hodnota = 0.2727 > 0.05). Zde tedy nastává rozkol mezi výsledky testů normality. Před stanovením závěru o rozdělení náhodného výběru se podívejme na vizualizaci datového souboru (viz obrázek 38). sirka hlavy (v mm) Obrázek 38: 3D graf a tečkový diagram s 95% elipsou spolehlivosti pro šířku hlavy a šířku tváře žen (v mm) 3D graf nám nabízí pohled na dvourozměrné normální rozdělení tvořené pospolitým kopcem s alespoň jednou výrazně odlehlou hodnotou a několika méně odlehlými hodnotami. Aby data pocházela z dvourozměrného normálního rozdělení, je potřeba, aby alespoň 81 bodů leželo uvnitř 95% elipsy spolehlivosti. Zbylé 4 body mohou ležet mimo elipsu spolehlivosti. Z tečkového grafu je zřejmé, že mimo elipsu spolehlivosti však leží 5 bodů. V tomto případě tedy grafická vizualizace podporuje výsledek Roystonova a Mardiova testu. Na základě výsledků testování a grafické vizualizace tedy zamítáme hypotézu o dvourozměrném normálním rozdělení náhodného výběru šířek hlavy a šířek tváře žen. Protože předpoklad dvourozměné normality obou náhodných výběrů je porušen, nemůžeme k testování hypotézy ze zadání použít parametrický test o rozdílu korelačních koeficientů p\ — pí. 62 Příklad 10.15. Test o rozdílu dvou korelačních koeficientů p\ — p2 Mějme datový soubor 05-one-sample-correlation-skull-mf.txt, proměnnou skuli.pH popisující největší výšku mozkovny a proměnnou face.H popisující morfologickou výšku tváře (viz sekce ??). Na hladině významnosti a = 0.10 zjistěte, zda korelační koeficient největší výšky mozkovny a morfologické výšky tváře mužů je statisticky významně větší než korelační koeficient největší výšky mozkovny a morfologické výšky tváře žen. Řešení příkladu 10.15 Datový soubor načteme příkazem read.delim() a odstraníme z něj chybějící hodnoty příkazem na.omit(). Operátorem [] vybereme z tabulky nejprve údaje o největší výšce mozkovny (skuli.pH), resp. morfologické výšce tváře (face.H) u mužů (sex == 'm') a následně údaje o největší výšce mozkovny (skuli.pH), resp. morfologické výšce tváře (face.H) u žen (sex == 'f'). Dále zjistíme rozsahy obou náhodných výběrů (length()) rozsahy naměřených hodnot (range()) každé proměnné. Datový soubor obsahuje údaje o největší výšce mozkovny a morfologické výšce tváře u 164 mužů, 370 data <- read.delim('00-Data//05-one-sample-correlation - skuli-mf.txt') 371 data <- na.omit(data) 372 head(data) 373 skuli.pHM <- data [data$sex == 'm', 'skuli.pH ' ] 374 face.HM <- data[data$sex == 'm', 'face.H'] 375 skuli.pHF <- data [data$sex == 'f, 'skuli.pH'] 376 face.HF <- data[data$sex == 'f', 'face.H'] 377 378 nl <- length(skuli.pHM) # i64 379 n2 <- length(skuli.pHF) # 78 380 range(skuli.pHM) # 127-149 381 range(face.HM) # 100-136 382 range(skuli.pHF) # 120-142 383 range(face.HF) # 95-119 přičemž naměřené výšky mozkovny nabývají hodnot v rozmezí 127-149 mm, naměřené morfologické výšky tváře nabývají hodnot v rozmezí 100-136 mm. Soubor dále obsahuje údaje o největší výšce mozkovny a morfologické výšce tváře u 78 žen, přičemž naměřené výšky mozkovny nabývají hodnot v rozmezí 120-142 mm a naměřené hodnoty morfologické výšky tváře nabývají hodnot v rozmezí 95-119 mm. Řešení příkladu vede na test o rozdílu dvou korelačních koeficientů. Aby bylo možné provést parametrický test o rozdílu dvou korelačních koeficientů, musí oba datové soubory splňovat předpoklad dvourozměrné normality. Před samotným testováním je tedy potřeba ověřit, zda oba výběry tento předpoklad splňují. Závěr o dvourozměrné normalitě obou náhodných výběrů stanovíme na základě Henze-Zirklerova testu (a = 0.05) v kombinaci s grafickou vizualizací dat. 384 MVN::mvn(cbind(skuli.pHM, face.HM), mvnTest = 'hz')$multivariateNormality # 0.719099 Protože p-hodnota Henze-Zirklerova testu, tj. 0.7191 je větší než 0.05, hypotézu o dvourozměrné normalitě náhodného výběru největší výšky mozkovny a morfologické výšky tváře u mužů nezamítáme na hladině významnosti a = 0.05. Ke stejnému závěru bychom došli také použitím Mardiova testu (p-hodnota pro koeficient šikmosti = 0.3853 > 0.05, p-hodnota pro koeficient špičatosti = 0.4620 > 0.05) i Roystonova testu (p-hodnota = 19998 > 0.05). Nyní se podíváme na grafickou vizualizaci náhodného výběru (viz obrázek 39). 3D graf nám nabízí pohled na dvourozměrné normální rozdělení, tvořené dvěma vzájemně se prolínajícími kopci s několika odlehlými hodnotami. Zda je množství odlehlých hodnot v normě, zjistíme pohledem na tečkový diagram. Aby data pocházela z dvourozměrného normálního rozdělení stačí, aby alespoň 95% bodů, tj. 156 bodů, leželo uvnitř elipsy spolehlivosti. Zbylých 8 bodů může ležet mimo elipsu spolehlivosti. Z grafu vidíme, že mimo elipsu spolehlivosti leží právě 8 bodů, což je v pořádku. O náhodném výběru největší výšky mozkovny a morfologické výšky tváře u mužů tedy na základě testování i grafické vizualizace předpokládáme, že pochází z dvourozměrného normálního rozdělení. Nyní prozkoumáme předpoklad normality největší výšky mozkovny a morfologické výšky tváře u žen. Jelikož p-hodnota Henze-Zirklerova testu, tj. 0.9997 je větší než 0.05, hypotézu o dvourozměrné normalitě náhodného výběru největší výšky mozkovny a morfologické výšky tváře žen nezamítáme na hladině významnosti 63 Obrázek 39: 3D graf a tečkový diagram s 95% elipsou spolehlivosti pro největší výšku mozkovny a morfologickou výšku tváře mužů (v mm) 385 MVN::mvn(cbind(skuli.pHF, face.HF), mvnTest = 'hz')$multivariateNormality # 0.9997521 a = 0.05. Ke stejnému závěru bychom došli také Mardiovým testem (p-hodnota pro koeficient šikmosti = 0.9674 > 0.05, p-hodnota pro koeficient špičatosti = 0.6223 > 0.05) i Roystonovým testem (p-hodnota = 0.4054 > 0.05). Obrázek 40: 3D graf a tečkový diagram s 95% elipsou spolehlivosti pro největší výšku mozkovny a morfologickou výšku tváře žen (v mm) 3D graf nabízí pohled na pospolité dvourozměrné normální rozdělení tvořené jedním souvislým kopcem (viz obrázek 40). Aby data pocházela z dvourozměnrého normálního rozdělení, je potřeba, aby 95% elipsa spolehlivosti pokrývala alespoň 74 bodů. Zbylé čtyři body se mohou realizovat mimo elipsu spolehlivosti. Pohledem na tečkový diagram zjišťujeme, že mimo elipsu spolehlivosti leží pouze tři body, což podporuje náš závěr o dvourozměrné normalitě náhodného výběru. O náhodném výběru největší výšky mozkovny a morfologické výšky tváře u žen tedy na základě testování i grafické vizualizace předpokládáme, že pochází z dvourozměrného normálního rozdělení. Protože oba výběry pochází z dvourozměrného normálního rozdělení, můžeme k testování použít parametrický test o rozdílu korelačních koeficientů p\ —p2- Zadaným úkolem je zjistit, zda korelační koeficient největší výšky mozkovny a morfologické výšky tváře mužů je statisticky významně větší než korelační koeficient největší výšky mozkovny a morfologické výšky tváře žen. Toto tvrzení je zněním alternativní hypotézy a nulové hypotéza tvoří doplněk k alternativní hypotéze. 64 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Korelační koeficient nejvčtší výšky mozkovny a morfologické výšky tváře mužů je menší nebo roven korelačnímu koeficientu nejvčtší výšky mozkovny a morfologické výšky tváře žen. Hi : Korelační koeficient nejvčtší výšky mozkovny a morfologické výšky tváře mužů je větší než korelační koeficient nejvčtší výšky mozkovny a morfologické výšky tváře žen. • matematická formulace nulové a alternativní hypotézy H0 ■ pi < P2 Pi ~ P2 < Po, kde p0 = 0 Hi : pi > p2 pi - p2 > Po, kde p0 = 0 (pravostranná alternativa) 2. Volba hladiny významnosti • Hladinu významnosti volíme podle zadání a = 0.10. 3. Testování kritickým oborem • Fisherovy Z-transformace výběrových korelačních koeficientů K výpočtu Fisherových Z-transformací je potřeba nejprve stanovit hodnoty výběrových korelačních koeficientů, a to pomocí příkazu cor(). Výběrový korelační koeficient výšky morkovny a morfologické výšky tváře mužů r\ = 0.33064, výběrový korelační koeficient výšky a šířky nosu peruánské populace vyšel r2 = 0.06417. Nyní vypočítáme Fisherovy Z-transformace obou výběrových korelačních koeficientů a konstanty po = 0 z nulové hypotézy. 7 1, 1 + Zl = 2lnT^ = llnl + 0.3306431 = lln 2 1 - 0.3306431 2 = 0.5 x 1.62842 = 0.81421 = 0.8142 Zo = - ln- 2 1-R 2 1. 1 + 0.06417166 1, - ln-= - ln 1.137144 2 1 - 0.06417166 2 0.5 x 0.1285199 = 0.06425995 1^1+Po 2 1 1 In 2 1 Po 1 + 0 0 ■lni 0.5 x 0 = 0 386 alpha <- 0.10 387 ksiO <- 0 388 rl <- cor(skull.pHM, face.HM) # 0.3306431 389 r2 <- cor(skull.pHF, face.HF) # 0.06417166 390 zl <- 1 / 2 * log((l + rl) / (1 - rl)) # 0.3435501 391 z2 <- 1 / 2 * log((l + r2) / (1 - r2)) # 0.06425996 65 • Testovací statistika n\— 3 712—3 0.3435501 - 0.06425996 164-3 1 78-3 0.2792901 ^0.00621118 + 0.01333333 0.2792901 ^0.01954451 0.2792901 0.1398017 1.997759 392 zw <- (zl - z2 - ksiO) / sqrt(l / (nl - 3) + 1 / (n2 - 3)) # 1.997759 • Kritický obor W = («i_„ ; oo) = (mi-o.io ; oo) = (m0.9o ; oo) = (1.281552; oo) 393 qnormCl - alpha) # 1.281552 • Závěr testování Protože realizace testovací statistiky zyy = 1.9978 náleží do kritického oboru, tj. zyy e W, Hq zamítáme na hladině významnosti a = 0.10. 4. Testování intervalem spolehlivosti • Interval spolehlivosti Pro výpočet intervalu spolehlivosti si nejdříve vypočítáme hodnotu sg. 1 1 + Tli — 3 «2—3 + ——- = V0.00621118 + 0.01333333 164 - 3 78 -3 ^0.01954451 = 0.1398017 = 0.1398 (d, 2) = (tanh (zx - z2 - ssui-a) 5 2) = (tanh (0.3435501 - 0.06425996 - 0.1398017 x «i_0.io) ; 2) = (tanh (0.2792901 - 0.1398017 x 1.281552) ; 2) = (tanh (0.100127) ; 2) = (0.09979368; 2) 66 394 sg <- sqrtCl / (nl - 3) + 1 / (n2 - 3)) # 0.1398017 395 dh <- tanh(zl - z2 - sg * qnorm(l - alpha)) # 0.0997938 • Závěr testování Protože po = 0 nenáleží do Waldova 90% empirického levostranného intervalu spolehlivosti, tj. po = 0 ^ IS, Hq zamítáme na hladině významnosti a = 0.10. 5. Testování p-hodnotou • p-hodnota p-hodnota = ~Pr(Zw > zw) = 1 - vt(zw < zw) = 1 - Pr(Zw < 1.997759) = 0.02287137 396 1 - pnorm(zw) # 0.02287137 • Závěr testování Protože p-hodnota = 0.02287 je menší než a = 0.10, Hq zamítáme na hladině významnosti a = 0.10. 6. Grafická vizualizace výsledků testování Rozdíl v míře závislosti mezi oběma náhodnými výběry vizualizujeme pomocí tečkového diagramu s lineárními regresními přímkami prokládajícími zobrazené body (viz obrázek 44). Koeficienty lineární regresní přímky pro největší výšku mozkovny a morfologickou výšku tváře mužů získáme pomocí funkce lm(), jehož argumentem bude vztah face. H M skul I. pH M vyjadřující závislost mezi proměnnou face. H M na ose y a proměnnou skull.pHM na ose x. Do grafu dále dokreslíme body největší výšky mozkovny a morfologické výšky tváře žen, které superponujeme lineární regresní přímkou získanou funkcí lm() s argumentem face.HF skuli.pHF. Nakonec do grafu doplníme popisek obsahující hodnoty výběrových korelačních koeficientů a legendu. 397 par(mar = c(5, 4, 1, 1), family = 'Times') 398 plot(skull.pHM, face.HM, xlab = '', ylab = 'morfologická výska tvare (v mm)', 399 xlim = c(120, 150), ylim = c(95, 135), las = 1, 400 pch = 21, col = 'blue', bg = 'cornflowerblue') 401 k <- lm(face.HM ~ skull.pHM)$coef 402 x <- seq(min(skull.pHM), max(skull.pHM), length = 1000) 403 y <- k[l] + x * k [2] 404 lines (x, y, col = 'darkblue' , lwd = 2) 405 406 points(skull.pHF, face.HF, pch = 21, col = 'red', bg = 'salmon') 407 k <- lm(face.HF ~ skull.pHF)$coef 408 x <- seq(min(skull.pHF) , max(skull.pHF) , length = 1000) 409 y <- k[l] + x * k [2] 410 lines(x, y, col = 'darkred', lwd = 2) 411 412 rl <- round(rl, digit = 4) 413 r2 <- round(r2, digit = 4) 414 mtext('nej vet si vyska mozkovny 415 mtext(bquote(paste(rho [1] == . 416 legend('topleft', pch =21, pt 417 col = c('blue', 'red'), 418 bty = 'n') (v mm)', side = 1, line = 2.3) (rl), ', ', rho [2] == .(r2))), side = 1, line = 3.7) bg = c('cornflowerblue', 'salmon'), legend = c('muzi', 'zeny'), 67 130 > 120 - st no - 100 T 125 130 135 140 145 nejvetsi vyska mozkovny (v mm) p, = 0.3306, p2 = 0.0642 Obrázek 41: Krabicový diagram největší výšky mozkovny a morfologické výšky tváře žen (v mm) 7. Interpretace výsledků: Na základě všech tří způsobů testování zamítáme hypotézu Hq. Korelační koeficient největší výšky mozkovny a morfologické výšky tváře u mužů je statisticky významně větší než korelační koeficient největší výšky mozkovny a morfologické výšky tváře u žen. Mezi největší výškou mozkovny a morfologickou výškou tváře u mužů existuje mírný stupeň přímé lineární závislosti (pi = 0.3306). Mezi největší výškou mozkovny a morfologickou výškou tváře u žen existuje velmi nízký stupeň přímé lineární závislosti (pi = 0.0642). ★ 68 Příklad 10.16. Test o rozdílu dvou korelačních koeficientů p\ — p2 Mějme datový soubor 13-two-samples-correlations-trunk.txt, proměnnou lowex.L popisující délku dolní končetiny v mm a proměnnou tru.L popisující délku trupu v mm (viz sekce ??). Na hladině významnosti a = 0.01 zjistěte, zda je korelační koeficient délky dolní končetiny a délky trupu u mužů menší než korelační koeficient délky dolní končetiny a délky trupu žen. Řešení příkladu 10.16 Datový soubor načteme příkazem read.delim() a odstraníme z něj chybějící hodnoty. Operátorem [] vybereme z tabulky nejprve údaje o délce dolní končetiny (lowex.L) a délce trupu (tru.L) u mužů (sex == 'm'), resp. u žen (sex == 'f). Dále zjistíme rozsahy obou náhodných výběrů a rozsahy naměřených hodnot každé proměnné. Datový 419 data <- read.delim('00-Data//13-two-samples-correlations-trunk.txt') 420 data <- na.omit(data) 421 lowex.LM <- data[data$sex == 'm', 422 tru.LM <- data[data$sex == 'm', 423 lowex.LF <- data[data$sex == 'i', 424 tru.LF <- data [data$sex == 'i', 425 426 nl <- length(lowex.LM) # 75 427 n2 <- lengthdowex . LF) # 100 428 range(lowex.LM) # 915-1114 429 range(tru.LM) # 4oo-548 430 range(lowex.LF) # 836-1076 431 range(tru.LF) # 323~492 soubor obsahuje údaje o délce dolní končetiny a délce trupu u 75 mužů, přičemž naměřené délky dolní končetiny nabývají hodnot v rozmezí 915-1114mm, naměřené délky trupu nabývají hodnot v rozmezí 400-548mm. Soubor dále obsahuje údaje o největší délce dolní končetiny a délce trupu u 100 žen, přičemž naměřené délky dolní končetiny nabývají hodnot v rozmezí 836-1076 mm a naměřené délky trupu nabývají hodnot v rozmezí 323-492 mm. Řešení příkladu vede na test o rozdílu dvou korelačních koeficientů. Před provedením parametrického testu o rozdílu korelačních koeficientů je potřeba ověřit předpoklad dvourozměrné normality obou náhodných výběrů. Závěr o dvourozměrné normalitě obou výběrů stanovíme na základě Henze-Zirklerova testu (a = 0.05) v kombinaci s grafickou vizualizací dat. 432 MVN::mvn(cbind(lowex.LM, tru.LM), mvnTest = 'hz')$multivariateNormality # 0.543489 Protože p-hodnota Henze-Zirklerova testu, tj. 0.5435 je větší než 0.05, hypotézu o dvourozměrné normalitě náhodného výběru délky dolní končetiny a délky trupu mužů nezamítáme na hladině významnosti a = 0.05. Ke stejnému závěru bychom došli také použitím Mardiova testu (p-hodnota pro koeficient šikmosti = 0.5599 > 0.05, p-hodnota pro koeficient špičatosti = 0.4298 > 0.05) i Roystonova testu (p-hodnota = 0.1876 > 0.05). Na 3D grafu vidíme dvourozměné normální rozdělení tvoření jedním kopcem (viz obrázek 42). Aby data pocházela z dvourozměrného normálního rozdělení stačí, aby alespoň 71 bodů leželo uvnitř 95% elipsy spolehlivosti. Zbylé 4 body mohou ležet mimo elipsu spolehlivosti. Z tečkového diagramu vidíme, že mimo elipsu spolehlivosti leží pouze dva body a jeden bod leží na hranici elipsy spolehlivosti. O náhodném výběru největší délky dolní končetiny a délky trupu mužů tedy předpokládáme, že pochází z dvourozměrného normálního rozdělení. Nyní prozkoumáme předpoklad normality délky dolní končetiny a délky trupu u žen. 433 MVN::mvn(cbind(lowex.LF, tru.LF), mvnTest = 'hz')$multivariateNormality # 0.587225 Jelikož p-hodnota Henze-Zirklerova testu, tj. 0.5872 je větší než 0.05, hypotézu o dvourozměrné normalitě náhodného výběru největší výšky mozkovny a morfologické výšky tváře žen nezamítáme na hladině významnosti a = 0.05. Ke stejnému závěru bychom došli také Mardiovým testem (p-hodnota pro koeficient šikmosti = 0.1769 > 0.05, p-hodnota pro koeficient špičatosti = 0.8360 > 0.05) i Roystonovým testem (p-hodnota = 0.6339 > 0.05). 'lowex.L ' ] 'tru.L'] 'lowex.L ' ] 'tru.L ' ] 69 Obrázek 42: 3D graf a tečkový diagram s 95% elipsou spolehlivosti pro délku dolní končetiny a délku trupu u mužů (v mm) Obrázek 43: 3D graf a tečkový diagram s 95% elipsou spolehlivosti pro délku dolní končetiny a délku trupu u žen (v mm) 3D graf nabízí pohled na dvourozměrné normální rozdělení tvořené jedním souvislým kopcem a několik odlehlými hodnotami (viz obrázek 43). Aby data pocházela z dvourozměnrého normálního rozdělení, je potřeba, aby odlehlých hodnot nebylo více než pět. Pohledem na tečkový diagram zjišťujeme, že mimo elipsu spolehlivosti leží pouze čtyři body. O náhodném výběru délky dolní končetiny a délky trupu žen tedy předpokládáme, že pochází z dvourozměrného normálního rozdělení. Protože oba výběry pochází z dvourozměrného normálního rozdělení, můžeme k testování použít parametrický test o rozdílu korelačních koeficientů p\ — pí. Zadaným úkolem je zjistit, zda korelační koeficient délky dolní končetiny a délky trupu mužů je statisticky významně menší než korelační koeficient délky dolní končetiny a délky trupu žen. Toto tvrzení je zněním alternativní hypotézy a nulové hypotéza tvoří doplněk k alternativní hypotéze. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Korelační koeficient délky dolní končetiny a délky trupu mužů je větší nebo rovek korelačnímu koeficientu délky dolní končetiny a délky trupu žen. Hi : Korelační koeficient délky dolní končetiny a délky trupu mužů je menší než korelační koeficient délky dolní končetiny a délky trupu žen. 70 • matematická formulace nulové a alternativní hypotézy H0 ■ Pi > P2 Pi ~ P2 > Po, kde p0 = 0 Hi : pi < p2 pi - p2 < Po, kde p0 = 0 (levostranná alternativa) 2. Volba hladiny významnosti • Hladinu významnosti volíme podle zadání a = 0.05. 3. Testování kritickým oborem • Fisherovy Z-transformace výběrových korelačních koeficientů K výpočtu Fisherových Z-transformací je potřeba nejprve stanovit hodnoty výběrových korelačních koeficientů. Výběrový korelační koeficient délky dolní končetiny a délky trupu mužů r\ = 0.05976, výběrový korelační koeficient délky dolní končetiny a délky trupu žen vyšel r2 = 0.285256. Nyní vypočítáme Fisherovy Z-transformace obou výběrových korelačních koeficientů a konstanty po = 0 z nulové hypotézy. 7 1, 1 + Zi = — ln- 2 1 - i?i 1 1 + 0.05976 1 = - ln —- = - ln 1.127116 2 1 - 0.05976 2 = 0.5 x 0.1196626 = 0.0598313 = 0.05983 Z2 = 2lnT^ = llnl + 0.285256 =lln 2 1 - 0.285256 2 = 0.5 x 0.5867888 = 0.2933944 1, 1+Po « ln - 2 1 - po - ln ——— = - ln 1 2 1-0 2 0.5 x 0 = 0 434 alpha <- 0.01 435 ksiO <- 0 436 rl <- cor(lowex.LM, tru.LM) # 0.05975781 437 r2 <- cor(lowex.LF, tru.LF) # 0.285256 438 zl <- 1 / 2 * log((l + rl) / (1 - rl)) # 0.0598291 439 z2 <- 1 / 2 * log((l + r2) / (1 - r2)) # 0.2933943 • Testovací statistika 71 _ Zi - Z2 - Co _ 0.0598291 - 0.2933944 - 0 ~ / i i i y 75-3 ~ 100-3 _ -0.2335653 ~ ^0.01388889 + 0.01030921 _ -0.2335653 ~ ^0.02419817 _ -0.2335653 ~ 0.1555576 = -1.501471 440 zw <- (zl - z2 - ksiO) / sqrt(l / (nl - 3) + 1 / (n2 - 3)) # -1.501471 • Kritický obor W = (—oo ; ua) = (-oo; m0.oi) = (-oo; -1.644854) 441 qnorm(alpha) # -2.326348 • Závěr testování Protože realizace testovací statistiky zyy = —1.5015 nenáleží do kritického oboru, tj. zyy ^ W, Hq nezamítáme na hladině významnosti a = 0.05. 4. Testování intervalem spolehlivosti • Interval spolehlivosti Pro výpočet intervalu spolehlivosti si nejdříve vypočítáme hodnotu sg. 1 1 + Tli — 3 Tl2 — 3 i i + ™—r: = ^0.01388889 + 0.01030928 75 - 3 100 -3 ^0.02419817 = 0.1555576 = 0.1556 (—2, h) = (—2 ; tanh (z\ — z2 — sgua)) (-2, h) = (-2 ; tanh (0.0598291 - 0.2933944 - 0.1555576m0.oi)) = (-2 ; tanh (-0.2335653 - 0.1555576 x (-2.326348))) = (-2; tanh (0.1283158)) = (-2; 0.1276162) 72 442 sg <- sqrtQ / (nl - 3) + 1 / (n2 - 3) ) # 0.1555576 443 hh <- tanh(zl - z2 - sg * qnorm(alpha)) # 0.1276162 • Závěr testování Protože po = 0 náleží do Waldova 99% empirického pravostranného intervalu spolehlivosti, tj. po = 0 G IS, Hq nezamítáme na hladině významnosti a = 0.01. 5. Testování p-hodnotou • p-hodnota p-hodnota = Pr(Zyy < zw) = Pľ(Zw < -1.501471) = 0.06661688 444 pnorm(zw) # 0.06661688 • Závěr testování Protože p-hodnota = 0.06662 je větší než a = 0.01, Hq nezamítáme na hladině významnosti a = 0.01. 6. Grafická vizualizace výsledků testování Míru závislosti obou náhodných výběrů vizualizujeme pomocí tečkového diagramu s lineárními křivkami prokládajícími naměřené body (viz obrázek 44). Koeficienty lineární přímky pro délku dolní končetiny a délku trupu mužů, resp. žen, získáme pomocí funkce lm(), analogicky jako v předchozích příkladech. 445 par(mar = c(5, 4, 1, 1), family = 'Times') 446 plot(lowex.LM, tru.LM, xlab = '', ylab = 'délka trupu (v mm)', 447 xlim = c(820, 1100), ylim = c(350, 550), las = 1, 448 pch = 21, col = 'blue', bg = 'cornflowerblue ' ) 449 k <- lm(tru.LF " lowex.LF)$ coef 450 x <- seq(min(lowex.LF) , max(lowex.LF) , length = 1000) 451 y <- k[l] + x * k[2] 452 lines(x, y, col = 'darkred', lwd = 2) 453 454 points(lowex.LF , tru.LF, pch = 21, 455 k <- lm(tru.LM " lowex.LM)$coef 456 x <- seq(min(lowex.LM) , max(lowex. 457 y <- k[l] + x * k [2] 458 lines(x, y, col = 'darkblue' , lwd 459 460 rl <- round(rl, digit = 4) 461 r2 <- round(r2, digit = 4) 462 mtext('délka dolni končetiny (v mm)', 463 mtext(bquote(paste(rho [1] == . (rl), 464 legend('topleft', pch = 21, pt.bg = c 465 col = c('blue ' , 'red'), legend 466 bty = 'n') col = 'red', bg = 'salmon') LM), length = 1000) = 2) side = 1, line = 2.3) ', ', rho [2] == .(r2))), side = 1, line = 3.7) ('cornflowerblue', 'salmon'), = c('muzi' , 'zeny') , 7. Interpretace výsledků: Na základě všech tří způsobů testování nezamítáme hypotézu Hq. Korelační koeficient délky dolní končetiny a délky trupu mužů není statisticky významně menší než korelační koeficient délky dolní končetiny a délky trupu žen. Mezi délkou dolní končetiny a délkou trupu mužů existuje velmi nízký stupeň přímé lineární závislosti 73 550 500 - 450 400 - 350 • muzi • zeny 0° CD© o" On o „ o ^o„ a o o o 0 c&o qgP9cP0° -1-1-1-1-1— 850 900 950 1000 1050 délka dolni končetiny (v mm) p, = 0.0598, p2 = 0.2853 1100 Obrázek 44: Krabicový diagram s 95% elipsou spolehlivosti pro délku dolní končetiny a délku trupu u žen (v mm) (pi = 0.0598). Mezi délkou dolní končetiny a délkou trupu žen existuje nízký stupeň přímé lineární závislosti (pi = 0.2853). 74 10.5 Test o rozdílu parametrů pi — P2 dvou alternativních rozdělení Nechť Xn,... Xitvj je náhodný výběr z alternativního rozdělení Alt(pi) a X21, ■ ■ ■ X2n2 Je na něm nezávislý náhodný výběr s alternativního rozdělení Alt(p2). Nechť po je konstanta. Na hladině významnosti a testujeme jednu z následujících tří hypotéz oproti příslušné alternativní hypotéze. #01 '■ Pi — P'ž = Po oproti Hn : P\ — P2 ^ Po (oboustranná alt.) #02 : Pi — P'ž < Po oproti H12 '■ Pi — P2 > Po (pravostranná alt.) Hqs : pi — P2 > po oproti h13 : Pi — P2 < Po (levostranná alt.) Test nazýváme dvouvýběrovým Z-testem o rozdílu pravděpodobností pi — P2. Testovací statistika má tvar M1-M2-po Zw = (10.6) /Mi(l-Mi) M2(l-M2) V íí! + jv2 kde Mi = -^r- Y^i^i Xu, kde Xu nabývá hodnot 0 nebo 1, je výběrový průměr prvního náhodného výběru, M2 = T^T ^i=i X-2i, kde X2Í nabývá hodnot 0 nebo 1, je výběrový průměr druhého náhodného výběru, N\ je rozsah prvního náhodného výběru, N2 je rozsah druhého náhodného výběru a po je konstanta z nulové hypotézy. Testovací statistika Zyy pochází asymptoticky ze standardizovaného normálního rozdělení, tj. zw= ; Mi-Ms-n, A /Mi(l-Mi) , M2(l-M2) V ífl + JÍ! To znamená, že pro malé rozsahy náhodných výběrů Ni a N2 nemusí mít rozdělení statistiky Zyy charakter standardizovaného normálního rozdělení. To je ale problém, neboť testování kritickým oborem, intervalem spolehlivosti i p-hodnotou je založeno na předpokladu, že Zyy má charakter standardizovaného normálního rozdělení. Proto pro malé rozsahy Ni či N2 náhodných výběrů nemůžeme statistiku Zyy k testování nulové hypotézy použít, neboť závěry takového testování by mohly být mylné. S rostoucími rozsahy náhodných výběrů Ni a N2 se však rozdělení statistiky Zyy čím dál více blíží ke standardizovanému normálnímu rozdělení, získává tak všechny jeho vlastnosti a závěry testování se stávají spolehlivě správnými. Zda máme dostatečný počet pozorování k provedení testu o rozdílu pravděpodobností pi — P2 prověříme podmínkami dobré aproximace, které mají obecně tvar iVipi(l-pi) >9 a N2p2{l-P2) > 9. (10.7) V případě, že skutečné hodnoty parametrů pi a p2 neznáme, nahrazujeme je jejich nestrannými odhady, tj. výběrovými průměry M\ a Mi. Podmínky dobré aproximace potom vypadají následujícím způsobem 7YiMi(l - Mi) > 9 a N2M2(1 - M2) > 9. (10.8) Pokud je podmínka 10.7 (resp. 10.8) dobré aproximace pro oba parametry pi i p2 splněna, můžeme test o rozdílu pravděpodobností pi — P2 použít, aniž bychom se vystavili riziku zavádějících výsledků testování. Pokud by však podmínka dobré aproximace pro jeden z parametrů nebo pro oba parametry pi a P2 nebyla splněna, nemůžeme test použít, dokud nerozšíříme jeden nebo oba datové soubory o další pozorování a nezvýšíme tak dostatečně rozsahy Ni i A^2 obou náhodných výběrů. Poznámka: Problému s nedostatkem pozorování při použití testu o rozdílu pravděpodobností pi — P2 se můžeme vyvarovat, pokud si ve fázi plánování experimentu spočítáme minimální potřebné rozsahy obou náhodných výběrů.Ve fázi sběru dat si potom již snadno ohlídáme, aby oba datové soubory obsahovaly potřebný počet pozorování, ideálně s nějakou rezervou. K výpočtu minimálního rozsahu náhodného výběru potřebujeme znát pouze předpokládanou hodnotu pravděpodobností pi a P2. 75 Kritický obor podle zvolené alternativní hypotézy má tvar Hu -Pi-P'ž^Po W = (-00 ; ua/2) U («i_a/2 ; 00) H12 : Pi - P2 > Po W = ; 00) Hiz ■ Pi - P2 < Po W = (-oo;ua) kde ua/2, Ui-a/2, ua, «i-a jsou kvantily standardizovaného normálního rozdělení, jejichž hodnoty získáme pomocí Ol a implementované funkce qnorm(). Interval spolehlivosti má podle zvolené alternativní hypotézy jeden z následujících tvarů Hn ■ Pi -P2^Po (d, h) = H12 ■ Pí ~P'2>P0 (d, 1) = H13 ■ Pí ~P'2 zw)} H12 : Pi - p2 > Po p-hodnota = Pr(Zw > zw) = 1 - Pr(Zw < zw) H13 : Pi - p2 < Po p-hodnota = Pr(Zw < zw) kde Zyy je náhodná veličina, zyy je realizace testovací statistiky Zyy (viz vzorec 10.6), tedy konkrétní číslo, a Pr(Zw < zw) je distribuční funkce standardizovaného normálního rozdělení, jejíž hodnotu získáme pomocí a implementované funkce pnormQ. 76 Příklad 10.17. Test o rozdílu parametrů P1-P2 Mějme k dispozici údaje o frekvenci výskytu epigenetického znaku sutura metopica (binomické proměnná) na lebkách jedinců z Anatolské populace (viz tabulka 1; (Eroglu, 2008)). Tabulka 1: Výskyt epigenetického znaku sutura metopica u jedinců Anatolské populace sutura metopica E muži 26 334 ženy 15 153 Na hladině významnosti a = 0.01 zjistěte, zda je frekvence výskytu epigenetického znaku sutura metopica u mužů Anatolské populace nižší než u žen Anatolské populace. Řešení příkladu 10.17 Nejprve se zorientujme v zadání příkladu. Ze zadání víme, že proměnná popisující frekvenci výskytu epigenetického znaku sutura metopica je binárního typu, tj. nabývá pouze hodnoty 1 (zástupný symbol pro úspěch, což je v našem případě výskyt epigenetického znaku sutura metopica) a hodnoty 0 (neúspěch, což je v našem případě absence epigenetického znaku sutura metopica). Náhodná veličina X popisující frekvenci výskytu epigenetického znaku sutura metopica u mužů tedy pochází z alternativního rozdělení, tj. X ~ Alt(pi), kde p\ je pravděpodobnost nastání úspěchu, tedy pravděpodobnost výskytu epigenetického znaku sutura metopica u mužů. Analogicky náhodná veličina Y popisující frekvenci výskytu epigenetického znaku sutura metopica u žen také pochází z alternativního rozdělení, tj. X ~ Alt(p2), kde p2 je pravděpodobnost výskytu epigenetického znaku sutura metopica u žen. Nejprve vytvoříme datový soubor pomocí příkazu data.frame(). 467 (data <- data.frame(sutura.metopica = c (26 , 15), suma = c (334 , 153), 468 row.names = c('muži', 'zeny'))) sutura.metopica suma 469 muzi 26 334 470 zeny 15 153 471 Z tabulky vidíme, že z celkového počtu 334 lebek mužů byl zaznamenán výskyt epigenetického znaku sutura metopica) na 26 lebkách. Dále z celkového počtu 153 lebek žen byl zaznamenán výskyt epigenetického znaku sutura metopica) na 15 lebkách. Rozsah náhodného výběru lebek mužů JVi = 334, rozsah náhodného výběru lebek žen A^2 = 153. Protože oba výběry pochází z alternativních rozdělení, použijeme pro ověření hypotézy ze zadání test o rozdílu parametrů p\ —p2- K použití tohoto testu je nejprve potřeba ověřit podmínku dobré aproximace zvlášť pro každý náhodný výběr. Jelikož přesné hodnoty parametrů p\ a p2 neznáme, musíme je v podmínce dobré aproximace nahradit jejich bodovými odhady, tedy hodnotami výběrových průměrů M\ a M2. iVi mi = -^j>i = +1 + 1 +1 + 1 +1 + 0 +■■■ + ()) Nx ^ 334 2=1 — = 0.07784431. 334 m2 = ^-É^ = 7^(1 + ■■■ + 1 + 0+ ■■■ + 0) N2 N2^~" 153" 2=. 0.09803922. 15 153 = 77 472 NI <- 334 473 N2 <- 153 474 xl <- 26 475 x2 <- 15 476 477 ml <- xl / NI # 0.07784431 478 m2 <- x2 / N2 # 0.09803922 Podmínky dobré aproximace potom dopočítáme následujícím způsobem iViMi(l - Mi) = 334 x 0.07784431 x (1 - 0.07784431) = 334 x 0.07784431 x 0.9221557 = 23.97605 479 NI * ml * (1 - ml) # 23.97605 Jelikož číslo 23.9761 > 5, je podmínka dobré aproximace pro muže splněna. N2M2(1 - M2) = 153 x 0.09803922 x (1 - 0.09803922) = 153 x 0.09803922 x 0.9019608 = 13.52941 480 N2 * m2 * (1 - m2) # 13.52941 Jelikož číslo 13.5294 > 5, je podmínka dobré aproximace pro ženy splněna. Obě podmínky dobré aproximace jsou splněny, můžeme tedy přistoupit k testování otázky ze zadání. Zde poznamenejme, že v zadání příkladu není zmínka o znění nulové hypotézy, pouze o záměru zjistit, zda je frekvence výskytu znaku sutura metopica u mužů nižší než u žen. Toto tedy bude znění alternativní hypotézy, zatímco nulovou hypotézu musíme vhodně dodefinovat. Testování nulové hypotézy potom provedeme v posloupnosti šesti kroků. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Frekvence výskytu epigenetického znaku sutura metopica u mužů je vyšší nebo rovna frekvenci znaku sutura metopica u žen. Hi : Frekvence výskytu epigenetického znaku sutura metopica u mužů je nižší než u žen. • matematická formulace nulové a alternativní hypotézy H0-Pi>P2 Pi -P2> Po, kde po =0 Hi : pi < p2 pi - p2 < po, kde p0 = 0 (levostranná alternativa) 2. Volba hladiny významnosti • Hladina významnosti a = 0.01 (viz zadání příkladu). 3. Testování kritickým oborem 78 • Testovací statistika (Mi - M2) - po Z, w Mi(l-Mi) , M2(l-M2) iVi + N2 (0.07784431 - 0.09803922) - 0 0.07784431(1-0.07784431) , 0.09803922(1-0.09803922) 334 + 153 -0.02019491 ^0.0002149239 + 0.0005779577 -0.02019491 0.02815815 -0.7171959 = -0.7172 481 alpha <- 0.01 482 pO <- 0 483 zw <- ((ml - m2) - pO) / (sqrt(ml * (1 - ml) / Nl + m2 * (1 - m2) / N2)) # -0.7171956 • Kritický obor W = (—oo ; ua) = (-oo; m0.oi) = (-oo; -2.3263) 484 qnorm(alpha) # -2.326348 • Závěr testování Protože realizace testovací statistiky zyy = —0.7172 nenáleží do kritického oboru, tj. zyy ^ W, Hq nezamítáme na hladině významnosti a = 0.01. 4. Testování intervalem spolehlivosti • Interval spolehlivosti mi(l-mi) m2(l - m2) (-1, h) = | -1; mi - m2 - ^/---+-—-uD ^^n^„^ /0.07784431(1 - 0.07784431) 0.09803922(1 - 0.09803922) -oo ; (0.07784431 - 0.09803922) - \j-K—-'- +-K—-%.oi ("I (-1 (-1 -0.02019491 - ^0.0002149239 + 0.0005779577 x (-2.326348)) -0.02019491 - 0.02815815 x (-2.326348)) 0.04531075) = (-1; 0.04531) 485 HH <- (ml - m2) - sqrt (ml * (1 - ml) / Nl + m2 * (1 - m2) / N2) * qnorm(alpha) # 0. 04531075 • Závěr testování Protože po = 0 náleží do Waldova 99% empirického pravostranného intervalu spolehlivosti, tj. p0 = 0 G IS, Hq nezamítáme na hladině významnosti a = 0.01. 5. Testování p-hodnotou 79 • p-hodnota p-hodnota = ~Pr(Zw < zw) = Pr(Zw < -0.7171956) = 0.2366267 = 0.2366 486 p.val <- pnorm(zw) # 0.2366267 • Závěr testování Protože p-hodnota = 0.2366 je větší než a = 0.01, Hq nezamítáme na hladině významnosti a = 0.01. 6. Grafická vizualizace výsledků testování Výsledek testování vizualizujeme pomocí sloupcového diagramu relativních četností (viz obrázek 45), který vygenerujeme příkazem rel.barplot.two(). Funkce rel.barplot.two() je implementována v RSkriptu Sbirka-AS-l-2018-funkce.R. Před použitím samotné funkce je tedy potřeba načíst RSkript příkazem sourceQ. 0.0 -1 1-1- muzi zeny pohlaví Obrázek 45: Sloupcový diagram relativních četností výskytu epigenetického znaku sutura metopica mužů a žen Anatolské populace 7. Interpretace výsledků: Na základě všech tří způsobů testování nezamítáme hypotézu Hq . Frekvence výskytu epigenetického znaku sutura metopica u mužů není statisticky významně nižší než u žen. ★ 80 Příklad 10.18. Test o rozdílu parametrů pi-p2 Mějme datový soubor 20-more-samples-probabilities-pubis.txt obsahující údaje o frekvenci výskytu tří stupňů změn kostního reliéfu (nepřítomné změny, malé změny, střední až výrazné změny) na vnitřní straně stydké kosti (os pubis) v blízkosti stydké spony sedací (symphysis publica) u žen z tří kosterních souborů (evropského původu, afrického původu a Inuitů). Pro více informací o datovém souboru viz sekci ??. Na hladině významnosti a = 0.05 testujte hypotézu, že frekvence výskytu středních změn na stydké kosti u žen afrického původu je větší nebo rovna frekvenci výskytu středních změn na stydké kosti u inuitských žen. Řešení příkladu 10.18 Datový soubor načteme příkazem read.delim() a odstraníme z něj NA hodnoty příkazem na.omit(). 487 data <- read.delim('00-Data/20-more-samples-probabilities-pubis.txt') 488 data <- na.omit(data) 489 head(data) origin absence trace .to.small moderate .to.large number.of cases 1 European 30 20 10 60 2 African 56 37 17 110 3 Inuits 16 6 13 35 490 491 492 493 Načtená tabulka obsahuje tři řádky, z nichž nás zajímají pouze řádky týkající se africké a inuitské populace (tj. řádek 2 a 3), a pět sloupců, z nichž nás zajímají pouze sloupce obsahující obsahující údaje o výskytu středních změn na vnitřní straně stydké kosti (trace.to.small) a celkový počet skeletů (number.of.cases). Pouze tyto údaje si vložíme do nové tabulky data.AI, a to pomocí operátoru [] a množinového operátoru %in%. 494 data. AI <- data [data$origin 7,in°/, c('African', 'Inuits'), c ('trace . to . small ' , ' number . of . cases')] 495 row.names(data.AI) <- c('African', 'Inuits') 496 data.AI trace .to.small number.of cases African 37 110 Inuits 6 35 497 498 499 Mějme nyní náhodnou veličinu X popisující frekvenci výskytu středních změn na vnitřní straně stydké kosti v populaci Afrických žen. Tato veličina nabývá hodnoty 1 (úspěch; výskyt středních změn na vnitřní straně stydké kosti) nebo hodnoty 0 (neúspěch; výskyt jiného typu změn (výrazná nebo žádná změna) na vnitřní straně stydké kosti). Předpokládáme tedy, že náhodná veličina X pochází z alternativního rozdělení, tj. X ~ Alt(pi), kde pi je výskytu středních změn na vnitřní straně stydké kosti v populaci afrických žen. Analogicky mějme náhodnou veličinu Y popisující frekvenci výskytu středních změn na vnitřní straně stydké kosti v populaci Inuitek. Tato náhodná veličina pochází rovněž z alternativního rozdělení, tj. Y ~ Alt(p2), kde p2 je pravděpodobnost výskytu středních změn na vnitřní straně stydké kosti v populaci Inuitek. Z tabulky data.AI vidíme, že z celkového počtu 110 skeletů afrických žen byl zaznamenán výskyt středních změn na vnitřní straně stydké kosti u 37 skeletů. Z celkového počtu 35 skeletů Inuitek byl zaznamenán výskyt středních změn na vnitřní straně stydké kosti u 6 skeletů. Rozsah náhodného výběru skeletů afrických žen JVi = 110, rozsah náhodného výběru skeletů Inuitek N2 = 35. Protože oba výběry pochází z alternativních rozdělení, použijeme pro ověření hypotézy ze zadání test o rozdílu parametrů pi — p2. Nejprve však ověříme podmínku dobré aproximace zvlášť pro každý náhodný výběr. Přesné hodnoty parametrů pi a p2 neznáme, nahradíme je proto v podmínce dobré aproximace jejich bodovými odhady, tedy hodnotami výběrových průměrů M\ a M2. 1 Nl 37 mi = — Y2Xi = Tiö = °-3363636- 81 1 ^2 m2 = — ^Xl = ^ = 0.1714286. _6_ N2^f"° 35 =i Podmínky dobré aproximace potom dopočítáme následujícím způsobem iViMi(l - Mi) = 110 x 0.3363636 x (1 - 0.3363636) = 110 x 0.3363636 x 0.6636364 = 24.55454 = 24.5545 500 NI * ml * (1 - ml) # 24.55455 Jelikož číslo 24.5545 > 5, je podmínka dobré aproximace pro populaci afrických žen splněna. N2M2(1 - M2) = 35 x 0.1714286 x (1 - 0.1714286) = 35 x 0.1714286 x 0.8285714 = 4.971429 = 4.9714 501 N2 * m2 * (1 - m2) # 4.971429 Jelikož číslo 4.9714 < 5, není podmínka dobré aproximace pro populaci Inuitských žen splněna. Jelikož podmínka dobré aproximace pro populaci Inuitských žen není splněna, nemůžeme otestovat nulovu hypotézu zadanou v zadání příkladu pomocí test o rozdílu parametrů p\—p2. Připomeňme, že test o rozdílu parametrů pi —p2 je asymptotickým testem, tj. testem, jehož výsledky jsou spolehlivé pouze pro dostatečné rozsahy obou náhodných výběrů, které nám zajišťuje právě splnění obou podmínek dobré aproximace. ★ Poznámka: Při analýze reálných dat bychom měli několik možností, jak se se situací podobnou situaci nastalé v příkladu 10.18 vypořádat. První možností je zamyslet se, zda by bylo možné doplnit datový soubor o další údaje. Doplnění údajů není v příkladě 10.18 preferováno, neboť data byla nasbírána v rámci studie (Stewart, 1970) a nasbírání dalších hodnot v podobné kvalitě na podobných skeletech by nebylo snadné. Druhou možností by bylo zvolit na otestování hypotézy ze zadání jiný test, který není tak striktní na splnění podmínek dobré aproximace. O takových testech si povíme více v kapitole ??. V praxi bychom však pravděpodobně volili postup, v rámci kterého bychom do striktně analytického postupu vložili špetku lidskosti, zamhouřili obě oči nad mírným porušením podmínky dobré aproximace a po řádném okomentování a zdůvodnění bychom test o rozdílu parametrů pi — p2 na analýzu nulové hypotézy ze zadání příkladu 10.18 použili. Je však nezbytné pracovat s tímto postupem velmi opatrně, a při závažnějším porušení byť jen jedné podmínky dobré aproximace volit raději jiné, vhodnější metody, nebo zauvažovat o rozšíření datového souboru. V případě příkladu 10.18 jsme zůstali z důvodu větší názornosti striktní, test nulové hypotézy ze zadání jsme neprovedli a příklad jsme hned po zjištění porušení podmínky dobré aproximace ukončili. Závěrem poznamenejme, že této nepříjemné situaci bychom se vyvarovali, kdybychom v rámci fáze plánování experimentu na základě podmínek dobré aproximace stanovili minimální rozsahy obou náhodných výběrů, a to analogicky, jako tomu bylo v příkladech ?? a ?? v kapitole ??. 82 Příklad 10.19. Test o rozdílu parametrů p\-p2 Mějme datový soubor 14-two-samples-probabilities-sexratio.txt. Nechť binární proměnná sex obsahuje údaje o pohlaví novorozenců a binární proměnná o.sib.N obsahuje údaje o počtu starších sourozenců novorozence (viz sekce ??). Na hladině významnosti a = 0.01 testujte hypotézu, že pravděpodobnost narození chlapce je stejná u novorozenců s žádným sourozencem a u novorozenců s jedním sourozencem. Řešení příkladu 10.19 Datový soubor načteme příkazem read.delimQ a odstraníme z něj NA hodnoty příkazem na.omitQ. 502 data <- read.delim('00-Data/14-two-samples-probabilities-sexratio.txt') 503 data <- na.omit(data) 504 head(data) sex o.sib.N 1 m 0 2 m 0 3 f 0 4 m 0 5 m 0 6 m 1 505 506 507 508 509 510 511 Načtená tabulka obsahuje dva sloupce. Sloupec sex, popisující pohlaví novorozence, disponuje dvěma variantami (m = mužské pohlaví a f = ženské pohlaví). Sloupec o.sib.N, popisující počet starších sourozenců, disponuje též dvěma variantami (0 = žádný satarší sourozenec, 1 = jeden starší sourozenec). Datový soubor nyní rozdělíme pomocí operátoru [] na dva subsoubory. První bude obsahovat pouze údaje týkající se novorozenců s žádným starším sourozencem (dataO) a druhý bude obsahovat údaje týkající se pouze novorozenců s jedním starším sourozencem (datal). Mějme nyní náhodnou veličinu X popisující frekvenci narození novorozenců mužského pohlaví s žádným 512 dataO <- data[data$o.sib.N == 0, ] 513 datal <- data[data$o.sib.N == 1, ] starším sourozencem. Tato veličina nabývá hodnoty 1 (úspěch; narození chlapce) nebo hodnoty 0 (neúspěch; narození děvčete).Předpokládáme tedy, že náhodná veličina X pochází z alternativního rozdělení, tj. X ~ Alt(pi), kde pi je pravděpodobnost nastání úspěchu, neboli pravděpodobnost narození chlapce u populace novorozenců s žádným starším sourozencem. Analogicky mějme náhodnou veličinu Y popisující frekvenci narození novorozenců mužského pohlaví s žádným starším sourozencem. Tato náhodná veličina pochází rovněž z alternativního rozdělení, tj. Y ~ Alt(p2), kde p2 je pravděpodobnost narození chlapce v populaci novorozenců s jedním starším sourozencem. Nyní zjistíme počet narozených chlapců v populaci novorozenců s žádným a s jedním starším sourozencem a rozsahy obou náhodných výběrů. Z tabulky vidíme, že z celkového počtu 595 novorozenců s žádným starším sourozencem 514 xl <- sum(data0$sex == 'm') 515 x2 <- sum (dat al $ sex == 'm') 516 NI <- length(data0$sex) 517 N2 <- length(datal$sex) 518 tab <- data.frame(x = c(xl, x2), N = c(Nl, N2)) bylo 310 novorozenců mužského pohlaví a že z celkového počtu 518 novorozenců s jedním starším sourozencem bylo 277 novorozenců mužského pohlaví. Rozsah náhodného novorozenců s žádným sourozencem Ni = 595, rozsah náhodného novorozenců s jedním starším sourozencem = 518. Protože oba výběry pochází z alternativních rozdělení, použijeme pro ověření hypotézy ze zadání test o rozdílu parametrů pi —p2- Nejprve však ověříme podmínku dobré aproximace zvlášť pro každý náhodný výběr. Protože přesné hodnoty parametrů pi a p2 neznáme, nahradíme je v podmínce dobré aproximace jejich bodovými odhady, tedy hodnotami výběrových průměrů M\ a Mi. 83 1 Nl 310 mi = — V x, =-= 0.5210084. 1 277 m2 = — V x, =-= 0.534749. =i 519 ml <- xl / NI # 0.5210084 520 m2 <- x2 / N2 # 0.534749 Podmínky dobré aproximace potom dopočítáme následujícím způsobem 7YiMi(l - Mi) = 595 x 0.5210084 x (1 - 0.5210084) = 595 x 0.5210084 x 0.4789916 = 148.4874 521 NI * ml * (1 - ml) # 5.8o4348 Jelikož číslo 148.4874 > 5, je podmínka dobré aproximace pro populaci novorozenců s žádným starším sourozencem splněna. N2M2(1 - M2) = 518 x 0.534749 x (1 - 0.534749) = 518 x 0.534749 x 0.465251 = 128.8745 522 N2 * m2 * (1 - m2) # 128.8745 Jelikož číslo 128.8745 > 5, je podmínka dobré aproximace populaci novorozenců s jendím starším sourozencem splněna. Obě podmínky dobré aproximace jsou splněny, můžeme tedy přistoupit k testování otázky ze zadání. V zadání máme uvedeno přesné znění nulové hypotézy, zbývá stanovit alternativní hypotézu tak, aby byla doplňkem k hypotéze nulové. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : Pravděpodobnost narození chlapce je shodná pro populaci novorozenců s žádným starším sourozencem a pro populaci s jedním starším sourozencem. Hi : Pravděpodobnost narození chlapce není shodná pro populaci novorozenců s žádným starším sourozencem a pro populaci s jedním starším sourozencem. • matematická formulace nulové a alternativní hypotézy H0 : pi = p2 pi-p2= p0, kde p0 = 0 Hi : pi ^ P2 Pi - P2 7^ Po, kde p0 = 0 (oboustranná alternativa) 2. Volba hladiny významnosti • Hladina významnosti a = 0.01 (viz zadání příkladu). 84 3. Testování kritickým oborem • Testovací statistika Zv (Mi - M2) - po Mi(l-Mi) , M2(l-m) iVi + N2 (0.5210084 - 0.534749) - 0 0.5210084(1-0.5210084) 0.534749(1-0.534749 595 -0.0137406 + 518 ^0.0004194263 + 0.0004802944 -0.0137406 0.02999534 -0.4580912 -0.45809 523 alpha <- 0.01 524 pO <- 0 525 zw <- ((ml - m2) - pO) / (sqrt(ml * (1 - ml) / NI + m2 * (1 - m2) / N2)) # -0.4580921 • Kritický obor W = (-00 ; ua/2) U zw)} = 2mm{Pľ(Zw < zw), 1 - Pr(Zw < zw)} = 2mm{Pľ(Zw < -0.4580921), 1 - Pľ(Zw < -0.4580921)} = 2 min{0.3234431, 0.6765569} = 2 x 0.32344310 = 0.6468862 = 0.6469 530 p.val <- 2 * min(pnorm(zw), 1 - pnorm(zw)) # 0.6^68863 • Závěr testování Protože p-hodnota = 0.6469 je větší než a = 0.01, Hq nezamítáme na hladině významnosti a = 0.01. 6. Grafická vizualizace výsledků testování , Výsledek testování doložíme sloupcovým diagramem relativních četností, jenž vykreslíme pomocí funkce rel.barplot.twoQ implementované v RSkriptu Sbirka-AS-l-2018-funkce.R (viz obrázek 46). 1.0 -i 0.8 1 0.6 > 'M 0.4 0.2 -0.0 - 52.1 % 53.47 % AI.9 % 46.53 % ■ muzske pohlaví ■ zenske pohlaví zadny jeden počet starších sourozenců Obrázek 46: Sloupcový diagram relativních četností narození chlapce v populaci novorozenců s žádným starším sourozencem a s jedním starším sourozencem 7. Interpretace výsledků: Na základě všech tří způsobů testování nezamítáme hypotézu Hq. Mezi pravděpodobností narození chlapce v populaci novorozenců s žádným starším sourozencem a s jedním starším sourozencem neexistuje statisticky významný rozdíl. 86 ★ 87 Příklad 10.20. Test o rozdílu parametrů P1-P2 Mějme datový soubor 26-two-samples-probabilities-palmar.txt obsahující údaje o vysokém, středním, nízkém a jiném zakončení dlaňových linií na pravé a levé straně 50 žen z mechské populace a 87 žen z populace Rajbanshi (viz sekce ?? a tabulky 2 a 3). Na hladině významnosti a = 0.05 zjistěte, zdaje frekvence výskytu nízkého zakončení tří dlaňových linií na levé straně u žen z mechské populace vyšší než frekvence výskytu nízkého zakončení tří dlaňových linií na levé straně u žen z populace Rajbanshi. Tabulka 2: Frekvence zakončení dlaňových linií u žen z mechské populace f-mech 1 r n Hi 5 6 11 Mi 16 15 31 Lo 19 22 41 others 10 7 17 Tabulka 3: Frekvence zakončení dlaňových linií u žen z populace Rajbanshi f-raj 1 r n Hi 27 37 64 Mi 17 26 43 Lo 21 10 31 others 22 14 36 Řešení příkladu 10.20 Z obou tabulek si vytáhneme pouze relevantní údaje, a sice počet výskytů nízkého zakončení tří dlaňových linií na levé straně u žen mechské populace (x\ = 19) a počet výskytů nízkého zakončení tří dlaňových linií na levé straně u žen z populace Rajbanshi (x2 = 21). Ze zadání příkladu dále víme, že rozsah náhodného výběru žen z mechské populace JVi = 50 a rozsah náhodného výběru žen z populace Rajbanshi = 87. 531 xl <- 19 532 x2 <- 21 533 NI <- 50 534 N2 <- 87 535 tab <- data.frame (x = c (xl , x2) , N = c (NI , N2) ) Naším úkolem je porovnat frekvenci výskytu nízkého zakončení tří dlaňových linií na levé straně u dvou indických populací. U mechské populace máme k dispozici počet úspěchů, tj. počet výskytů nízkého zakončení tří dlaňových linií na levé straně. Náhodná veličina X popisující výskyt nízkého zakončení tří dlaňových linií na levé straně u žen z mechské populace pochází z alternativního rozdělení, tj. X ~ Alt(pi), kde pi je pravděpodobnost nastání úspěchu, tedy pravděpodobnost výskytu nízkého zakončení tří dlaňových linií na levé straně. Podobně u populace žen z Rajbanshi máme k dispozici počet výskytů nízkého zakončení tří dlaňových linií na levé straně. Náhodná veličina Y popisující výskyt nízkého zakončení tří dlaňových linií na levé straně u žen z populace Rajbanshi pochází z alternativního rozdělení, tj. X ~ Alt(p2), kde pi je pravděpodobnost výskytu nízkého zakončení tří dlaňových linií na levé straně. Řešení příkladu vede na situaci, kdy pravděpodobnost pl porovnáváme s pravděpodobností p2, tedy na dvouvýběrový test o rozdílu pi —pi-Před použitím tetsu je potřeba ověřit pro každý náhodný výběr zvlášť podmínku dobré aproximace. Protože přesné hodnoty parametrů pi a p2 neznáme, nahradíme je v podmínce dobré aproximace jejich bodovými odhady, tedy hodnotami výběrových průměrů M\ a M2. mi 1 Nl — Y 19 50 0.38. 88 1 N'2 21 ^a;2 = — = 0.2413793. 7 = 1 536 ml <- xl / NI # 0.38 537 m2 <- x2 / N2 # 0.2413793 Podmínky dobré aproximace potom dopočítáme následujícím způsobem /ViMi(l - Mi) = 50 x 0.38 x (1 - 0.38) = 50 x 0.38 x 0.62 = 11.78 538 NI * ml * (1 - ml) # 11.78 Jelikož číslo 11.78 > 5, je podmínka dobré aproximace pro ženy z mechské populace splněna. N2M2(1 - M2) = 87 x 0.2413793 x (1 - 0.2413793) = 87 x 0.2413793 x 0.7586207 = 15.93103 539 N2 * m2 * (1 - m2) # 15.93103 Jelikož číslo 15.93103 > 5, je podmínka dobré aproximace pro ženy z populace Rajbanshi splněna. Obě podmínky dobré aproximace jsou splněny, můžeme se tedy zaměřit na otázku v zadání. Naším úkolem je zjistit, zda je frekvence výskytu nízkého zakončení tří dlaňových linií na levé straně u žen z mechské populace vyšší než frekvence výskytu nízkého zakončení tří dlaňových linií na levé straně u žen z populace Rajbanshi. Protože v zadání není zmínka o nulové hypotéze, je tvrzení z předchozí věty zněním alternativní hypotézy. Zbývá stanovit nulovou hypotézu tak, aby byla doplňkem k alternativní hypotéze. 1. Stanovení hypotéz • slovní formulace nulové a alternativní hypotézy Hq : frekvence výskytu nízkého zakončení tří dlaňových linií na levé straně u žen z mechské populace je menší nebo rovna frekveni výskytu nízkého zakončení tří dlaňových linií na levé straně u žen z populace Rajbanshi. Hi : Frekvence výskytu nízkého zakončení tří dlaňových linií na levé straně u žen z mechské populace je vyšší než frekvence výskytu nízkého zakončení tří dlaňových linií na levé straně u žen z populace Rajbanshi. • matematická formulace nulové a alternativní hypotézy H0 ■■ Pi < P2 Pi ~ Vi < Po, kde p0 = 0 Hi : pi > p2 pi - p2 > po, kde p0 = 0 (pravostranná alternativa) 2. Volba hladiny významnosti • Hladinu významnosti volíme v souladu se zadáním a = 0.05. 3. Testování kritickým oborem 89 • Testovací statistika = (Mi -M2) -po /Mi(l-Mi) , M2(l-M2) V JVi ~r iV2 _ (0.38 - 0.2413793) - 0 /0.38(1-0.38) 0.2413793(1-0.2413793) V 50 + 87 _ 0.1386207 ~ V0.004712 + 0.002104774: _ 0.1386207 ~ ^0.006816774 _ 0.1386207 ~ 0.08256376 = 1.678953 = 1.6790 540 alpha <- 0.05 541 pO <- 0 542 zw <- ((ml - m2) - pO) / (sqrt(ml * (1 - ml) / Nl + m2 * (1 - m2) / N2)) # 1.678953 • Kritický obor W = (tíl-a ; oo) = («1-0.05; oo) = («0.95 ; oo) = (1.6449; oo) 543 qnorm(l - alpha) # 1.644854 • Závěr testování Protože realizace testovací statistiky zyy = 1.6790 náleží do kritického oboru, tj. zyy e W, Hq zamítáme na hladině významnosti a = 0.05. 4. Testování intervalem spolehlivosti • Interval spolehlivosti , /mi(l -mi) m2(l - m2) (d, h) = [ mi - m2 - \ j-—--1--—-«i-a/2 ; 1 0.38 - 0.2413793 - ,/ + 0-2413793(1 - 0.2413793)^^ . 1 (0.1386207 - ^0.004712 + 0.002104774 x 1.644854; l) (0.1386207 - 0.08256376 x 1.644854; 1) (0.0028153; 1) 90 544 dh <- (ml - m2) - sqrt (ml * (1 - ml) / Nl + m2 * (1 - m2) / N2) * qnorm(l - alpha) # 0. 002815394 • Závěr testování Protože po = 0 nenáleží do Waldova 95% empirického levostranného intervalu spolehlivosti, tj. po = 0 ^ IS, Hq zamítáme na hladině významnosti a = 0.05. 5. Testování p-hodnotou • p-hodnota p-hodnota = 2Vt{Zw > zw) = 1 - Pr(Zw < zw) = Pr(Zw < 1.678953) = 0.04658058 = 0.04658 545 p.val <- 1 - pnorm(zw) # 0.04658058 • Závěr testování Protože p-hodnota = 0.04658 je menší než 0.05, Hq zamítáme na hladině významnosti a = 0.05. 6. Grafická vizualizace výsledků testování Rozdíly mezi oběma populacemi vizualizujeme sloupcovým diagramem relativních četností (viz obrázek 47). 546 source('Sbirka-AS-I-2018 -funkce . R ' ) 547 data <- cbind(c(19, Nl - 19), c(21, N2 - 21)) 548 par(mar = c(5, 4, 1, 1), family = 'Times') 549 rel.barplot.two(data, col = c('orange2', 'darkseagreen3'), 550 border = 'grey40 ' , 551 names = c('Mech', 'Rajbanshi ' ) , 552 legend = c('nizke zakončeni', ' jine zakončeni'), 553 density = 60, col.number = 'black', main = '', 554 xlab = 'populace', ylab = 'relativní četnost') 7. Interpretace výsledků: Na základě všech tří způsobů testování zamítáme hypotézu Hq. Frekvence výskytu nízkého zakončení dlaňových linií na levé straně u žen z mechské populace je statisticky významně vyšší než frekvence výskytu nízkého zakončení dlaňových linií na levé straně u žen z populace Raj banshi. ★ 91 1.0 -i 0.8 0.6 0.4 0.2 - 62% 75.86 % ■ nizke zakončeni □ jine zakončeni 0.0 -1 1-1- Mech Rajbanshi populace Obrázek 47: Sloupcový diagram relativních četností zakončení dlaňových linií z levé strany u žen z mechské populace a žen z populace v Rajbanshi 92