Zadání příkladů - cvičení č.l - 15-9-23 Příklad č.l (porovnání dvou typů modelů) (přednáška) Model rozděleni pravděpodobností je modelem náhodné proměnné X, např. (1) model rozdělení pravděpodobnosti náhodné proměnné X šířka dolní čelisti, nebo (2) model rozdělení pravděpodobnosti náhodné proměnné X hrubost kožních řas u dospělých zdravých žen. Statistický model je modelem náhodné proměnné Y\X (Y kauzálně závisí na X), např. (1) model závislosti náhodné proměnné Y šířka dolní čelisti na proměnné X pohlaví, nebo (2) model závislosti náhodné proměnné Y hrubost kožních řas u dospělých zdravých žen na proměnné X BMI. Všimněme si, že náhodné proměnné označujeme X anebo Y podle toho, jaký model je charakterizuje. Příklad č.2 (jednoduchý náhodný výběr) V jednoduchém náhodném výběru o rozsahu n z populace s konečným rozsahem N má každý prvek stejnou pravděpodobnost vybrání. Pokud vybíráme bez vracení (opakování), mluvíme o jednoduchém náhodném výběru bez vraceni (Dalgaard, 2008). Pokud vybíráme s vracením, mluvíme o jednoduchém náhodném výběru s vracením. Mějme množinu M. s N = 10 prvky a chceme z ní vybrat n = 3 prvky (a) bez vracení, (b) s vracením. Kolik máme možností? Jak vypadá jedna takováto možnost, pokud M = {1,2,..., 10}? Zopakujte to samé pro N = 100, n = 30 a množinu M = {1,2,..., 100}. Příklad č.3 (jednoduchý náhodný výběr) Mějme skupinu lidí označených identifikačními čísly (ID) od 1 do 30. Vyberte (a) náhodně 5 lidí z 30-ti bez návratu, (b) náhodně 5 lidí ze 30-ti s návratem a nakonec (c) náhodně 5 lidí ze 30-ti bez návratu, přičemž lidé s ID od 28-mi do 30-ti mají pravděpodobnost vybrání 4x vyšší než lidé s ID od 1 do 27. Příklad č.4 (normální rozdělení) Mějme náhodnou proměnnou X (může to být např. výška postavy desetiletých dívek) a předpokládejme, že tato náhodná proměnná má normální rozdělení s parametry fi (střední hodnota) a a2 (rozptyl), což zapisujeme jako X ~ N(fi,a2), fi = 140.83, a2 = 33.79. Normální rozdělení představuje model rozdělení pravděpodobnosti pro tuto náhodnou proměnnou. Vypočítejte pravděpodobnost Pr(a < X < b) = Pľ(X < b) - PľX < a) = Fx(b) - Fx(a), kde a = - ka, b = + ka, k = 1, 2, 3. Nakreslete hustotu rozdělení pravděpodobnosti, vybarvěte oblast mezi body a a b a popište osy x a y tak, jako je uvedeno na obrázku 1. 120 130 140 150 160 170 120 130 140 150 160 170 120 130 140 150 160 170 vyska (cm) vyska (cm) vyska (cm) Pr(135.02 O, j = 1,2, p G (—1,1) jsou parametry. Potom 6 = (/ii,/12, o"2, o"|, p). Výraz v exponentu můžeme zapsat jako _ 1 / x - /i! \T / of po-i0-2 \ 1 / x _ /ii \ 2 V Ž/ - P2 y V P^i^ of / \V-fJ-2 )' Marginální rozdělení1 jsou X ~ N (pi, o^) a "K ~ a?" (/í2, čt|), p je koeficient korelace2(Viz obrázek 5) Příklad č.15 (dvojrozměrné normální rozdělení) (1) Nakreslete hustotu dvojrozměrného normálního rozdělení a^a*, E) pomocí funkce image() a superponujte ho s konturovým grafem hustoty toho stejného rozdělení pomocí funkce contour(). (2) Nakreslete hustotu dvojrozměrného normálního rozdělení N2(n, S) pomocí funkce persp(). Hustotu rozsekejte na 12 intervalů, kde hodnoty v těchto intervalech budou odpovídat barvám terrain.colors(12). Použijte následující parametry: • Pl = 0, P2 = o, 01 = 1,02 = 1. »P = 0; • Pl = 0, P2 = o, 01 = 1,02 = 1. >P = 0.5; • Pl = 0, P2 = 0, 0"! = 1-2, 02 = 1,P = 0.5 Vzorové řešení je uvedeno na obrázku 5. 1 Marginální rozdělení je rozdělení náhodné proměnné, zde X nezávisle na y a naopak Y nezávisle na X. 2Z tohoto příkladu je zřejmé, že na dostatečný popis dvojrozměrného normálního rozdělení potřebujeme pět parametrů, t.j. střední hodnotu a rozptyl pro marginální rozdělení náhodných proměnných X a Y a korelační koeficient p = p(X,Y) popisující sílu lineárního vztahu X &Y. 6 H! = O, n2 = O, o, = 1, o2 = 1, p = O Hi = O, n2 = O, o, = 1, o2 = 1, p = 0.5 H, = 0, n2 = 0, o, = 1, o2 = 1.2, p = 0.5 Obrázek 5: Hustoty dvojrozměrného normálního rozdělení při různých parametrech (první řádek -konturový graf; druhý řádek - perspektivní trojrozměrný graf v podobě plochy); čím je p odlišnější od nuly, tím více se kontury liší od kruhů (mění se na elipsy); se zvyšujícím se rozdílem mez o\ a o"2 se zvětšuje rozdíl rozptýlení koncentrických kruhů ve směru jednotlivých os (říkáme, že rozdíl variability proměnných laľse zvětšuje.) Příklad č.17 (standardizované normální rozdělení) Náhodný vektor (X, Y)T má dvojrozměrné normální rozdělení N2 (0, S), kde 0 = (0, 0)T aS - ' 9 s hustotou 0 (x, y) = f (x, y) =-- exp p 1 x2 — 2pxy + y2 kde (x,y)T G IR2, p G (—1,1) jsou parametry, potom 0 = (0,0,1, l,p). Výraz v exponentu můžeme psát jako í ŕ x \ i 1 p \ (x 2\yJ \ p 1/ \y marginální rozdělení jsou obě N(0,1) a p je koeficient korelace. 7 Příklad č. 18 (standardizované normální rozdělení) Nechť náhodnou proměnnou X ~ N(pi, af) je největší výška mozkovny (skuli.pH; v mm) a náhodnou proměnnou Y ~ N(p2,a2) je morfologická výška tváře (face.H; v mm). Nechť X a Y mají dvojrozměrné normální rozdělení s parametry (pi, p2)T a a\ a P Jsou parametry kovarianční matice S. Když od náhodné proměnné X odpočítáme její střední hodnotu p\ a tento rozdíl podělíme odmocninou z rozptylu (o"i), dostaneme náhodnou proměnnou Zx, která má asymptoticky normální rozdělení se střední hodnotou p\ = 0 a rozptylem o\ = 1, což zapisujeme jako Zx ~ N(0,1). Pokud od náhodné proměnné Y odečteme její střední hodnotu p2 a tento rozdíl podělíme odmocninou z rozptylu (o"2), dostaneme náhodnou proměnnou Zy, která má asymptoticky normální rozdělení se střední hodnotou /i2 = 0a rozptylem o\ = 1, což zapisujeme jako Zy ~ ÍV(0,1). Potom (Zx, Zy)t má standardizované dvourozměrné normální rozdělení N2(fi, S) s parametry [i = (0,0)T a o\ = 1, o"2 = 1 a p jsou parametry kovarianční matice S. Příklad č.19 (dvourozměrné normální rozdělení) Simulaci pseudonáhodných čísel z N2(fi, S) můžeme v R vytvořit následujícími způsoby: 1. použitím funkce mvrnorm() z knihovny MASS; 2. použitím funkce rmvnorm() z knihovny mvtnorm 3. použitím funkce rnorm() a následujícího algoritmu: Nechť X1 ~ ÍV(0,1) a X2 ~ ÍV(0,1); potom (Yi, F2)T ~ iV2(/x, S), kde fi = (fiu p2)T je vektor středních hodnot a o\ a a| a p jsou parametry kovarianční matice S, přičemž síla lineárního vztahu Yi a Y2 je daná velikostí a znaménkem p; = aiXi + p± &Y2 = a2(pXi + a/1 — p2X2) + P2- Nasimulujte pseudonáhodná čísla Yi a Y2 z ^(a*, S). Vypočítejte dvourozměrný jádrový odhad hustoty (Yi, Y"2)T pomocí funkce kde2d(). Nakreslete jej také pomocí funkce image() a superponujte jej konturovým grafem hustoty dvourozměrného normálního rozdělení N2(n,H) pomocí funkce contour(). Hustotu rozsekejte na 12 intervalů, kde hodnoty v těchto intervalech budou odpovídat barvám terrain.colors(12). Při simulaci použijte následující parametry: (a) pi = 0, ^2 = 0, 01 = 1, 0-2 = 1, p = 0; (1) n = 50, (2) n = 500 (b) pi = 0, ^2 = 0, o"i = 1, 0-2 = 1, p = 0.5; (1) n = 50, (2) n = 500 (c) pi = 0, ^2 = 0, O"! = 1, a2 = 1.2, p = 0.5; (1) n = 50, (2) n = 500 Vzorové řešení viz obrázek 6. 8 Simulace pseudonahodnych cisel z N2(u., I) _funkce mvrnorm; N = 300_ i-1-1-1-1-1-1-r -3-2-10 1 2 3 (i! = 0, [i2 = 0, a-i = 1, o2 = 1, p = 0 Simulace pseudonahodnych cisel z N2(u., I) funkce rmvnorm; N = 300 0 0 o V \X° ° CM - o / y i 7 offff iV\r^\ 00 \\\\ \ \ ° ID olgl o A Jtí^ vil A jo J o o • o v O o Tr^o.o2—■ 0 ° 0 ^ - -3-2-10 1 2 3 Hi = 0, n2 = 0, o, = 1, o2 = 1, p = 0 Simulace pseudonahodnych cisel z N2(^i, E) funkce rnorm; N = 300 Hi =0, n2 = 0, o, = 1, o2 = 1, p = 0 Simulace pseudonahodnych cisel funkce mvrnorm; N = 300 Simulace pseudonahodnych cisel funkce rmvnorm; N = 300 Simulace pseudonahodnych cisel funkce rnorm; N = 300 Hi = 0, p-2 = 0, o, = 1, o2 = 1, p = 0 Hi = 0, p-2 = 0, o, = 1, o2 = 1, p = 0 -3-2-10 1 2 Hi =0, n2 = 0, o, = 1, o2 = 1, p = 0 Obrázek 6: Hustoty dvourozměrného normálního rozdělení Příklad č.23 (binomické rozdělení, binomický experiment) Experiment sestávající z fixního počtu Bernoulliho experimentů (ozn. N) se nazývá binomický experiment. Pravděpodobnost úspěchu označme p, pravděpodobnost neúspěchu q = 1 — p. Náhodná proměnná X je počet pozorovaných úspěchů po dobu experimentu. Pravděpodobnost X = x za podmínky, že X pochází z binomického rozdělení Bin(N,p), píšeme jako Pr(X = x) x px(l-p) n-x x = 0,l,...,N (1) (Ugarte a kol. 2008). Střední hodnota E[X] = Np a rozptyl Var[X] = Np(l — p). Naprogramujte a zobrazte v R pravděpodobnostní funkci a (kumulativní) distribuční funkci pro 5m(5,0.5). Řešení viz obrázek 7. 9 Pravdepodobnosti funkce binomického rozděleni Bin(5,0.5) Distribuční funkce binomického rozděleni Bin(5,0.5) ~i-r 2 3 Obrázek 7: Pravděpodobnostní a distribuční funkce binomického rozdělení Bin(5, 0.5) Příklad č.26 (Poissonovo rozdělení; počet havárií za týden) Pokud každý z 50 milionů lidí řídí v Itálii řídí auto následující týden nezávisle, potom pravděpodobnost smrti při autonehodě bude 0.000002, kde počet úmrtí má binomické rozdělení Bin(50rmZ, 0.000002) anebo limitní Poissonovo rozdělení s parametrem A = 50rraZ x 0.000002 = 100. Příklad č.27 (Poissonovo rozdělení; pruské armádní jednotky) Nechť početnosti úmrtí X jako následek kopnutí koněm v Pruských armádních jednotkách (Bort-kiewicz, 1898) mají Poissonovo rozdělení s parametrem A, tj. X ~ Poiss(X). Pravděpodobnost, že někdo bude smrtelně zraněný v daném dni, je extrémně malá. Mějme 10 vojenských jednotek za 20-letou periodu s rozsahem M = 200 (200 = 10 x 20), kde, při početnostech úmrtí n = 1, 2, 3,4, 5+ v dané jednotce a v daném roce, zaznamenáváme také početnosti vojenských jednotek mn při daném n, kde M = ^mn (viz tabulka). Vypočítejte očekávané početnosti, za předpokladu X ~ Poiss(X), kde A = (2) n 0 1 2 3 4 5+ mn 109 65 22 3 1 0 Příklad č.28 (podíl chlapců a dívek v rodinách) Nechť X představuje početnost chlapců mezi dětmi v rodinách. Zde můžeme předpokládat, že X Bin(N,p), tj. rodina může mít vychýlený poměr pohlaví dětí ve směru k chlapcům nebo k dívkám. V realitě tedy můžeme mít velmi mnoho rodin jen s chlapci nebo jen s děvčaty a nemáme dostatek rodin s poměrem pohlaví blízkým 51 : 49 (poměr chlapců ku dívkám). Z toho nám vyplývá, že rozptyl početnosti chlapců bude ve skutečnosti větší než rozptyl předpokládaný binomickým rozdělením Bin(n, P). Příklad č.29 (overdispersion v binomickém modelu) V klasické studii poměru pohlaví u lidí z roku 1889 na základě záznamů z nemocnic v Sasku (více informací viz Lindsey a Altham, (1998)) zaznamenal Geissler (1889) rozdělení počtu chlapců v rodinách. Mezi M = 6115 rodinami s N = 12 dětmi pozoroval následující početnosti chlapců (n jsou početnosti chlapců a mn početnosti rodin s n chlapci). 10 n 0 1 2 3 4 5 6 7 8 9 10 11 12 mn 3 24 104 286 670 1033 1343 1112 829 478 181 45 7 Vypočítejte mn za předpokladu, že početnosti chlapců X v rodinách mají binomické rozdělení s parametry 7T En=0 iVM 0.5192 (3) a N = 12, ozn. X ~ Bin(N,n). Příklad č.30 (overdispersion v Poissonově modelu) Mějme početnosti úrazů n mezi dělníky v továrně, kde početnosti dělníků mn při daném n (viz tabulka) (Greenwood a Yule (1920)). n 0 1 2 3 4 >5 mn I 447 132 42 21 3 2 Vypočítejte očekávané početnosti dělníků za předpokladu, že početnosti úrazů na dělníka X mají Poissonovo rozdělení s parametrem ^ nm„ A = 4^-- = 0.47. (4) Ozn. X ~ Poiss(A). Příklad č.31 (binomické rozdělení, simulační studie) Vygenerujte pseudonáhodná čísla X (početnosti úspěchů) opakovaná M-krát (M = 1000) z Bin(N, p), kde iV = 5 a p = 0.5. Vytvořte tabulku vygenerovaných (simulovaných) i teoretických relativních početností (pro n = 0,1,... ,5). Superponujte histogram vygenerovaných pseudonáhodných čísel s teoretickou pravděpodobnostní funkcí (viz obrázek 8). Histogram vygenerovaných pseudonáhodných čísel superponovaný teoretickou AT ^\ Obrázek í pravděpodobnostní funkcí Bin(N,p). 11 Příklad č.32 (binomické vs normální rozdělení) Nechť Xn ~ Bin(N,p), potom můžeme aproximovat binomické rozdělení normálním následovně: XN ~ N(Np, Np(l - p)), kde také platí ZN= fN~Np ~JV(Q,1). Ukažte, že CLV platí pro N = 100 a p = 0.5 na tři desetinná místa. Příklad č.33 (normální rozdělení, simulační studie) Na základě simulační studie prověřte, že pokud X ~ iV(150,6.25), potom [X]n ~ iV(150, 6.25/n). Použijte n = 30. Pro každou simulaci X vypočítejte aritmetické průměry xm, m = 1,2,..., M, kde M = 500 000. Superponujte je histogramem v relativní škále s teoretickou křivkou hustoty pro Xn. Vypočítejte Pr(Xn > 151) ze simulovaných dat a porovnejte tento výsledek s teoretickou (očekávanou) pravděpodobností. Řešení viz obrázek 9. i-1-1-1-1 i-1-1-1 i-1-1-1-1 146 148 150 152 154 148 149 150 151 149.0 149.5 150.0 150.5 151.0 průměry průměry průměry n-5 n-30 n-100 Obrázek 9: Histogram vygenerovaných průměrů superponovaný teoretickou křivkou hustoty Xn. Příklad č.34 (normální rozdělení, simulační studie) — — 2 2 Nechť X ~ N(fii,af) a Y ~ ÍV(/í2, o"!). Potom Xni — Yn2 ~ N(fii — fi2, ^ + ^J). Generujte pseu-donáhodná čásla X a Y rozdělení N(fij,<7j), j = 1,2, kde fii = 100, o\ = 10, fi2 = 50, 02 = 9 při (a) rii = 4, n2 = 5, (b) n\ = 100, n2 = 81. Pro každou simulaci XaY vypočítejte rozdíl xm — ym, m = 1,2,... M, kde M = 1000. Superponujte histogram těchto rozdílů v relativní škále s teoretickou křivkou hustoty rozdílu Xni — Ýri2. Pro případ (a) i (b) vypočítejte Pr(Vni — Ýri2 < 52) na základě empirického (vygenerovaného) a teoretického rozdělení Xni — Yn2. 12 o co _. I-1-1-1-1 ° I-1-1-1-1 30 40 50 60 70 46 48 50 52 54 rozdil průměru rozdil průměru ni=4, n2 = 5 ni = 100, n2 = 81 Obrázek 10: Histogram vygenerovaných rozdílů průměrů superponovaný teoretickou křivkou hustoty rozdělení rozdílu výběrových aritmetických průměrů Příklad č.35 (statistika) Mějme náhodný výběr (X\, X2, ■ ■ ■, Xn)T, kde Xi G IR, i = 1, 2,..., n, potom příklady statistik jsou: • 7i = £ľ=i*i 9. Pro každou simulaci X vypočítejte z\y,m, m = 1,2,..., M, kde M = 1000. Superponujte histogram vygenerovaných testovacích statistik v relativní škále s teoretickou křivkou hustoty Zyy- 13 realizace statistiky Zw N= 10 , p= 0.1 , Hp= 0.9 2 3 4 5 -3 -2 realizace statistiky Zw N= 10 , p=0.5 , Hp= 2.5 -5 -A realizace statistiky Zw N= 10 , p= 0.9 , Hp=0.9 -3 -2 realizace statistiky Zw N= 100 , p=0.1 , Hp= 9 2 3 N= 100 , p=0.5 , Hp= 25 realizace statistiky Zw N= 100 , p= 0.9 , Hp=9 Obrázek 11: Histogram vygenerovaných testovacích statistik v relativní škále superponovaný s teoretickými křivkami hustoty. Příklad č.36 mluví o použití jednovýběrové testovací statistiky pro parametr binomického rozdělení (pravděpodobnost) pro různé pravděpodobnosti a různé početnosti. Pokud není Haldova odmínka splněná,není možné testovací statistiku použít. Příklad č.37 (testovací statistika, simulační studie) Na základě simulační studie prověřte, že pokud (a) X ~ N(fi,a2), kde fi = 0, a2 = 1 a (b) X ~ [(1 — p)N(fi,a2) + pN(fi,af)], kde fi = 0, a2 = 1, p = 0.05, o\ = 2, potom testovací statistika F = ^n^S má asymptoticky Xn-i rozdělení o n — 1 stupních volnosti. Použijte rozsahy náhodných výběrů n = 15 a n = 100. Pro každou simulaci X vypočítejte FpoZ)m, m = 1,2,..., M, kde M = 1000. Superponujte histogram vygenerovaných testovacích statistik v relativní škále s teoretickou křivkou hustoty F. 14 X~N(0,1) X~N(0,1) I-1-1-1-1 o |-1-1-1-1-1 O 10 20 30 40 60 80 100 120 140 160 realizace testovací statistiky F realizace testovací statistiky F n=15 n=100 Obrázek 12: Histogram vygenerovaných testovacích statistik v relativní škále superponovaný s teoretickými křivkami hustoty ÍV(0,1). X~(1 -p)N(0,1 )+p*N(0,2) X~(1 -p)N(0,1 )+p*N(0,2) 200 realizace testovací statistiky F realizace testovací statistiky F n=15 n=100 Obrázek 13: Histogram vygenerovaných testovacích statistik v relativní škále superponovaný s teoretickými křivkami hustoty (1 — p)N(0,1) + pN(0,2). 15