Ústav matematiky a statistiky Přírodovědecká fakulta Masarykova univerzita
Analýza prežívania
Vybrané kapitoly z neparametrických metód pre cenzurované dáta, zadania príkladov a domácich úloh a
niektoré riešenia
Stanislav Katina katina@math.muni.cz
10. decembra 2014
Obsah
1 Udalosti a cenzurovanie 1
1.1 Cenzurovanie I. typu.................................... 2
1.2 Cenzurovanie II. typu.................................... 2
1.3 Progresívne cenzurovanie.................................. 3
1.4 Ľubovoľné a náhodné cenzurovanie............................ 3
1.5 Intervalové cenzurovanie I. typu.............................. 3
1.6 Intervalové cenzurovanie II typu.............................. 4
2 Základné charakteristiky prežívania 5
3 Odhady základných charakteristík prežívania 5
3.1 Empirická funkcia prežívania................................ 6
3.2 Odhady funkcie prežívania................................. 7
3.3 Odhady rizika a kumulatívneho rizika........................... 11
3.4 Odhady rozptylu kumulatívneho rizika a funkcie prežívania............... 12
3.5 Empirické intervaly a pásy spoľahlivosti.......................... 14
3.5.1 Intervaly spoľahlivosti ............................... 14
3.5.2 Intervaly spoľahlivosti pre kumulatívne riziko .................. 15
3.5.3 Intervaly spoľahlivosti pre funkciu prežívania................... 15
3.6 Odhady strednej hodnoty a mediánu prežívania..................... 16
3.7 Odhad strednej hodnoty a mediánu zostatkového života a ich rozptyl ......... 18
4 Testy na porovnanie kriviek prežívania 21
5 Charakteristiky definované sčítacím procesom 26
6 Softvérová implementácia charakteristík prežívania 28
7 Štruktúra funkcií na cvičenia a do DÚ 30
Katina, S., 2014: Analýza prežívania
1
1   Udalosti a cenzurovanie
Udalosťou nazývame ukončenie pozorovania z dôvodu zlyhania alebo smrti pacienta - do konca sledovaného obdobia.
Príklady udalostí:
• overall survival - smrť z akéhokoľvek dôvodu,
• progression-free survival - prvé znaky progresie choroby alebo smrť,
• disease-free survival - prvé znovuobjavenie sa choroby alebo smrť,
• event-free survival - prvé znovuobjavenie sa choroby, objavenie sa inej špecifikovanej choroby alebo smrť,
• disease-specific survival (cause-specific survival) - smrť ako dôsledok špecifikovanej choroby,
• relapse-free survival (recurrence-free survival) - prvé znaky recidívy (opakovania sa) chodoby,
• time-to-progression - prvé znaky progresie choroby.
V mnohých prípadoch v danom pokuse musíme pozorovania ukončiť z dôvodu iného ako je zlyhanie alebo smrť pacienta. Teda do konca sledovaného obdobia dôjde k úmrtiu len niektorých subjektov (pacientov), zatiaľ čo u ostatných k úmrtiu do konca sledovaného obdobia buď nedôjde alebo sa tieto subjekty z pozorovania stratia. Cenzurované pozorovanie obsahuje v sebe len čiastočnú informáciu o čase úmrtia.
Predpokladajme, že n subjektov s rovnakými znakmi začneme pozorovať v čase t = 0. Pre ich časy do zlyhania Tí, T2,..., Tn platí, že sú nezávislé, rovnako rozdelené náhodné premenné, kde náhodná veličina Ti (i = 1,2,... ,n) má distribučnú funkciu F (t), o ktorej obyčajne predpokladáme, že je absolútne spojitá. Teoreticky by sme mali pozorovať subjekty tak dlho, pokiaľ nezískame úplný náhodný výber časov do zlyhania pre všetkých n pozorovaní. Ide o veľmi dlhý časový interval, a preto nemôžeme túto štúdiu alebo experiment z rôznych dôvodov urobiť až do konca. Musíme preto vhodným spôsobom rozhodnúť o ukončení pozorovaní, k čomu nás možu donútiť napr. technické alebo ekonomické dôvody, a teda získame neúplné cenzurované dáta. V medicíne sa môžu vyskytnúť z nasledujúcich príčin:
• ukončenie štúdie (termination of the study): pacient prežije časový interval experimentu,
• konkurenčné riziko (competing risk): pacient zomrie z iného dôvodu, ako v dôsledku sledovanej choroby,
• prerušenie/vy sadenie liečby (drop-out): pacient preruší liečbu a odíde z kliniky predčasne, napr. z dôvodu zlých vedľajších účinkov liečby, pacient sa sám rozhodne nepokračovať v liečbe,
• strata z ďalšieho sledovania (loss to follow-up): pacient sa rozhodne presťahovať a nemáme o ňom už žiadne informácie.
Za predpokladu, že experiment prebehne bez cenzurovania, i-ty jedinec z populácie o rozsahu n zomrie za čas Tt od začiatku sledovania. Čas, ktorý mienime sledovať i-teho jedinca označíme C,t
(a > 0).
Podľa spôsobu ukončenia experimentu môžeme rozdeliť cenzurovanie do nasledovných typov:
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
2
• cenzurovanie I. typu,
• cenzurovanie II. typu,
• progresívne (zrýchlené) cenzurovanie,
• ľubovoľné a náhodné cenzurovanie,
• intervalové cenzurovanie I. typu,
• intervalové cenzurovanie II. typu.
1.1    Cenzurovanie I. typu
Predpokladáme, že všetkých n jedincov vstupuje do experimentu súčasne. Ide o cenzurovanie časom, pretože si zvolíme pevné číslo tc, ktoré nazveme fixovaný cenzurujúci čas, teda jediná príčina cenzurovania je plánované ukončenie experimentu. Takto získame informáciu o prvých d pozorovaniach j^1) < jH2) < ... < T^d\ kde T^ < tc < T^d+1^ a o ostatných n — d nepozorovaných zlyhaniach vieme iba to, že nastali až po cenzurujúcom čase tc. Počet skutočne pozorovaných zlyhaní d je náhodná veličina, ktorá môže nadobúdať hodnoty 0,1,... ,n.
Namiesto pozorovania Ti pozorujeme Xi,X2,... ,Xn, kde
Tí, Ti < tc, pre necenzurované X,t
Xi = min (Tí,tCJ - ,
1 tc,li > tc, pre cenzurované Xi
Pre distribučnú funkciu X platí Pr(T > tc) > 0 v bode x = tc, pretože Pr(X = tc) = Pr(T > t c) = 1 — F (t c) > 0. To znamená, že X je zmiešaná náhodná premenná so spojitým a diskrétnym komponentom. Kumulatívna distribučná funkcia M (X) náhodnej premennej X bude teda spojitá do času tc (v čase < tc) a v čase tc skočí do jednotky (diskrétny komponent). Majme binárnu premennú 5, ktorá indikuje, či bol čas zlyhania pozorovaný alebo cenzurovaný,
^      í 1, Ti < tc, pre necenzurované Xi \ 0,Tj > tc, pre cenzurované X í
Skutočným pozorovaniam potom zodpovedá náhodný vektor (Xi, ô i), i = 1,2,..., n.
Pozn.: Treba si uvedomiť, že ak (ô = 0 a T < tc) implikuje, že čas zlyhania bol presne T = tc, čo nastane s nulovou pravdepodobnosťou, ak je T spojitá premenná. Pre diskrétne premenné môžeme definovať tc rovné poslednému dostupnému času zlyhania, ktorý môže byť pozorovaný. Potom Pľ(ô = 0 A T < tc) Ý 0.
1.2   Cenzurovanie II. typu
Aj pri tomto type cenzurovania vstupuje do experimentu súčasne všetkých n jedincov. Zvolíme pevné číslo d, ktoré nazveme fixovaný počet zlyhaní. Ukončenie potom nastáva po vopred zvolenom počte d zlyhaní, tj. ide o cenzurovanie zlyhaním, kde d = [np] + l,p E (0,1) ,T^ <       <....< T^d\
Odtiaľ vidíme, že d < n a všetkých zostávajúcich n — d pozorovaní má čas zlyhania väčší ako T^d\ Dostávame tak celistvú informáciu o prvých d pozorovaniach. Náhodnou veličinou je čas trvania experimentu T^d\
Namiesto pozorovania Tj pozorujeme X1}X2,... ,Xn, kde
X = miníT- T(d)) = { Ti,Ti ~ ^'' pľe necenzúrované Xí 1 \      , Tí >      , pre cenzurované X,-t
Skutočnému pozorovaniu potom zodpovedá náhodný vektor (Xí,ô,j), kde
) l,Tj < T(d\ pre necenzúrované X,-t 0,Tí > T^d\ pre cenzurované X,-t
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
3
1.3 Progresívne cenzurovanie
Cieľom tohto typu cenzurovania je čo najväčšie skrátenie času experimentu. Znovu vstupuje do experimentu súčasne všetkých n jedincov. Teraz máme dve možnosti postupu.
Zrýchlené cenzurovanie I. typu. Ide o cenzurovanie zlyhaním, kedy zvolíme čísla tci,i = 1,2,..., k, ktoré nazveme fixované cenzurujúce časy. Vieme, že tc\ < tC2 < ... < tck. V čase tCi vyradíme rrii subjektov, t.j. v čase tc\ vyradíme ni\ subjektov, v čase tC2 vyradíme ni2 subjektov, ..., v čase tck vyradíme m k subjektov. Po fc-tom kroku máme vyradených ni\ + rri2 + ... + rrik subjektov. Náhodnou veličinou je počet skutočne pozorovaných zlyhaní d E {0,1,..., n}.
Zrýchlené cenzurovanie II. typu. Ide o cenzurovanie časom, kedy zvolíme čísla dÍ7 ktoré nazveme fixované počty zlyhaní. Vyradenie teda nastáva po vopred zvolenom počte d,i zlyhaní, kde di = [npi] + l,pi G (0,1), t.j. po di zlyhaniach vyradíme ni\ subjektov, po d<i zlyhaniach vyradíme ni2 subjektov, po dk zlyhaniach vyradíme m k subjektov. Po fc-tom kroku máme vyradených ni\ + m2 + ... + nik subjektov. Náhodnou veličinou je čas trvania experimentu.
1.4 Ľubovolné a náhodné cenzurovanie
V predchádzajúcich typoch cenzurovania vstupuje do experimentu súčasne všetkých n jedincov a časy zlyhania sú v takom poradí, v akom nastali. Ale takýto ideálny prípad nemusí nastať.
Čas do cenzurovania i-teho jedinca označme d, potom C\, C2, ■ ■ ■, Cn sú nezávislé, rovnako rozdelené náhodné premenné, kde náhodná veličina d (i = 1,..., n) má distribučnú funkciu G (t).
Pozorujeme X\, X2,..., Xn, kde
^       .  /rn  „ .     í Ti, T,j_ < C i, pre necenzurované X,t
Xi = mm (Ti, d) = S n rp ^ n ,        , v
Ci,li > Ci, pre cenzurované Aj
Skutočnému pozorovaniu potom zodpovedá náhodný vektor (Xi, ô i), i = 1,2,... ,n, kde Xi — min (T, C i) ,
^ _ í ^,Ti < Ci, pre necenzurované Xi \0,T > Ci, pre cenzurované X i
Ak Ci považujeme za ľubovoľné konštanty a nie za náhodné veličiny, kde napr. Ci = c,i = 1,2,... ,n, potom hovoríme o ľubovolnom cenzurováni. Cenzurovanie I. typu je jeho špeciálnym prípadom.
Ak považujeme d za náhodné veličiny, potom hovoríme o náhodnom cenzurovaní a predpokladáme nezávislosť náhodných veličín d a, T. Tento typ cenzurovania vzniká často v medicínskych štúdiách.
Pozn.: Je to však dosť silný predpoklad. Treba si uvedomiť súvislosť s typom cenzúry drop-out (z dôvodu vedľajších účinkov liečby). Ale práve za tohoto predpokladu je funkcia vierohodnosti ľahko spočítateľná.
1.5 Intervalové cenzurovanie I. typu
Majme n subjektov. Označme Ti,i = 1,2,... ,n, nepozorovatelné časy zlyhania. Skutočnému pozorovaniu potom zodpovedá náhodný vektor (Xi} ó~i), kde Xi = Ci sú časy cenzúr a ôi = I(T < Ci), t.j.
£ _ í 1, 7i < Cj, pre necenzurované Xi \0,T > Ci, pre cenzurované X i
Príklad 1 (nádor pľúc, animálny model) Laboratórne myši sú injektované látkou, ktorá spôsobuje nádor. Kedže tento druh nádoru nie je smrteľný, je potrebné myš najprv zabiť, aby sme zistili, či bol nádor indukovaný, t.j. po časovom úseku náhodnej dĺžky C je myš zabitá, aby sme zistili, či sa nádor vyvinul alebo nie. Endpoint záujmu je čas T do objavenia sa nádoru. pred
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
4
1.6   Intervalové cenzurovanie II typu
Majme opäť n subjektov. Označme T,i = 1,2,... ,n, nepozorovatelné časy zlyhania. Vieme len, že T nastalo buď vnútri nejakého náhodného časového intervalu, pred jeho ľavou hranicou alebo po jeho pravej hranici. Označme Cu a C u časy dvoch vyšetrení a indikačné funkcie definujeme nasledovne ôu = I(Ti < Cu), ô2i = I(Cu < T < C2í) a ô3i = I(Ti > C2í), t.j.
^       í ^,Tí < Cu, pre necenzurované X;t \ 0, Tj > Cu, pre cenzurované X;t '
^       íl, Cu < T < C2i, pre necenzurované X,t \0,T > C2i, pre cenzurované X,t
a nakoniec 53i = 0.
Potom máme nasledovné tri možnosti:
1. udalosť mohla nastať niekedy pred prvým vyšetrením Cu, kde 5U = 1 a 82i = ó^í = 0,
2. udalosť mohla nastať niekedy medzi prvým a druhým vyšetrením, t.j. v intervale (Cu,C2í), kde ôu = 0, 52i = 1 a ó3i = 0,
3. udalosť sa do druhého vyšetrenia nevyskytla, t.j. mohla nastať niekedy po C2Í (ale nevieme kedy), kde ôu = 0, 82% = 0 a 53i = 0.
Nech Xu = Cu a X2Í = Cii- Skutočnému pozorovaniu potom zopovedá náhodný vektor
{Xii, X21, Ôu, Ô2í).
Všimnime si, že ôsi nie je potrebné použiť, pretože nemáme ďalšie vyšetrenie po Cii- Keby sme mali Csí alebo aj ďalšie (po ňom nasledujúce) vyšetrenia, hovorili by sme zovšeobecnenom intervalovom cenzurovaní.
Príklad 2 (nádor pľúc, pacienti) Pacienti navštevovali kliniku opakovane každých A až 6 mesiacov, kde pozorovania sú buď intervaly (Cu,C2í) ak sa retrakcia prsníka vyskytla medzi poslednými dvoma návštevami alebo {C2%, 00), ak sa do C21 retrakcia nevyskytla.
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
5
2   Základné charakteristiky prežívania
Príklad 3 (k-ty moment času prežívania) Nech nezáporná náhodná veličina T je charakterizovaná funkciou prežívania S (T). Nech je k-ty moment, E[Tk], konečný, E[Tk] < oo,k G N.
a) Ukážte, že platí E [T] = ^íeNo Pr(T > t) = ^íeNo 1 — F (T) = ^2te^0S(T). Použite pri tom definíciu strednej hodnoty E [T] = ^2te^0tPr(t) a pomocné tvrdenie 1 + 1 + ... + 1 = J^=1 1 =
t-krát
b) Ukážte, že platí E[T] = jQ  S (t) dt. Použite pri tom definíciu strednej hodnoty E [T] = fQ tf(t)dt,
aplikujte vlastnosti súm z (a) ako aj J0°° s (í) dt = J^°(Jq ldx)s (í) dt. Výpočet Vám uľahčí metóda per-partes.
c) Pomocou metódy per-partes ukážte, že E[Tk] = k J0°° tk~1S(t)dt. DÚ
Príklad 4 (stredná hodnota zostatkového života) Pomocou metódy per partes ukážte, že platí posledná rovnosť v nasledovnom vzťahu DÚ
f°° (u — t) f(u)du     f°° S(u)du
">•*>=Ev -ť|T ž *i =' S(t)   = Lwr-
Príklad 5 (funkcia prežívania)  Ukážte, že platí posledná rovnosť v nasledovnom vzťahu DU
/■°° ^ m     íí\,m\     , í ŕ du
o (t) = /    xj(x)dx = exp I  /  \{x)dx   = exp (—A(t)) =-—- exp
J t \Jo / mrl{ť)       \   J q mrl{u)
Príklad 6 (riziko) Ukážte, že platí posledná rovnosť v nasledovnom vzťahu DÚ
\{t) = -?L\nS(t) = ^j\ = (^-mrl{t)
dt      K '     S (t)     \dt     w     J mrl{t)
Príklad 7 (hustota) Ukážte, že platí posledná rovnosť v nasledovnom vzťahu DÚ
,/ n       9     .     w N„, N     f d       .      \   mrl(0)        (    ŕ du f (t) = ~^S(t) = X(t)S(t) = ( -mrlit) + 1 ) /    ^'  exp '
dt  w      w  w     \dt     w      J (mrl(t)Y       V   Jo rnrl(u) Táto rovnosť vyplýva zo vzťahu f (t) = \(t)S(t) a z výsledku predchádzajúcich dvoch príkladov.
3   Odhady základných charakteristík prežívania
Príklad 8 (Akútna myelogénna leukémia, AML) AML pacienti boli po absolvovaní chemoterapie a zmiernení príznakov náhodne rozdelení do dvoch skupín. Prvá skupina (skupina A) dostala udržujúcu chemoterapiu a druhá (kontrolná; skupina B) nie (pozri tabuľku). Cieľom bolo zistiť, či udržujúca chemoterapia predlžuje čas remisie (opätovného zhoršenia stavu). cvič.
Tabuľka 1: Dáta AML
skupina	čas po kompletný relaps (v týždňoch)	n	udalostí	cenzúr
A	9,13,13+, 18, 23, 28+, 31, 34,45+, 48,161+	11	7	4
B	5,5, 8, 8,12,16+, 23, 27, 30, 33,43,45	12	11	1
Tri náhľady na problém analýzy AML dáta:
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
6
problém 1: po odstránení cenzurovaných pozorovaní:
problém 2: po ošetrení cenzurovaných pozorovaní, ktoré zoberieme do úvahy akoby boli udalosťami (zlyhaniami) a
problém 3: berúc do úvahy cenzurované pozorovania.
Tabulka 2: Základné charakteristiky AML podľa skupín (v týždňoch)
	problém 1	problém 2	problém 3
skupina	A B	A B	A B
x	25.1 21.7	38.5 21.3	52.6 22.7
x	23.0 23.0	28.0 19.5	31.0 23.0
Príklad 9 (časy do zlyhania alebo cenzúry, AML) Nakreslíte časy do zlyhania a časy do cenzúry pre AML dáta ako horizontálne úsečky ukončené bodom typu • (pre zlyhania) a bodom typu + (pre cenzúry). Na y-ovej osi označte pacientov skratkami Al-All a B1-B12 a na x-ovej osi zobrazte časy po desiatich týždňoch od 0 po 160.
A11 A10
B12 B11 B10
I-1-1-1-1-1-1-1-1-1-1-1-1-1-1-1-1
0     10    20    30    40    50    60    70    80    90   100   110   120   130   140   150 160
cas(v týždňoch)
Obr. 1: Časy do zlyhania (ozn. •) a časy do cenzúry (ozn. +) pre AML dáta
3.1   Empirická funkcia prežívania
Majme diskrétny prípad. Predpokladajme, že časy do zlyhania sú už zoradené, t.j. označenie t^ < t(2) < • • • < — ni přeznačíme na t\ < t2 < ... < tj. Potom empirickú funkciu prežívania
vypočítame nasledovne
Sn(t)
^pozorovaní > t _ #{U > t} _ Yn=i HU > t)
n
n
n
kde n je rozsah súboru. Ide teda o podiel tých pacientov, ktorí sú stále v riziku. Sn(t) je sprava spojitá schodovitá funkcia so schodmi smerom nadol v každom ti. Ide o konzistentný odhad skutočnej funkcie prežívania S (t). Exaktné rozdelenie nSn(t) pre každé fixované t je Bi(n,p), kde n je počet pozorovaní
a p = Pr(T > t). Na základe centrálnej limitnej vety pre každé fixované t platí asymptoticky Sn(t) ~ N(p,p(l — p)/n).
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
7
Príklad 10 (programovanie; empirická funkcia prežívania) Naprogramujte v algoritmus
na výpočet empirickej funkcie prežívania. cvič.
Príklad 11 (empirická funkcia prežívania, AML) Vypočítajte empirickú funkciu prežívania vo všetkých časoch t pre skupinu A (1) pre problém 1, (2) problém 2 a (3) pre problém 3. Zobrazte tieto tri funkcie do jedného obrázka ako schodovité funkcie. Použite tri rôzne typy čiar alebo tri rôzne farby a pridajte legendu a popis osí. cvič.
Riešenie (čiastkové - pre problém 3; pozri tabulku)
t	0	9	13	18	23	28	31	34	45	48	161
Sn(t)	11 11	10 íi	8 11	7 11	6 11	5 11	4 11	3 11	2 11	i íi	0 íi
problém 1 problém 2 a 3
I
50
I—
100
150
cas (v týždňoch)
Obr. 2: Empirická funkcia prežívania pre AML dáta (problém 1, 2 a 3)
Príklad 12 (programovanie; početnosti di a tíj) Naprogramujte v <® algoritmus na výpočet d,t
a rii pre všetky časy zlyhania ti. Výstupom bude tabulka so stĺpcami ti, di a n,i. cvič.
3.2   Odhady funkcie prežívania
Klasický Kaplan-Meierov (KM) odhad funkcie prežívania je definovaný nasledovne
i:ti<t
kde A(£j) = Aj = ^ (ide o maximálne vierohodný odhad rizika), di je počet zlyhaní (úmrtí) v čase ti < t & rii počet ľudí v riziku tesne pred časom t. KM odhad kumulatívneho rizika A-km (t) = — ln Skm (t) = — Y^i-t <tm — • Ak je najväčší čas cenzurovaný čas cmax > í/, potom podľa definície S (í/) = S (cmax) a keďže funkcia vierohodnosti nezávisí na S (t) pre t > cmax, nie je možný žiadny odhad za cmax.
Príklad 13 (programovanie; KM odhad funkcia prežívania) Naprogramujte v Qt algoritmus
na výpočet KM odhadu funkcie prežívania Skm (t) vo všetkých časoch zlyhania. cvič.
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
8
Príklad 14 (KM odhad funkcia prežívania, AML) Vypočítajte S km (t) vo všetkých časoch zlyhania a porovnajte ho s Sn (t) pre skupinu A. Zobrazte tieto dve funkcie do jedného obrázka ako schodovité funkcie. Použite dve rôzne typy čiar alebo dve rôzne farby a pridajte legendu a popis osi. cvič.
Riešenie (čiastkové; pozri tabulku)
t	0	9	13 13+	18	23	28+	31	34	45+	48	161+
Jn (t) Skm (t)	íi íi 1	10 íi 0.91	8 8 11 11 0.82 0.82	7 11 0.72	6 11 0.61	5 11 0.61	4 11 0.49	3 11 0.37	2 11 0.37	i íi 0.18	0 íi 0.18
Skm (0) -Skm (9) = Skm (13) Skm (13+ Skm (18) Skm (23) Skm (28+ Skm (31) Skm (34) Skm (45+ Skm (48) S'km (161-
Skm (0) x = Skm (9) x = Si
íi-i íi
10-1 10
km (13) x -g-
Skm (13) x Skm (18) x ^ = Skm (23) x SKM (23) x ^
Q 4-1
(31) x ^
'km
= Skm (34) x ^2 = Skm (34) x ^ l = Skm (48) x ^
Skm (13)
Skm (23)
Skm (34) = Skm (48)
-1—
100
cas (v týždňoch)
KM EFP
I—
150
Obr. 3: Dva odhady funkcie prežívania - Skm (t) a Sn - pre AML dáta (skupina A) Modifikovaný KM odhad funkcie prežívania je definovaný nasledovne
SKMmod(í) = _ [ (1 - Aj), kde 1 - Aj
i:ti<t
Lj=0   < n4-j
Breslowov (B) odhad funkcie prežívania definujeme nasledovne
SB (í) = exp (-Ana (í)) = _ [ exp (-A,
ak d{ = 1 ak dj > 1
i:ti<t
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
9
kde Ana (t) = ^2i,i.t.<t \ = ^2i-t<t ^~ Je Nelson-Aalenov (NA) odhad kumulatívneho rizika.
Flemingom a Harringtonom (FH) modifikovaný B odhad funkcie prežívania definujeme nasledovne
SFHmodB (t) = exp (— ApHmodNA (t)
' je FH modifikovaný NA odhad kumulatívneho rizika.
kde ApHmodNA - Y,í:ti<t [Ľj=0
Príklad 15 (programovanie; odhady funkcie prežívania) Naprogramujte v ^tt algoritmus na výpočet SjcMmod (t), S b (t) a SpHmodB (t) v každom čase zlyhania. cvič.
Príklad 16 (porovnanie; odhady funkcie prežívania, AML)  Vypočítajte nasledovné odhady funkcii prežívania - S km (t), SKMmod(t), S b (t) a SpHmodB(t) pre AML dáta (skupina A). Zobrazte tieto štyri funkcie do jedného obrázka ako schodovité funkcie. Použite štyri rôzne typy čiar alebo štyri rôzne farby a pridajte legendu a popis osí. cvič.
t  d n
KM
B
FHmodB
KMmod
#[1,]     9  1   11   0.9090909 0.9090909 0.9131007 0.9131007
#[2,]   13  1   10  0.8181818 0.8181818 0.8262077 0.8262077
8  0.7159091 0.7159091 0.7291257 0.7291257
7  0.6136364 0.6136364 0.6320630 0.6320630
5  0.4909091 0.4909091 0.5174894 0.5174894
4  0.3681818 0.3681818 0.4030212 0.4030212
2  0.1840909 0.1840909 0.2444447 0.2444447
#[3,] 18 1
#[4,] 23 1
#[5,] 31 1
#[6,] 34 1
#[7,] 48 1
KM a KMmod B a FHmodB
-1-
100
cas (v týždňoch)
I—
150
Obr. 4: Štyri odhady funkcie prežívania - Skm (t), SxMmod (t), Sb (t) a SpHmodB (t) - pre AML dáta (skupina A)
Príklad 17 (odhady funkcie prežívania)  Vypočítajte odhady SKMmod{t), S b (t) a S FHmodB {t) v
čase 12 a 88 (pozri tabuľku). pred
t	di	rii
4.5	1	70
11.5	2	68
16.0	1	65
20.7	2	55
20.8	1	53
31.0	1	47
34.5	1	45
46.0	1	34
61.0	1	25
87.5	5	15
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
10
Pre 2 zhody v čase 12 platí:
SKM (12) = (69/70) (66/68) = 0.9567
5KMmod (12) = (69/70) (67/68) (66/67) = (69/70) (66/68) = 0.9567
Druhý prípad predstavuje úpravu Skm (t) pri zlome zhôd rozdelením času 11.5 na 11.48 a 11.52, z
čoho je zrejmá invariantnosť Skm (t) na zhody.
SB (12) = exp [- (1/70 + 2/68)] = 0.9572
^FHmodB (12) = exp [- (1/70 + 1/68 + 1/67)] = 0.9570
Pre 5 zhôd v čase 88 platí:
Škm (88) = 0.5294
5FHmodB (88) = 0.5395
SB (88) = 0.5709
SFHmodB (t) dáva vo všeobecnosti odhad bližšie ku Skm (t) a je menší ako Sb (t) pri zhodách.
Príklad 18 (porovnanie; odhady funkcií prežívania)  Vypočítajte nasledovné odhady funkcii prežívania - Skm (t), SxMmod (t), S b (t) a Spu-modB (t) pre dáta (pozri tabuľku). Zobrazte tieto štyri funkcie do jedného obrázka ako schodovité funkcie. Použite štyri rôzne typy čiar alebo štyri rôzne farby a pridajte legendu a popis osí. Použite dáta z predchádzajúceho príkladu.
20
"T"
40
"T"
60
I
80
E x
KM a KMmod B
FHmodB
I
20
60
"T"
80
Obr. 5: Štyri odhady funkcie prežívania - Skm (t), SKMmod (t), Sb (t) a SFHmodB (t) ~ Pre AML dáta (skupina A)
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
11
9	#			t	d	n		KM		KMmod		B		FHmodB
10	#	[U	4	5	1	70	0	9857143	0	9857143	0	9858158	0	9858158
11	#	[2,]	11	5	2	68	0	9567227	0	9567227	0	9572435	0	9570334
12	#	[3,]	16	0	1	65	0	9420039	0	9420039	0	9426294	0	9424225
13	#	[4,]	20	7	2	55	0	9077492	0	9077492	0	9089677	0	9084623
14	#	[5,]	20	8	1	53	0	8906218	0	8906218	0	8919781	0	8914822
15	#	[6,]	31	0	1	47	0	8716724	0	8716724	0	8732004	0	8727148
16	#	[7,]	34	5	1	45	0	8523019	0	8523019	0	8540099	0	8535351
17	#	[8,]	46	0	1	34	0	8272342	0	8272342	0	8292578	0	8287967
18	#	[9,]	61	0	1	25	0	7941449	0	7941449	0	7967421	0	7962991
19	#[10 ,]		87	5	5	15	0	5294299	0	5294299	0	5708907	0	5395385
Rozdiely medzi odhadmi funkcií prežívania sú malé, s výnimkou veľmi malého n alebo pri veľkom počte zhôd. Odhad Skm (t) je invariantný na zhody. Odhad Sb (t) môžeme vylepšiť použitím odhadu ^FHmodB (t) v prípade velkého počtu zhôd. Ak e~x « 1 — x, pre malé x, potom sú oba odhady Skm (t) a Sb (t) podobné v prípade, že nárasty dA sú malé, čo znamená, že v riziku je veľa osôb. Oba odhady sú v skutočnosti asymptoticky ekvivalentné (n —► oo). Avšak e~x > 1 — x =>- Sb (t) > Skm (t) v prípade, že ide o konečné náhodné výbery (ale ak posledný príspevok je smrť, potom Skm (t) = 0, Sb (t) > 0).
Príklad 19 (porovnanie odhadov funkcie prežívania) Použitím Taylorovho rozvoja (prvého rádu) mínus logaritmu A^a (t) ukážte, že Akm^) je prvým členom tohoto rozvoja.
3.3   Odhady rizika a kumulatívneho rizika
Príklad 20 (riziko a časové jednotky) Závislosť hodnoty rizika na jednotkách času.
Pr(t < T < t + At\T >t) = \
Ak At = \ dňa, potom \{t) = f = 0.75 na deň.
3
i.
Ak At =     týždňa, potom \(t) = -jr = 5.25 na týždeň.
Odhad rizika v časovom intervale t,L <t < ti+1 je definovaný nasledovne
rií(ti+1 - U)''
Príklad 21 (riziko, AML) Vypočítajte odhad rizika \(t) v intervale U < t < ti+1 v čase 26 (pre skupinu A).
A(23) = \ = 0.143
Ä(26) = Ä(23) = 7{^_2i) = 0.018
Príklad 22 (programovanie; riziko) Naprogramujte v ^St algoritmus na výpočet odhadov rizika \(t) a \(t) pre všetky íj.
Príklad 23 (riziko, AML) Vypočítajte odhad rizika \(t) a porovnajte ho s odhadom \(t) v intervale ti < t < ti+i pre všetky ti.
20	#	cas	n . i d	i	lambda.KM		lambda	. INT
21	#1	9	11	1	0	0909	0 .	0227
22	#2	13	10	1	0	1000	0 .	0200
23	#3	18	8	1	0	1250	0 .	0250
24	#4	23	7	1	0	1429	0 .	0179
25	#5	31	5	1	0	2000	0 .	0667
26	#6	34	4	1	0	2500	0 .	0179
27	#7	48	2	1	0	5000		NA
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
12
Príklad 24 (programovanie, kumulatívne riziko) Naprogramujte v Qt algoritmus na výpočet odhadov kumulatívneho rizika Akm a A-na Vre všetky t. cvič.
Príklad 25 (kumulatívne riziko, AML) Vypočítajte Akm a A^a Vre všetky t (pre skupinu A). Zobrazte tieto dve funkcie do jedného obrázka ako schodovité funkcie. Použite dve rôzne typy čiar alebo dve rôzne farby a pridajte legendu a popis osi. cvič.
Obr. 6: Dva odhady kumulatívneho rizika - Akm a Ana _ pre AML dáta (skupina A)
3.4   Odhady rozptylu kumulatívneho rizika a funkcie prežívania
Greenwoodov (G) odhad rozptylu kumulatívneho rizika definujeme nasledovne
d2 (t) = Var g ÄKM (t)  = Y" ——
d,
^tni(rii - di)'
plug-in odhad rozptylu kumulatívneho rizika (nazýva sa aj Kleinov odhad; maximálne vierohodný odhad)
d2K(t) = VarK [A(í)l = ]T
di(rii - d,
i:ti<t
binomický odhad rozptylu kumulatívneho rizika
dl(t) = VarB [Ä (í)j = Y,
^ n2 (n,-, — 1):
i:U<t     1 V   * '
NA odhad rozptylu kumulatívneho rizika (nazýva sa aj Aalenov, Poissonov alebo Tsiatisov odhad)
ťxT (t) = Var
i:ti<t
n
2 '
FH modifikovaný NA odhad rozptylu kumulatívneho rizika
A
FHmodNA
fdi-\
«)| = E E
i:ti<t \j=0
i^í - j)'
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
13
Odhad d2G(t) je o niečo väčší ako cr\(t) v prípade vyšších hodnôt odhadu kumulatívneho rizika. Vo všeobecnosti platí crK(t) < cr|(£) < a\(ť). Ak sa nevyskytujú zhody, potom cr|(£) = a\(ť). Všetky odhady sú menšie ako <TG(t). Odhad <x|(£) môžeme v prípade zhôd vylepšiť použitím odhadu crFH(£). Odhad d* (t) má menšiu výchylku ako (xK(t). V praktických aplikáciách sú odhady rozptylu rôzne v závislosti od množstva zhôd v analyzovaných dátach. Ignorovaním sumy ^i.í<í dostaneme odhady rozptylov rizika v čase t, ktoré označujeme ako Oq(íj), 3K(ti), cr|(£j), a\(ti) a crFH(£j).
Príklad 26 (rozptyl kumulatívneho rizika, AML)  Vypočítajte odhady rozptylu kumulatívneho rizika (y2G{t), (y2K{t), (y2B{t), o2T{t) a crFjj{ť) v čase 26 pre AML dáta (skupina A). pred
^(26) = Var G [Arm (26)] ?K(26) = VarK [Ä (26) £|(26) = VarB í A (26)
11(11-1)	4- 1 T 10(10-	-1)	+ 8(8-1)
1(11-1) 113 -1-	1(10-1) 103	+	1(8-1) 1 83 1
1(11-1)	1(10-	1)	1 1(8-1)
112(11-1)	r 102(10	-1)	T 82(8-l)
ll2 ^ 102	+ j£r +	1 72	= 0.0543
7(7-1)
0.0619
73
^ = 0.0478
= 0.0543
Príklad 27 (programovanie, rozptyl kumulatívneho rizika) Naprogramujte v * na výpočet odhadov rozptylu kumulatívneho rizika aG(t), (y2K{t), (y2B{t), cr2T(t) a a2FH(t).
algoritmus
cvič.
Príklad 28 (rozptyl kumulatívneho rizika, AML) Vypočítajte odhady rozptylu kumulatívneho rizika o2G{t), (y2K{t), (y2B{t), o2T{t) a cŕFjj{ť) pre AML dáta (skupina A).
cvič.
28	#	cas	ni	di	rozptyl . G	rozptyl.K	rozptyl.B	rozptyl . T	rozptyl.FH
29	#1	9	11	1	0.0091	0.0075	0.0083	0.0083	0.0083
30	#2	13	10	1	0.0202	0.0165	0.0183	0.0183	0.0183
31	#3	18	8	1	0.0381	0.0302	0.0339	0.0339	0.0339
32	#4	23	7	1	0.0619	0.0477	0.0543	0.0543	0.0543
33	#5	31	5	1	0. 1119	0.0797	0.0943	0.0943	0.0943
34	#6	34	4	1	0.1952	0.1266	0.1568	0.1568	0.1568
35	#7	48	2	1	0.6952	0.2516	0.4068	0.4068	0.4068
Greenwoodov odhad rozptylu KM odhadu funkcie prežívania je definovaný ako
Var,
G
Srm (t) — Skm (t)
Var
ln Skm (ŕ)  = Skm (ŕ) ^o(ŕ)
Aalen-Johansenov odhad rozptylu KM odhadu funkcie prežívania píšeme ako
VarAj
Skm (t)  —  Skm (t)
Var
Ana (í)  = Skm (í) 32T(ť).
Oba odhady rozptylu funkcie prežívania odhadnú rozptyl menší ako v skutočnosti (underestimation) pre malé a stredne velké výbery. V priemere je Greenwoodov odhad bližšie k skutočnému rozptylu. Podobne ako pre KM odhad môžeme zapísať aj rozptyl B odhadu funkcie prežívania (Klein 2003)
Vare
SB(í)  = SB(Í)
Var
Arm (í)  =  SB(í) a2G{t).
Var^a
SB(í)  = SB(Í)
Var
Ana (í)  = SB(í)
t-
Pokiaľ máme v dátach zhody, použijeme nasledovné rozptyly FHmodB odhadu funkcie prežívania
Var q
S
FHmodB
(t)
S
FHmodB
(t)
Var
Arm (t)
S
FHmodB
(t) O&t).
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
14
Varn a
S
FHmodB
Si
FHmodB
Var
A-NA (t)    —   SpHmodB (t) 0"|H(í).
Vo všeobecnosti potom odhad rozptylu funkcie prežívania definujeme ako
Var
S (t)  = S (t)
Var
In S (t)
Príklad 29 (rozptyl funkcie prežívania, AML)  Vypočítajte odhad rozptylu KM a B odhadu funkcie prežívania v čase 13 pre skupinu A pomocou Greenwoodovej formule. pred
V^gJSkm(13)] = 0.822(TI(I^TI + -KarNA[SB( 13)] = 0.832(^ + 751) = 0.0125
0.0136
Príklad 30 (programovanie, rozptyl funkcie prežívania) Naprogramujte v       algoritmus na výpočet odhadov rozptylu KM a B odhadu funkcie prežívania. cvič.
Príklad 31 (rozptyl funkcie prežívania, AML)  Vypočítajte odhad rozptylu KM a B odhadu funkcie prežívania vo všetkých časoch pre skupinu A pomocou Greenwoodovej formule. cvič.
36	rozptyl		.funkcie.			prezivania(survobjA.KM)	
37	#	kontrola		--	- Gl	je  rozptyl vypocitany	funkc iou
38	#	cas	ni	di	r oz	ptyl.Gl  rozptyl.G2	
39	#1	9	11	1		0.0075 0.0075	
40	#2	13	10	1		0.0135 0.0135	
41	#3	18	8	1		0.0195 0.0195	
42	#4	23	7	1		0.0233 0.0233	
43	#5	31	5	1		0.0270 0.0270	
44	#6	34	4	1		0.0265 0.0265	
45	#7	48	2	1		0.0236 0.0236	
46	B .	rozptyl .		funkcie.prezivania(survobjA.B)			
47	#	kontrola		--	- BI	je  rozptyl vypocitany	funkc iou
48	#	cas	ni	di	r oz	ptyl.BI  rozptyl.B2	
49	#1	9	11	1		0.0069 0.0069	
50	#2	13	10	1		0.0125 0.0125	
51	#3	18	8	1		0.0180 0.0180	
52	#4	23	7	1		0.0217 0.0217	
53	#5	31	5	1		0.0253 0.0253	
54	#6	34	4	1		0.0255 0.0255	
55	#7	48	2	1		0.0243 0.0243	
3.5   Empirické intervaly a pásy spoľahlivosti 3.5.1    Intervaly spoľahlivosti
Pri výpočte hraníc empirických intervalov spoľahlivosti (ISs) sa používajú nasledovné princípy
1. Waldov princíp, skóre princíp a vierohodnostný princíp (všetky tri vychádzajúce z funkcie vierohodnosti),
2. princíp transformácie funkcie vierohodnosti v S (t) pomocou nejakej funkcie g (S (t)) aplikovaný na Waldov princíp (hranice IS vypočítané pomocou g (S (t)) sa spätne transformujú do škály S (t)), kde podľa g (S (t)) rozoznávame nasledovné škály
• g (S (t)) = S (t) - škála funkcie prežívania
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
15
g (S (t)) = ln S (t) - škála logaritmu funkcie prežívania (log škála funkcie prežívania: škála kumulatívneho rizika) - spätná transformácia pre funkciu prežívania exp (ln S (t)),
g (S (t)) = ln(— lnS'(í)) log-log škála funkcie prežívania (log škála kumulatívneho rizika; škála logaritmu kumulatívneho rizika)
spätná transformácia pre funkciu prežívania exp (exp (— ln S (t))), spätná transformácia pre kumulatívne riziko exp (— ln S (t)).
g (S (t)) = arcsin ^(S(t))1^2^    arcus sínus odmocninová škála funkcie prežívania
spätná transformácia pre funkciu prežívania (sin (arcsin ((S(t)Ý^2
g (S (t)) = arcsin ^exp ^_ln^W^ _
arcus sínus odmocninová škála kumulatívneho
rizika
spätná transformácia pre kumulatívne riziko —2 ln sin ( (arcsin (exp
lnS(ť)
Pre extrémne malé alebo extrémne velké hodnoty t sa môže stať, že empirický interval spoľahlivosti pre S (t) má hranice mimo intervalu (0,1). Ak budeme predpokladať asymptoticky normálne rozdelenie nejakej transformácie S (t), ktorej hodnoty nie sú ohraničené, môžeme sa tomuto problému vyhnúť.
3. princíp úpravy hraníc IS pomocou efektívneho rozsahu súboru v čase t (n*(t) alebo n** (t)) z dôvodu výskytu cenzúr v dátach (resp. rozdielu medzi rozptylom S (t) v čase t vypočítaným pre cenzurované dáta a rozptylom za predpokladu, že cenzúry v dátach nie sú)
4. princíp korekcie hraníc IS z dôvodu zlých štatistických vlastností (konzervatívny alebo liberálny IS, t.j. predpokladáme, že nominálny koeficient spoľahlivosti 1 — a je iný ako teoretický)
3.5.2 Intervaly spoľahlivosti pre kumulatívne riziko
Pri tvorbe IS pre A(í) máme tri možnosti.
1. Akm(í) a af;(t),
2. Ana(č) a ^t(^) (namiesto môžeme použiť aj cr|(£) alebo SkW);
3-  ApjjmodNA(t) a °fhW;
3.5.3 Intervaly spoľahlivosti pre funkciu prežívania
Pri tvorbe IS pre S (t) máme päť možností podľa použitých odhadov funkcií prežívania a ich rozptylov
1- Skm (t) a Var g SKM (t)
2. SB (í) a VarNA SB (í)
3. SB (t) a VarAj SB (t)
4. S-
FHmodB
(í) a Var^a
S-
FHmodB
(t)
5. S-
FHmodB
(í) a VarQ
Si
FHmodB
(t)
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
16
Príklad 32 (intervaly spoľahlivosti pre kumulatívne riziko) Naprogramujte v Qt 100 x (1 — a)% intervaly spoľahlivosti pre kumulatívne riziko
(a) v škále kumulatívneho rizika,
(b) v log škále kumulatívneho rizika.
Pri výpočte použite AkmÍí) a Ap^modNA (t) a ich rozptyly. cvič.
Príklad 33 (intervaly spoľahlivosti pre kumulatívne riziko) Vypočítajte v 95% intervaly spoľahlivosti pre kumulatívne riziko v oboch škálach pre AML dáta (skupina A). Zobrazte do okna 2x1 odhady kumulatívneho rizika a IS pre kumulatívne riziko v tvare písmena „I" v každom čase zlyhania v oboch škálach. Použite (a) Akm^) a (b) ApHmodNA (t) a ich rozptyly. cvič.
Príklad 34 (intervaly spoľahlivosti pre funkciu prežívania) Naprogramujte v      100 x (1 — a)% intervaly spoľahlivosti pre funkciu prežívania
(a) v škále funkcie prežívania,
(b) v log škále funkcie prežívania a
(c) v log-log škále funkcie prežívania.
Pri výpočte použite Skm (t) a SFHmodB (t) a ich rozptyly. cvič.
Príklad 35 (intervaly spoľahlivosti pre funkciu prežívania) Vypočítajte v ^íl 95% intervaly spoľahlivosti pre funkciu prežívania vo všetkých troch škálach pre AML dáta (skupina A). Zobrazte do okna 3x1 odhady funkcie prežívania a IS pre funkciu prežívania v tvare písmena „I" v každom čase zlyhania vo všetkých troch škálach. Použite (a) Skm (t) a (b) SpHmodB (t) a ich rozptyly. cvič.
3.6 Odhady strednej hodnoty a mediánu prežívania Pre odhad strednej hodnoty času zlyhania E [T] platí
(j,
S (ŕ) dt,
kde S (t) je KM odhad funkcie prežívania a £max je maximum pozorovaných časov. Tomuto odhadu hovoríme aj urezaný.
V diskrétnom prípade odhad strednej hodnoty vypočítame nasledovne
max(ímax,cmax) 1
í+i
U)S{U
i=0
kde í0 = 0a/<nje počet rôznych zlyhaní a tj = £max. Ak ímax < cmax, kde cmax je maximálnym časom do cenzúry, potom odhad strednej hodnoty počítame na intervale (0,cmax). Pre KM odhad funkcie prežívania je odhad strednej hodnoty prežívania nedefinovateľný, ak je posledné pozorovanie cenzúra. Je možné definíciu modifikovať tak, že odhad je nulový za posledným pozorovaním. Takýto odhad strednej hodnoty je síce vychýlený smerom k nule, ale doteraz nebola objavená lepšia alternatíva. Podľa použitých odhadov funkcií prežívania rozlišujeme nasledovné odhady strednej hodnoty
A^KM = /ÍKMmod, /iB a /^FHmodB■
Potom pre odhad rozptylu strednej hodnoty času prežívania platí
-, 2
Var
o
S (u) du
d2(t)dt
a v diskrétnom prípade
Var [fi] =
í:íj<max(ŕmax,cmax) 1
E
(tj+1 - tj)S(tj)
j'.ti ^tj ^Iliax(ŕmaxiCmax) 1
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
17
56
57
58
59
60
61
62
63
64
Za S (t) dosadíme Skm (t), Sb (t) alebo SFHmodB (t). Podobne ako pri rozptyle funkcie prežívania dostaneme päť nasledovných rozptylov
VarG [/íkm] = VarG [//RMmod], V"arNA [//b], V"arAj [//km], V"arNA [/tFHmodB] a VarG [/iFHmodB]-
Príklad 36 (priemerný čas prežívania, AML)  Vypočítajte priemerný čas prežívania fi, jeho rozptyl Var\p\ a 100x(l — a)% intervalov spoľahlivosti pre /i pre AML dáta (skupina A). pred
/2 = 52.6 týždňa, Var^i] = 19.82,95%IS = (13.792,91.408) týždňa
Príklad 37 (programovanie, priemerný čas prežívania) Naprogramujte v <® algoritmus na výpočet obsahu pod krivkou prežívania. cvič.
Príklad 38 (priemerný čas prežívania, AML) Vypočítajte priemerný čas prežívania pre AML dáta (skupina A). Pri výpočte použite Skm (t) a S b (t). Okomentujte výsledky. Nakreslite krivku odhadu funkcie prežívania a pomocou funkcie polygonO vyfarbite obsah pod ňou.
# výsledok 60.30259
# kontrola
print(survobj A.B,rmean="individual " )   # 60.3
Keďže Skm (t) < SB (t), potom /xKM < //KM.
Príklad 39 (programovanie, rozptyl priemerného času prežívania) Naprogramujte v Q$ algoritmus na výpočet rozptylu priemerného času prežívania. cvič.
Príklad 40 (rozptyl priemerného času prežívania, AML)  Vypočítajte rozptyl a smerodajnú odchýlku priemerného času prežívania pre AML dáta (skupina A). Pri výpočte použite Skm (t) a S b (t). Okomentujte výsledky.
rozptyl.mu(survobjA.KM,error="G") # 393.1735 rozptyl.mu(survobjA.B,error="T")     # 414.8907
rozptyl.mu(survobjA.B,error="G")     #  655.1461   (ako  přednastavené  v R) sqrt(rozptyl.mu(survobj A.KM,error = "G"))   # 19.8286 sqrt(rozptyl.mu(survobj A.B,error = "T"))   # 20.36887 sqrt(rozptyl.mu(survobj A.B,error = "G"))   # 25.59582
Medián času prežívania je 50-ty percentil íq.5- Medián funkcie prežívania je potom S^o.s) = 0.5. Výberový medián je definovaný ako prvý čas, v ktorom S (to.5) = 0.5. Teda (Jl = to.5 < S-1 (0.5). Niekedy je potrebné použiť lineárnu interpoláciu pre to.5 v podobe
A   S (U) - 0.5 S (U) - 0.5
/w - +   -     -š(u+i) ~tl + m 'tl to'5 tl+1"
kde f (U) = s^^fa+i) ^ čo je pravdepodobnosť zlyhania v intervale (U+i — ti) škálovaná jeho dĺžkou.
Horná a dolná hranica intervalu spoľahlivosti (IS) pre medián je definovaná na základe IS pre S (t) v danom čase, t.j. horná hranica IS pre medián je prvý čas, v ktorom je horná hranica IS pre S (t) väčšia alebo rovná 0.5 (rovnaký vzťah platí pre dolnú hranicu IS pre medián). To korešponduje s narysovaním horizontálnej úsečky na grafe krivky prežívania, tj. přetnutím tejto úsečky s krivkou prežívania, dolnou a hornou hranicou IS pre S (t) v danom čase.
Nech tp je p-ty kvantil rozdelenia T (100 x p-ty percentil), potom
S(tp) = Pr(T >tp) = l-p,tp< S-^l - p).
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
18
Keďže KM krivka prežívania je schodovitá funkcia, inverzia S,_1(íp) nie je jednoznačne definovaná. Odhad kvantilu bude potom
tp = min{ti : S (U) < 1 - p}. Aplikovaním delta metódy na VarG[S(tp)] dostaneme
Var [t%
VarG[S(tp)] IfítpW
S(up) - S(lp)
lp Up
kde up = max{íj : S (ti) > 1—p + e}a/p = min{£j : S (ti) < 1 — p — e} pre i = 1,2,..., J < n. kde u korešponduje s hornou hranicou, / s dolnou, J je počet rozdielnych časov zlyhania, e je veľmi
malé číslo. Najčastejšie e = 0.05 akceptovateľné, ale musí byť velké, ak \lp dosadíme Skm(£); 5b (t) alebo SpHmodB (t)- Potom dostaneme t
Ur,
0. Za S (t)
-jxkm) _ £(KMmod)   J(B) ^ ^(FHmodB)
, Lp     a, Lp
Var
Í(KM) tp
Var
r(KMmod) tp
Var
ífB)
Lp
a Var
7(FHmodB) tp
Príklad 41 (medián, AML) Vypočítajte medián času prežívania fi, jeho rozptyl Var■[//] a 100x(l —
a)% intervalov spoľahlivosti pre í0.5 P^e ^4ML daŕa (skupina A). pred
A4 = ŕo.5, V"ar[
= V«r-o[5(F0.5)l [/(Í0.5)]2
^0.5 = 31 týždňov, Mo.5 = max{£j : S (U) > 0.55} = 23 a /o.5 = min{£j : S (U) < 0.45} = 34
7(21) = S(u0.5)-Š(T0.5) _ S(23)-5(34) _ 0.6136364-0.3681818 _ g
í).5-«0.5
yCľíTIl - í 0-16419327^ 1/ ar [ól\ — k 0.02231405 ^
34-23
11
54.144, JVar[31] = 7.358,95% IS = (16.578,45.422) týždňa
Príklad 42 (programovanie, kvantily) Naprogramujte v     funkciu na výpočet odhadu kvantilov času prežívania tp, Var[tp] a 100 x (1 — a)% intervalov spoľahlivosti pre tp. DU
Príklad 43 (medián, AML) Vypočítajte medián ťo.5 a jeho 100 x (1 — a)% intervalov spoľahlivosti pre AML dáta (skupina A).
65
66
67
68
69
70
kvantily.km(survobj A.KM,p  =  0.5,   eps  =  0.05,   conf.coef   = 0.975)
# kvp sd(kvp) DHIS HHIS #[1,]     31     7.3583   16.578 45.422
kvantily.km(survobj A.B,p =  0.5,   eps  =  0.05,   conf.coef  = 0.975)
# kvp sd(kvp) DHIS HHIS #[1,]     34    8.5516   17.2392 50.7608
3.7   Odhad strednej hodnoty a mediánu zostatkového života a ich rozptyl
Odhad strednej hodnoty zostatkového života v čase t nazývame priemerný zostatkový život v
čase t a vypočítame ho v spojitom prípade ako
mri (í)
t b(u)du
a v diskrétnom prípade ako
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
19
kde U < t < U+i.
Za S(-) dosadíme Skm('), Sb(') alebo í>FHmodB (') • Rozptyl priemerného zostatkového života odhadneme nasledovne
Var [mrl(í)]
S2 (t)
S(x)dx
crG(u) du
1 2 rt
S(u)du
u t
ď^iiŕjdu
pre spojitý prípad, kde u G (t,tmax), a ako
1
Var [mrl(í)]
S'2(í)
E
^t^ti^maxftmaxjCmax)-!
(ť,+1 - ť)S(ťi) +
E
(tj+1 - t j) S (t j)
a2G(tt
E
{tj+i - t j) S (t j)
J-U+i
^
i:ti<t
pre diskrétny prípad, kde ti < t < ti+1. Za S (t) dosadíme Skm (t), Sb (t) alebo SpHmodB (t). Potom dostaneme mrlxM (t), mrie (t) a mrlpHmodB (t) päť nasledovných rozptylov
Vare
mri
KM
Vare
mri
KMmod
Var^a
mrh
VarAJ
mri
KM
Var^a
mri
FHmodB
a V^arc
mri
FHmodB
Medián zostatkového života (median residual life, median remaining lifetime) definujeme ako
—r, ,    ,       ,    ,          .Síti) - 0.5S(t)     .       .    5(íi) - 0.5S(t) mrl(í) = (U -t) + (ti+1 - ti)J^-^—^ = (U-t) + -K
S(U) - S(ti+1)
f (U
kde U < t < U+i, S (U) > 0.5S(t) a S(U+i) < 0.5S(t). To znamená, že ak sme mali pravdepodobnosť prežitia S (t) v čase t, potom budeme mať nažive polovicu z tohoto množstva, t.j. 0.5S(t), v čase t + mrl(í). Rozptyl odhadneme nasledovne
Var
mrl(í)
— ÍS^
Za S (t) dosadíme Skm (t), Sb (t) alebo SFHmodB (t)- Potom dostaneme mrlxM (t), mrie (t) a mrlpHmodB (t) päť nasledovných rozptylov
Varc
mrl
KM
Var q
mri
KMmod
Var
na
mrh
Var
aj
mri
KM
VarNA
mri
FHmodB
a Var q
mri
FHmodB
Príklad 44 (programovanie, priemerný zostatkový život) Naprogramujte v     funkciu na výpočet odhadu strednej hodnoty zostatkového života v čase t. cvič.
Príklad 45 (priemerný zostatkový život, AML) Vypočítajte priemerný zostatkový život v čase t = 0 a t = 30 pre AML dáta (skupina A). Nakreslite krivku odhadu funkcie prežívania a pomocou funkcie polygonO vyfarbite obsah pod ňou za bodom t. cvič.
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
20
71
72
mri(survobjA.KM,   tcko  =  0)$mrl   #  52.64545   52.64545 52.64545 mri(survobjA.KM,   tcko  =  30)$mrl   #  45.70000   52.73077 52.57692
Príklad 46 (programovanie, rozptyl priemerného zostatkového života) Naprogramujte v funkciu na výpočet odhadu rozptylu priemerného zostatkového života v čase t. DÚ
Príklad 47 (rozptyl priemerného zostatkového života, AML)  Vypočítajte odhadu rozptylu priemerného zostatkového života v čase t = 0 a t = 30 pre AML dáta (skupina A). DU
Príklad 48 (programovanie, medián zostatkového života) Naprogramujte v 1í funkciu na výpočet odhadu mediánu zostatkového života v čase t. DU
Príklad 49 (medián zostatkového života, AML)  Vypočítajte medián zostatkového života v čase
t = 0 a t = 30 pre AML dáta (skupina A). DÚ
Príklad 50 (programovanie, rozptyl mediánu zostatkového života) Naprogramujte v funkciu na výpočet odhadu rozptylu mediánu zostatkového života v čase t. DÚ
Príklad 51 (rozptyl mediánu zostatkového života, AML)  Vypočítajte odhadu rozptylu mediánu zostatkového života v čase t = 0 a t = 30 pre AML dáta (skupina A). DÚ
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
21
4   Testy na porovnanie kriviek prežívania
Príklad 52 (nádor pľúc) Nech tij,i = 1,... ,rij,j = 1,2 sú časy do zlyhania (úmrtia) od diagnostiky nádoru pľúc v mesiacoch, kde j = 1 predstavuje I. typ terapie a j = 2 zasa II. typ terapie (pozri tabuľku). Otestujte H0 : S± (t) = S2 (t) oproti alternatíve H\ : S± (t) 7^ S2 (t) pomocou Sw a Smw nesledovne (1) Sw = Wy a Sw = Wx, (2) Smw = Uy a Smw = Ux- Vždy presne naformulujte H\. Okomentujte výsledky. cvič.
tij	52 240	19 53	15 43	340	133   111 231	378 49
skup	1 2	2 1	1 2	2	1      1 2	1 1
Určenie poradí
t (i)	15	19	43	49	52	53	111	133 231	240	340 378
skup	1	2	2	1	1	1	1	1 2	2	2 1
r(i) í r(2) í	1	-	-	4	5	6	7	8	-	- 12
	-	2	3	-	-	-	-	9	10	11
Príklad 53 (asymptotická normalita Smw) Pre ííi = 7 a íí2 = 5 porovnajte v <5lt asymptotické rozdelenie Smw s JeJ exaktným rozdelením. Na výpočet asymptotickej hustoty použite funkciu dnorm() a na výpočet asymptotickej distribučnej funkcie použite funkcie dnorm() a cumsum(). Na výpočet exaktnej hustoty použite funkciu dwilcoxO a na výpočet exaktnej distribučnej funkcie použite funkciu pwilcoxO. Teoretické a exaktné rozdelenie superponujte v podobe (t) hustoty, (2) distribučnej funkcie a (3) qq-diagramu s qq-priamkou (na x-ovej osi bude sekvencia x od teoreticky možného min(SMw) P° teoreticky možné max(SMw) a na y-°veJ osí teoretické kvantily y vypočítané pomocou funkcie qnorm(); qq-priamka bude prechádzať bodmi (s0.25,ío.75) a (^0.25,^0.75)^- cvič.
Príklad 54 (asymptotická normalita Smw) Porovnajte v^ň asymptotické rozdelenie Smw s JeJ exaktným rozdelením pre (t) n\ = 5 a n>2 = 50, (2) n\ = 50 a n>2 = 50, (3) n\ = 50 a n>2 = 100 a (4) ri\ = 100 a ri2 = 100. cvič.
Príklad 55 (nádor pľúc pokrač.) Nech tij,i = 1,... ,rij,j = 1,2 sú časy do zlyhania (úmrtia) od diagnostiky nádoru pľúc v mesiacoch, kde j = 1 predstavuje I. typ terapie a j = 2 zasa II. typ terapie (pozri tabulku). Otestujte v Hq : Var\S\ (t)] = Var[S2 (t)] oproti alternatíve H\ : Var\S\ (t)] 7^ Var[S2 (t)] pomocou S$t a ^st- Okomentujte výsledky. cvič.
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
22
Riešenie (aj v ^B$):
Siegel-Tukey test
tij skup	52 1	240 2	19 2	53 1	15 1	43 2	340 2	133 1	111 1	231 2	378 49 1 1
Ä<" R?	9	6	3	11	1	5	4	10	12	8	2 7
Príklad 56 Naprogramujte v*© test Hq : Var[S± (t)] = Var[S2 (t)] oproti alternatíve H\ : Var[S± (t)] ^ Var[S2 (t)] pomocou S$t a S^ použitím algoritmu 1. Okomentujte výsledky. cvič.
Príklad 57 (dvojčatá; Wilcoxonov znamienkový test) Pri skúmaní agresívnosti dvojčat psychológovia zaznamenali nasledovné skóre (pozri tabulku). Otestujte Hq na a = 0.05, že prvorodené dvojča je agresívnejšie. Použite Wilcoxonov znamienkový test, t.j. štatistiky (í) Sw, (2) Sw a (3) Sw s korkeciou rozptylu na zhody. Porovnajte výsledky. pred.
i	1	2	3	4	5	6	7	8	9	10	11
Xi dvojča 1	85	70	76	68	90	72	77	90	70	71	87
Yi dvojča 2	87	76	75	64	95	72	65	89	65	80	80
Zi = Xi — Yi	-2	-6	1	4	-5	0	12	1	5	-9	7
\Zi\	2	6	1	4	5	-	12	1	5	9	7
Ri	3	7	1.5	4	5.5	-	10	1.5	5.5	9	8
Príklad 58 (asymptotická normalita S kw) Pre ni = n2 = n3 = 5 porovnajte v asymptotické rozdelenie Skw s JeJ exaktným rozdelením. Na výpočet asymptotickej hustoty použite funkciu dchisqO a na výpočet asymptotickej distribučnej funkcie použite funkcie dchisqO a cumsum(). Na výpočet exaktnej hustoty použite funkciu dwilcoxO a na výpočet exaktnej distribučnej funkcie použite funkciu dKruskalWallis() v knižnici SuppDists. Teoretické a exaktné rozdelenie superponujte v podobe (í) hustoty, (2) distribučnej funkcie a (3) qq-diagramu s qq-priamkou (na x-ovej osi bude sekvencia x od teoreticky možného min(SKw) P° teoreticky možné max(S'xty) a na y-°veJ osi teoretické kvantily y vypočítané pomocou funkcie qchisqO; qq-priamka bude prechádzať bodmi (2?o.25) s0.75) a (^0.25)^0.75) alebo alternatívne bude qq-priamku reprezentovať os prvého a tretieho kvadrantu). D U
Príklad 59 (asymptotická normalita Skw) Porovnajte v asymptotické rozdelenie Skw s JeJ exaktným rozdelením pre (í) n\ = n>2 = 5 a n% = 50, (2) n\ = 122 = n% = 50, (3) n\ = n>2 = 50 a n3 = 100 a (4) m = n2 = n3 = 100. DÚ
Príklad 60 (WBC pokrač.) Majme pacientov s akútnou myeloidnou leukémiou a v rámci nich skupinu AG-pozitívnych (výskyt určitých špecifických indikátorov choroby v kostnej dreni). Pre chorobu je charakteristické, že s počtom bielych krviniek (white blood cells counts, WBC) vzrastná závažnosť choroby. Nech ti,i = 1,2,..., 17 sú časy do zlyhania v týždňoch (pozri tab.). (a) Otestujte Hq : Si (ŕ) = S2 (í) = S3 (í), alternatíva a) H\ : Si (t) 7^ S j (t), i < j; i, j = 1,2,3, b)
Si(í) y S2(í) y S3(í), c) Si(í) í S2 (í) í S3(í). Použite (í) SKW, (2) Sj, (3) Sc a (4) SL. (b) Vypočítajte Kendalov korelačný koeficient t medzi časom do zlyhania a príslušnosťou do skupiny 1,2, a 3 pomocou Sj. (c) Porovnajte Sc a S;v- (d) Otestujte odklon od trendu, (e) Graficky znázornite príslušnosť do skupiny voči času do zlyhania spolu s krivkou spájajúcou priemerný čas do zlyhania v každej skupine. Okomentujte adekvátnosť linárneho trendu. pred.
Riešenie:
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
23
Tabuľk
a:
WBC	U
2300	65
750	156
4300	100
2600	134
6000	16;;
10500	108
10000	121
17000	4
5400	39
7000	143
9400	56
32000	26
35000	22
100000	1
100000	1
52000	5
100000	65
Rozdeľme AG-pozitívnych pacientov do troch skupín nasledovne
• Skupina 1: WBC > 100000, nľ = 3, (1, 1, 65)
• Skupina 2: WBC G (10000,100000), n2 = 6, (108, 121, 4, 26, 22, 5),
• Skupina 3: WBC < 10000, n3 = 8, (65, 156, 100, 134, 16, 39, 143, 56)
Tabuľka:
sk.l	1	1					
sk.2	-	-	4 5	-	22	26 -	
sk.3	-	-	-	16	-	- 39	
poradie	1.5	1.5	3 4	5	6	7 8	
sk.l	-	65	-	-	-	-	-
sk.2	-	-	-	108	121	-	-
sk.3	56	65	100	-	-	134 143	156
poradie	9	10.5	12	13	14	15 16	17
i?i = 4.50, R2 = 7.833, R3 = 11.5625
skw = ^Jiy EjU 7^-3 (n + !) = TŤ^iš p-hodnota= 0.0924
3 x 4.52 + 6 x 7.833
x 11.56252
-3xlž
4.762662.
SL = Y?j=1 n j (L j - M j) Rj = 3 x(0- 14) x 4.5 + 6 x (3 - 8) x 7.833 + 8 x (9 - 0) x 11.5625 = 408.5
Var [SL] = ZU ni (Li ~ Mrf = t3 x (° " 14)2 + 6 x (3 - 8)2
+8 x (9 - O)2] = 187.99732
ZL = 408.5/187.9973 = 2.172903, p-hodnota=0.0298
SKW ~ (ZL)2 = 4.762662 - 2.1729032 = 0.0412, p-hodnota=0.8392
n-1-1-r
1.5       2.0       2.5 3.0
skupina
Obr. 8: Lineárny trend pre WBC dáta
Príklad 61 (nádor pľúc pokrač.) Nech tij,i = 1,... ,rij,j = 1,2 sú časy do zlyhania (úmrtia) od diagnostiky nádoru pľúc v mesiacoch, kde j = 1 predstavuje I. typ terapie a j = 2 zasa II. typ terapie (pozri tabuľku). Otestujte Hq : S± (t) = S2 (t), alternatíva H\ : S± (t) ^ S2 (t). Použite (í) Skw, (2) Sj, (3) Sc a (4) Sl- Vždy presne naformulujte H\.
Príklad 62 (pokrač. WBC) Majme pacientov s akútnou myeloidnou leukémiou a v rámci nich skupinu AG-pozitívnych (výskyt určitých špecifických indikátorov choroby v kostnej dreni). Pre chorobu je charakteristické, že s počtom bielych krviniek (white blood cells counts, WBC) vzrastná
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
24
závažnost choroby. Nech ti,i = 1,2,..., 17 sú časy do zlyhania v týždňoch prislúchajúce zoradeným WBC (pozri tab.). Vypočítajte Kendallov korelačný koeficient r. Otestujte nezávislosť medzi počtom bielych krviniek a časmi do zlyhania na hladine významnosti a = 0.05. DTJ
WBC	U	Cij	díj
750	156	0	16
2300	65	5	9
2600	134	1	13
4300	100	3	10
5400	39	5	7
6000	16	7	4
7000	143	0	10
9400	56	3	6
10000	121	0	8
10500	108	0	7
17000	4	4	2
32000	26	1	4
35000	22	1	3
52000	5	1	2
100000	1	1	0
100000	1	1	0
100000	65	0	0
Riešenie:
cíj = #(í tí, T WBCi) pod i, c = Y, Cij, dtJ = #(| U, t WBd) pod i, d = Y,dij, c = 33, d = 101, = t —   c~d   — —0 5
í^ = á» = 0.032,
z~ = -0.5/^0.032 = -2.801 a p-hodnota=0.005.
Príklad 63 (WBC pokrač.) Vypočítajte Spearmanov korelačný koeficient r g. Otestujte nezávislost medzi počtom bielych krviniek a časmi do zlyhania pomocou Zg na hladine významnosti a = 0.05.
Príklad 64 (testy pre cenzurované dáta) Majme dáta z klinickej štúdie zhrnuté v nasledovnej tabuľke (pozri tabuľku), (a) Vytvorte kontingenčné tabuľky v každom čase zlyhania ti, i = 1,2,... ,7 použitím celkového počtu subjektov v riziku rii v čase ti, celkového počtu zlyhanídi v čase ti, celkového počtu subjektov prvej skupiny v riziku riu v čase ti a celkového počtu zlyhaní du subjektov prvej skupiny v čase ti. (b) Vypočítajte stredné hodnoty Ea[du\, rozdiely empirických a očakávaných početností du — Ea[dii], ako aj rozptyly Vara[du\. (c) Otestujte H0 : Ai (t) = X2 (t) oproti H\ : Ai (t) = 9X2 (t) pomocou testovacích štatistík Qgw, Qcm, Qtw a Qpp- (d) Nakreslite Kaplan-Meierove odhady funkcie prežívania pre obe skupiny do jedného obrázka, (e) Vypočítajte (t) 9p, Var[0MH] a 95%IS pre 9p, (2) Omh a 95%IS pre QMH a (3) 0*MH. cvič.
Riešenie:
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
25
U	Tli	di	nu	dá	E0[du]	dá	-Eolditl	Var0[dií]
3	10	1	5	1	0.50		0.50	0.2500
5	9	1	4	1	0.44		0.56	0.2469
7	8	1	3	1	0.38		0.62	0.2344
12	6	1	1	0	0.17		-0.17	0.1389
18	5	1	1	1	0.20		0.80	0.1600
19	4	1	0	0	0.00		0	0
20	3	1	0	0	0.00		0	0
suma				4	1.69	2.31		1.0302
Q = 2.312/1.0302 = 5.179674, p-hodnota=0.02285261. zq = ^2.3171.0302 = 2.275890, p-hodnota=2 x 0.01142630 = 0.02285259 KT v každom čase t f.
	di	a,i	sp	d2	a2	sp	d3	a3	sp	d i	04	sp
skup. 1	1	4	5	1	3	4	1	2	3	0	1	1
skup. 2	0	5	5	0	5	5	0	5	5	1	4	5
sp	1	9	10	1 8		9	1	7	8	1	5	6
	db	a5	sp	de	a6	sp		07	sp			
skup. 1	1	0	1	0	0	0	0	0	0			
skup. 2	0	4	4	1	3	4	1	2	3			
sp	1	4	5	1	3	4	1	2	3			
testovacie kritérium	Q	p-hodnot a
GW	4.695652	0.03023902
CM	5.197242	0.02262275
TW	4.970637	0.02578116
PP	4.732935	0.02959035
>
N Q.
o
00 o
CO
o
o
cg o
o o
10
15
20
-r
25
cas
Obr. 9: Kaplan Meierove odhady funkcie prežívania
Príklad 65 (AML, pokrač.) Nakreslíte (a) kumulatívne riziko AxM,j(t), j = 1,2 pre obe skupiny, (b) kumulatívne riziko ANA,j(t),j = 1,2, (c) Kaplan-Meierove krivky prežívania SKM,j(t),j = 1,2. (d) Otestujte H0 : \i(t) = X2 (t) oproti H\ : \i(t) = 9\2(t) pomocou testovacích štatistík Qmh a Qpp. Použite funkcie survdiff () s argumentami rho=0 (Qmh) a rho=l (Qpp). (e) Otestujte Hq : Ai (t) = X2 (t) oproti H\ : Ai (t) = 9X2 (t) pomocou testovacích štatistík Qgw , Qcm, Qtw a Qpp- (f) Vypočítajte (í) 0P, Var[0MH] a95%IS pre 9p, (2) Qmh a 95%IS pre Qmh a (3) Q*mh-
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
26
KM kumulativně riziko pre AML data
NA kumulativně riziko pre AML data
>
£
ZJ xl
skup M skup NM
i
i      i      i      i r
10 20 30 40 50 cas do relapsu {v týždňoch)
C T-
i «
T"
skup M skup NM
r
T"
T"
T"
T"
0      10     20     30     40 50 cas do relapsu (v týždňoch)
Obr. 10: Odhady kumulatívneho rizika pre AML dáta
Obr. 11: Kaplán Meierove odhady funkcie prežívania pre AML dáta
5   Charakteristiky definované sčítacím procesom
Vo formuláciách sčítacím procesom (Xi,ôi) nahradíme (Ni (t) ,Yi (t)), kde N i (ť) je počet pozorovaných udalostí v intervale (0,t) v jednotke i,
y , , _ Í 1 jednotka i je v riziku v čase t (pozorujem ju) ^ ' ~ \ 0 inak
Táto formulácia obsahuje dáta s pravým typom cenzúr ako špeciálny prípad, teda Y,j_ (t) = I ({Ti > t}) a Ni (t) = I ({Ti <t, ô = 1}). Všimnime si, že N (t) je zprava spojitá a Y (t) zľava spojitá. Y (t) je príkladom predikovateľného procesu, ktorého hodnoty v čase t sú známe nekonečne krátko pred t, v čase t~, ak nie skôr. S čítací proces je stochastický proces začínajúci v čase 0, ktorého trajektória je zprava spojitá funkcia so skokmi veľkosti 1. Pre N (t) potom platí {N (t) : t > 0} , N (0) = 0 [7].
Odhad kumulatívneho rizika je definovaný na základe agregovaného procesu Y (t) = Y,j_ (t). Ň (t) = Y,i Ní (t), dŇ (t) = AŇ (t) = N (t) - Ň (t~), kde N (t) je suma udalostí do času t vrátane,
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
27
Y (t) je počet jednotiek v riziku v čase t (formálne ide o počet jednotiek v riziku v časovom intervale (t — e, t) pre malé e).
Príklad 66 Majme náhodný vektor (Xj, 5J), definovaný nasledovne (pre nejakú fiktívnu i-tu štatistickú jednotku, t.j. subjekt)
1) (Xí,5í) = (3,0), t.j. v čase X{ = 3 je cenzúra, Ni(t) = N (3) = 0, y (3) = y (3) = 1 ^(iVi(3),yi(3)) = (0,l),
2) (Xí,5í) = (4,1), t.j. v čase Xt = 4 je udalosí (zlyhanie), N (A) = 1, Y* (4) = 1, t.j. (Wť(4),yť(4)) = (l,l),
3) Ak máme viac udalostí: (N (0.5), Y,t (0.5)) = (1,1), (Ni (2) , Y,t (2)) = (2,1).
Klasický Kaplan-Meierov (KM) odhad funkcie prežívania je definovaný nasledovne
SkmÍj) = Yl [i-AAkmÍU
i:ti<t
kde AArm (ti
je odhad rizika a KM odhad kumulatívneho rizika Arm (t) = — ln S km (t)
y (u)
Ei:ti<tln (X - A^Km(U
Breslowov (B) odhad funkcie prežívania definujeme nasledovne
SB (t) = exp Í-ANA (í)) = * [ e
i:ti<t
kde Ana (t) = Z/ŕt <t^i = Eŕt <t    3e Nelson-Aalenov (NA) odhad kumulatívneho rizika. Greenwoodov (G) odhad rozptylu kumulatívneho rizika definujeme nasledovne
dl(t) = Var g \AKM(t) \ = Y,
AN (s)
^tY (s) [F (s) - AN (s
plug-in odhad rozptylu kumulatívneho rizika (nazýva sa aj Kleinov odhad; maximálne vierohodný odhad)
AŇ(s) (Y (s) - AN (s))
dl(t) = Var p A (t) \ = Y
i:ti<t
Y3 (s)
binomický odhad rozptylu kumulatívneho rizika
~2 m    v v- AiV(g) (Y(s)-AŇ(s))
aB(t) = Var p A (t)  = > -_2---,
L      J     t?<t      Y (s)(Y(s)-l)
NA odhad rozptylu kumulatívneho rizika (nazýva sa aj Aalenov, Poissonov alebo Tsiatisov odhad)
AŇ(s)
d2T(t) = Var ANA (t) \ = Y
tít Y (s)
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
6   Softvérová implementácia charakteristík prežívania
Q$ knižnica library(survival) ponúka nasledovné možnosti. Označme
surv.obj <- survfit(Surv(cas, status)~ 1,  type="...", error="...",conf. type
= "...")). '
Potom
1. odhady funkcie prežívania vypočítame ako
(a) S km (t)-' type="kaplan-meier" (přednastavené);
(b) S b (t): type="fleming-harrington";
(c) SFHmodB(t): type="fh2";
2. odhady rozptylu ako
(a) Var g S km (t) ■ error="greenwood" (přednastavené);
(b) Var g S b (t)
: error="tsiatis";
3. typy intervalov spoľahlivosti (IS) ako
(a) žiadny: conf. type="none";
(b) škála funkcie prežívania: conf. type="plain";
(c) škála kumulatívneho rizika: conf. type ="log" (přednastavené);
(d) škála logaritmu kumulatívneho rizika: conf. type="log-log".
Ďalšími argumentami sú
4- koeficient spoľahlivosti conf. int=0. 95 (přednastavené); 5. úprava spodnej hranice IS, kde argument
(a) conf. lower="usual" (nemodifikovaná dolná hranica IS);
(b) conf. lower="peto" (používa Petov efektívny rozsah súboru n**(t));
(c) conf. lower= "modified" (používa Dorey-Korn modifikáciu).
Priemerný vek prežívania a jeho smerodajná odchýlka ako aj medián a jeho smerodajná odchýlka sa vypočítajú ako
1. print(surv.obj,print.rmean=TRUE) alebo
2. print(surv.obj, rmean="individual").
Na rozlíšenie typu cenzurovania je dôležitú počet argumentov funkcie Surv(). Ak sú dva, t.j. Surv (cas, status), ide o pravý typ cenzurovania. Ak sú tri, t.j. Surv (cas, cas 1, status), potom ide o intervalové cenzurovanie. Pomocným argumentom je type=". . . ", kde rozlišujeme
1. type="right" (pravý typ), type="interval" (intervalový typ cenzurovania I. typu, kde interval (—oo, ti) označujeme (NA, ti));
2. type="interval2" (intervalový typ cenzurovania II. typu; kde interval je typu (tu,t2^ alebo
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
29
3. interval (ti, oo), ktorý označujeme (ti, NA)). Dolnou hranicou intervalu môže byť aj 0 a hornou hranicou ímax. Výstupmi objektu surv.obj a summary(surv.obj) sú nasledovné
1. rozsah - n,
2. počet jedincov v riziku v jednotlivých časoch - n. risk,
3. počet udalostí (zlyhaní) v jednotlivých časoch - n. event, 4- počet cenzur v jednotlivých časoch - n. censor,
5. odhad funkcie prežívania v jednotlivých časoch - surv,
6. odhad odmocniny z rozptylu funkcie prežívania v jednotlivých časoch - std. err,
7. dolná hranica IS pre funkciu prežívania v jednotlivých časoch - lower,
8. horná hranica IS pre funkciu prežívania v jednotlivých časoch - upper,
9. časy, v ktorých nastalo
• zlyhanie alebo cenzúra - time (pre surv.obj).
• zlyhanie - time (pre summary(surv.obj)).
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
30
7   Struktura funkcií na cvičenia a do DU
Sruktúra funkcie pre odhad kumulatívneho rizika, jeho rozptylu, intervaly spoľahlivosti a pás spoľahlivosti pre kumulatívne riziko a graf, tj. funkcie, ktorá bude „všetko " počítat a zobrazovat v jednotlivých časoch zlyhania to, t\,..., max(£max, cmax); je nasledovná (číselný výstup je vo forme tabuľky TAB. kum .riz. IS, pozri nižšie)
kumulativně.riziko.all(cas, status, type="KM", error="G",
conf.coef=0.95,conf.type="plain", plot=FALSE,plot.conf.int=FALSE,plot.conf.band=FALSE),
kde sú možné nasledovné
1. kombinácie odhadov
(a) type="KM" pre A^-^f(í) a error="G" a2G(t); přednastavené;
(b) type="NA" pre K^^if) a error="T" a2T(t); namiesto v\{t) môžeme použii aj v%{f) alebo o2K(t), potom error="B", resp. error="K"');
(c) type="FHmodNA"prekFHmodNA{t) a error="FH" d2FH{t);
2. typy IS
(a) škála kumulatívneho rizika: conf .type ="plain" (přednastavené);
(b) škála logaritmu kumulatívneho rizika: conf. type="log";
(c) škála arcus-sínusov á odmocninová: conf .typ e=" ar c-sin".
3. votby kreslenia obrázkov (číselné výpočty sú vždy výstupom funkcie)
(a) plot=FALSE, plot. conf. int=FALSE, plot. conf. band=FALSE - nekreslíme žiaden obrázok;
(b) plot=TRUE, plot. conf .int=FALŠE, plot. conf .band=FALSE - kreslíme obrázok len s odhadom kumulatívneho rizika;
(c) plot=TRUE, plot. conf. int=TRUE, plot. conf. band=TRUE - kreslíme obrázok s odhadom kumulatívneho rizika, a intervalmi a pásmi spoľahlivosti pre kumulatívne riziko.
Štruktúra funkcie pre funkciu prežívania, jej rozptylu, intervaly spoľahlivosti a pás spoľahlivosti pre funkciu prežívania a graf, tj. funkcie, ktorá bude „všetko " počítať a zobrazovať v jednotlivých časoch zlyhania t0, t\,..., max(tmax, cmax); je nasledovná (číselný výstup je vo forme tabuľky TAB. fcia.prez. IS, pozri nižšie)
funkcia.prezivani a.all(cas,status,type="KM",error="G",
conf.coef=0.95,conf.type="plain", plot=FALSE,plot.conf.int=FALSE,plot.conf.band=FALSE),
kde sú možné nasledovné 1. kombinácie odhadov
S km {t) ,
SB(t) ,
(a) type="KM" pre SKM(t) a error="GKM" VarG
(b) type="B" pre SB(t) a error="NAB" VarNA
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
31
(c) type="B" pre SB(t) a error="AJB" VarAJ SB(t)
(d) type="FHmodB" pre SFHm0dB(t) a error="NAFH" VarNA SFHmodB(t)
(e) type="FHmodB" pre SFHmodB{t) a error="GFH" VarG SFHmodB{t)
2. typy IS
(a) škála funkcie prežívania: conf. type ="plain" (přednastavené);
(b) škála logaritmu funkcie prežívania: conf. typ e ="log";
(c) škála logaritmu kumulatívneho rizika: conf. type="log-log";
(d) škála arcus-sínusov á odmocninová: conf. type="arc-sin".
3. voľby kreslenia obrázkov (číselné výpočty sú vždy výstupom funkcie)
(a) plot=FALSE, plot. conf. int=FALSE, plot. conf. band=FALSE - nekreslíme žiaden obrázok;
(b) plot=TRUE, plot. conf. int=FALSE, plot. conf .band=FALSE - kreslíme obrázok len s odhadom kumulatívneho rizika;
(c) plot=TRUE, plot. conf. int=TRUE, plot. conf. band=TRUE - kreslíme obrázok s odhadom kumulatívneho rizika, intervalmi a pásmi spoľahlivosti pre kumulatívne riziko.
Vyššie uvedené funkcie budú načítavať vo svojom tele ďalšie funkcie
1. početnosti (time, status) s výstupom v podobe tabulky obsahujúcej časy do zlyhania ti (ozn. ti), počty jedincov v riziku n,i tesne pred časom ti (ozn. ni) a počty zlyhaní dí v čase ti (ozn. di);
2. kumulativně, riziko (ti,ni, di, type="KM", error="G"), ktorej výstupom je tabuľka
TAB. kum. riz obsahujúca stĺpce ti, di, ni, Lambda a Var [Lambda] (nastavenia argumentov type a error kopírujú nastavenia funkcie kumulativně, riziko. all ());
3. kumulativně, riziko .IS(TAB. kum. riz, conf .type ="plain"), ktorej výstupom je tabuľka TAB.kum. riz. IS obsahujúca stĺpce ti, di, ni, Lambda, Var [Lambda], DH.IS a HH.IS (nastavenie argumentov conf .type a conf. coef kopírujú nastavenie funkcie kumulativně, riziko .all () );
4- funkcia.prezivania(ti,ni, di, type="KM", error="G"), ktorej výstupom je tabuľka
TAB. fcia.prez obsahujúca stĺpce ti, di, ni, Sa Var[S] (nastavenia argumentov type a error kopírujú nastavenia funkcie funkcia.prezivania. all ());
5. funkcia.prezivania. IS(TAB. fcia.prez, conf. coef=0. 95, conf. type="plain "), ktorej výstupom je tabuľka TAB. fcia.prez. IS obsahujúca stĺpce ti, di, ni, S, Var [S], DH.IS a HH.IS (nastavenia argumentu conf. type a conf. coef kopírujú nastavenie funkcie funkcia, prezivania. all ()).
Ďalšou funkciou bude funkcia
zakl. char.prezivania(TAB.fcia.prez, conf. coef=0.95,t,epsilon=0.05, plot.mu=FALSE, plot.mrl=FALSE),
(10. decembra 2014)
Katina, S., 2014: Analýza prežívania
32
ktorá bude načítavať vo svojom tele ďalšie funkcie (každá z nich bude mať výstup v podobe vektora obsahujúceho parameter, Var [parameter], DH. IS a HH. IS a celkový výstup je tabuľka so štyrmi riadkami (ozn. ako priemer, medián, priemerný.zostatkový.život a medián.zostatkové- ho.života)
a štyrmi stĺpcami; vstupom je okrem iného samotný odhad funkcie prežívania, preto nie je nutné mať ako vstup typ odhadu a jeho rozptyl)
1. priemerný. cas.prezivania(TAB. fcia.prez, conf. coef=0. 95,plot.mu=FALSE), ktorej výstupom bude vektor obsahujúci parameter fi, Var [parameter] Var[fľ\, DH.IS dolná a HH.IS horná hranica Waldovho empirického 100 x (1 — a)% IS pre strednú hodnotu času prežívania (argument plo t. mu znamená, že funkcia nakreslí obsah pod krivkou prežívania pomocou funkcie polygonO);
2. medián, casu. prezivania (TAB. fcia.prez, conf. coef=0. 95,plot.mrl=FALSE), ktorej
výstupom bude vektor obsahujúci parameter fi, Var [parameter] Var\pľ\, DH. IS dolná a HH. IS horná hranica Waldovho empirického 100 x (1 — a)% IS pre medián času prežívania (argument plot. mri znamená, že funkcia nakreslí za bodom t obsah pod krivkou prežívania pomocou funkcie polygonO );
3. priemerný. zostatkový. život (TAB. fcia.prez, t, conf. coef=0. 95), ktorej vstupom je okrem iného čas t a výstupom vektor obsahujúci parameter mrl(t), Var [parameter] Var[mrl(t)], DH. IS dolná a HH. IS horná hranica Waldovho empirického 100 x (1 — a)% IS pre strednú hodnotu zostatkového života;
4- medián, zostatkového. zivota(TAB. fcia.prez, t, conf. coef=0. 95), ktorej vstupom je okrem iného čas t a výstupom vektor obsahujúci parameter mrl(t), Var [parameter] Var[mrl(t)], DH.IS dolná a HH.IS horná hranica Waldovho empirického 100 x (1 — a)% IS pre medián zostatkového života.
(10. decembra 2014)