M5VM05 Statistické modelování 9. Zobecněné lineární modely Jan Koláček (kolacek@math.muni.cz) Ustav matematiky a statistiky, Přírodovědecká fakulta, Masarykova univerzita, Brno podzim 2013 Jan Koláček (PřF MU) M5VM05 Statistické modelování podzim 2013 1/57 Motivace V reálném světě má mnoho procesů jiný, než lineární vztah závislosti. Např. v ekonomii se ukazuje, že mnoho vztahů má logaritmickou závislost, k vysvětlení procesů v přírodních vědách se užívají reciproké, mocninné i další vztahy. Vysvětlovaná veličina popisující pravděpodobnost přežití člověka, v případě určité nemoci a určitého způsobu léčby, může z definice pravděpodobnosti nabývat hodnot pouze z intervalu [0,1], což by v případě klasického lineárního modelu bylo možné zajistit jen za přijetí určitých omezení na parametry modelu. Také normalita chyb je často nesplněným předpokladem klasického lineárního regresního modelu. Připomeňme, že normalita se vyznačuje nezávislosti střední hodnoty a rozptylu. Typicky např. u ekonomických veličin s rostoucí střední hodnotou obvykle roste rozptyl náhodné veličiny, přičemž náhodné chyby mají v těchto případech často nesymetrická, kladně sešikmená rozdělení. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 2/57 Základní pojmy a definice I Definice 1 Mějme parametrický prostor G C K"1. Řekneme, že systém m-parametrických hustot TZg = {f(y;9):0=(01.....6mf e 0} je regulární, jestliže platí (1) O C IR"1 je otevřená borelovská množina. (2) Množina M = {y G R" - f (y',6) > 0} nezávisí na parametru 6. (3) Pro každé y E M existuje konečná parciální derivace /ľ(y;0) = ä^ (í = i.....«)• Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 3/57 Základní pojmy a definice II Definice 1 (4) Pro všechny d = (6lr..., 0m)T E 0 platí kde F(y;0) je odpovídající distribuční funkce. (5) Pro všechny 6 = (8i,... ,8m)T E O je integrál ain/(y;0)ain/(y;0) M 30; de. dF(y;0) i,j=l,...,m konečný a matice J = J(0) = (/!y(0))^-_1 je pozitivně definitní. Matice J se nazývá Fisherova informační matice o parametru 6. Jan Koláček (PřF MU) M5VM05 Statistické modelováni podzim 2013 4 / 57 Základní pojmy a definice Definice 2 Nechť / E ^řeg- Pak náhodný vektor U = U(0) = (U1(6),...,Um(0))T se složkami Ui = Ui(6) = ain/(Y;0) d6i se nazývá skórový vektor příslušný hustotě/. Jan Koláček (PřF MU) M5VM05 Statistické modelováni podzim 2013 5 / 57 Základní pojmy a definice I Věta 3 (1) Je-li f E F™eg 3 pro i,j = 1,..., m existují pak fí!(r,e) = EV (d) = 0 a d2f(y;Q) DU(0)=J(0) Jan Koláček (PřF MU) M5VM05 Statistické modelováni podzim 2013 6 / 57 Základní pojmy a definice II Věta 3 (2) Platí-li navíc pro i,j = 1,... ,m = o, pak kde J(0) = -E(U'(0)), U'(fl) = au,-(e) 30; y=i Jan KoláCek (PřF MU) M5VM05 Statistické modelováni podzim 2013 7 / 57 Základní pojmy a definice Uvažujme náhodný výběr Y„ = (Yi,..., Yn)T z rozdělení/ e J7^. Označme M = {y £ E :f(y;d) > 0}. Pak sdružená hustota Značení: funkce: n. vektory: U* u* matic, fce: J /v„(y;0) = flf(yi, ô), y = (yi.....y«)' e f* = K = E(0;y) u*(0) U*(0) J(0) J«(0) ln/(y*;0) in/yB(y;0) f31n/(Yt;fl) i,-39;—" ^31njV,(Y;fl) V 331 ' 31n/(Yt;fl)\T de„, J 31n/Y„(Y;8)^T * ' 3 ŕ? m Jan KoláCek (PřF MU) ■••ÍMáhy)á'ny^v(y^));; ,V1 vl5VM05 Statistické modelovaní i 2013 Základní pojmy a definice I Věta 4 Uvažujme náhodný výběr Y„ = (Yi,..., Yn)T z rozdělení s hustotou f E F™eg. (1) Pokud pro i, j = 1,... ,tn existují fíi(y;e) = pak EU*n(0) =0 d2f(y;Q) DU;(0) = n](d) . Jan Koláček (PřF MU) M5VM05 Statistické modelováni podzim 2013 9 / 57 Základní pojmy a definice II Věta 4 (2) Platí-li navíc pro i,j = 1,... ,m o, (tj. f je regulární i v 2. derivacích), pak E(U£'(0)) = -nj(fl), /rete dU*{9)\ U*(0) 30; Jan KoláCek (PřF MU) M5VM05 Statistické modelováni podzim 2013 10 / 57 Základní pojmy a definice I Věta 5 Mějme náhodný výběr\n = (Yi,... ,Yn)T z rozdělení s regulárni hustotou f £ ?7eg- Označme M = {y 6 R :/(y; 0) > 0}. Necht pro všechna y E M, 0 E 0 a i,j = 1,..., m existují druhé parciální derivace hustoty f [y; 6). (1) Pak platí A nm(oj(0)). Dále platí -u;(0)Tj(0)-1u;(0) A X2(m). Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 11 / 57 Základní pojmy a definice II Věta 5 (2) Platí-li navíc, že f je regulární i v 2.derivacích, tj. J (y; e) o, pak matice náhodných veličin Vm- 1 au'm\ 1 (d2ln{d;Y) n y=i s.; -J(0)- Jan Koláček (PřF MU) M5VM05 Statistické modelováni podzim 2013 12 / 57 Základní pojmy a definice Definice 6 (a) Věrohodnostní funkcí rozumíme funkci vektorového parametru 6 Wy)=f(r,o) (b) logaritmickou věrohodnostní funkcí nazýváme funkci Z(0;y) =lnL(0;y) (c) Řekneme, že odhad 0mle = 0MLE(Y) je maximálně věrohodný odhad (MLE) vektorového parametru 6, pokud platí L(6MLE;\) > L(6;\) pro všechna 6 G G. Jan KoláCek (PřF MU) M5VM05 Statistické modeloval podzim 2013 13 / 57 Základní pojmy a definice Věta 7 Mějme náhodný výběr\n = (Yi,...,Yn)T z rozdělení s regulární hustotou f e ?Ťeg- Označme M = {y E~R :f(y;6) > 0}. Nechi pro všechna y e M, 6 E © a i,j = 1,..., m existují druhé parciální derivace hustoty f [y; 6) a platí (1) ^(0mle-0) ~ NmfrW)-1) (2) W = (?mle - 0)TnJ(0)(0mle - 0) ~ x2{m) , tzv. Waldova statistika. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 14 / 57 Základní pojmy a definice Definice 8 Řekneme, že pozorovaní pochází z rozdělení exponenciálního typu, pokud jeho pravděpodobnostní funkce (v případě diskrétních rozdělení) či hustota (v případě spojitých rozdělení) je tvaru f(y)=exp{a(y)b(8)+c(8)+d(y)}, kde 8 je (neznámý) tzv. přirozený parametr a a{y),b{8),c{8),d{y) jsou známé funkce. Pokud a a(y) = y, říkáme že pravděpodobnostní funkce, popř. hustota je v kanonické formě. • v konkrétním rozdělení figurují další neznámé parametry, nazveme je tzv. rušivými parametry. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 15 / 57 Základní pojmy a definice V dalším budeme uvažovat pouze regulární a kanonické formy spolu s podmínkou b(6) = 8 a přitom zavedeme do označení jeden rušivý parametr

0, d (y) jsou známé funkce, a pokud i/>(<^>) = > 0, <^> > 0 je tzv. faktor měřítka (scale factor) oj > 0 je známá apriorní váha. Tato forma se také nazývá škálovou formou hustoty exponenciálního typu Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 Základní pojmy a definice Věta 9 Mějme náhodnou veličinu Y z rozdělení s regulárni hustotou f exponenciálního typu: /(3/) = exp{^P+%,<ř)}. (1) Pak EY = i\Q) Necht navíc platí kdef"(Y;6) = &&?±, pak de2 DY = i'(6)ip( 0. f(y) = V2 : exp 1 fy-f 7(8) = exp < 2^ ii/2 i, ,^ ( = 1. Jan Koláček (PřF MU) M5VM05 Statistické modelováni podzim 2013 19 / 57 Základní pojmy a definice Příklad 11 (Binomické rozdělení) Mějme Z~Bi(n,n), n 6 N,7r 6 (0,1). pak fz(z) = ©7Tz(l-7T)"-z = exp{zln(T^)+nln(l-7r)+ln©} pro z = 0,... ,n, přičemž EZ = ]i = nn a DZ = nn(í — n). Pravděpodobnostní funkce není ve škálové formě, proveďme reparametrizaci = ln 71 1 - 71 = ln nn n — nn = ln n — ]i n = 1+é a 1 — n = 1 1 + e9' Jan KoláCek (PřF MU) M5VM05 Statistické modelováni podzim 2013 20 / 57 Základní pojmy a definice Tedy d(y,(f>) ^ 7(0) =nln(l + ee) m = í = 1 ĺú = 1 \,f>2. jsou neznámé parametry a x,- jsou známé kovariáty. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 29 / 57 Příklad Jan Koláček (PřF MU) M5VM05 Statistické modelování podzim 2013 30 / 57 Příklad Jestliže Y,- ~ G(a,/3(- = jsou pro i= í,...,n nezávislé náhodné veličiny (EY(- = Hi = a/3,), g(}ii) = \tí}Ií = + fÍ2xi Je logaritmická linkovací funkce, pi,p2 a a = ^ Jsou neznámé parametry (a je rušivý parametr) a x,- jsou známé kovariáty. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 31 / 57 Příklad -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 X Obrázek : Ukázka GLM modelu s linkovací funkcí g(ji) = \nfi pro náhodnou veličinu Y s gamma rozdělením. Jan Koláček (PřF MU) M5VM05 Statistické modelování podzim 2013 32 / 57 Příklad Příklad 15 Poissonovská regrese. Yi ~ Po(m) jsou pro i = 1,... ,n nezávislé náhodné veličiny (EY{ = g(Hi) = \nni = fii+frxi je logaritmická linkovací funkce, f>\,f>2. jsou neznámé parametry a x,- jsou známé kovariáty. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 33 / 57 Příklad Příklad Příklad 16 Binomická regrese: Y i ~ Bi(rii, ni) jsou pro i = 1,... ,n nezávislé náhodné veličiny, kde je logistická li n kovací funkce, fii, f$2 Jsou neznámé parametry a x,- jsou známé kovariáty. Například ve farmaceutickém experimentu může být počet pacientů, kterým byla podána dávka x,- nového léku a Y i počet pacientů dávající pozitivní odpověď na danou dávku x,- nového léku. y. Jestliže pozorujeme, že roste spolu s x,-, hledáme model, ve kterém 7i(- je funkcí x,-, hodnot 0 < 7T(- < 1. Proto model TCj = fÍ2xi nem vhodný, avšak /3i + /32X,- = ln í j obvykle pracuje dobře. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 35 / 57 Příklad Jan Koláček (PřF MU) M5VM05 Statistické modelování podzim 2013 36 / 57 Praktický příklad Příklad 17 V souboru „motak.Rdata" jsou uložena data o lovu tetřeva dravcem jménem Moták pilich (Circus cyaneus) v závislosti na výskytu tetřeva. Označme Y, procento zkonzumovaných tetřevů a x,- počet tetřevů v dané oblasti. Teorie zabývající se chováním těchto dravců navrhují k modelování použít vztahu e(Y) = m = txxf ô + xf kde Y i má Gamma rozdělení. Je tedy třeba odhadnout neznámé parametry cc a S. Užitím linkovací funkce inverse dostáváme 1 _ 1 3 Definování nových parametrů fio = l/tx.afii=S/tx. dostáváme lineární vztah 1 „ „ 1 .3 ' Jan Koláček (PřF MU) M5VM05 Statistické modelováni Praktický příklad Konzumace tetřeva motákem 40 60 80 100 120 počty tetřeva Obrázek : Aplikace Gamma regrese s linkovací funkcí g(]i) — na data motak. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 38 / 57 Odhady neznámych parametrů v GLM Všimněme si, že rozdělení náhodných veličin Y,- jsou stejného typu a logaritmus sdružené věrohodnostní funkce má tvar f (0;y) = tm-,Vi) = t ■ Odhad neznámych parametru metodou maximální věrohodnosti dostaneme řešením rovnic typu Podle věty 5 konverguje matice druhých parciálních derivací skoro jistě k matici —]n, která je při regularitě systému hustot negativně definitní. Jan KoláCek (PřF MU) M5VM05 Statistické modeloval podzim 2013 39 / 57 Řešení věrohodnostních rovnic I Věta 18 Mějme náhodný výběr Y = (Y\,..., Y„ )T, který se řídí zobecněným lineárním modelem s lín kovací funkcí g(m) = *i P = Vi i=l,...,n. Předpokládejme, že pro i = 1,... ,n existují příslušné derivace 7'(0;)'7"(^i) 3 platí EYt = m = 7'(6i) DYj = 7"(6i)Uf)- Pak u; = u;^) = Íx-^-^d^ (6) Jan KoláCek (PřF MU) M5VM05 Statistické modelování podzim 2013 40 / 57 Řešení věrohodnostních rovnic II Věta 18 což lze zapsat maticově U* = U£(j8) = (Lij.....U*k)T = XTWQr (8) Jn = J„(j6)= (/í)*.=i=XTWX, (9) kde r = (ri(j8).....r„(/3))T r,- = Y,- - # = Y,--g~L{x]p) 1 / 3^;" W = í/Íflg{H7i(j8),...,H7„(j8)} H7," m. ydfJ. 1 /„T , 2 Jan Koláček (PřF MU) M5VM05 Statistické modelování podzim 2013 41 / 57 Řešení věrohodnostních rovnic Řešíme tedy věrohodnostní rovnice dl* _ _ w = - = r^- = yx^Yi-1li)^ = o /=i k Ty nejsou lineární vzhledem k neznámým parametrům, musí se řešit numerickou iterací. Jan KoláCek (PřF MU) M5VM05 Statistické modelováni podzim 2013 42 / 57 Newton-Raphsonova metoda Chceme-li najít řešení systému nelineárních rovnic U£(/3) = 0, lze použít následující iterativní postup: O Nejprve provedeme linearizaci pomocí Taylorova rozvoje v okolí bodu /30, kde /S0 je nějaký počáteční odhad: U*(/S) » U*(/S0) + U*'(/S0)(/S -/S0). Protože U£(/3) = 0, pak po jednoduchých úpravách dostaneme /3«/30-[u*'(/30)]_1U*(/30). O Odhady parametrů v s-tém kroku jsou získány ze vztahu Iterační proces popsaný v předchozím bodě pokračuje tak dlouho, dokud -ís+l) ~(s) /T - £ Rí 0. Jan KoláCek (PřF MU) M5VM05 Statistické modelováni podzim 2013 43 / 57 Metoda skórování Alternativní procedurou k Newton-Raphsonově metodě je tzv. metoda skórování, kdy se matice druhých parciálních derivací (/S) nahradí její střední hodnotou, tj. maticí — Jn(/S), kde Jn(/S), je informační matice. Druhý iterační krok pak upravíme takto: )8 = )8 + j.. Jan KoláCek (PřF MU) M5VM05 Statistické modeloval podzim 2013 45 / 57 Důsledek Hypotézu H0 : CTß = 0 zamítáme na hladině významnosti a, pokud platí W >*?-«( Ui_a, /TI 1 2' kde opět Fisherovou informační matici Jn(/3) aproximujeme maticí Jn(/3MLE). Jan Koláček (PřF MU) M5VM05 Statistické modelováni mle^ podzim 2013 47 / 57 Ověřování vhodnosti modelu Definice 20 Maximální GLM, který označíme GLMmax, splňuje následující podmínky (1) Maximální model je zobecněný lineární model se stejným typem rozdělení jako zkoumaný GLM model. (2) Maximální model a zkoumaný mají stejnou linkovací funkci. (3) Počet parametrů maximálního modelu je roven počtu vysvětlovaných veličin n, maximálně věrohodný odhad parametru ^max je n-rozměrný vektor fimax. Definice 21 Minimální GLM, který označíme GLMm(n, splňuje následující podmínky (1) Minimální model je zobecněný lineární model se stejným typem rozdělení jako zkoumaný GLM model. (2) Minimální model a zkoumaný mají stejnou linkovací funkci. (3) Počet parametrů minimálního modelu je roven 1, maximálně věrohodný odhad parametru fimin je skalár j}min. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 48 / 57 Submodel Definice 22 Mějme zobecněný lineární model s maticí plánu X„xjt a vektorem neznámých parametrů Submodel, který označíme GLMsub, splňuje následující podmínky (1) Submodel je zobecněný lineární model se stejným typem rozdělení jako zkoumaný GLM model. (2) Submodel a zkoumaný model mají stejnou linkovací funkci. (3) Vektor neznámých parametrů físub E IR^a matice plánu Qnxq, pro kterou platí Qnxq = "X-nxk^kxq- Aby GLMsui, byl submodelem modelu GLM, musí každý sloupec matice Q patřit do obalu sloupců matice X. To bude splněno právě tehdy, bude-li Q typu Qnxq = ^-nxk^kxq- Je třeba si uvědomit, že GLMsuj, je speciálním případem modelu GLM. Platí-li tudíž pro náhodný výběr Y model GLMSU},, platí pro Y také model GLM. Jan Koláček (PřF MU) M5VM05 Statistické modeloval podzim 2013 49 / 57 Deviace Deviace v zobecněných lineárních modelech je obdobou rozptylu u klasických lineárních regresních modelů. Deviace je tedy kritériem vhodnosti zobecněného lineárního modelu. Jak bude patrné z definice, metoda maximální věrohodnosti totiž odpovídá hledání minima deviace modelu. Definice 23 Mějme modely GLM a GLMmax. Nechf náhodný výběr Y se řídí modelem GLMmax- Skálovou deviací modelu GLM (scaled deviance) rozumíme statistiku D = 2 fQ3m„;Y)-fQ3;Y) kde f>max,f> jsou odpovídající maximálně věrohodné odhady. Jan KoláCek (PřF MU) M5VM05 Statistické modelováni Ověřování vhodnosti submodelu Věta 24 Mějme základní model GLM s /}