šansside suhet (OR, odds ratio), Coxi võrdeliste riskide mudeliga aga riskimäärade suhet (HR, hazard ratio). Erinevates analüüsides on leitud, et Coxi mudelitest saadud hinnangud on täpsemad, Coxi regressioon on teatud juhtudel võimsam, Coxi mudeli hindamine võtab rohkem aega (Staley et al., 2017) ning et šansside suhte ja riskimäärade suhte omavaheline erinevus on seda suurem, mida pikem on jälgimisaeg ning mida suurem on haigestumusrisk või haiguse ja riskiteguri vahelise seose tugevus (Symons ja Moore, 2002). 11 Märgime veel, et jälgimiseelsete ja jälgimisaegsete juhtude koos analüüsimiseks on küll välja töötatud spetsiaalseid meetodeid, näiteks eksponentsiaalse kalde meetod (Maziarz et al., 2019), kuid siiani ei ole need skaleeritavad ülegenoomsete analüüside jaoks.

12 2 Analüüsimetoodika 2.1 Elukestusanalüüs

Elukestusanalüüsi (survival analysis) kasutatakse selliste andmete puhul, kus huvi pakub aeg algmomendist mingi kindla sündmuse toimumiseni. Algmomendiks võib olla näiteks inimese sünd või uuringusse kaasamise aeg, huvipakkuvaks sündmuseks näiteks haigestumine, mingi sümptomi esinemine või surm. Ajavahemikku algmomendist lõppsündmuseni nimetatakse elu kestuseks (survival time). Tsenseerimine Sageli ei ole kõikide subjektide elukestus teada, sellisel juhul nimetatakse nende elukestust tsenseerituks. Tavaliselt on tsenseerimine tingitud sellest, et uuritaval ei ole sündmus enne katse lõppu või analüüsi tegemist toimunud. Inimese elukestus tsenseeritakse ka siis, kui inimene sureb (kui surm ei ole huvipakkuvaks sündmuseks) või lahkub uuringust enne selle lõppemist. Tsenseerimise korral on peamiselt tegemist paremalt poolt tsenseerimisega, mis tähendab seda, et uuritav sündmus ei toimunud enne uuringu lõppu. Tähtsaks eelduseks on see, et tsenseerimine oleks mitteinformatiivne ehk toimuks huvi pakkuvast sündmusest sõltumatult. Paremalt tsenseeritud andmete korral tähendab see, et sub jekti tsenseerimine mingil ajahetkel ei anna sündmuse toimumise aja kohta muud infot peale selle, et tsenseerimisaeg on väiksem sündmuse toimumise ajast. Haigestumisandmete uurimisel ei tohi seega inimese tsenseerimine oleneda tema tervislikust seisundist – kui näiteks kontroll gruppi kuuluvad inimesed lahkuvad uuringust seetõttu, et nad jäid väga haigeks, ja ravigruppi kuuluvad inimesed seetõttu, et nad tunnevad end juba tervelt, siis on tegemist informatiivse tsenseerimisega – ning kõiki klassikalisi tsenseeritud andmete uurimiseks mõeldud analüüsi meetodeid rakendada ei saa (Kalbfleisch ja Prentice, 2002: 13). Selles töös kasutatavate geenivaramu andmete puhul on tegemist mitteinformatiivse paremalt tsenseerimisega – teada on viimane ajahetk, mil subjektil ei olnud

huvipakkuvat sündmust toi 13 munud – ning tsenseerimine on tingitud sellest, et geenidoonor ei ole analüüsi tegemise ajaks haigust saanud. Sel juhul võetakse tsenseerimisajaks kuupäev, mil toimus viimane ühendamine (linking) diagnoosiandmetega. Väike osa tsenseerimistest toimub aga ka seetõttu, et geeni doonor sureb enne huvipakkuva haiguse saamist. Kuigi me sel juhul teame, et inimene ei saa kunagi huvipakkuvat haigust, ja seega ei ole tegelikult tegemist mitteinformatiivse tsenseerimi sega, siis eeldusel, et suremus ei mõjuta seost riskiteguri ja haiguse vahel, on Coxi võrdeliste riskide mudeli parameetreid ka sel juhul võimalik õigesti hinnata (Therneau et al., 2021). Ajaskaala valik Kestusandmete analüüsimisel on oluline valida sobiv algmoment ja ajaskaala. Sageli on algmo mendi valikuks kas inimese sünniaeg või uuringuga liitumise aeg. Kliiniliste uuringute puhul hakatakse inimesi tihti jälgima alates näiteks haigestumisest või mingi ravi alustamisest ja huvi pakub näiteks inimese tervenemine, surm või sümptomite taaste ke. Sel juhul on loomulikuks algmomendiks uuringusse kaasamise aeg ja ajaskaalaks uuringus oldud aeg (time-on-study). Siis võrreldakse uuritavat teiste uuritavatega, kes on uuringus olnud sama kaua. Inimese vanuse arvesse võtmiseks lisatakse sel juhul kovariaadina mudelisse ka uuringuga liitumise vanus. (Thiébaut ja Bénichou, 2004) Epidemioloogiliste uuringute puhul, kus huvi pakub aeg mingi haiguse tekkeni, kasutatakse ajaskaalana enamasti inimese vanust ja algmomendiks on seega sünniaeg (Thiébaut ja Bénichou, 2004). Sellisel juhul ei ole uuringuga liitumise aeg seotud huvipakkuva haigusega ja seega ei ole aeg uuringuga liitumisest haigestumiseni väga oluline. Arvestades, et enamiku haiguste puhul on vanus olulisim riskitegur, siis on tähtis just inimese vanuse arvesse võtmine, ja kui ajaskaala on

inimese vanus, siis võrreldaksegi uuritavaid teiste samaealistega. Arvestada tuleb aga vasakult tõkestatusega (left truncation), see tähendab, et uuringusse saavad jõuda vaid need inimesed, kes on elus uuringu alguses. Vasakult tõkestatuse mittearvestamine võib viia nihkega hinnanguteni (Klein ja Moeschberger, 2003). Samal ajal on näidatud, et kui epidemioloogiliste uuringute puhul on huvipakkuv argument 14 tunnus ajas muutumatu ega ole seotud liitumisaegse vanusega, siis annab nende erinevate aja skaalade kasutamine sarnaseid tulemusi (Korn et al., 1997; Ingram et al., 1997; Canchola et al., 2003; Thiébaut ja Bénichou, 2004). Selles töös kasutame elukestusandmete analüüsimisel aja skaalana aega alates geenivaramuga liitumisest ja kovariaadina lisame vanuse liitumise ajal. Eespool tutvustatud mõisteid on illustreeritud joonisel 1. Sellel on numbritega 1-6 tähista tud erinevad vaatlusalused subjektid, kes on kas paremalt tsenseeritud, vasakult tõkestatud või jälgimiseelsed või jälgimisaegsed juhud. Aeg 1 2 3 4 5 6 Geenivaramuga liitumine Analüüsi tegemine

Joonis 1. Subjektid 1-4 on geenivaramuga liitunud. Subjektid 1 ja 4 on paremalt tsenseeritud, subjekt 2 on jälgimisaegne juht ja subjekt 3 on jälgimiseelne juht. Subjektid 5 ja 6 on vasakult tõkestatud. Sellel joonisel on beeži joonega tähistatud aeg, mil indiviid ei olnud geenidoonor, ja sinise joonega on märgitud geenidoonoriks olemise aeg. Pidevjoonega on tähistatud see aeg, mil toi munud sündmused on analüüsi läbiviimise ajaks teada, kriipsjoonega see aeg, mille kohta info analüüsi tegemise ajal puudub. Risti ja ringjoonega on tähistatud vastavalt surm ja haigestu mine. Subjektid 5 ja 6 kujutavad vasakult tõkestatust, nemad ei ole geenivaramuga liitunud 15 ja nende haigestumise kohta meil info puudub. Subjekt 3 on saanud huvipakkuva haiguse en ne geenivaramuga liitumist ja on seega uuritava haiguse suhtes jälgimiseelne juht. Subjekt 2 on haigestunud pärast liitumist ja on seega uuritava haiguse suhtes jälgimisaegne juht. Subjekt 4 on paremalt tsenseeritud surma tõttu ja subjekt 1 on paremalt tsenseeritud seetõttu, et huvipakkuv sündmus ei olnud toimunud enne analüüsi tegemist.

2.2 Võrdeliste riskide mudel See peatükk koos alapeatükkidega põhineb raamatul (Collett, 2015), kui ei ole viidatud teisiti. Anname enne võrdeliste riskide mudeli tutvustamist ülevaate elukestusanalüüsi tähtsamatest mõistetest ja definitsioonidest. Olgu T mittenegatiivne elukestust tähistav juhuslik suurus tihedusfunktsiooniga f(t) ja jaotus funktsiooniga F(t). Üleelamisfunktsioon S(t) on tõenäosus, et huvipakkuv sündmus toimub pärast ajamomenti t > 0: S(t) = P(T > t) = 1 − F(t). Riskifunktsioon h(t) iseloomustab hetkelist sündmuse toimumise riski ajamomendil t > 0, kui on teada, et see ei toimunud enne ajamomenti t. Pideva aja korral on riskifunktsioon kujul P(t 6 T < t + ∆t| T > t) h(t) = lim ∆t→0 ja selle väärtused on vahemikus [0,∞). ∆t=f(t) S(t) Eelmise võrduse põhjal saame, et t > 0 korral kehtivad järgmised seosed: h(t) = − d dtlog S(t), 16 millest

S(t) = exp(−H(t)) (1) ja H(t) = −log S(t), (2) kus H(t) = R t 0h(u)du on kumulatiivne riskifunktsioon, mis kirjeldab kumulatiivset sündmuse toimumise riski ajaks t. Selleks, et leida elukestust prognoosivaid tunnuseid ja nende mõju riskifunktsioonile, kasuta takse võrdeliste riskide mudeleid (proportional hazards models). Olgu meil n vaatlust, p kirjeldavat tunnust ja olgu i-nda (i = 1, . . . , n) vaatluse riskifunktsioon hi(t) ning xi = (x1i, . . . , xpi)Ttema argumenttunnuste väärtuste vektor. Olgu h0(t) baasriski funktsioon ehk riskifunktsioon juhul, kui kõigi seletavate tunnuste väärtusteks on 0. Siis võrde liste riskide mudeli korral saab i-ndale vaatlusele vastava riskifunktsiooni kirjutada kujul hi(t) = h0(t)ψ(xi), kus ψ(xi) on ajast t mittesõltuv xi funktsioon, mis kirjeldab riskifunktsioonide suhet ehk riski määrade suhet (HR, hazard ratio). Et riskimäärade suhe ei saa olla negatiivne, kasutatakse tihti ψ(xi) = exp(βTxi), kus β = (β1, β2, . . . , βp)T on argumenttunnustele vastavate parameetrite vektor. Eelnevat kokku võttes on võrdeliste riskide mudeli kujuks hi(t) = h0(t) exp(βTxi). (3)

Seega kahe indiviidi ajabkorral, kelle argumenttunnuste vektorid ning riskifunktsioonid on 17 vastavaltxa jaxb ningha (t) jahb (t), kehtib võrdelisteriskide eeldus:hb (t)= h0 (t) exp(β T xa ) ha (t)h0 (t) exp(β T xb )

exp{β T (xa − xb )} ehk riskimäärade suhe ei sõltu ajastt. Viies mudeli (3) kujule log hi (t)h0 (t) =β T xi , näeme, et tegemist on lineaarse mudeliga riskimäärade suhte logaritmile. Parameetriliste võrdeliste riskide mudelite korral on baasriskifunktsioon h0 (t)määratud mingi kindla jaotuse parameetritega. Poolparameetrilisel juhul hinnatakseainult parameetridβja baasriskifunktsiooni ei määrata. 2.2.1 Coxi võrdeliste riskide mudel Coxi võrdeliste riskide mudel, mis on üks populaarsemaid meetodeidelukestusanalüüsis, on poolparameetriline mudel, sest baasriskifunktsiooni h0 (t) mudelis (3) ei määrata, küll aga lei takse hinnangud regressioonikordajatele βj (j = 1, . . . , p) ja huvipakkuvatele riskimäärade suhetele. Parameetrite hindamine Coxi mudeli parameetrid saab hinnata suurima tõepära meetodil. Kasutusel on niinimetatud osaline tõepärafunktsioon, mis ei sõltu baasriskist h0(t) ega kasuta täpseid vaadeldud sündmuse toimumise aegu, vaid ainult nende aegade omavahelist järjestust. Olgu meil vaatluse all n subjekti, kellest r-il toimub sündmus, ja ülejäänud n −r subjekti on paremalt tsenseeritud. Järjestame sündmuse toimumiste ajad t(1) < . . . < t(r). Tähistagu iga j = 1, . . . , r korral R(j) nende inimeste hulka, kellel ei ole toimunud sündmust ja kes 18 ei ole tsenseeritud enne hetke t(j), ning olgu x(j) argumenttunnuste vektor indiviidil, kellel esineb sündmus hetkel t(j). Selliste tähistuste korral on Coxi (1972; 1975) tuletatud mudelile (3) vastava osalise tõepärafunktsiooni kujuks L(β) = Yr j=1 P exp(βTx(j) ) l∈R(j)exp(β Txl). Olgu U(β) logaritmitud osalise tõepärafunktsiooni esimene tuletis β järgi: U(β) = ∂log L(β) ∂β. Suurima osalise tõepära hinnangud leitakse siis võrrandist U(β) = 0, mille lahendamiseks kasutatakse tavaliselt iteratiivseid meetodeid, näiteks Newton-Raphsoni meetodit (Therneau ja Grambsch, 2000: 41). Coxi võrdeliste riskide mudeli eelduseks on see, et riskimäärade suhe oleks ajas

muutumatu, ehk et riskimäärad oleksid võrdelised. Selle eelduse täidetust saab kontrollida näiteks Schoenfeldi jääkide abil (vt lisa 1). Võrdsed sündmuste toimumiste ajad Siiani oleme teinud eelduse, et kõik sündmuste toimumiste ajad on erinevad. Praktikas võib aga juhtuda, et näiteks mõõtmiste ebatäpsuse tõttu esineb andmetes ka võrdseid sündmuse toimu miste aegu. Väga tihti on sündmuse toimumise aeg teada päeva täpsusega, kuid samal päeval saab sündmus toimuda mitmel indiviidil; seda näeme eelkõige suurtes valimites, nagu ka geeni varamu andmebaasis. Kui andmetes esineb võrdseid sündmuse toimumiste aegu, siis muutub vastav osaline tõepärafunktsioon arvutuslikult keeruliseks ja kasutatakse erinevaid lähendeid. Enamasti kasutatakse lihtsat Breslow’ lähendit (Breslow, 1972; Peto, 1972), kus samadel aega del toimuvate sündmuste puhul eeldatakse, et need on erinevad ja toimuvad üksteise järel: 19 L(β) ≈Yr j=1 nP exp(βTsj ) l∈R(j)exp( βTxl) odj, kus iga j = 1, . . . , r korral on D(j) nende indiviidide hulk, kellel toimus sündmus hetkel t(j), sj =Pi∈D(j)xi on nende indiviidide argumenttunnuste vektorite summa, kellel toimub sünd mus hetkel t(j), ja dj on ajahetkel t(j)toimunud sündmuste arv. Veidi keerulisema ja täpsema lähendi pakkus Efron (1977), mis on vaikimisi kasutusel ka R-i paketi survival (Therneau, 2020) funktsioonides: L(β) ≈Yr j=1 Qdj k=1 nP exp(βTsj ) l∈R(j)exp(βT xl) −(k − 1)d−1 j P l∈D(j)exp(βT xl) o.

Paneme tähele, et juhul kui kõik sündmuste toimumiste ajad on erinevad ehk dj = 1 iga j = 1, . . . , r korral, siis on mõlemad lähendid täpselt samad ja võrdsed esialgse Coxi osalise tõepärafunktsiooniga. Kui võrdseid sündmuste toimumiste aegu on vähe, siis annavad mõlemad lähendid sarnaseid tulemusi. Selles töös kasutame võrdsete haigestumisaegade puhul Efroni tõepärafunktsiooni lähendit. Martingaalijäägid Coxi mudeli jaoks on defineeritud erinevaid jääke: näiteks Cox-Snelli jäägid mudeli üleüldise sobivuse hindamiseks, martingaalijäägid kovariaatide sobiva funktsionaalse kuju määramiseks, hälbimuse jäägid erindite leidmiseks, Schoenfeldi jäägid ning skoorijäägid võrdeliste riskide eelduse kontrollimiseks ja erindlike vaatluste leidmiseks. Siin peatükis kirjeldame lähemalt martingaalijääke. Martingaalijääkide nimetus tuleb sellest, et need on võimalik tuletada, kasutades juhuslike prot sesside, täpsemalt martingaalide teooriat. Neid jääke kasutatakse näiteks selgitavate tunnuste funktsionaalse kuju määramiseks. Martingaalijääkide arvutamiseks on vaja hinnangut kumula 20 tiivsele riskifunktsioonile. Riski- ja üleelamisfunktsiooni hindamisest on kirjutatud lisas 2. Olgu meil nagu ennegi vaatluse all n indiviidi, kellest r-il toimub sündmus ja ülejäänud n−r on paremalt tsenseeritud, ning olgu meil hinnatud p argumenttunnusega Coxi mudel. Olgu δi = 0, kui i-s (i = 1, . . . , n) vaatlus on tsenseeritud, ja δi = 1 muidu. Olgu i-nda indiviidi vaatlusajaks tija olgu talle hinnatud Nelson-Aaleni kumulatiivne riskifunktsioon Hˆi(t) (vt lisa 2 valem (8)). Martingaalijääk i-nda indiviidi jaoks on defineeritud järgmiselt:

mˆ i = δi − Hˆi(ti) = δi − Hˆ0(ti) exp(βˆTxi). Martingaalijääkide väärtused on −∞ja 1 vahel, kusjuures kõik tsenseeritud vaatluste martingaalijäägid on negatiivsed. Paneme tähele, et suurus δi on i-nda indiviidi vaadeldud sündmuste arv ajavahemikus alg momendist hetkeni ti (kas 0 või 1 vastavalt sellele, kas indiviid on tsenseeritud või mitte). Hinnatud kumulatiivset riskifunktsiooni Hˆi(ti) aga saab interpreteerida kui oodatavat sünd muste arvu indiviidil i ajavahemikus algmomendist hetkeni ti. See tähendab, et i-nda indiviidi martingaalijääki saab interpreteerida kui tema ajavahemikus algmomendist hetkeni ti vaadel dud ja oodatud sündmuste arvu vahet. Seega saab martingaalijääkide abil leida need subjektid, kellel hinnatud mudeli põhjal toimus sündmus liiga vara või liiga hilja – absoluutväärtuselt suu red negatiivsed jäägid vastavad nendele subjektidele, kelle elukestus on pikk, kuid kellel mudeli põhjal pidanuks sündmus toimuma varem, ühelähedased jäägid vastavad neile, kelle elukestus oli lühem, kui võiks arvata mudeli põhjal. Martingaalijääke saab kasutada Coxi mudeli lähendamiseks lineaarse regressiooniga. Selles töös kasutatakse martingaalijääke, hindamaks SNP-de mõju haigestumisele. Koos negu i-nda indiviidi mudelis kasutatavate argumenttunnuste vektor kahte tüüpi tunnustest xi = (gi, z1i, . . . , zpi)T, kus g tähistab SNP-d ja zj (j = 1, . . . , p) on ülejäänud kovariaadid. Olgu nende argumenttunnuste parameetriteks vastavalt βg ja β1, . . . , βp. SNP mõju hindamiseks sobitatakse esmalt Coxi mudel nii, et kaasatakse ainult mittegeneetilised argumenttunnused ehk 21 hinnatakse nii-öelda nullmudel SNP suhtes: hi(t) = h0(t) exp(β1z1i + . . . + βpzpi).