Bayesovské metody (seminár) September 13, 2016 Doporučená literatúra: Anděl, J., Matematická statistika, SNTL/ALFA, Praha, 1985. Grendár, M., Bayesovská statistika (http://www.savbb.sk/ grendar/pdf) Hušková, M., BAYESOVSKÉ METODY, UNIVERZITA KARLOVA, Praha, 1985. (http://www.karlin.mif.cuni.cz/ huskova) Pázman, A., BAYESOVSKÁ STATISTIKA, Univerzita Komenského, Bratislava, 2003. Tento text v zásade sleduje knižku A.Pázmana. 1. Úvod a motivačný príklad Definícia 1.1. Majme pravdepodobnostný priestor (íl, A, P). Náhodné javy A\, A2,... G A tvoria úplný systém javov, ak platí 00 (1.1) Ai n A3■ = 0, i ŕ j, a (J At = íi. i=l Poznámka. Úplný systém javov môže byt aj konečný. Veta 1.1. (Vzorec pre úplnú pravděpodobnost). Nech A\, A2,... je úplný systém javov v pravdepodob-nostnom priestore (íl, A, P) taký, že (1.2) P(AA>0, i = 1,2,.... Potom platí 00 (1.3) P(B) = YJP{B\AAP(AX). i=l Dôkaz: P(B) = P(Bílíl) = P (B íl U=! Ai) = P ({JZi(B n A)) = £~i P(B n Ai) = £~i P(B\Al)P(Al). * Veta 1.2. (1. Bayesov vzorec). Nech A\,A'2,... je úplný systém javov v pravdepodobnostnom priestore (íl, A, P) taký, že P(Ai)>0, i =1,2,.... 1 2 Ak P (B) > O, tak platí (14) P(A \E) PWAi)P(Ai) 7-l2 ( ' ( 3' ' ~ T,r=i P(B\AAP(AX) ' 3 ~ lj 2" Dôkaz: Pre lubovolné j je Veta 1.3. (2. Bayesov vzorec). Nech Ai,A'2,... je úplný systém javov v pravdepodobnostnom priestore (íl, A, P) taký, že P(At) > 0, i = 1,2,... . ďalej A e A, že P (Ä) > 0 a B e A. Platí „,DI ^ P(AnA<)>0}P(Ai)P(A|Ai)P(g|AnAi) ( } ( 1 }" EŠľ^W^W • Dôkaz: spravte si sami. Poznámka. Vety 1.1, 1.2 a 1.3 platia aj v prípade, že úplný systém javov je konečný. Poznámka. P (Aj) v Bayesových vzorcoch sú tzv. apriórne pravdepodobnosti a P(Aj\B) aposteriórne pravdepodobnosti (po vykonaní pokusu s výsledkom B). Poznámka. V prípade 1. Bayesovho vzorca ide o riešenie situácie, ked máme hypotézy A\,..., ktoré sa navzájom vylučujú, ale vyčerpávajú všetky možnosti. Poznáme ich (apriórne) pravdepodobnosti P(Ai). Nastal jav A a poznáme pravdepodobnosti P(A\Ai). Pýtame sa na (aposteriórne; nové, ktoré berú do úvahy skutočnost, že mastal A) pravdepodobnosti P(Ai\A) V prípade 2. Bayesovho vzorca ak nastal jav A, pýtame sa na pravděpodobnost javu B. Poznámka. Nieje vždy jednoduché volit správny pravdepodobnostný model pre výpočet podmienených pravdepodobností. Príklad 1.1. (Lekárska diagnostika). Vieme, že určitou (konkrétnou) chorobou Ch trpí 1% populácie. Choroba je diagnostikovaná na základe vyšetrenia, ktorého spolahlivost je (i) 95% ak vyšetrovaná osoba trpí chorobou Ch (ii) 70 % ak vyšetrovaná osoba netrpí chorobou Ch. Vyšetrujeme náhodne zvolenú osobu. Určte pravděpodobnost správnej diagnózy, ak výsledok vyšetrenia je (a) pozitívny (podlá výsledku vyšetrenia je osoba chorá) (b) negatívny (podlá výsledku vyšetrenia je osoba zdravá). Riešenie: Označme jav A - vyšetrovaná osoba trpí chorobou Ch (je chorá) B - výsledok vyšetrovania je pozitívny Zo zadania vieme P (A) = 0.01 (pravděpodobnost, že vybraná osoba je chorá) Táto pravděpodobnost sa volá prevalencia alebo tiež apriórna pravděpodobnost choroby Vyšetrenie (spolahlivost vyšetrenia) sa charakterizuje dvomi charakteristikami, a síce pravdepodobnostou P(B\A) = 0.95 tzv. citlivost testu alebo aj senzitivita testu pravdepodobnostou P(B\A) = 0.7 tzv. špecificita testu, (a) Máme určit vlastne P(A\B) (lebo v tomto prípade výsledok testu bol pozitívny, teda test hovorí, že vyšetrovaná osoba je chorá (diagnóza je, že pacient je chorý) a my máme určit pravděpodobnost správnej dignózy). 3 Zo zadania vieme, žeP(A) = 0.01, P (A) = 0.99, P(B\A) = 0.95 aP(B\A) = 1-P(B\A) = 1-0.7 = 0.3. Podia Bayesovho vzorca (A, A sú hypotézy) P(A\B) =_P(A)P(B\A)_ =_0.01-0.95_= V 1 ' P(A)P(B\A) + P(A)P(B\A) 0.01-0.95 + 0.99-0.3 Je to aj aposteriórna pravděpodobnost, že pacient je chorý, ak výsledok testu bol pozitívny. Je to prekvapivý výsledok, čakali by sme "omnoho lepší" výsledok. Celkom máme 29 700 + 950 = 30 650 pozitívnych výsledkov, z toho správne pozitívnych je 950, čiže 950 P(A\B) =-= 0.030995. v 1 ' 30650 _ (b) Analogicky (zase A, A sú hypotézy) P(Ä\B) = _ P(Ä)P(B\Ä) _ =_0^7_= V 1 ' P(A)P(B\A) + P(A)P(B\A) 0.99-0.7+ 0.01-0.05 Je to aposteriórna pravděpodobnost, že pacient nie je chorý, ak výsledok testu bol negatívny. Naozaj celkovo máme 69 300 + 50 = 69 350 negatívnych výsledkov, z toho správne negatívnych je 69 300 a teda __ 69300 pravděpodobnost správnej diagnózy u negatívnych výsledkov testu je P(A\B) = - = 0.99928. 69350 2. Bayesova veta Uveďme si základné tvrdenia. Veta 2.1. (O podmienenej hustote). Nech A a v sú a—konečné miery a nech združená hustota náhodných vektorov Y € lZr a Z € lZn~r vzhladom k jj, = A x v je p(y, z). Nech (2.1) q(y)= / p(y,z)du(z) je marginálna hustota náhodného vektora Y. Potom podmienená hustota vektora Z pri pevnom Y = y je rovná ,00, , , , Jp(y|z)/9(y), ak n) = 1 (Bm je systém borelovských podmnožin priestoru lZm). Poznáme a) podmienenú hustotu pravdepodobnosti vektora Y vzhladom k a—konečnej miere A(y), za predpokladu, že je daná hodnota 0 náhodného vektora 0 (2-4) /(y|0), b) hustotu pravdepodobnosti náhodného vektora 0 vzhladom k a—konečnej miere v{ff) (2.5) vr(0). Podlá Vety 2.1 združená hustota pravdepodobnosti (©', Y')' vzhladom k /i = v x A sa rovná (2.6) fl(0,y) = T(0)/(y|0), marginálna hustota pravdepodobnosti veličiny Y (predikčná hustota) je /*(y)= / ir(o)f(y\o)dv(o) Ju a podmienená hustota 0, za predpokladu, že je známa hodnota Y = y, sa rovná g(0,y) vr(0)/(y|0) (2.7) w(9\y) /* (y) In^o)f(y\e)du(ey Pre dáta y máme pravdepodobnostný model, t.j. predpokladáme, "že sa riadia" hustotou f(y\0) (v prípade spojitých údajov) alebo pravdepodobnostnou funkciou (v prípade diskrétnych údajov). Obyčajne hustotu zapisujeme ako f(y;0), ale v tomto prípade hustotu zapisujeme ako f(y\0), aby sme zdôraznili, že pravdepodobnostný model je hustotou pre dáta za predpokladu, že hodnota parametra je 0. Ďalej predpokladajme, že sme schopní vyjádřit "náš názor" (naše presvedčenie) týkajúce sa parametra 0 vo forme apriórnej hustoty tt(O). Túto hustotu musíme vediet "skonštruovat" separátne (inou cestou, inak) než z nameraných dát. Poznámka. Ak Y má (absolútne) spojité rozdelenie (A je Lebesgueova miera) a 0 tiež absolútne spojité rozdelenie, tak tt(9), f(y\0) a ir(0\y) sú "obyčajné" hustoty. Ak Y má spojité rozdelenie a 0 diskrétne (ľ(9) je sčítacia miera), 0 G {0\, O2,...}, tak tt(0) je pravde-podobnostná funkcia, f(y\9) je hustota a Efc>i7r(é'fc)/(ylé'fc) je pravdepodobnostná funkcia. Ak Y má diskrétne rozdelenie a 0 spojité, tak tt(O) je hustota, f(y\0) je pravdepodobnstná funkcia a 7r(0|y) je hustota. Ak Y je diskrétny náhodný vektor, Y G {yi, 1/2,...} a 0 tiež diskrétny náhodný vektor, 0 G {0\, Oi, •••}, tak vr(0), /(y|0) aj A^\y,) = ^-fa ' lfl v i = 1,2, sú pravděpodobnostně funkcie. 5 Ešte iný pohlad na Bayesov vzorec. V prípade diskrétnych náhodných premenných 8 a ľ nech A\,Am je úplný systém borelovských množín na reálnej osi a B\,Bn nejaké borelovské množiny. Nech sú známe podmienené pravdepodobnosti P(Bi\Ak) = P (Y eBj|G6 Ak) a tzv. apriórne pravdepodobnosti P (Ak) = P (Q G Ak). P (Ak) vyjadrujú "očakávania" javov Ak ešte pred pozorovaním javov Bi. Vplyvom pozorovania javu Bi sa pravděpodobnost "očakávania" P (Ak) zmení na podmienenú pravděpodobnost P(Ak\Bi) určenú Bayesovým vzorcom (1.4) (alebo (2.7)) P(Bl\Ak)P(Ak) Vzorce (2.4)-(2.7) samozrejme platia ak 0 a Y sú náhodné veličiny. 3. Príklad (podia Anděl, str. 279) Výroba prebieha každý deň. pravděpodobnost vyrobenia chybnej súčiastky je p. Predpokladajme, že táto pravděpodobnost sa zo dňa na deň (mierne) mení a rozdelenie pravdepodobnosti veličiny p sa dá aproximovat beta rozdelením s hustotou (3.1) ff(p) = ^_p»-i(i_p)"-i) 0 0, b > 0 sú parametre, B (a, b) = xa~1(l — x)b~1dx je beta funkcia. Keď výroba prebieha "ustálené", môžeme "dostatočne presne" odhadnut hodnoty parametrov a tieto odhady považovat za skutočné hodnoty a a b. Ak teraz spravíme pokus - náhodne vyberieme n výrobkov (s vrátením), pričom m z nich bude chybných, pravděpodobnost vybratia práve m chybných súčiastok (ak je pravděpodobnost vyrobenia chybnej súčiastky rovná p) je (3.2) S(m\p) = {^jpm(l-p)n~m, m = 0,1,..., n. Podlá (2.7) je aposteriórna hustota 7r(p|m) náhodného parametra p pri zistenom počte m chybných súčiastok z n vybraných súčiastok rovná Ap)f(y\p) _ TOjP^d-p^^a-p)"-7" (3.3) 7r(p|m) J-00^(p)f(y\p)dp ^1^-1(1 -p)b-1(^)pm(l - p)"-™dp 1 pa+m-l^_p^b+n-m-l^ 0 0+, b —> 0+, tak n*(p) n**(p) = p_1(l - p)-1, 0 0+, b —> 0+ dostávame (4.1) Tr(plm) -> 7r**(p|m) = —----p^1™"1 (1 - p)^™"1™"1. B (m, n — m) Tento výsledok by sme dostali aj priamo dosadením funkcie tt** za apriórnu hustotu. Ak teda platí 0 < m < n, je tt**(p\m) "obyčajná" hustota, ak ale m = 0 alebo m = n, nedostaneme konečnú limitu funkcie 7r(p|m). Nezápornú meratelnú funkciu voláme nevlastná hustota. Ďalšie vlastnosti Bayesovho vzorca. Ak pre každé dve rôzne yi,y2 v (2.4) platí, že podiel j^1] J, nezávisí od 0, tak ir(0\yi) = i"(0|y2) pre každé 0. Z bayesovského hladiska sú výsledky experimentov y1 a y2 "informačne ekvivalentné". V bayesovskom prístupe dôležitú úlohu má funkcia vierohodnosti (4.2) ly[B) cc f(y\0). Je to trieda všetkých reálnych funkcií premennej 0, líšiacich sa od f(y\0) len o multiplikatívny člen nezávislý od 0. Z Bayesovho vzorca vyplýva, že štatistické postupy založené na aposteroórnej pravdepodobnosti spĺňajú tzv. princíp vierohodnosti: Ak y, y* sú dva výsledky rôznych experimentov s tým istým parametrickým priestorom í~ž a ak ly(9) = ly'(9), tak výsledky y a y* vedú k tým istým štatistickým záverom o parametri 0. 5. Pravidlo retazenia pre nezávislé pozorovania Uvažujme dva nezávislé experimenty s tým istým parametrickým priestorom í~ž. Experiment I je daný systémom podmienených hustôt (5.1) {/(y|0): eefi}, v ktorom pozorujeme y G 3^ a experiment II nezávislý na experimente I, daný systémom podmienených hustôt {g(z\o): eefi}, v ktorom pozorujeme z e Z. V kombinovanom experimente II, v ktorom pozorujeme (nezávislú) dvojicu (y, z) GjxZ, daným systémom podmienených hustôt {/i(y,z|0) = /(y|ŕ%(z|0): 9eí]} pri apriórnej hustote tt(O) dostávame z (2.7) aposteriórnu hustotu 7r(0)/(y|ŕ%(z|0) (5.2) Jnn(0)f(y\0)g(z\0)du(0Y 7 Teraz realizujme experiment II, ale za apriórnu hustotu neuvažujme 7r(0), ale aposteriórnu hustotu z experimentu I, teda AO)f(y\e) <0\y) Podlá (2.7) dostávame aposteriórnu hustotu br(0|y)]fl(z|0) fnn(6)f(y\6)dv(6y (5.3) Jn[ir(0\y)]g(z\0)dv(0) lfnn(0)f(y\0)dv(0) 1(1 AO)f(y\0) 1 g(z\6)dv(6) fnn(0)f(y\0)dv(0)\ 7r(0)/(y|0)g(z|0) Ja*(0)f(y\0)g(z\0)dv(Oy Hustoty (5.2) a (5.3) sú zhodné. Toto pravidlo môžeme rozšířit na lubovolnú postupnost nezávislých experimentov. Aposteriórnu hustotu získanú z prvých k experimentov možno použit ako apriórnu hustotu v ďalšom (k + 1) — vom experimente. 6. Bayesov vzorec a postačujúce štatistiky Nech Y 6 J G S„, 8 € íí € Bm a r je meratelné zobrazenie r : lZn —> lZk. Takéto zobrazenie voláme štatistika. (Poznamenávame len, že ak m = k, štatistika je odhadom.) Štatistiku voláme postačujúcou ak podmienené hustoty pre y, za podmienky, že t (y) = t, nezávisia od 0. Z učebníc matematickej štatistiky (pozri napr. Anděl, str. 262) je známe, že štatistika r je postačujúca práve vtedy, ak existujú reálne (nezáporné a meratelné) funkcie h(y), q(t,0) také, že (6.1) S(y\e) = h(y)q(T(y),6) pre každé y, 0 (veta o faktorizácii). Pripomíname len, že ak Y je spojitá náhodná veličina (alebo náhodný vektor), tak f{y\0) je hustota a ak Y je diskrétna náhodná veličina (alebo náhodný vektor), tak f{y\0) je pravdepodobnostná funkcia. Označme 7r*(0|t)t=1-(y) aposteriórnu hustotu, ak bola pozorovaná hodnota t štatistiky t. Veta 6.1. (Veta o postačujúcej štatistike). Ak t (y) je postačujúca štatistika, tak platí vr(0|y) = 7r*(0|t)t=x(y). Dôkaz: Ak 0 aj Y sú diskrétne náhodné veličiny (alebo náhodné vektory), a r(y) je postačujúca štatistika, tak aposteriórna hustota (vlastne pravdepodobnostná funkcia) ((i 2] (0M = h(y)q(T(y),0)n(0) = q(r(y), 0)vr(0) [ ' J J2lh(y)q(T(y),9M9) E* sMy), ©M©)" Presnejšie v tomto prípade realizácie náhodnej veličiny (alebo náhodného vektora) 0 sú {0i,02, •••} a realizácie náhodnej veličiny (alebo náhodného vektora) Y sú {yi,y2, •••} Platí niO,) = P(& = Oj), f(yi\Oj) = P(Y = y2 0 = 9 j), HyMriy^OjMOj) _ g(t(yi), 0X0,-) E»i MyOXyO, 0X0,0 E,->i sMyJ, 0X0,0' 8 Pre diskrétne náhodné veličiny (resp. náhodné vektory) dostaneme analogicky ako vo vete 2.1, že združená pravdepodobnostná funkcia g(0j,yi) = 7r(0j)/(yi|0J). Aposteriórna hustota náhodnej veličiny 0 (vlastne pravdepodobnostná funkcia), ak bola pozorovaná hodnota t štatistiky r je ^(ej-|t)t=T(y)=P{e = ej-|Y = yi: T(y2) = t}, j = 1,2,.... Ak označíme náhodné udalosti Aj-{® = ej}, j = i,2,..., Bl-{y = yl: r(yj = t}, i =1,2,..., tak |t)t=T(y) = P{0 = 0,1 Y = y, : r(yi) = t} = P(^| U 2*) = = = = Ey<; T(y<)=t -P(Q = ^ Y = yj) = Ey<; T(y<)=t[-P(e = = ^ ® = ^)] Ey,: T(yi)=t p(y = Vi) Ey,: T(yi)=t{Es>i n® = os)p(y = Yl\& = es)} Ey<; T(y<)=t WO/fol^)] E^^JE^i^)/^)}' Ak t je postačujúca štatistika, tak podlá (6.1) vr*(6»J|t)t=x(y) Ey,: T(yi)=t{Es>i AOsMvM^), es)} Es>i Aos)q(t, o„)' Zo vztahu (6.2) dostávame, že ak r je postačujúca štatistika a t je pozorovaná hodnota tejto štatistiky, tak aposteriórna pravdepodobnostná funkcia je <0j\yi- T(yJ = t) teda vetu sme v prípade diskrétnych náhodných veličín (alebo vektorov) 0 a Y dokázali. Vo všeobecnom prípade nájdeme dôkaz napr. v knihe A.Pázmana. Jft V nasledujúcich kapitolách si zavedieme mieru množstva informácie a ukážeme, že vo všeobecnosti hustota 7r*(0|t) obsahuje menej informácie o neznámom parametri 0 než hustota 7r(0|y). 7. I-divergencia a informácia získaná z experimentu Bayesovské ponímanie znamená, že pred experimentom je všetka informácia o parametri 0 zhrnutá v apriórnej hustote tt(0) a úplná informácia o parametri 0 po experimente je zhrnutá v aposteriórnej hustote 7r(0|y). Niekedy sa vyžaduje charakterizovat množstvo tejto informácie jediným číslom (podobne ako charakterizovat rozdelenie pravdepodobnosti jednou číselnou charakteristikou). Táto úloha samozrejme je diskutabilná a nie je jednoznačná. Jedna z možností je vyjádřit "ekonomicky užitočnú" informáciu podobne ako v kapitole 23. Cielom môže byt aj vyjadrenie "presnosti", s ktorou hustota tt(0) alebo ir(0\y) určujú hodnotu parametra 0 v danom experimente. Tu ide o miery variability rozdelenia pravdepodobnosti náhodnej premennej 0. Môžeme merat aj nedostatok informácie obsiahnutej v tt(0) nejakou "mierou neurčitosti" apriórneho rozdelenia. Jedna z takýchto možných mier vznikla vo fyzike a neskôr sa osvedčila v rôznych vedných oblastiach. 9 Je to entropia. Budeme sa v ďalšom zaoberat aj vyjadrením stredného množstva informácie, ktoré možno získat z experimentu. Vyjadríme ho pomocou strednej I-divergencii. I-divergencia meria "odlišnost" dvoch distribúcií. 8. Entropia Definícia 8.1. Nech 0 € íl £ Bm je diskrétna náhodná veličina, ktorá nadobúda hodnoty {61,62, •••} s pravděpodobnostami {pi,P2, •••}, teda P(& = 6i) = pi i = 1, 2,.... Pravdepodobnosná funkcia {Oí,Pí}í, i = 1, 2,... určuje rozdelenie pravdepodobnosti P. Entropia takto daného rozdelenia pravdepodobnosti P sa rovná Ent(P) = - y^pilnpj. i Ak množina hodnôt {61,62, ■■■} je konečná (obsahuje n hodnôt), tak entropia je z uzavretého intervalu (0,lnn). Minimálnu hodnotu (rovnú 0) nadobúda, ked pre jeden index Íq je pio = 1 a pre ostatné indexy i 7^ io Je Pi = 0- V tomto prípade s istotou nastáva jediná z hodnôt, a síce 8io, a nie je neurčitost. Maximálnu hodnotu (rovnú ln n) nadobúda entropia, ak pi = pre všetky i. Vtedy je neurčitost, ktorá z hodnôt {6\, 62,6n} nastane, najväčšia. Kladieme OlnO = 0, lebo limz^0+ (zmz) = 0- Definícia 8.2. Nech 0 € í~ž € Bm je spojitá náhodná veličina, ktorej hustota (voči Lebesgueovej miere) je f (6). Potom entropia (spojitého) rozdelenia pravdepodobnosti P určeného hustotou /(•) je Ent(P) = - f ln f (6) f (6)d6 = -£ (ln/(©)). Jsi Podotýkame len, že ent(&) má iné vlastnosti ako entropia diskrétneho rozdelenia Ent(&) (napr. nadobúda hodnoty z intervalu (—00,00)). 9. I-divergencia Definícia 9.1. Nech 0 G lZm je diskrétna náhodná veličina, ktorá nadobúda hodnoty {6\, 62, ■■■} s pravde-podobnostami {p\,p2, teda P(& = 6i) = pi i = 1,2,... a 0* nadobúda tie isté hodnoty {61,62, ■■■} s pravděpodobnost ami {qi, í/2, •••}, teda Q(&* = 6i) = qi i = 1, 2,... . P a Q sú teda dve rozdelenia pravdepodobnosti, pričom pre tieto dve rozdelenia platí, že ak qi > 0, tak aj pi > 0 (miera Q je absolútne spojitá vzhladom k miere P). I-divergencia I(P, Q) týchto mier je (9.1) /(P;Q) = ^Píin^ ak pi ln < 00, inak položíme I(P, Q) = 00. Kladieme 0 ln ^ = 0. Definícia 9.2. Nech 8 £ í] € Bm je spojitá náhodná veličina s hustotami fp(6) a Jq(6) (voči Lebesgueovej miere), pričom pravdepodobnostná miera Q je absolútne spojitá voči miere P (t.j. ak pre nejakú borelovskú množinu A platí P(A) = 0, tak aj Q(A) = 0). I-divergencia I(P, Q) týchto mier je (9.2) I(P,Q) = Jjnj^fP(6)d6 = £P ak Jn ln Íf\q\ fp(6)dd existuje, inak I(P, Q) = 00. ln f p (Q) /q(®) 10 Poznámka. Niekedy sa namiesto I(P, Q) píše I(fp, f q). Platí tvrdenie, že vždy I(P, Q) > 0 a I(P, Q) = 0 práve vtedy ak P = Q (pozri knihu A.Pázmana, str. 24). I(P, Q) nie je metrika (vzdialenost medzi rozdeleniami), lebo I(P,Q) ^ I(Q,P). Miera P má vo výraze pre I(P,Q) dominantné postavenie, lebo vzhladom na ňu sa berie stredná hodnota. Príklad 9.1. Majme n—rozmerný náhodný vektor X a jeho dve rôzne normálne rozdelené rozdelenia pravdepodobnosti P a Q dané hustotami /i(x) = -^-e4(x " Mj'S-^x - Mi) t {p Q} V ' (2ir)ídeíÍ(Si) Spočítajte I(fP, fQ). Riešenie. /p(x) 1 det(So) 1 1/ w/v-,-1/ 7^(x) = 2 det(Sp) ~2^~^p'p <-x ~ + 2 (x ~ ~ vMq - Mp)) £q (x - /iP - - /ip)) = = 2 det(Sp) " 2(X " Mp) p (X " Mp) + 2(X " Mp) « (X " Mp)+ + (X - /ip/E^/ip - //q) + Í(/lp - ^g/E^/lp - //q). Počítajme /(/p(x)>/0(x)) = fp[ln^]=/wnIn^/P(x)dK = = ^Jln^||^/p(x)dx-^ i(x-Mp)'S-1(x-Mp)/p(x)dx+^ i(x-Mp)'S51(x-Mp)/p(x)dx+ + / (x-Mpľ^ÍMp -Mg)/p(x)dx+ / i(/iP -/^'E^/ip -/iQ)/p(x)eŕx = = fp {^ln SIIpT} " fp {^(X " ^)'SpX(X - Mp)} +fp {i(X - Mp)'S^(X - Mp)} + +5p {(X-zip/Eg^/ip -mq)| +5p |^(Mp -/íq/Eq^Mp - Mq)| = 1, det(So) ti lr„ i„ , 1. i. = 2 ln ďellšpT " 2 + ír 2{S^} + ä(Mp " EQ " V prípade, že Sp = Sq = S (rozdelenia sa líšia len v strednej hodnote), dostávame 2I(P, Q) = (Mp - Mq/E-Vp - /xQ). Príklad 9.2. Majme diskrétnu náhodnú premennú X € {0,1, ...,n} a jej dve rôzne binomické rozdelenia pravdepodobnosti P a Q dané hustotami (pravdepodobnostnými funkciami) Spočítajte I(P,Q). Riešenie. fj(x)=r)e?(l-6j)n-x, x = 0,l,...,n, je{P,Q} m-—— = a;ln--h (n — x m-. Sqíp) eQ v ' i-eQ 11 preto podia (9.1) n I(P,Q) = Y, , d p 1 - 8 p x ln---h (n — x) ln-— "p{l - 6P)n-x = n9p ln + n(l - 6P) ln--^. 10. Súvis I-divergencie s Fisherovou informačnou maticou (Spojitý prípad.) Majme štatistický experiment s výberovým priestorom Y G y £ B„ a hustotami {/(•|0): eeneBm} (voči Lebesgueovej miere). Ak pre malé A G í] sa hustoty f(-\9) a f(-\9 + A) značne odlišujú, tak na základe experimentu budeme môct dobre rozlišovat 9 od 9 + A. Preto budeme môct v tomto prípade "presne" odhadovat parameter 9. Odlišnost f{-\9) od f(-\9 + A) je vhodné merat pomocou I-divergencie (ío.i) /e,e+A = J]n /(y|g^A)/(ylg)rfA(y)- (Diskrétny prípad.) Ak máme štatistický experiment s výberovým priestorom Y € {yi,y2, •••} a hustotami (pravdepodobnostnými funkciami) {/(•|0): eeneBm}, tak odlišnost f(-\9) od f(-\9 + A) je vhodné merat pomocou I-divergencie (io.2) w = E^/(y^fA)/(y^)- Veta 10.1. (Súvis I-divergencie a Fisherovej informačnej matice.) Ak pre skoro všetky y G }" G B„ je f(y\9) dva razy spojité diferencovatelnou funkciou a ak možno zaměnit poradie derivovania a integrovania vo výraze í-f(y\0)d\(y), y i potom platí 2/e,e+A = A'm(0)A + o(||A||2), kde m(0) je Fisherova informačná matica s prvkami =)2 i -,Am) a o(||A||2) je taká funkcia ||A||2, že liniA^o —M A , 2 = 0- o(\\A\\2) |A| Dôkaz: Dôkaz spravíme pre spojitý prípad. Pre diskrétnu náhodnú veličinu Y sa spraví analogicky. Pomocou Taylorovej vety dostávame (10.3) ln/(y|0 + A) = ln/(y|0) + ""^'"'A + ^A'" ľľ^'"yA + ô(| |A| |2). 91n/(y|e) 1 y92ln/(y|e) 89' 2 ÔBidOj Pretože (10.4) / AJ m9)dX(y) = í ± d-^l^f{y\9)dX{y) 12 z (10.1) pomocou (10.3) a (10.4) dostávame 2W = 2^1n-^L_/(y|^A(y) = = "2 Jyd^^Af(y\e)dX(y)-JyA'd2^ /(y|e)dA(y) = = A'm(6>)A + o(||A||2). * 11. Stredné množstvo informácie z experimentu Bayesovská odpoveď na otázku, kolko informácie poskytuje daný experiment, vyplýva z porovnania apriórnej a aposteriórnej hustoty. V experimente danom výberovým priestorom 3^ a meraniami Y, teda náhodnou veličinou (náhodným vektorom) Y a s hustotami (alebo pravdepodobnostnými funkciami) {/(•|0): eeneBm} (voči Lebesgueovej miere alebo voči sčítacej miere), množstvo informácie o parametri 0, (keď výsledok experimentu bol y), možno vyjádřit pomocou výrazu (11.1) I[n(-\y),n(-)}= f [ln ak 0 je (absolútne) spojitý (voči Lebesgueovej miere) náhodný vektor, alebo T(0i|y)" 7r(0|y)d0, (11.2) ln- vr(02) T(0i|y), ak 0 je diskrétny náhodný vektor s hodnotami {0\, O2,...}. Výraz (11.1) resp. (11.2) je vždy nezáporný, je nulový práve ak tt(0) = ir(0\y). Z Bayesovho vzorca (2.7) vyplýva, že v takomto prípade T(»)/(y|0) w(9\y) Jnn(0)f(y\0)dv(0y f(y\0) = / 7r(0)/(y|0)^(0). Pretože lavá strana predchádzajúcej rovnosti závisí od 0 a pravá strana nie, musí byt funkcia vierohodnosti pre dané y konštantná. Vektor y neobsahuje žiadnu informáciu o parametri 0. Pri jednom výsledku y máme množstvo získanej informácie dané vztahom (11.1) resp. (11.2). Prediktívna hustota náhodného vektora Y je f* (y) a preto stredné množstvo informácie, ktoré možno získat z experimentu sa rovná (11.3) v spojitom prípade, alebo W = fY{/[vr(-|y),7r(-)]}= / 7[7r(.|y),7r(.)]/*(y)dA(y) Jy (11.4) 13 ak Y je diskrétny náhodný vektor. Pozrime sa na iné interpretácie Iexp. Uvažujme spojitý prípad a predpokladajme, že môžeme zamieňat poradie integrácie. Pomocou (2.6) a (2.7) upravíme (11.3) (11.5) Iexp=£Y{I[7T(-\y),7T(-)}} = /[7r(.|y),7r(.)]r(y)dA(y)= / / [ln7r(0|y)] w(9\y)du(9) ) f*(y)d\(y) y Jy \Ju \nw(9\y)]g(9,y)du(9)d\(y)- / [\nw(9)]g(9,y)du(9)d\(y). y Ju JyJu Spočítajme teraz strednú I-divergenciu f(-\0) voči /*(•) I[f(-\9),r(-)H0)du(9) = u 9(0, y) \nf{y\9)]f{y\9)d\{y)\^9)du{9)-u \Jy J JuJy ln <0\y) f(y\9)n(9)dv(9)d\(y) \lnw(9\y)]g(9,y)du(9)d\(y)- / [\nw(9)] g(9, y)dv(0)d\(y) = Iexp. lyJU Jy Ju Stredné množstvo informácie, ktoré možno získat z experimentu sa rovná strednej I-divergencii f(-\0) voči /*(•)• Spočítajme I-divergenciu združenej hustoty g(0,y) voči súčinu marginálnych hustôt f*(y) a it(9). Pomocou (2.7) a (11.5) dostávame Uxy m ^ 9(0,y) y Ju ln 7r(0|y) g{6,y)d{vx\){6,y)= / [\Rg{9,y)]g{9,y)dv{9)d\{y)-Jy Ju g(9,y)du(9)d\(y) - / [ln7r(0)] g(9, y)dv(0)d\(y) = Iexp. Jy Ju Stredné množstvo informácie, ktoré možno získat z experimentu sa rovná I-divergencii združenej hustoty g(0,y) voči súčinu marginálnych hustôt f* (y) a it(9). Preto dostávame (pozri Poznámku pod Definíciou 9.2) Veta 11.1. Stredné množstvo informácie získanej z experimentu Iexp je vždy nezáporné, pričom Iexp = 0 vtedy a len vtedy, keď pozorovaný vektor Y a parametre 0 sú nezávislé náhodné vektory, t.j. keď v Y nie je obsiahnutá žiadna informácia o 0. Poznámka. Na str. 29 v knižke A.Pázmana je dokázané tvrdenie, že stredná informácia získaná z pozorovania r(y) (t(-) je nejaká štatistika) je menšia alebo rovná strednej informácii získanej z pozorovania y. Rovnost nastáva práve vtedy, ak t (y) je postačujúca štatistika. "Náhrada" pozorovaného vektora y hodnotou t (y) nemôže zväčšit množstvo informácie o 0, nanajvýš ho zachovat. Toto nastane práve vtedy, ak je t(-) postačujúca štatistika. 12. Využitie I-divergencie v asymptotike Bayesovho vzorca Zopakujme si Zákon velkých čísel - Chinčinovu vetu (pozri Anděl, str. 183). 14 Veta 12.1. (Chinčinova.) Nech X\,X2,... je postupnost nezávislých náhodných veličín, ktoré majú rovnaké rozdelenie s konečnou strednou hodnotou jj,. Potom pre n —> oo 1 ^ n ^—^ podlá pravdepodobnosti. V kapitolke 5. "Pravidlo retazenia pre nezávislé pozorovania" sme ukázali, že postupne pridávané nezávislé experimenty "akumulujú" informáciu. Majme zvolené nejaké apriórne rozdelenie (môže byt aj subjektívne). Nezávisle opakujeme ten istý experiment tak, aby skutočná hodnota parametra 0 bola tá istá (ale neznáma). Označme ju 9*. Ukážeme, že pri velkom počte nezávislých opakovaní experimentu, výsledné aposteriórne rozdelenie sa v limte koncentruje do bodu 9*. Veta 12.2. Nech 0 je diskrétny náhodný vektor s hodnotami {9\, 9Ž, ■■■} ■ Apriórne pravdepodobnosti tt(Oí) > 0 pre každé i = 1, 2,... . Teda apriórna pravdepodobnostná funkcia diskrétnho náhodného vektora 0 je {7r(0i)}i>1. Nech yi,y2, ■■■,y„ sú pozorované hodnoty v n nezávislých experimentoch s tou istou hustotou f(-\9), čiže yi,y2, ■■■,y„ sú realizácie spojitých nezávislých rovnako rozdelených náhodných vektorov Yi, Y2,Y„. Nech pre každé Q i 7^ 9 j platí f(-\9i) 7^ f(-\9j). Potom pre aposteriórnu pravdepodobnostnú funkciu platí 'l, ak 9X = 9*, o, ak e, ^ e*, lim 7r(02|y") pričom y™ = (y1; y2,yj. Dôkaz: Platí 71 f(yn\Ok) = l[f(yi\ok), k = 1,2,.... i=l Teda podlá Poznámky pod Vetou 2.2 (Bayesovou) platí Ľfc^(W(y»l*fc) E^^f^j Ľ^^exp^}' kde Náhodné veličiny ln fc) ]) , k = 1,2,n majú strednú hodnotu So- f(Yk\e*) -7(/(-|»*),/(-|»*)) <0, Ä; = 1,2,..,n, pričom táto stredná hodnota je rovná 0 práve ak f(-\0i) = f(-\9*), čo nastane práve vtedy ak 9i = 9*. Teda ak 9i = 9*, tak podlá pravdepodobnosti 1 y, /(Yfclft) n tí f(Yk\e*) 0. 15 Ak Oj ^ O*, tak podlá pravdepodobnosti iy,nMl,a 0, pričom Y" = (Yi,Y2,...,Y„)'. * Dokázaná veta je len jednoduchším prípadom asymptotických tvrdení o aposteriórnej hustote. V knihe J.Anděla v kapitole o bayesovských metódach nájdete dôkaz, že aposteriórna hustota asymptoticky (pri velkom počte pozorovaní) málo závisí od apriórnej hustoty. Podobný zmysel má aj veta 2.1 v skriptách M.Huškovej. 13. Princíp neurčitosti Od vzniku bayesovskej štatistiky bola snaha stanovit, aké apriórne rozdelenie treba zvolit, ak nemáme žiadnu informáciu o 0. Ak 0 je absolútne spojitý (vzhladom k Lebesgueovej miere), pričom P(& € íí) = 1, íí € Bm a Lebesgueova miera /lí(íí) > 0, podlá princípu neurčitosti sa volí apriórne rozdelenie rovnomerné na í~ž, teda tt(9) = k > 0, 0 G í~ž (obyčajne sa volí k = 1). V prípade, že 0 je diskrétny, jeho apriórne rozdelenie pravdepodobnosti volíme tiež rovnomerné. Často sa stáva, že /lí(íí) = 00 (alebo diskrétna náhodná veličina 0 nadobúda spočítatelne vela hodnôt). Potom apriórna hustota je nevlastná. Príklad 13.1. Nech náhodná veličina Y má binomické rozdelenie pravdepodobnosti s pravdepodobnost-nou fumkciou f(y\0) f(y\6) = Q &y(í - 9)n~y, y = 0,1, 2,n (n je známe celé nezáporné číslo). O parametri O G (0,1) nemáme žiadnu apriórnu informáciu. Podlá princípu neurčitosti volíme apriórnu hustotu n(6) = 1, 0 < 6 < 1. Podlá (2.7) je aposteriórna hustota parametra O rovná qv(1 _ 0)n~y tt(%) = , [ '-- O<0<1. V |y; B{y+l,n-y+\y Príklad 13.2. Nech Y\, Y2,Yn je náhodný výber z normálneho rozdelenia JV(0, a2), kde a2 je známe kladné číslo. O parametri 0 G íl = 1Z nemáme žiadne apriórne informácie. V tomto prípade za apriórnu hustotu parametra 0 (podlá princípu neurčitosti) volíme n(9) = 1, e e n 16 " 2íj2 Si=l(ž/i (nevlastná hustota). Pretože f(yi,y2,-,Vn\6) = ,0 L n aposteriórna hustota parametra O je podia (2.7) t(%i,ž/2, -,y„) = vr(6»|y) In aV2Tr e 2, ■^Eľ=i(^-y)2e-^(y-o)2 M0-y)2 Aposteriórna hustota parametra 6 je JV t/, — \ n Princíp neurčitosti (alebo princíp rovnakej pravdepodobnosti) vedie v niektorých prípadoch k protirečeniam. Napríklad ak máme experiment s hustotami f(y\9), v ktorom P (Q € (0,1)) = 1. Zaujíma nás ale "nový" parameter j3 = 82. Tento tiež nadobúda hodnoty z (0,1). Podlá princípu neurčitosti by oba parametre mali byt rozdelené rovnomerne, t.j. s apriórnymi hustotami (13.1) 7r(0) = 1, 7t(/3) = 1. Platí veta Veta 13.1. (Anděl, str. 46.) Nech X má spojitú distribučnú funkciu F (x). Predpokladajme, že F'(x) = f (x) existuje všade s výnimkou najviac konečne vela bodov. Nech t je rýdzo monotónna funkcia, ktorá má všade deriváciu. Položme Y = t(X). Označme r inverznú funkciu k t. Potom Y má hustotu g(y) = f(T(y))\T'(y)\. Označme t(6) = 62, čiže r(ß) = = 6. Teda platí n(ß) = ir(y/ß) (13.1). dß 2Vß Toto je v spore s 14. Jeffreysova apriórna hustota Predchádzajúce úvahy viedli štatistikov k záveru, že bolo by vhodné namiesto princípu neurčitosti volit apriórne rozdelenie tak, aby aposteriórne rozdelenie nezáviselo od parametrizácie modelu. V prípade, že 0 G lZm je (absolútne) spojitý (vzhladom na Lebesgueovu mieru /i), p(& € íí) = 1, pričom /j,(íl) > 0, je riešenie tohto problému vo vete 14.2. Ešte pred jej sformulovaním si zopakujme niekolko definícií a viet. Definícia 14.1. Nech náhodný vektor Y = (Y\, Yi, Yn)' má hustotu f(y\9) (vzhladom k a—konečnej miere /i), pričom 0 = (6\,9m)'. Predpokladajme, že platí: (A) 0 G í~ž, kde í~ž je neprázdna otvorená množina v lZm. (B) Množina M = {y : /(y, 9) > 0} nezávisí od Q. df(v 0} (C) Pre skoro všetky y e M (vzhladom k fi) existujú parciálne derivácie f[(y,9) = —^—, i = o6i 1,2,...,m. (D) Pre každé i a pre všetky 0 G í~ž platí JM //(y, 9)d/j,(y) = 0. (E) Pre každú dvojicu (i, j) existuje konečný integrál m ŕm ľ fí&Mk^ff flw f \ MlÁ9) = L ŕ(y,o) /fr'Wy)- (F) Matica M(0) s prvkami {Wl(0)}ij = Mij(o) je pozitívne defmitná pre každé 0 G í~ž. 17 Ak sú splnené predpoklady (A) až (E), tak M(0) sa nazýva Fisherova informačná matica. Ak sú splnené predpoklady (A) až (F), tak tak hovoríme, že systém hustôt je regulárny. Poznámka. Definícia 14.1 je platná aj ked náhodný vektor Y = (Y\,Y'2, ...,Yn)' má diskrétne rozdelenie pravdepodobnosti. Vtedy jeho "hustota" je jeho pravdepodobnostná funkcia. V bode (D) podmienka Im ň(y^d)d^(y) = 0 Je vlastne £ r//(Y|0)/j(Y|0) /í (Y, e) /(y, e) 0 a Fisherova informačná matica má prvky Mi j (9) £9 - f2SY^ Poznámka. V učebnici J. Anděla, str. 261 je dokázaná veta, ktorá tvrdí, že v prípade regulárneho systému hustôt {/(y|0), 9 g SI}, ak existujú derivácie f"j = q/1 q/1—, i, j = 1,2, ...,m a pre všetky 9 g SI platí £g /"(Y|0) /(Y|0) 0, i,j = 1,2,...,m, tak Mij(e) -£9 89,89j d2 ln/(Y|0) 89,89, Definícia 14.2. Zobrazenie f z množiny A do množiny B nazývame prostým, ak platí {xi g A, x2 e A, xi ^ x2} {f(Xl) ^ f(x2)}. Definícia 14.3. Nech f je zobrazenie z 7?.r do TU. Ak je y = /(x), kde y = (1/1, ...,yr) a x = (x\, ...,xr), položme t/i = fi(xi,xr), i = 1,2,r. Povieme, že zobrazenie f je regulárne v množine M c lZr, ak platí: 1. M je otvorená. 2. Funkcie /1, ...,/r majú parciálne derivácie prvého rádu spojité v M. 3. Pre každé x g m platí, že jakobián Df(x) 7^ 0, pričom / Žíl 8x\ Df (x) = det 8h_ \ 8xr 8fr 8fr \ 8x\ 8xr ) Veta 14.1. (O substitúcii, Anděl, str. 47.) Nech f je zobrazenie otvorenej množiny P c TU na Q c TU. Nech f je regulárne a prosté v P s jakobiánom Df(x). Nech m c Q je borelovská množina a F meratelná reálna funkcia. Potom platí / F(x)dx= / F [f(u)] |Df(u)|eŕu, J m J í-1 (m) akonáhle jeden z integrálov existuje. Veta 14.2. (Anděl, str. 291.) Nech náhodný vektor Y má pri danom parametri 0 g SI g Bm hustotu /(y|0) (alebo pravdepodobnostnú funkciu {/(yi|0)}i>i) • Predpokladajme, že systém hustôt {/(y|0) : 9 g SI} je regulárny a má Fisherovu informačnú maticu M(0), ľ(9) je Lebesgueova miera. Nech platí (14.1) Položme (14.2) 0 < / f(y\0)\M(0)\ždO < 00. 1 Jnf(y\9)\M(9)\U9 Nech H je regulárne a prosté zobrazenie SI na SI* g Bm. Označme 77 = H(0) a f*(y\r)) = /(y|H-1(í7)). Potom je {/*(y|í?), f] g SI*} regulárny systém hustôt. Ak označíme M*(77) jeho Fisherovu informačnú maticu, potom pre lubovolnú množinu B g Bm, H C SI platí (14.3) cf(y\0)\M(0)\*dO h(b) c/*(y|T/)|M*(T/)|2eŕT7. 18 Dôkaz je uvedený v knižke J.Anděla, str.291. Z Vety 14.2 vyplýva, že pri apriórnej hustote parametra 0 rovnej funkcii |M(0)|^ (alebo íubovolnému kladnému násobku tejto funkcie) je aposteriórna hustota parametra 0 rovná c/(y|0)|M(0)|5 a je to hustota. Vzorec (14.3) zaručuje, že aposteriórna pravděpodobnost lubovolnej "rozumnej" množiny B nezávisí od toho, či sa vychádza od parametra 0 alebo od nejakej jeho hladkej transformácie r\ = 11(0). Na lavej strane (14.3) je totiž výraz p(9 G B|y) a na pravej strane p(r] G H(B)|y) = p(0 G B|y). Pri Jeffreysovej voľbe apriórnej hustoty parametrov 0 a r\ sú obe tieto aposteriórne pravdepodobnosti rovnaké a nemožno dôjsť k paradoxnému výsledku ako v prípade princípu neurčitosti. Jeffreys v r. 1946 [An invariant form for the prior probability in estimation problems. Proc. Roy. Soc. A, 186, 453-461] navrhol "univerzálnu" apriórnu hustotu (vzhladom na Lebesgueovu mieru) ttj(0) oc [deí(M(6>))]Í Jeffreys navrhol svoju hustotu len pre jednorozmerný parameter. V takomto prípade sa jeho hustota osvedčila, lebo vtedy vždy možno zmenit parametrizáciu modelu tak, aby informačná matica (vlastne Fisherova miera informácie - číslo) bola konštantná, t.j. nezávislá od hodnoty parametra. Toto je možné len výnimočne, ak parameter 0 je mnohorozmerný. Jeffreysova apriórna hustota býva často nevlastná. Príklad 14.1. Nech náhodná veličina Y má binomické rozdelenie pravdepodobnosti s pravdepodobnost-nou funkciou f(y\9) f(y\6) = Q &y(í - 9)n~y, y = 0,1, 2,n (n je známe celé nezáporné číslo). Nájdite Jeffreysovu apriórnu hustotu. Riešenie: Binomické rozdelenie spĺňa predpoklady Vety 14.2. Dostávame dln f (y\9) y n —y y — n9 89 = 9 ~ 1 - 9 = 9(1-9)' Fisherova miera informácie je M (9) = £g dln f (Y\6) 89 £{Y-ni 92(l-9)2 ~ 9(1 lebo £(Y) = n9 a disperzia "D(Y) = £(Y — £(Y))2 = n9(l — 9). Jeffreysova apriórna hustota je rovná k , =, kde k > 0. vW^W) Ked zvolíme apriórnu hustotu pre parameter 9 V n y/6(l-6) dostávame z (2.7) aposteriórnu hustotu parametra 9 Tr(9)f(v\9) 9-i(l-9)-i(n)9y(l-9)n-y \ J^(9)f(y\9)d9 />"i(l - 9)~i (ny)9v(l - 9)-vd9 B{y+\,n-y+\) Majme teraz záujem o parameter j3 = 92. Podlá Vety 13.1 je aposteriórna hustota parametra j3 rovná (14.4) ■n(P\v) = •••> ž/n|0> = -~-e 2o-2 2^i=iAw» tiež, že (2ir)-(jn ln /(y|0, cr) = In 2vr - nln o- - —2 ^(y, - S)2. Platí 2cr2 2=1 a) Ak parameter cr poznáme, tak Fisherova miera informácie je rovná (podlá Poznámky nad Definíciou 14.2) d2 ln/(Y|0,o-)" dff2 n 7žž ■ M{9) = -£g Pri známom cr to je konštanta a preto Jeffreysova apriórna hustota parametra 9 je n(9) = 1, 9 eTZ. b) Ak parameter 9 poznáme, ale nepoznáme parameter cr, tak Fisherova miera informácie je rovná (podlá Poznámky nad Definíciou 14.2) M(a) = -£„ d2 ln/(Y|9, cr) do 2n 7žž ■ Preto Jeffreysova apriórna hustota parametra cr je 1 7r(cr) = -, cr > 0. cr c) Ak nepoznáme ani parameter 9 ani parameter cr, tak prvky Fisherovej informačnej matice sú (podlá Poznámky nad Definíciou 14.2) 'O2 \nf(Y\9,a)' Mn(9,a) = -£g dff2 n 20 M12(6,a) M22(9) = -Sg d2 \nf(Y\9,a) deda d2 ln/(Y|6»,cr) 2n da2 M21(6,a) = M12(6). Jeffreysova apriórna hustota je ir(6,a) = -Ar, 6»e^, cr>0. Poznámka. J. Anděl v učebnici na str. 294 poznamenáva, že väčšinou v prípade c) predchádzajúceho príkladu sa za apriórnu hustotu volí 1 ir*(9,a) en, cr > o. Zdôvodňuje sa to tým, že sa dá mnohokrát dopredu předpokládat, že 9 a a sú nezávislé a apriórna hustota tt*(9, cr) sa rovná súčinu apriórnych hustôt 7r(0)7r(cr) (ako keby tt(9) a it( o. Platí Preto M{9) \nf{y\9) = -\nyl-9 + y\n9 'Ô2 In f (y\9) 892 Sa Preto Jeffreysova apriórna hustota je itj(9) oc —=. v 9 15. Hierarchizácia apriórneho rozdelenia V niektorých prípadoch možno apriórnu hustotu zapísat v tvare (15.1) vr(0) = ^71-1(017)11-2(7)^(7), kde 7ri(0|7) je hustota, ktorej tvar je známy, ale neznáme sú parametre 7 G Q (realizácie náhodného vektora T), ako aj neznáma je hustota 7^(7) vzhladom na a—konečnú mieru jj, (my uvažujeme alebo Lebesgueovu mieru alebo sčítaciu mieru). Poznámka. Upozorňujeme len, že "všeobecný" zápis (15.1) chápeme: (i) Ak je 0 spojitá aj T spojitá, tak 7ri(0|7) je hustota známeho tvaru a 7^(7) je neznáma (apriórna) hustota náhodného vektora T. Táto je nenulová na Q. Potom (15.2) vr(0) = / 7n(0|7)7r2(7)d7. Jg (ii) Ak je 0 spojitá a T diskrétna s hodnotami {71,72, •••}, tak {7ri(0|7i)}i>i sú hustoty známeho tvaru a {i"2(7j)}j>i Je neznáma (apriórna) pravdepodobnostná funkcia parametra 7 G {71, 72,...}. Potom (15.3) vr(0) = ^ 711(017^2(7,). i>l 21 (iii) Ak je 0 diskrétna s hodnotami {#i,#2, •••} a T spojitá, tak {íti(0j |t)}j>i je pravdepodobnostná funkcia známeho tvaru a 7^(7) je neznáma (apriórna) hustota náhodného vektora T. Táto je nenulová na Q. Potom (15.4) MOj)}j>i = { í ^i(0j\l>2(l)dl) • Ug ) j>i (iv) Ak je 0 diskrétna s hodnotami {0\, O2,...} a T diskrétna s hodnotami {7-,, 72,...}, tak {tti (Oj |Tí)}j>i sú pravděpodobnostně funkcie známeho tvaru (pre každé i G {1, 2,...}) a {tt2(~/í)}í>i je neznáma (apriórna) pravdepodobnostná funkcia parametra 7 G {7i,72, •••}• Potom (15.5) MOj)h>i = \ e^(0,17^(7,) [ i>l Vztah (15.1) znamená hierarchickú štrukturalizáciu. Hustota 7^(7) je primárna a hustota tt(O) je sekundárna v tejto chierarchizácii. Pôvodne sme mali experiment {/(y|0): 6 G í~ž} a teraz uvažujeme experiment {g(yh): 7 e G}, kde ff(y|7)= / f(y\0)iri(0\~i)dp(0). V tomto experimente je apriórna hustota 7^(7). Táto sa volí rovnomerná (t.j. podlá princípu neurčitosti), alebo akákolvek iná. Jaj volba nieje tak kritická, ako volba apriórnej hustoty tt(0) v pôvodnom experimente {/(y|0): 6 G íž}. Príklad 15.1. (Podlá knižky A. Pázmana, str. 9.) Novovyvinutý prístroj (alebo metóda) umožňuje detekovat chorobu, resp. nepřítomnost choroby s 95% pravdepodobnostou. To znamená, že zo 100 pacientov, ktorým prístroj indikoval chorobu, približne 5 pacientov touto chorobou netrpí. Naopak, zo 100 pacientov, ktorým prístroj indikoval, že sú zdraví, v skutočnosti približne 5 je chorých. Náhodná veličina Y—indikácia choroby prístrojom má realizácie yi — prístroj chorobu indikoval alebo 1/2— prístroj chorobu neindikoval. Parameter 9 nadobúda hodnoty 9\— vyšetrovaná osoba je chorá alebo 62— vyšetrovaná osoba je zdravá. Teda 7r(#) je pravdepodobnostná funkcia popisujúca apriórne rozdelenie pravdepodobnosti ochorenia v experimente {f(y\9) : 0 G {#i,#2j-}- Zo zadania platí /(íft|0i) = /(íte|02) = O,95 a f(y2\91) = f (yi\92) = 0,05. Ľudia sú kategorizovaní podlá určitých predispozícií k chorobe na s kategórií 7i,...,7s a podlá štatistík pravdepodobnostná funkcia ochorenia v jednotlivých kategóriách je {ti(^-|7í)}Li, i =1,2,.., s, 22 Toto je apriórne rozdelenie pravdepodobnosti ochorenia v niektorej z kategórií. Nepoznáme však apriórnu primárnu hustotu 7^(7) na T = {1, 2,s}. Môže sa zvolit rovnomerná. Namiesto experimentu {f(y\9) : O G {#1, #2}} s (neznámou) apriórnou pravdepodobnostou tt(9) budeme teraz uvažovat experiment {g(y\j) : 7 G {1, 2,s}}, pričom {S(žfc|7)}?=i = |e^I^)ti(^-|7)| • Apriórna hustota v tomto experimente je 7T2(A) (rovnomerná na {1, 2,s}. Predchádzajúce úvahy zhrnieme vo vete Veta 15.1. (Pázman, str. 47.) Ak platí (15.1), tak aj aposteriórna hustota má hierarchickú štruktúru t^ly) = / 7ri(é*l7,y)7i"2(7|y)rfM(7)I ■h kde m , = /(y|0)7ri(fl|7) ni[ l7'yj Jnf(y\9)M0h)du(e) je aposteriórna hustota pre 0 pri daných parametroch 7 a kde fl'(y|7)i"2(7) ""2(7! y) 16. Empirické metódy určenia apriórnej hustoty Pre predikčnú hustotu platí (pozri časí Inferencia v 2. kapitole) (16.1) /*(y)= / f(y\0)ir(0)dv(0). Keby sme poznali /*(y), riešením integrálnej rovnice (16.1) by sme mohli dostat apriórnu hustotu tt(0). Z náhodného výberu Yi, Y2,Y„ vieme len odhadnut hustotu (resp. pravdepodobnostnú funkciu) f* (y) (napríklad jadrovými odhadmi hustoty). Potom môžeme riešit integrálnu rovnicu (16.1) a takto získat odhad 7r(0) apriórnej hustoty tt(O). Pri tomto postupe odhadovanie predikčnej hustoty sa obyčajne nerobí bayesovskými metódami, preto tento postup nie je (obyčajne) čisto bayesovský. V prípade, že v bayesovských metódach použijeme namiesto apriórnej hustoty tt(O) jej odhad tť(9), hovoríme o empirických bayesovských metódach. Treba tu zdôraznit dva fakty: - Odhad 7r(0) treba získat z pozorovaní, ktoré nezávisia od pozorovaní v experimente, v ktorom používame bayesovské metódy. - Odhad 7r(0) možno získat len z takých pozorovaní, v ktorých sa parameter 0 náhodne mení podlá tt(O). Nie je možné preto odhadovat subjektívnu apriórnu pravděpodobnost. Príklad 16.1. (Podlá knižky A. Pázmana, str. 48.) Nech Y\,Y'2, ...,Yn je náhodný výber z Poissonovho rozdelenia s pravdepodobnostnou funkciou (16.2) f(y\0) = e-9^r, y = 0,1,2,..., y]- pričom parameter O sa náhodne mení podlá hustoty tt(9) (vzhladom k Lebesgueovej miere) a P(0 G (0,oo)) = 1. Ak aj pozorovanie Yn+\ je poissonovsky rozdelené s parametrom 9 (realizáciou náhodnej 23 premennej 0) a nezávislé od Y\,Y2, ...,Yn, aposteriórne rozdelenie (hustota) pri získanom pozorovaní yn+i, teda 7r(0|y„+i) je podlá (2.7) rovné ir(6)f(yn+1\6) H^Vn+1)- ^K{6)f{yn+1\6)d6 a jeho stredná hodnota £ [it(6\yn+i)\ je (16.3) £ br(%„+i)] J™6n(6)f(yn+1\6)d0 J0°°iT(6)f(yn+1\6)d6 • Táto stredná hodnota sa niekedy považuje za bayesovský odhad parametra 9. Po dosadení podlá (16.2) do (16.3) dostávame dá A\ C í (0\ m (ž/n+1 + 1)! J (ž/n+1 + 1) , 1A (!6-4) £ 7r(%„+i) =--=-———- (yn+1 + l). f°°n{e)e-°e—d6 f ^ Výraz (16.4) vieme spočítat len ked tt(9) poznáme. Ak tt(9) nepoznáme a n je dost velké číslo, postupujeme nasledovne: Tl(v) Označme n(y) počet výskytov hodnoty y v n nezávislých pozorovaniach Y\,Y2, ...,Yn. -je relatívna n početnost výskytu hodnoty y v týchto n nezávislých pozorovaniach a je to frekvenčný odhad pravdepodobnosti f* (y). Teda podlá (16.4) empiricky získaná aposteriórna stredná hodnota £ [it(6\yn+ij\ (odhad tejto strednej hodnoty) je £ vr(6 yn+i)\ = —-.-—t- (yn+i + 1). n{y„ +1) 17. Konjugované systémy apriórnych hustôt Konjugovaným systémom apriórnych hustôt, priradeným k danému štatistickému modelu experimentu {/(y\9) : 0 g íí} nazývame taký parametrizovaný systém V apriórnych hustôt na íl, že a) ak 7r(-) g V, tak aj aposteriórna hustota vr(-|y) g V pre každé y g J-', b) prechod od 7r(-) g "P k vr(-|y) sa deje jednoduchou, lahko vypočítatelnou zmenou nejakých parametrov. 17.1. Konjugované systémy pre hustoty exponenciálneho typu V užšom zmysle konjugovaný systém je priradovaný hustotám, ktoré sú exponenciálneho typu, t.j. ktoré sú tvaru (17.1) /(y|0) = exp {-^(y) + ť(y)7(0) - n(6)} , kde Y g 3^ c TU1, 0cííc Km, t : y TZk, 7 : íl -> TU, ip : y -> TI, k : íl K sú dané funkcie. Konjugované apriórne hustoty majú tvar (17.1 a) 7Tc,;(0) = K(c,l)exp{c'j(6) - Ik(0)} , kde c g TZk, l g TZ sú také, aby " ec'7(0) - H0)dv{e) = [K(C; < ^ n 24 [K(c + t(y),Z + l)]-1 < oo pre každé y £ J. Po dosadení do Bayesovho vzorca sa přesvědčte, že ak apriórna hustota je vrcí(-), tak aposteriórna hustota je 7rc+t(yy+1(-), teda v skutočnosti ide o velmi jednoduchý prepočet parametrov (17.2) c->c + t(y), l^l+l. Problémom môže byt výpočet K{c, l). Metoda má význam hlavne ak rozdelenia ttc,i(-) patria k rozdeleniam známym z teórie pravdepodobnosti, kde normovači výraz K(c, l) už bol stanovený. Systém konjugovaných apriórnych hustôt {ttcj(0)} priradený k hustotám exponenciálneho typu (17.1) možno podstatne rozšířit, ak uvažujeme aj konvexné zmesi hustôt v tvare k (17-3) ^Wí7rCiiíi(0), i=l kde ují > 0, Yli=i wi = 1- Platí veta Veta 17.1. (Pázman, str. 51.) Ak tt(0) má tvar (17.3), tak aposteriórna hustota je tiež konvexnou zmesou tvaru k i=l kde koeficienty v tejto konvexnej kombinácii sú K(Cl,h) ^ (y) = ^ Vfc K(cj,lj) ^(Cj+t(y),ZJ + l) Dôkaz: Ak dosadíme (17.3) do Bayesovho vzorca (2.7), dostávame ^7rCiiíi(0)/(y|0) 0, & > 0 a 5(a, 6) = TTr^T- Aby 7rCi;(6») bola hustota, musí platit [k(c, O]"1 = ľ ec{i-e)lm-cde = ľ e^-1 [1-6)^-^-^6 = b(c+i, im-c+i) = r(c + 1)r(žm~c+1). Jo Jo T(lm + 2) Pravidlá "prepočítavania" parametrov sú dané vztahom (17.2), preto v tomto prípade c —> c + y, l —> l + 1. a aposteriórna hustota je n(6\y) = nc+yM1(6) = K(c + y,l + 1)0C+»(1 - 0)«+V™-°-v = =_r((ž + i)TO + 2)_QC+ _ {l+1)m_c_y T(c + y + l)T((l + l)m-c-y+l) K ' Z tohoto vztahu vieme napríklad spočítat aposteriórnu strednú hodnotu rr ia\ w C + V+ 1 = (ž + i)m + 2' čo je bayesovský odhad (vlastne jeho realizácia) parametra 6 (pozri nižšie). V praxi je postup taký, že v súlade s kapitolou 16 nájdeme odhady neznámych parametrov apriórnej hustoty tt(6) z konjugovaného systému hustôt (obyčajne nie bayesovskými postupmi) a neznáme parametre nahradíme ich odhadmi. V tomto príklade potom dostávame (Z + 1)to + 2 Príklad 17.2. (Poissonovský experiment, A. Pázman, str. 53.) Nech Y má Poissonovo rozdelenie pravdepodobnosti s pravdepodobnostnou funkciou (hustotou vzhladom k sčítacej miere) av f(y\6) = e-e— oce-e+yln6, y = 0,1,2,... , y]- teda 3^ = {0,1, 2,}, 0 G íl = (0, oo). Hustota f(y\6) patrí do exponenciálnej triedy hustôt. Konjugovaná apriórna hustota (voči Lebesgueovej miere) na íl má podlá (17.1 a) tvar ttcJ(6) = K(c, l)e~w6c. Pre gama funkciu platí T (a) = J0°° xa~1e~xdx, a > 0. Preto [tf(c,/)]-!= H e-w6cd6=-±- H e-uu^-1du~V{C+l) ]c+l / íc+1 o 1 Jo 1 26 pre c > —1, l > 0. Pravidlá "prepočítavania" parametrov sú dané vztahom (17.2), preto v tomto prípade c —> c + y, l —> l + 1. a aposteriórna hustota je vr(%) = vrc+y,í+1(e) = Ä"(c + y,l + l)e-W>0c+v = ^,+ 'J ^, e"^1)^. T{c + y + 1) Z tohoto vztahu vieme napríklad spočítat aposteriórnu strednú hodnotu rr ia\ w c + V+l swm\ = -j^i-, čo je bayesovský odhad (vlastne jeho realizácia) parametra 9 (pozri nižšie). V praxu opát neznáme parametre c a l nahradíme ich (nie bayesovskými) odhadmi získanými z "iného" experimentu. Teda dostávame c + y+1 £[tt(%)] 1+ 1 Príklad 17.3. (A. Pázman, str. 54.) Nech Yi,Y2,...,Yn je náhodný výber z rovnomerného rozdelenia na (0, 9), teda 0 G íl = (0, oo). 7, £ J= (0, oo) má hustotu (vzhladom k Lebesgueovej miere) f(Vi\8) e-1, ak Vi e (o,6»), [0, ak yi e (9, oo). Hustota náhodného vektora Y = (1^, Y2, ...,Yn)' sa nedá zapísat v tvare (17.1), ale keď definujeme í (y) = max í/í, môžeme písat kde f(y\9) = l[f(yt\9) = 9-nl{t{y)^)(9), 2ľ(a,6)(él) 1, ak 9 e (a,b), 0, ak 9 ^ (a, 6). Len podotýkame, že í (Y) je postačujúca štatistika. Konjugovanú apriórnu hustotu zvolíme (17.4) ncj(9) = (l - V^K^é-H-1!^)^) (přesvědčte sa, že ttcj(9) je vskutku hustota). Preto (17.5) /(y|0Ki(0) oc č-("+í)i(max{Cit(y)}i0o)(č). Zo vztahov (17.4) a (17.5) vidíme pravidlá pre "prepočítavanie" parametrov c —> max{c, í(y)}, l —> l + n. Príklad 17.4. Majme náhodný výber Y\, Y2,Yn z rozdelenia N(9, a2), pričom a2 poznáme. Parameter 0 G íl = 1Z. Podlá Príkladu 13.2 je hustota náhodného vektora Y = (Y\, I2,Yn)' rovná \2 f(yi,y2,--,yn\0) = /(y|i -\- e 2a2 Ei=l(ž/i 27 1 e 2, (2vr) ~an Ako konjugovaný systém sa v tomto prípade v praxi berie množina všetkých rozdelení N(a,b2), kde a G K, b>0. Teda 1 2vr b -e 262 ôk(č-«)2 Preto hustota aposteriórneho rozdelenia parametra 0 je vr(%) /(y|flK,b(fl) Úpravou dostaneme -2M^-I/)2e-2^-«)2 n&2y + acr^ 2a- :(0-yf 2b2 (9 - a)2 nb2 + a2 b2 a2 vr(%) ~ N nb2y + aa2 v nb2 + a2 b2a2 nb2 + a2 ' nb2 + a2 Stredná hodnota aposteriórneho rozdelenia je váženým priemerom parametra a (strednej hodnoty apriórneho rozdelenia, ktorú poznáme pred experimentom) a odhadu y založeného len na výberových hodnotách (na hodnotách získaných v experimente). Disperzia aposteriórneho rozdelenia je menšia než disperzia apriórneho rozdelenia b2 (experimentom zmenšujeme rozptyl). Disperzia aposteriórneho rozdelenia je ale tiež menšia ako —, čo je disperzia Y. Keby sme čerpali informácie o 8 len z experimentu, tak najlepší nestranný odhad n _ parametra 9 by bol Y. Bayesovský odhad využíva apriórnu informáciu a je "presnejší" ako ten najlepší frekventistický odhad, ktorý apriórnu informáciu nevyužíva. Ukážme si ešte príklady viacrozmerných konjugovaných systémov. Príklad 17.5. (Normálne rozdelenie s neznámou strednou hodnotou a disperziou, A. Pázman, str. 55.) Nech Y\,Yn je náhodný výber z rozdelenia N(8\, ^-), teda neznáma stredná hodnota je 9\ G 1Z a neznáma inverzná disperzia je 92 G (0, 00), 0 = (0i, 02)' G 1Z x (0, 00). Hustota Yi je a preto f(yí\0i,e2) = f(yl\e) f(yi,m,-,yn\o) = f(y\0) fexp{-|(,-,l)2} ^rTT e 2 2 (2tty 1 {(Eľ=i^ Er=i%2)(!i§)-f(-in^2+^ —— exp *• 2 Í02)} Pri označení (17.1) je t (y) = (Eľ=i Eľ=i%2)'> 7(0) = (M2, -f)', k(0) = §(-ln02+ 0?02) a pravidlá pre "prepočítavanie" parametrov sú podlá (17.2) ci c'2 c2 Ešte spočítame K(c,l). Podlá formuly pod (17.1 a) [K(c,l)] •2 Z ďž TeClW2 ~ -VlV2 -c2-ddidd2 (— 00,00) x (0,oo) 28 oo r />oo nl ni (li L d92 oo nl 6.? e 2 o oo rd _ 6 2 e 2 C'2 2eC1^2-f^-C^d91 nl) nl J\d6i 02nl in ci \ 2 e 2 nl) dOx d92 o de* oo nl 6,? e 2 v^ nl> 1 27r 2Tll oo nl—1 j — oo d6o C2'nl — c1 2 2 e_t72 2«; d6»2. Z podmienok na parametre gama-hustoty dostávame ci e c2 > 0, l > 0, alebo Substitúciou c2 1 cieK,c2>r, — < l < 0. In n dostávame z posledného integrálu C2nl — cf 2nl Mel)]-1 0o=U 2vr 1» aP kde a C2 a p ■■ nl + 1 2 " 2 Príklad 17.6. (Multinomický experiment, A. Pázman, str. 56.) Majme diskrétny náhodný vektor Y = (Yi, Y'2,!>)', ktorého zložky sú diskrétne náhodné veličiny. Môžu nadobúdat hodnoty 0,1,..., n. Vektor Y má multinomické rozdelenie pravdepodobnosti s pravdepodobnostnou funkciou /(ž/i, yr\eu 9r) = f(y\0) = P(Y1 = yi,Yr = yr) 9f. ž/l!.-ž/r! kde ž/í £ {0,1,n}, Y^i=i Vi = n> ^ (0> -0> Í = 1> 2,r, Xlj=i @j = 1- Takouto hustotou (pravdepodobnostnou funkciou) je matematicky popísaný napríklad experiment, v ktorom robíme n nezávislých pokusov. Každý pokus može rezultovat v jednom z r výsledkov, pravděpodobnost nastatia j—teho výsledku je 6j. Náhodná veličina Yi je počet nastatí i—teho výsledku v týchto n pokusoch. Konjugovaný systém hladáme v tvare r i=l (c = (ci,...,cr)', y = {(ž/i,ž/2,-..,yn) : Vi e {0,l,...,n}, J2l=1yt = n}, fl = {0 = (9u...,6r) : 6j G (0,1), j = 1, 2,r, Y^íj=i &j = !}• Pravidlá pre "prepočítavanie" parametrov sú Cj Cj + ž/í i = 1, 2, ...,r. Rozdelenie tohoto tvaru je známe. Je to Dirichletovo rozdelenie s hustotou r(ai + ... + ar) aa%-\ 1^=1^) l\ 1 29 (pozri napr. http://en.wikipedia.org/wiki/Dirichlet_distribution). Teraz už lahko dostaneme vyjadrenie Kc 18. Aproximácie integrálov Pri bayesovskej inferencii je klučové určenie aposteriórnej hustoty ir(0\y). Túto dostávame zo vzorca (2.7), pričom musíme spočítat integrál j(l'K(0)f(y\0)dv(0) (normovači faktor). V ďalšom budeme potřebovat spočítat strednú hodnotu aposteriórneho rozdelenia vektora 0, čiže opát len integrály Jn 9i it(0)f(y\0)dv(0), i = 1,2,.., m. Pokial je dimenzia dim(0) malá, môžu sa použit "bežné" numerické metódy približného výpočtu týchto integrálov. Načrtneme dve metódy vhodné pre dim(0) = 1. Obe metódy využívajú poznatok, že velké počty pozorovaní je tvar funkcie vierohodnosti f(y\9) (ako funkcie O pri daných, nameraných y) podobá tvaru hustoty normálneho rozdelenia. Teda uvažované integrály možno " rozumne" pretransformovať do tvaru $(í)e~~dí -oo s danou funkciou $(í). 18.1. Metóda kvadratury (Hermiteova-Gaussova kvadratura) Táto metóda využíva ako pomocný nástroj systém Hermiteových polynómov {Hk(t)}^_0, kde Hfc(.) = (-i)va|k*a. Integrál (18.1) sa podlá tejto metódy rovná <5>{ť)e-2dt = Y,HJ<5>{tJ) + E, 3 = 1 pričom tj, j = 1, 2,k sú korene Hermitovho polynomu H k (t) a Hj sú určené nasledovne 2fc+1fc!v/7? 12 i =-- H Chyba E sa rovná Hj = —rr-,— j = 1,2,...,k. 2fc(2fc)! [l) Podrobnejšie pozri napr. v knižke Ralston, A, A first course in numerical analysis, MCGRAW-HILL, INC, 1965. Tam sú uvedené aj hodnoty t j a H j pre j = 1,2,3, 4, 5. 18.2. Metóda Laplaceovej transformácie Metóda Laplaceovej transformácie využíva poznatok, že funkcia vierohodnosti f(y\6) (ako funkcia 9 pri daných, nameraných y) je pre veké počty pozorovaní značne koncentrovaná okolo bodu 9 = argmax/(y|#) 30 (9 je odhad získaný metodou maximálnej vierohodnosti). Touto metodou sa počíta integrál tvaru b(9)e-h^d9, -oo kde h(9) = — ln f(y\9). Funkcie b(9) a h(9) aproximujeme kvadratickými Taylorovými rozvojmi v okolí bodu b{9) = b{9) + db{9) d9 + 1 d2b(9) 2 dff2 )2 = b0 + b1(9-9) + ^b2(9-9)2, 2 dff2 (9-9f =h0 + h2( Lineárny člen sme v druhom vztahu vynechali, lebo dh{9) d9 0. Integrál (18.2) aproximujeme integrálom I = e + ^2(e-š)2 e 2 h2_e-!f(9 - 9)2 d9 d9 b0 + 2h9 Metodu možno zovšeobecnit aproximáciou funkcie b(9) rozvojom vyššieho rádu. 19. Niektoré simulačně metody generovania nezávislých realizácií náhodnej veličiny V mnohých prípadoch použitia bayesovských metód ide o určenie aposteriórnej hustoty, alebo o určenie strednej hodnoty aposteriórneho rozdelenia (teda o integrál). Základný princíp simulačných metód je generovanie postupnosti realizácií nezávislých náhodných veličín z nejakého (daného) rozdelenia. Hustotu potom aproximujeme histogramom, momenty aproximujeme výberovými priemermi. Tieto aproximácie sú tým lepšie, čím je väčší rozsah simulovaného súboru. Balíky počítačových programov obyčajne umožňujú generovat náhodný výber (resp. pseudonáhodný výber) z rovnomerného rozdelenia R(0,1). Bližšie sa s mechanizmom činnosti generátora (pseudo)náhodných čísel a so spôsobmi testovania vytvorenej postupnosti môžeme oboznámit napr. v skriptách Kalas, J., Pekár, J., Simulačné metódy, Bratislava: MFF UK, 1991. Tu uvedieme tri metódy. Viac o simulovaní náhodných veličín sa dočítate v knižke Antoch, J., Vorlíčková D., Vybrané metody statistické analýzy dat, ACADEMIA, Praha, 1992. 19.1 Metóda inverznej transformácie Táto metóda umožňuje generovanie realizácií nezávislých (jednorozmerných) náhodných veličín Yi,Y2,... z rozdelenia aké má náhodná veličina Y s danou distribuňou funkciou Fyiy) = PiY < y) a kvantilovu funkciou Fy1^) = iní {y : Fy(y) > u}, 0 < u < 1 . Toto rozdelenie môže byt diskrétne alebo aj spojité. Metóda je založená na nasledujúcom tvrdení: Veta 19.1. Nech náhodná veličina U má rovnomerné rozdelenie na (0,1). Nech Fyiy) je lubovolná distribučná funkcia. Potom náhodná veličina Y = F^iU) má rozdelenie s distribučnou funkciou Fy(y). Dôkaz nájdeme v knižke J. Anděla, str. 6. Postup pri generovaní je jednoduchý. Ak x\,X2, ■■■ je náhodný výber z rozdelenia R(0,1), tak F~1(xi), F~1(x2),... je náhodný výber z rozdelenia s (požadovanou) distribučnou funkciou Fy(y). 31 19.2 Zamietacia metoda Metoda je založená na tvrdení: Veta 19.2. Nech a > 0. Náhodný vektor U G lZn má hustotu (vzhladom k Lebesgueovej miere) /u(u) a súčasne podmienená hustota (inej) náhodnej veličiny (jednorozmernej) V, teda fv\u(v\u) Je rovnomerná na intervale (0, a f (u)) 1 ak v G (0, a/(u)), /V|uMu)= { «/u(u)' I 0, ak v i (0,a/(u)) práve vtedy, keď náhodný vektor (U7, K)' G lZn+1 je rovnomerne rozdelený na množine {(u',«)' e : ueTZn, /u(u) > 0, v G (0,a/(u))} . Dôkaz: nájdeme v knihe Devroye, L., Non-uniform random variate generation, Berlin, Springer, 1986. Algoritmus generovania postupnosti vektorov y1,y2, ••• G Ti" ktorá je náhodným výberom z rozdelenia s hustotou /y(y): - Zvolme nejakú hustotu gu(u), u G 7?™, z ktorej vieme lahko generovat náhodný výber, pričom táto hustota je taká, že existuje a > 0, že pre každé u G lZn platí /y (u) < agu (u). - Generujme náhodný výber ui,U2,... zodpovedajúci hustote gu(u). - Na i—tom kroku generujeme Xi ako (jednorozmerný) jednobodový náhodný výber z rovnomerného rozdelenia R(0, agu(ui)). Podlá Vety 19.2 to ale znamená, že body (uí,Xí) sú nezávisle generované z rovnomerného rozdelenia na množine G = {(u', x)' e lZn+1 : 0i, ale ak Xi > /y(uj), tak bod Ui vyradíme z ďalšieho uvažovania. - Je zrejmé, že takto vybrané vektory (y^o^)' sú rovnomerne rozdelené na množine F = {(ý,x)' enn+1: /Y(y)>o, o 0}, ak stacionárne rozdelenie markovovského retazca je spojte a /(•) je jeho hustota, alebo ku každému z bodov množiny {x : p(x) > 0}, ak stacionárne rozdelenie markovovského retazca je diskrétne a p(-) je jeho pravdepodobnostná funkcia. Táto skutoňost zaručí, že integrály tvaru Jx $(x)/(x)eíA(x) (resp. $(x) Xlxex P(x)) možno aproximovat sumou 33 — Y^ii=i ^(x^). Rozbor tohto problému patrí do teórie markovovských procesov (pozri napr. Robert, Ch., P. Méthodes de Monte Carlo par Chaínnes de Markov. Paris: Economica, 1996). MCMC metódy možno rozdělit do dvoch (hlavných) skupín: a) metódy odvodené od algoritmu Hastinga a Metropolisa, b) metódy odvodené od algoritmu Gibsa. 21.1. Algoritmus Hastinga a Metropolisa Opis algoritmu. Zvolíme lubovolnú triedu podmienených hustôt (resp. podmienených pravdepodob-nostných funkcií) na X {q(-\y): y e X} tak, aby platilo: a) Ak pre A C X je JA /(x)dA(x) > 0 (£xeAp(x) > 0), tak JA g(x|y)dA(x) > 0 (£xeA g(x|y) > 0) pre každé y G X. b) Ľahko možno simulovat náhodné výbery zodpovedajúce hustote (pravdepodobnostnej funkcii) q(-\y). c) existuje symetria v zmysle g(x|y) = g(y|x); x, y e X, alebo je známy analytický výraz pre funkciu (x, y) € X x X —> t/(x|y). Hustoty (pravděpodobnostně funkcie) g(- |y) sú pomocné a nazývajú sa inštrumentálne. Lubovôla pri výbere inštrumentálnej hustoty je obmedzená požiadavkou, aby príslušný markovovský proces bol ergodický. Príkladom vhodnej hustoty je t/(x|y) = g (y) nezávisle na x (pozri knižku A. Pázmana). Nech x^") G X je vektor vygenerovaný algoritmom v n—tom kroku. Vektor x(™+1) G X dostaneme v (n + 1)—vom kroku takto: 1. Pri danom x^™) generujeme vektor y(") G X zodpovedajúci hustote (pravdepodobnostnej funkcii) g(-|x(")). 2. Vypočítame výraz , , , , f /"(y^Wx^lyH) 1 fl(x(") y(")) =min J 1 j y ^ . ■ ■ l /(x("))g(y(")|x(")) J 3. S pravdepodobnosťou iž(x("),y(")) zvolíme x(™+1) = y^ a s pravdepodobnosťou 1 — R(yó-n\ y*-™-1) zvolíme x(™+1) = x^™). Poznámka. Na stanovenie iž(x("),y(")) nepotrebujeme vediet normovači koeficient v hustote /(x) (v pravdepodobnostnej funkcii p(x)). Toto využijeme ak požadovaná hustota (pravdepodobnostná funkcia) je aposteriórna, kde stačí použit súčin it{ff)f(y\ff) namiesto 7r(0|y) (resp. 7r(0i)/(y|0i), i = 1,2,... namiesto n(0i\y),i = 1,2,...) . Výber bodu x^™+1) v algoritme je náhodný a pravděpodobnost, s ktorou je vyberaný, nezávisí od prdchádzajúcich x^, i < n. Dostaneme teda markovovský proces. Vzhladom na značnú lubovôlu pri výbere inštrumentálnych hustôt sa javí prekvapujúce, že práve hustota /(x) zodpovedá stacionárnemu stavu tohto procesu. Práve táto lubovôla umožňuje vybrat inštrumentálne hustoty tak, aby sa stacionárny stav dosiahol čo najrýchlejšie, resp. výpočtovo najjednoduchšie. Označme P(A\x) pravděpodobnost, že proces sa (na lubovolnom) kroku dostane zo stavu x do stavu, ktorý patrí do množiny A. P(A\x.) voláme aj prechodové jadro (transition probability kernel). Platí veta (Pázman, A., Bayesovská štatistika, str. 65) 34 Veta 21.1. Pri akejkolvek volbě inštrumentálnej hustoty (pravdepodobnostnej funkcie) spĺňajúcej podmienku a), algoritmus Hastinga a Metropolisa generuje realizáciu markovovského procesu a pre každú me-ratelnú množinu A platí (i) ak /(x) je hustota / P(A\x)f(x)dX(x) = / /(x)cřA(x) pre každú meratelnú A C X ■JX j A (A(-) je Lebesgueova miera), (ii) ak p(xi), i = 1,2,... je pravdepodobnostná funkcia P(A\xi)p(xi) = p(x) pre každú meratelnú A C X. i>l xěA Teda /(x) je hustota stacionárneho rozdelenia (p(x) je pravdepodobnostná funkcia stacionárneho rozdelenia) tohto procesu. 21.2. Algoritmus Gibsa Tento algoritmus sa používa v mnohorozmernom prípade, ked: - X = X\ x X'2 x ... x Xm, pričom každé x e X možno rozložit na m komponent x = (x\,xm)' Xi € X{, i = 1,2,...,m. Predpokladáme, že \(Ä) = Xi(Ai) x ... x \m(Am), ak A = Aľ x ... x Am, kde A{ C X{. (Toto platí napríklad pre Lebesueovu mieru.) - pre každé i = 1,2,m sa dajú generovat náhodné výbery z množiny Xi zodpovedajúce podmieneným hustotám týchto komponentov, teda /iMí^jlj^ij^i)- Algoritmus predpisuje prechod od x^™) = (x^, ...,Xm^)' ku x(™+1) = (x\n+1\ ...,Xm+1^)' po jednotlivých komponentoch takto: - generujeme výber x^+1^ z hustoty fi(-\x^2n\ ...,x^m), - generujeme výber x;2 z hustoty j2\-\%\ ixz , ■■■,xrn ), - generujeme výber xa(n + 1) z hustoty fa(-\x\n+1\ x2n+1\ x^Xm"1), atd, - generujeme výber xm(n+ 1) z hustoty fm(-\x\n+1\x2n+1\ ...,x^^). Platí veta (Pázman, A, Bayesovská štatistika, str. 67) Veta 21.1. Postupnost x^^x^2),... generovaná algoritmom Gibsa je realizáciou markovovského procesu a /(x) = f(xi,...,xm) je hustota pravdepodobnosti stacionárneho rozdelenia tohto procesu. Dôkaz ergodičnosti príslušného markovovského procesu nájdeme napr. v knihe Robert, Ch, P. Méthodes de Monte Carlo par Chainnes de Markov. Paris: Economica, 1996. Poznámka. Poznamenávame len, že ak potrebujeme integrovat nejakú funkciu h(xi), ktorá závisí len od i—teho komponentu vektora x x xex / = j h(xi)f(x)d\(x) (I=Y1 M^)p(x)) tak túto aproximujeme sumou n k=l 35 teda z postupnosti x^1), x^2\ ... používame iba príslušnú postupnost i—tých komponentov: x\2\... . Odtial vyplýva zovšeobecnený algoritmus Gibsa: Vytvoríme nejaký priestor Z a najdeme takú hustotu (pravdepodobostnú funkciu) g(x, z) na X x Z, že - hustota /(x) je marginálna ku g(x,z), - pomocou algoritmu Gibsa aplikovaného na g(x, z) generujeme postupnost vektorov (x^1), z^1)), (x^2\ z^2)), - na aproximáciu integrálov typu / = / $(x)/(x)dA(x) (7=^ $(x)p(x)) použijeme sumu TI j = ^£*(*ífc)) *:=i a teda používame postupnost ktorá zodpovedá hustote /(x) (pravdepdobnostnej funkcii p(x)). Takýto postup použijeme, ked generovanie podlá hustoty g(x, z) je jednoduchšie, než podlá hustoty /(x) (pravdepodobnostnej funkcie p(x)). Je tu možnost rôznej volby g(x, z) a tým aj rôznych variantov algoritmu Gibsa. 22. Bayesovské bodové a intervalové odhady a testy pre jednorozmerný parameter 22.1. Bayesovské bodové a intervalové odhady pre jednorozmerný parameter (Dve nasledujúce kapitoly sledujú text z knižky J.Anděla.) Nech dimŕ? = 1. Za bodový bayesovský odhad parametra 9 sa v niektorých prípadoch berie stredná hodnota aposteriórneho rozdelenia (presnejšie za realizáciu bayesovského odhadu sa berie stredná hodnota aposteriórneho rozdelenia). V prípade príkladu z 3.kapitoly je hustota aposteriórneho rozdelenia "(Pl™) = W—Ĺ-lhZ--,Pa+m-\l ~p)h+n~m-\ 0 < p < 1 B(a + m,b + n — m) a stredná hodnota aposteriórneho rozdelenia 1 pa+m-1{í_p)b+n-m-lpdp_ B(a + m + l,b + n-m) _ o B(a + m,b + n — m) B(a + m,b + n — m) _ T(a + b + n) T(a + m + l)T(b + n - m) _ a + m ~ Y(a + m)Y(b + n-m) T(a + b + n+l) ~ a + b + n' tyl Ak je n velké (v zrovnaní s a, b), tak tento odhad je blízko hodnoty —. Toto je odhad parametra p bez ohladu na apriórne rozdelenie). Je to v zhode so skúsenostou: Ak n je velké "vela vieme" o p bez ohladu na predchádzajúce vedomosti. Ak naopak nevieme nič o p, teda nerobíme žiaden výber - m = 0, n = 0, odhad je len z predchádzajúcich vedomostí, teda —(stredná hodnota B (a, b) rozdelenia). Iná metóda je za bayesovský odhad zobrat modus aposteriórneho rozdelenia (presnejšie za realizáciu bayesovského odhadu zobrat modus aposteriórneho rozdelenia). V tomto prípade p k}, kde k je konštanta určená podmienkou n(p\m)dp = 1 — a. 22.2. Bayesovské testy v prípade jednorozmernáho parametra Budeme pokračovat v úvahách z predchádzajúcej kapitoly. Ak už máme čísla D a H zvolené tak, aby splňovali podmienku (22.2), môžeme pristúpit aj k bayesovskému testovaniu hypotéz o parametri p. Predpokladajme, že chceme testovat hypotézu H0 : p e A, kde A je nejaký interval obsiahnutý v (0,1). V prípade, že celý interval A padne mimo (D, H), hypotézu Hq zamietame. V tomto prípade je zrejme aposteriórna pravděpodobnost menšia alebo rovná a. Nebolo by však rozumné založit test len na hodnote P(A\m) (- pravděpodobnost, že náhodná veličina s hustotou aposteriórneho rozdelenia ir(p\m) sa realizuje v intervale A) bez konfrontácie s intervalom (D, H). Keby bol totiž interval A velmi krátky, mohla by byt pravděpodobnost P(A\m) velmi malá dokonca aj v prípade, že 37 by A obsahoval modus alebo strednú hodnotu aposteriorného rozdelenia. To by viedlo k zamietnutiu Hq, hoci by poloha intervalu A vzhladom k aposteriórnej hustote nebola nijako "podozrivá". Niekedy sa hovorí o teste hypotézy H0 ■ P = Po, kde po € (Oj 1) Je dané číslo. V tomto prípade A je uzavretý interval obsahujúci jediný bod pg. Je jasné, že v tomto prípade je P(A\m) = 0, takže podmienka P(A\m) < a je splnená triviálne. O zamietnutí Hq teda rozhoduje len poloha bodu po vzhladom k (D, H). Keby sme nemali nejakú podmienku vztahujúcu sa k intervalu (D, H), hypotézu H0 by sme vždy zamietli. Tu by sa dalo argumentovat, že hypotéza H0 nie je rozumná, lebo dopredu vieme, že jej pravděpodobnost je rovná nule. Ked však padne bod po mimo (D, H), existuje také jeho okolie Aq, že je disjungtné s (D, H) a P(Ao\m) = jAo ir(p\m)dp > 0. Môže sa povedat, že ide vlastne o test hypotézy, že p je z nejakého okolia bodu pg. 23. Štatistické rozhodovanie (Kapitola sleduje knižku A. Pázmana: Bayesovská štatistika.) Teória štatistického rozhodovania je založená na ekonomickom princípe, podlá ktorého má rozhodovanie maximalizovat priemerný výnos, resp. minimalizovat priemernú stratu rozhodovania. Teória sa dá formulovat aj bez použitia apriórneho rozdelenia, teda nebayesovsky. Bayesovský prístup má tiež svoje opodstatnenie. Základy teórie sa dajú sformulovat dokonca aj bez pozorovaní (t.j. bez experimentu). Takýto prístup je síce značne obmedzujúci, v teórii však umožňuje elementárne definovat základné pojmy, ktoré potom možno lahko preniest na prípad rozhodovania s experimentom. 23.1. Základné pojmy rozhodovania bez experimentu Základné pojmy si vysvetlime na elementárnom príklade: Treba rozhodnut o spôsobe výstavby objektu v prípade, že máme pochybnosti o stave podložia, pričom nemáme možnosti ho preskúmat. Z pozície stavebníka sú dva možné stavy podložia: 9\: podložie je pevné, vhodné na stavbu, 6'2: podložie je slabé a stavat sa môže len so zosilnenými základmi. Rozhodovatel (stavebník) má volit medzi tromi rozhodnutiami: cti: nestávat vôbec, 0,2'. stavat štandardným spôsobom, 0,3: stavat so zosilnenými základmi. Označme L(9, a) finančnú stratu, ktorú utrpí stavebník ak podložie je v stave 9 a je prijaté rozhodnutie a. Hodnoty L(9, a) sú v nasledujúcej tabulke (záporná strata = zisk) a\9 6»! 92 ai 5 000 0 a2 -30 000 50 000 a3 3 000 -20 000 Rozhodnutie 0,3 je lepšie ako rozhodnutie a\ pri akomkolvek stave podložia. Rozhodnutia a\ a ŕi2 sú neporovnatelné, kým nemáme apriórne (alebo nejaké experimentálne) informácie o stave podložia. Chápeme to v 38 takom zmysle, že keď je podložie pevné (9i) lepšie je rozhodnut sa pre ŕi2 (stavat štandardným spôsobom) a keď je podložie slabé (#2) lepšie je, keď sa rozhodneme pre a\ (nestávat). Ak ale vieme s dost velkou apriórnou pravdepodobnostou, že podložie môže byt v stave 62 (slabé), uprednostníme aj rozhodnutie pred ŕi2- V prípade, že nemáme apriórne informácie o stave podložia, rozhodnutie a\ možno z ďalšieho rozhodovania vylúčit. Zovšeobecnime tieto úvahy: Nech íl je množina stavov prírody alebo (prozaickejšie) íl je parametrický priestor. Nech 2Í je množina možných rozhodnutí. Nech L : íl x 2Í —> 1Z je stratová funkcia, o ktorej budeme předpokládat, že je zdola ohraničená. Rozhodnutie ai považujeme za ekvivalentné s a2 (píšeme ai ~ 32), ak pre každé 0 G íl platí L{0,3\) = L(0,3.2)- Podobne ai je rovnomerne nie horšie ako a2 (píšeme ai ^ sľž), ak pre každé 0 G íl platí L{0,3\) < L{6,3.2). Ak však naviac platí, že existuje 0* G íl také, že L{0*,&i) < L{6*, 32), potom ai je rovnomerne lepšie ako a2 (píšeme ai -< 32). Zrejme relácia ^ definuje čistočné usporiadanie množiny rozhodnutí 2Í. Nemusí existovat rovnomerne najlepšie rozhodnutie (pozri predchádzajúci príklad). Môžeme však použit dôležité náhradné pojmy, a síce pojem prípustného rozhodnutia a pojem úplnej množiny rozhodnutí. Definícia 23.1. Rozhodnutie a G 2Í je prípustné, ak neexistuje také b G 2Í, že b -< a. Množina rozhodnutí C C 2Í je úplná, ak ku každému rozhodnutiu a G 2Í — C existuje rozhodnutie b G C, že b -< a. Označme symbolom V množinu všetkých prípustných rozhodnutí. Rozhodnutia, ktoré nie sú prípustné, možno z rozhodovania vylúčit. Podobne, ak sa podarí nájst úplnú množinu C, možno vylúčit rozhodnutia patriace do 2Í — C. Množina Cq C 2Í sa nazýva minimálne úplná, ak Cq je úplná a je podmnožinou každej úplnej množiny. Veta 23.1. Ak V je úplná, tak sa rovná minimálnej úplnej množine. Dôkaz: Ak C je úplná a a ^ C, tak a ^ V, pretože existuje b G C, že b -< a. Teda V C C. Odtial vyplýva, že ak je V úplná, je aj minimálne úplná. Jft Ak teda "P je úplná možno množinu rozhodnutí 2Í redukovat na množinu V a ďalšia redukcia už nie je možná. Ak máme k dispozícii apriórnu hustotu tt(9), stredná strata pri rozhodnutí a je / L(0,3)ir(0)d\(0) resp. ^ L{6X, a)7r(02). Jn i>i Rozhodnutie a^ G 2Í možno považovat za optimálne (vzhladom na túto strednú stratu), ak a, = argmin/L(MWW) resp. a. = arg min £ L(0t, 3M0t). Takéto rozhodnutie sa volá bayesovské. V súlade s tým, čo sme povedali o prípustných rozhodnutiach, platí, že ak množina V je úplná, tak rozhodnutie a^ je alebo prípustné, alebo existuje prípustné rozhodnutie b^ -< a^, pričom platí í L{6,3lľ)-n{6)d\{6)= í L(0,b7r)7r(0)dA(0) Ja Ja (a analogicky aj v diskrétnom prípade). Návod na dôkaz tohto tvrdenia je v knižke A. Pázmana, str. 74. 23.2. Rozhodovanie na základe experimentu 39 V príklade v kapitole 23.1 je rozumné, aby si stavebník, prv než sa rozhodne stávat, dal preskúmat podložie. To znamená, že príslušní odborníci vykonajú merania. Výsledkom týchto meraní je vektor údajov y. Tento je samozrejme ovplyvnený náhodnými efektmi, ktoré v závislosti od stavu podložia 0 dajú vystihnut hustotami pravdepodobnosti f{y\0) resp. príslušnými pravdepodobnostnými funkciami. Prv než sa vykoná rozhodnutie v rozhodovacom probléme, je rozumné vykonat rôzne pozorovania, ktoré súhrnne nazývame experiment. Matematicky sme experiment charakterizovali triedou hustôt (alebo pravde-podobnostných funkcií) {/(y|0): e e íl}. Ulhou štatistika je nájst rozumné pravidlo, ako každému výsledku experimentu y G y(€ž £>„) priradí niektoré rozhodnutie a G 2Í, t.j. nájst vhodné zobrazenie A: y € ^ a € 21. Toto zobrazenie sa nazýva rozhodovacia funkcia. Ak máme danú rozhodovaciu funkciu A, tak rizikom, alebo rizikovou funkciou nazývame funkciu (0, A) G íl x D -> r(0, A)= f L(8, A(y))/(y|0)dA(y) ■Jy y ktorá vyjadruje strednú stratu pri stave 0 a pri rozhodovaní podlá rozhodovacej funkcie A. Symbolom D tu značíme množinu rozhodovacích funkcií dovolených v danom rozhodovacom probléme. Základnou požiadavkou na každú A G -D je, aby existoval uvedený integrál. Ak ten čo rozhoduje chce minimalizovat rizikovú funkciu, ocitá sa formálne v tej istej situácii, ako v prípade rozhodovania bez experimentu. Namiesto trojice (fi,2l,L(0,a)), ktorá sa uvažovala v kapitole 23.1, má teraz trojicu (íl, D, r(0, A)). Mnohé pojmy možno teda priamo preniest z kapitoly 23.1. Rozhodovacie funkcie sú vo vztahu Ai ^ A2, ak pre každé 0 G íl platí r(0, Ai) < r(0, A2), čím definujeme rovnomerné usporiadanie rozhodovacích funkcií. Rozhodovacia funkcia Ai je prípustná, ak neexistuje iná rozhodovacia funkcia A2, ktorá by bola rovnomerne lepšia, teda neexistuje A2 pre ktorú platí A2 -< A\. Podobne ako v 23.1 definujeme úplné množiny rozhodovacích funkcií a platia presne tie isté všeobecné súvislosti medzi prípustnostou a úplnostou. Pretože usporiadanie rozhodovacích funkcií reláciou ^ je čiastočné, vo všeobecnosti neexistuje rovnomerne najlepšia rozhodovacia funkcia. Ak však poznáme apriórnu hustotu tt(9), porovnávame rozhodovacie funkcie na základe stredného rizika (alebo strednej straty), teda na základe integrálu (23.1) R(A) = / r{0,A)-K{0)dv{0) L(0,A(y))/(y|0)dA(y) ir(0)dv(0). >a J n i-Jy Bayesovskou rozhodovacou funkciou voláme tú rozhodovaciu funkciu, pre ktorú platí Ajr = arg min R(A). 40 Tu však končí analógia s rozhodovaním bez experimentu. Pri rozhodovaní bez experimentu minimalizujeme integrál Jn L(0, a)ir(9)d\(9) resp. ^i>1 L(0i, a)7r(0j) vzhladom na prvky množiny 21, teraz potrebujeme minimalizovat dvojný integrál Jn L(0, A(y))/(y|0)cřA(y)J ir{0)dv{0) vzhladom na zobrazenia A, ktorých oborom hodnôt je mnžina 21 (množina možných rozhodnutí), teda vzhladom na podstatne komplikovanejšiu štruktúru. Ukazuje sa však, že za pomerne všeobecných podmienok táto minimalizačná úloha má elegantné riešenie, ktoré je dané v nasledujúcej vete. Veta 23.2. (Základná veta o bayesovských rozhodovacích funkciách.) Nech stratová funkcia L(0,a) je meratelná a zdola ohraničená. Nech D je množina všetkých meratelných zobrazení z J-' do 21. Nech pre skoro všetky y €ž y existuje riešenie úlohy (23.2) Aw(y) = arg min / L(0, a)it{6\y)dv{6), ktoré je meratelnou funkciou y. Potom toto riešenie je bayesovskou rozhodovacou funkciou pri apriórnej hustote tt(0). Dôkaz nájdete v knižke A. Pázmana na str. 76. Poznamenávame len, že v skutočnosti nieje ani potrebné určit celú bayesovskú rozhodovaciu funkciu AT. Ak zrealizujeme experiment a jeho výsledok je y, stačí riešit minimalizačnú úlohu (23.2) pre túto jedinú hodnotu y. Často sa popri bayesovských rozhodovacích funkciách uvažuje aj tzv. minimaxná rozhodovacia funkcia Am definovaná vztahom Am(y) = arg min < maxr(0, A) > , wj 6 ai lJy 'y Podlá Vety 23.2 bayesovská rozhodovacia funkcia pri danom y je Aw(y) = arg min | - ^ L(02, P(-)M02|y) j = arg min |-^(lnF^W^Iy) j = = arg mini -^(lnP(0í))7r(0í|y) + ^(ln7r(0í|y))7r(0í|y) l = arg min/(ir(-|y), P(-)). P('' [ i>i i>i J P^'') Minimalizujeme teda /—divergenciu. Podlá Poznámky pod Definíciou 9.2 toto minimum sa dosahuje ak P(-) a 7r(-|y) sa zhodujú. Teda bayesovská rozhodovacia funkcia sa v tomto prípade rovná A(y) =vr(-|y). Tento výsledok je zaujímavým informačným zdôvodnením používania Bayesovho vzorca. Žiadne iné rozdelenie pravdepodobnosti na í~ž nepostihuje tak dobre informáciu získanú z experimentu ako práve aposteriórne rozdelenie. 24. Bayesovské odhady a testy z hladiska teórie štatistického rozhodovania Teória štatistického rozhodovania, ktorá vznikla ako ekonomicky motivovaná teória rozhodovania využívajúca experiment (pozorovania, merania), ovplyvnila myslenie v teoretickej štatistike, hlavne v teórii odhadu. 24.1. Bayesovské bodové odhady z hladiska teórie štatistického rozhodovania Ak v experimente {/(y|0): een} zvolíme 2Í = íž a ak stratovú funkciu L(0,a) zvolíme tak, aby vyjadrovala odchýlku odhadu od skutočnej hodnoty parametra, tak bayesovské rozhodovacie funkcie sú vlastne (bayesovskými) odhadmi parametra 0 (teda odhadovacími štatistikami). Takými stratovými funkciami sú napríklad L(0,a) = ||a-0||2 42 alebo L(0,a) = ]TK Bayesovské odhady sa porovnávajú s "klasickými" odhadmi frekvenčnej štatistiky, hlavne s odhadom maxima vierohodnosti. Naopak, pri "klasických" odhadoch sa overujú vlastnosti formulované v teórii rozhodovania, ako sú přípustnost (admisibilita), invariantnost, optimálnost vzhladom na niektorú stratovú funkciu. Je vhodné specializovat základnú vetu o rozhodovacích funkciách (Vetu 23.2) na uvedené stratové funkcie. Veta 24.1. Nech S je symetrická pozitívne definitná matica a nech L(0,a) = (0-a)'S(0-a). Nech 7r(0) je apriórna hustota. Potom bayesovská rozhodovacia funkcia (t.j. bayesovský odhad) sa rovná podmienenej strednej hodnote 0 pri danom y (í In91f(9\y)n(9)dv(9) \ Ax(y) = áľ(0|y) /n/(0|yM0) £(0|y) je jednou z bayesovských rozhodovacích funkcií. Stredná strata pri tejto rozhodovacej funkcii sa rovná £ {(& - áľ(0|y))'S(0 - £(0|y))} = tr {S £ (cov(&\y))} , kde cov(&\y) je kovariančná matica náhodného vektora 0|y. Dôkaz: nájdete v knihe A. Pázmana, str. 82. Poznamenávame len, že ak vo vete 24.1 volíme S = I, tak L(0,a) = ||a-0||2. Veta 24.2. Nech u(-) je Lebesgueova miera na íž = lZm. Nech m L{9,Sí) = YJ\ar-0l\. Nech 7Tj(0j|y) je marginálna aposteriórna hustota parametra 9i *i(0i|y)= / 7r(0|y)d0i...d0i_id0i+i...d0n Ju™-1 a nech /ii(y) je medián tejto hustoty. Potom y e y -> (Mi(y),...,Mm(y))' je bayesovská rozhodovacia funkcia (bayesovský odhad). Dôkaz nájdete v knihe A. Pázmana na str. 83. Používaným odhadom v bayesovskej štatistike je aj odhad 0(y) = arg max 7r(0|y), 43 ktorého realizácia je 9(y) = arg max ir(0\y). Z bayesovho vzorca resp. z (2.7) vyplýva, že ak hustota tt(O) je konštantná, tak odhad maxima aposteriórnej hustoty je totožný s odhadom maxima vierohodnosti ktorého realizácia je 0(y) = arg max ln/(y|0), 9(y) = arg max ln/(y|0). Odhad 0(y) sa preto nazýva bayesovským odhadom maxima vierohodnosti. Podlá (2.7) možno písat 9(y) = arg max [ln/(y|0) +ln7r(0)] a preto sa tento odhad nazýva aj penalizovaným odhadom maxima vierohodnosti, pričom ln7r(0) je penalizácia. Dá sa ukázat, že odhad 0(y) nedostaneme zo žiadnej stratovej funkcie v zmysle kapitoly 23, môžemeho však vyjádřit ako limitu bayesovských rozhodovacích funkcií (pozri knihu A. Pázmana, str. 84). 24.2. Bayesovské intervalové odhady Vo frekvenčnej štatistike sa ako intervalové odhady parametra 0 používajú oblasti spolahlivosti. Sú to náhodné borelovské množiny, ktoré s predpísanou pravdepodobnosťou pokrývajú (pevný) parameter 0. Ich konštrukcia je spojená s konštrukciou optimálnych testov a môže byt velmi komplikovaná. V bayesovskej štatistike je situácia jednoduchšia. Používajú sa oblasti ohraničené krivkou konštantnej aposteriórnej hustoty 0(y) = {6 e í. : vr(0|y) > c(y)}, kde c(y) je zvolené tak, aby platilo (24.1) / 7r(6>|y)d^(6>) = 1 - a, j O (y) kde 1 — a je predpísaná spolahlivost. Teda aposteriórna pravděpodobnost toho, že 0 g C(y), je práve 1 — a. V prípade, že sa rovnost nedá realizovat, volíme za c(y) supremum z tých čísel, pre ktoré platí 7r(6>|y)d^(6>) > 1 - a. 0(y) Niekedy sa oblast 0(y) nazýva HDP oblast (highest posterior density region). Zrejme platí nasledujúca veta Veta 24.3. Nech y g 3^ je dané a 7r(0|y) je aposteriórna hustota vzhladom na mieru ľ(9) a nech pre W c ti platí / 7r(6>|y)d^(6>) = 1 - a. Jw Potom v(0(y)) < v(W). 24.3. Bayesovské testy 44 Vhodným výberom priestoru rozhodnutí 2Í a stratovej funkcie L(0,a) dostaneme rozhodovacie problémy podobné úlohám testovania hypotéz vo frekvenčnej štatistike. Podstatný rozdiel je v tom, že nemôžeme uvažovat hypotézy, ktorých apriórna pravděpodobnost je nulová (lebo ich aposteriórna pravděpodobnost je tiež nulová) a tiež v tom, že optimalita testov sa posudzuje úplne iným spôsobom než vo frekvenčnej štatistike. 24.3.1. Jednoduchá hypotéza a jednoduchá alternatíva Nech í~ž = {Oq, Oi}. Priestor rozhodnutí bude dvojprvkový 2Í = {a0,ai}, kde ao znamená, že prijímame hypotézu 0 = Oq a ai znamená, že prijímame alternatívnu hypotézu 0 = 0\. Stratová funkcia je L(ŕ?o,ai) = wq > 0, L(#i,ao) = wi > 0 L(#o,ao) = ^(^i,ai) = 0- Apriórne rozdelenie je diskrétne tt(Oq) > 0, ir(9i) = 1 — tt(Oq) > 0. Experiment je daný hustotami /(y|#o), /(y|#i) vzhíadom na (cr—konečnú) mieru A(-). Môžu to byt "obyčajné" hustoty vzhíadom na Lebesgueovu mieru A(-), alebo dve pravděpodobnostně funkcie {/(yi|0o)}i>i, {/(yil^i)}i>i " vtedy je miera A(-) sčítacia. Bayesovská rozhodovacia funkcia (teda bayesovský test) sa v súlade s Vetou 23.2 rovná A(y) = arg min [L(0O, a)7r(0o|y) + L(9U a)-7r(e1 |y)] = a€{a0,ai} Íai, ak wO7r(0o|y) < Wi7r(0i|y), a0, ak wivr(0i|y) < wO7r(0o|y). Nulovú hypotézu zamietame, ak prijímame alternatívnu hypotézu, teda ak A(y) = ai a to je práve vtedy ak (24.2) w0ir(00\y) < iunr(ei|y). Podía (2.7) je 7r(0|y) = —-2—^ —1-; teda (24.2) je splnené práve vtedy ak w17t(e1) c- f(y\o0y Porovnáme tento test s "klasickým" testom, preto napíšme oblast zamietania nulovej hypotézy 0 = Oq ako kde 24.3 c = 0 . Wl7r(6'l) Oblast zamietania nulovej hypotézy je určená podielom vierohodností, podobne ako v Neymannovej-Pearsonovej leme. Zopakujme si ju Veta 24.4. (Neymanova-Pearsonova lema.) Nech k danému a G (0,1) existuje také kladné číslo g, že pre množinu Wg = {y. /(y|0i)>s/(y|0o)} platí (24.4) / f(y\e0)d\(y) = a. ■Jwa 45 Potom pre lubovolnú množinu W G Bn splňujúcu podmienku (24.5) / f(y\e0)d\(y) = a Jw platí /(y|0i)dA(y) > / /(y|0i)dAy). Wg JW Dôkaz nájdeme v knižke J. Anděla na str. 239. Treba si uvědomit, že oblastou zamietania nulovej hypotézy Wg je určený test s hladinou významnosti (pravdepodobnostou chyby 1. druhu) a (t.j. pravděpodobnost že test zamieta nulovú hypotézu, ked ona platí, je a). Pričom tento test je najsilnejší spomedzi všetkých testov s hladinou významnosti a. Všetky testy s hladinou významnosti a sú určené oblastami zamietnutia W, pre ktoré platí (24.5). Chyba druhého druhu je taká, ked test nezamieta nulovú hypotézu, pričom platí alternatívna. Test určený oblastou zamietania W má pravděpodobnost chyby druhého druhu /(y|0i)dA(y) = l-/3 y-w a P je sila testu. Podia Neymanovej-Pearsonovej lemy je oblast zamietania pre najsilnejší test s hladinou významnosti a daný ako W° = V: TylT' 'u)- pričom konštantu g počítame zo vzorca (24.4). Pri bayesovskom testovaní naproti tomu namiesto chýb 1. a 2. druhu máme aposteriórne pravdepodobnosti oboch hypotéz 7r(ŕ?o|y), ""(^i|y) a strednú stratu (rizikovú funkciu, pozri (23.1)) rovnú R(A)=tt(00) / L(e0,A(y))f(y\e0)d\(y) + 7t(e1) / L(9U A(y))/(y|01)dA(y). Jy Jy Počítajme tt(0o) / L(0o,A(y))/(y|0o)dA(y) = Jy \(00)JyL(00,ai)f(y\00)d\(y), ak w0Tr(00\y) < w17r(e1\y) n(60)fyL(60,ao)f(y\60)d\(y), ak wMO^y) < w0n(00\y) _ (tt(00)w0 Jy f(y\e0)d\(y), ak y e y - W(c) I 0, inokedy = 7r(0oM, / f(y\00)d\(y) = tt(00)w0 Jy-w(c) kde a počítame podlá (24.5). Podobne f(y\e0)d\(y) - / f(y\e0)d\(y) y Jw(c) 7r(0o)wo(l - a), 71(0!) / L(01,A(y))f(y\01)d\(y) = Tr(01)w1f3, Jy kde p je sila testu, t.j. Jw^ /(y|^i)^A(y). Konštantu c počítame z (23.4). Pre strednú stratu vyššie uvedeného bayesovho testu dostávame R(A) = n(60)w0(l -a)+ 7r(0i)wi/3 = ir(0i)iui(c(l - a) + P). 46 24.3.2. Súčasné testovanie niekolkých jednoduchých hypotéz Na rozdiel od frekvenčnej štatistiky, v bayesovskej štatistike nerobí problémy rozhodovanie medzi niekolkými hypotézami. Nech 3^ = {0i, @k}, 2Í = {ai,a^}, kde a^ znamená, že prijmeme hypotézu 0 = Oj. Predpokladáme, že tt(Oj) > 0 pre každé j. Ďalej jwij > 0, ak i ŕ j, I 0, ak i = j. Potom A (y) = arg min V" L(0i; a)7r(0i|y) = a€{ai,...,afe} -f—* k k ak ^Wijir(0i\y) < ^wwir(0i\y) pre každé l. 24.3.3. Testovanie zloženej hypotézy proti zloženej alternatíve Nech íl = n0 U fii, kde íl0 n fii = 0. Nech 2Í = {a0,ai}. Nech /n ir{0)dv{0) > 0, i = 0,1. Stratovú funkciu volíme v princípe tak, aby {0, ak 0 g n,-, > 0, ak 0 £ flj. Obvyklé volby stratovej funkcie pre 0 g ííj sú L(0, aj) = > 0, ak i/ j, alebo L(0,a1) = K1 min 110-0*11. V zmysle Vety 23.2 hypotézu 0 g ííi prijmeme (t.j. hypotézu 0 g ÍIq zamietneme), keď í L(0,ai)7r(0|y)^(0) < í L(0,ao)7r(0|y)^(0). Takéto testovanie zloženej hypotézy proti zloženej alternatíve je typické pre diskriminačnú analýzu. (Podrobnejšie pozri v knižke J. Anděla str. 319-323.) Tvorba bayesovských testov môže byt jednoduchšia ako tvorba testov vo frekvenčnej štatistike.