Vybrané kapitoly z analýzy prežívania Testovanie hypotéz Stanislav Katina 1 Ústav matematiky a statistiky Prírodovedecká fakulta Masarykova univerzita v Brne ZS 2013 f evropský 1 sociální MINISTERSTVO ŠKOLSTVÍ, MLÁDEŽE A TĚLOVÝCHOVY fond V ČR EVROPSKÁ UNIE INVESTICE DO ROZVOJE VZDĚLÁVÁNI Stanislav Katina IT*1* OP Vzdělávání pro konkurenceschopnost Testy na porovnanie kriviek prežívania Prehfad testov Delenie testov na porovnanie kriviek prežívania • typ dát: necenzurované alebo cenzurované • počet porovnávaných kriviek prežívania: k = 2 alebo k>3 • typ pozorovaní: nepárové alebo párové • typ alternatívy: všeobecná alebo zoradená • crossing efekt (prekrižovanie sa kriviek prežívania): neexistuje alebo existuje • časový efekt liečby: neexistuje alebo existuje Stanislav Katina Poďakovanie Tento učební text vznikl za přispění Evropského sociálního fondu a státního rozpočtu ČR prostřednictvím Operačního programu Vzdělávání pro konkurenceschopnost v rámci projektu Univerzitní výuka matematiky v měnícím se světě (CZ. 1.07/2.2.00/15.0203) Stanislav Katina Testy na porovnanie kriviek prežívania Prehfad testov Neparametrické testy porovnania kriviek prežívania pre necenzurované dáta m testy porovnania dvoch kriviek prežívania (k = 2) • Wilcoxonov test (W) • Mann-Whitney test (MW) • Siegel-Tukey test (ST) • testy porovnania viac kriviek prežívania (k > 3) • Kruskal-Wallis test (KW) • Jonckheere test (J) • Cuzick test (C) • Letest(L) Stanislav Katina Testy na porovnanie kriviek prežívania Prehfad testov Testy na porovnanie kriviek prežívania Prehfad testov Neparametrické testy porovnania kriviek prežívania pre cenzurované dáta • testy porovnania dvoch kriviek prežívania (k = 2) • Gehan-Wilcoxon test, zovšeobecnený Wilcoxonov test (GB) • Cox-Mantel test, log-rank test (CM) • Tarone-Ware test (TW) • Peto-Peto test (PP) • testy porovnania viac kriviek prežívania (k > 3) • Gehan-Breslow test, zovšeobecnený Wilcoxonov test, zovšeobecnený Kruskal-Wallis test (GB) • Cox-Mantel test, log-rank test (CM) • Mantel-Haenszel test, log-rank test (MH) • Peto-Peto test (PP) Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Wilcoxonov test Predpoklady • Xi, ...,Xni je náhodný výber (NV) z nejakého spojitého rozdelenia • Y-\,Yn2 je NV z rovnakého spojitého rozdelenia a je oproti prvému rozdeleniu posunuté o nejakú konštantu S • veličiny Xi, ...,Xni a - S,Y„2 - S majú rovnaké rozdelenie • oba výbery sú nezávislé Hypotézy m Ho:5 = 0(S1(ř) = S2(ř),Vř) • Hi : ô ^ 0 (Si (ŕ) ŕ S2{t), pre aspoň jedno t) Stanislav Katina S2 (r) = S (r) Testované hypotézy m nulovná hypotéza H0 : Si (t) m alternatíva hypotéza H-\ : • Si (t) ŕ S2 (t) • Si (t) * S2 (t) • S-\ (t) y S2 (t) pre Vř S (t) je funkcia prežívania st st - 10 Sw - E0 [Sw] Zw = -/,/ ro n ~ N (" 1 J Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Mann-Whitney test Predpoklady • Xi, ...,Xni je NV z nejakého spojitého rozdelenia • Y-\,yn2 je NV z rovnakého spojitého rozdelenia a je oproti prvému rozdeleniu posunuté o nejakú konštantu S • oba výbery sú nezávislé • nech (X/, Y}) sú možné páry pozorovaní, pre ktoré môže nastat bud X, < Y} alebo X, > Y} Hypoŕézy • Ho:5 = 0(S1(ř) = S2(ř),Vř) sŕ • H-i : 5 > 0 (Si (ŕ) -< S2(ŕ), pre aspoň jedno ŕ) Stanislav Katina Stredná hodnota a rozptyl Sw Var0[Sw\t] Wilcoxonov test Ak A7i, n2 > 10 /71/72(/7+ 1) 12 Hi 1) 1 - , 21 ,x^f/fŕ/-1 n (n2 - 1) ^ 7 v-' ■M/,/ yy Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Mann-Whitney test Testy na porovnanie dvoch kriviek prežívania Mann-Whitney vs Wilcoxonov test Stredná hodnota a rozptyl SMW : Eo [Smw] : Var0 [SMW] = Var0 [Sw] Mann-Whitney test Ak n-!, n2 > 10 n^n2 n^n2 +a?2 + 1) 12 y Smw - Eq [Smw] », fn *^ zmw = —. ro--N (U, 1) V vdľQ [0mw\ Stanislav Katina Asymptotická normalita rozdelenia n2) nl <- 7; n2 <- 5; n <- nl+n2 SUMmax <- sum(6:12); SUMmin <- sum(l:5; Smw.max <- nl*n2-SUMmin+n2*(n2+l)/2 Smw.min <- nl*n2-SUMmax+nl*(nl+1)/2 x <- Smw.min:Smw.max fx <- dwilcox(x, nl, n2); Fx <- pwilcox(x, nl, E.Smw <- nl*n2/2; Var.Smw <- nl*n2*(n+1)/12 fx.norm <- dnorm(x,E.Smw,sqrt(Var.Smw)) Fx.norm <- cumsum(dnorm(x,E.Smw,sqrt(Var.Smw)))/sum(dnorm(x,E.Smw,sqrt(Var par(mfcol=c(1,3) ) plot (x, fx,type="h", go1="black", sub="nl = 7,n2 = 5", xlab="x", ylab="f(x)",main= "Hustota MW statistiky",ylim=range(c(fx.norm,fx) points(x, fx,pch=16,cex=0.1) lines(x,fx.norm,go1="red", plot(x. Fx,type="s" ylab="F(x)" ,lwd=5) , Gol="blaGk", sub="nl=7,n2=5", xlab="x", ,main= "Distribučná funkcia MW statistiky") points(x, Fx,pch=16,cex=0.7) lines(x,Fx.norm,co1="red",lwd=5) abline(h=0:1, col="gray20",lty=2) p.mw <- seq(0.001,0.999,length=length(x)) q.mw <- qnorm(p.mw,E.Smw,sqrt(Var.Smw)) plot(x,q.mw,main = "qq-diagram MW rozdelenia", ylab="kvantily", sub="nl=7,n2=5" pch=16,cex=0.7,asp=l) xl <- quantile(x,0.25); x2 <- quantile(x,0.75) yl <- quantile(q.mw,0.25); y2 <- quantile(q.mw,0.75) bl <- (y2 - yl)/(x2 - xl); b0 <- yl - bl * xl abline(a=b0,b=bl,lwd=2) Stanislav Katina Ux vyjadruje počet dvojíc X,, Vý, kde platí X, < Yý Uy vyjadruje počet dvojíc X, V}, kde platí X, > Vý Uy = n,n2+n^ni+')-WY Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Mann-Whitney vs Wilcoxonov test Hustota MW statistiky Distribučná funkcia MW statistiky qq-diagram MW rozdelenia 5 10 15 20 25 30 3í 5 10 15 20 25 30 3í 20 30 4C Obr.: Rozdelenie Mann-Whitneyho štatistiky SMw (^1 =7,n2 = 5) Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Siegel-Tukey a Levene test Testy na porovnanie dvoch kriviek prežívania Siegel-Tukey test Predpoklady m efekt nového typu liečby spôsobí rast alebo klesanie krivky prežívania oproti pôvodnému typu liečby • liečba nespôsobuje zmenu priemernej odpovede, ale výsledná odpoved • Xi, ...,Xni je náhodný výber (NV) z nejakého spojitého rozdelenia • Yi,Yn2 je náhodný výber (NV) z nejakého spojitého rozdelenia • oba výbery sú nezávislé Hypotézy • HQ : Var{S-\{t)) = Var{S2{t)), Vŕ • Hi : Var(Si(ŕ)) ŕ Var(S2{t)), pre aspoň jedno t Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Siegel-Tukey test Stredná hodnota a rozptyl Sst (resp. SfT): Eo [Sst] VarQ [SST] n 10 zst= sst-Eo [Sst]^n{q 1} Stanislav Katina Podstata Siegel-Tukey testu je nasledovná: m poradie priradíme najmenšiemu pozorovaniu • poradie R2 priradíme najväčšiemu pozorovaniu • poradie Rz priradíme druhému najmenšiemu pozorovaniu • poradie f?4 priradíme druhému najväčšiemu pozorovaniu atd. Siegel-Tukey štatistika "i /=1 kde daná suma prislúcha súčtu poradí pre prvý NV Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Levene test Podstata Leveneho alternatívy ST testu (Levene testu) je nasledovná: m odchýlky Dx = \X - /jx\ a Dy = IY - hy\ m odchýlky {cŕ, = |x, -x|}^ a {d; = \yj-y\}ji, • cř(1) < cŕ(2) < ... < cf(n), n = ni+n2 Levene štatistika S ĺ - S|'f - y^Rdiff,(j), /=1 kde ŕ?c//ff)(/) predstavujú poradia odchýlok od priemeru pre prvý NV Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Leven e test Testy na porovnanie dvoch kriviek prežívania Príklad Stredná hodnota a rozptyl SfT Eo [Sst] = £o Var0 [SST] = Var0 q alt °st oa/ŕ n^n2 {n<[+n2 + /{) 12 Levene test Ak n-,, n2 > 10 Z, =Z alt sl SfT - E0 [SfT] Var? [SST] N(0,1) Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Príklad 52 240 19 53 15 43 340 133 111 231 378 49 4>n 0 1 1 0 0 1 1 0 0 1 0 0 '(') 15 19 43 49 52 53 111 133 231 240 340 378 4>n 0 1 1 0 0 0 0 0 1 1 1 0 r1' 1 - - 4 5 6 7 8 - - - 12 - 2 3 - - - - - 9 10 11 - Wilcoxon test 35 /=1 Eo [SW] = 5-^± Var0 [Sw] = 7x5(^+5+1) Zw = (35-32.5)/6.157651 0.405999, p-hodnota=0.6847 Stanislav Katina Example Nech ŕ/,, / = 1,...,A7y,y = 1,2 sú časy do zlyhania (úmrtia) od diagnostiky nádoru pľúc v mesiacoch, kde j = 1 predstavuje I. typ terapie aj = 2 zasa II. typ terapie, nech ^ = 0 ak y = 1 a = 1 aky = 2 (pozri tab.). Otestujete H0 ■. (t) = S2 {t), alternatíva H, : (t) ^ S2 (t). Použite aj SMW, SST, SfT. Vždy presne naformulujte H^. 52 240 19 53 15 43 340 133 111 231 378 49 i'ij 0 1 1 0 0 1 1 0 0 1 0 0 Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Príklad Mann-Whitney test SMW = mn2-Y ŕ?í2) + "1 ("1 + 1 )/2 = 20 E0 [Smi/i/] = ^2/2 = 17.5, \/a/b[SML/i/] = n^n2(n + ^)/^2 = 37.91667 ZST = (2° - 17.5)/37.91667 = 0.405999, p-hodnota=0.6847 Siegel-Tukey test 52 240 19 53 15 43 340 133 111 231 378 49 0 1 1 0 0 1 1 0 0 1 0 0 r1' 9 - - 11 1 - - 10 12 - 2 7 r<2> - 6 3 - - 5 4 - - 8 - - (2) e0[Ssr] = Zsr = (26 '=1 5(7 26 -V,Var0[Sw]=7-^±ll 32.5)/6.157651 = -3.166792, p-hodnota=0.2912 Stanislav Katina Testy na porovnanie kriviek prežívania Prehfad testov Testy na porovnanie viacerých kriviek prežívania Kruskall-Wallis test Testované hypotézy • nulovná hypotéza H0 ■. S, (t) = S, {t) = S (t) • alternatíva hypotéza H-i : • Sj {t) ^ Sj (t) pre aspoň jedno i, j (KW test) • S,(ŕ)* Sy-(ŕ)(J,C,Ltesty) • S/ÍOČSyíOíJ.C.Ltesty) pre Vr,/' xs A y > y • diskordantná (neusporiadaná) pokiaf platí x, < xy A y > y alebo x,- < xs A y > y • ak platí x, = xj alebo y = y, nejde ani o konkordantný ani o diskordantný vzfah • C + D < n (n - 1), C je počet konkordantných dvojíc, D počet diskonkordantných dvojíc medzi (x-i,yi),(xn,yn) Stanislav Katina Alternatívna podoba Jonckheere štatistiky k-\ k k-\ k /=1 y'=1+/ /=1 y=1+/ kde E?=ľ ĽJLi+/"/"y = maxvs Sj Stredná hodnota a rozptyl Sj _ a72(2a7 + 3)-E"/2(2a7/ + 3) E[Sj] = 0, Var0 Alternatívna podoba Jonckheere testu 18 Sj-E0 Sj ^JvarQ Sj N(0,1) Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Zovšeobecnený Kendallov korelačný koeficient Kendallov korelačný koeficient r = 77(7^1) Ľľ=i E"=i s/gn (x, - x;) s/gn (y - y), čo je identické s f = 7í(7Pi) Eľ=ľ E"=/+i s/gn (x, - x;) sign (y - y), kde í 1 sign(xi-Xj) = Xy 1 ak Xj < Xj { 0 ak x, = Xy f 1 ak y > y sign{y, - y j) = = -0.5/VÖTÖ32 p-hodnota=0.002 0.032 -2.801 Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Cuzick test Testy na porovnanie viacerých kriviek prežívania Cuzick test Označenia • nech Rjj je poradie /y-teho pozorovania v združenom NV • nech Zy je s/fóre prislúchajúce NV, do ktorého /y-te pozorovanie patrí • n = EJLi nj Cuzickova štatistika n, k SC = ^2^2zijRij; /=1 ;'=1 Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Zovšeobecnený Spearmanov korelačný koeficient Označenia • nech (Xi, Vi )T,(Xn, Yn)T je výber z dvojrozmerného rozdelenia • nech R<\,..., Rn sú poradia veličín X1? ...,X„ • nech Q-i,Qn sú poradia veličín ,Yn Spearmanova štatistika S/v — ^ RjQj /=1 Stanislav Katina Stredná hodnota Eo[Sc]=(J2i^E[Z] = ln{n + '\)E[Z], kde £ [Z] = y!j=<\ zjPj, k je počet skupín, z/y = z,- =y, p, = n;/n Rozpŕy/ Var[Sc] nz (A7 + 1) 12 kde Var [Z] = z2p; - (E [Z])2 Cuzickov test 2 _ Sc - E0 [Sc] °~ ^VarQ [Sc] Var [Z], N (0,1) Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Zovšeobecnený Spearmanov korelačný koeficient Stredná hodnota E[SN] = n n + 1 Rozptyl Var[SN] 1 /A7(A72-1) n-1 V 12 Spearmanov test kde Rs y/{n-A)Vár[SN] Var [Rs] = ^ VřT^\Rs ~ W(0,1), (S/v-E[S/v]),E[f?s] = 01 Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Zovšeobecnený Spearmanov korelačný koeficient Testy na porovnanie viacerých kriviek prežívania Le test Vztah Spearmanovho a Pearsonovho korelačného koeficientu Ak (Xi, Yi)T,(Xn, Yn)T je výber z dvojrozmerného normálneho rozdelenia s korelačným koeficientom pXjy, potom platí 7T pX)Y = 2sin (^-) Rs Vztah Spearmanovho korelačného koeficientu a Cuzickovej štatistiky Cuzickova štatistika Sc je rovná Spearmanovej štatistike SN, kde jedna premenná predstavuje zoradenú (ordinálnu) premennú a druhá spojitú premennú Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Le test Stredná hodnota E[SL} = 0 Rozptyl Le test Var[SL] = n-^l^^-^ŕ /=1 Zl= \~£°iff;1 - a/(0,1) VVar0 [SL] Stanislav Katina Označenia • rij je rozsah /-teho NV • Lj = J2j<í nj = # pozorovaní vo všetkých výberoch naíavo od /-teho výberu, L-i = 0 • Mj = J2j>í nj = # pozorovaní vo všetkých výberoch napravo od /-teho výberu, Mk = 0 • Lj - Mj g {-n, n) m Rj je priemerné poradie pre /-ty výber Le štatistika SL = y£nj(Lj-Mj)Rj /=1 Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Zovšeobecnenia Test odklonu od trendu Skw SI Var [SL] Xk-2 Všeobecný tvar testovacej štatistiky SA = Y,n's'r'- /=1 kde rij sú rozsahy jednotlivých NV, s, skóre prislúchajúce jednotlivým NV a r, sú priemerné charakeristiky polohy prislúchajúce jednotlivým NV Potom (SA-E0[SA]y Var0 [SA] Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Príklad Testy na porovnanie viacerých kriviek prežívania Príklad Example (pokrač.) Nech tjj,i = '\,...,r)j,j = 1,2 sú časy do zlyhania (úmrtia) od diagnostiky nádoru pľúc v mesiacoch, kde j = 1 predstavuje I. typ terapie aj = 2 zasa II. typ terapie, nech ^ = 0 ak y = 1 a ty = 1 ak j = 2 (pozri tab.). Otestujte H0 ■. (t) = S2 (t), alternatíva ■. (t) ^ S2 (t). Použite aj Skw, Sj,SC- Vždy presne naformulujte H^. J 52 240 19 53 15 43 340 133 111 231 378 49 0 1 1 0 0 1 1 0 0 1 0 0 SKW = sw SMW = sc Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Príklad Rozdeľme AG-pozitívnych pacientov do troch skupín nasledovne nasledovne 1: WBC > 100000, n-, = 3, (1, 1, 65) • Skupina i. nuw / luuuuu, n-\ - o, v i, i, v^>) m 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) sk.1 1 1 sk.2 4 5 22 26 sk.3 16 39 poradie 1.5 1.5 3 4 5 6 7 8 sk.1 65 sk.2 108 121 sk.3 56 65 100 134 143 156 poradie 9 10.5 12 13 14 15 16 17 Stanislav Katina Example (pokrač.) Majme pacientov s akútnou myeloidnou leukémiou (AML) 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žnost choroby. Nech ti, i = 1,2,17 sú časy do zlyhania v týždňoch (pozri tab.). Otestujte H0 : Si (t) = S2 (t) = S3 (ŕ), alternatíva a) H, : S, (t) ŕ Sj (t), / < j; i J = 1,2,3, b) S, (t) ? S2 (t) ? S3 (t), c) S, (ŕ) ? S2 (ŕ) 5 S3 (0- Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Príklad Ri = 4.50, R2 = 7.833, R3 = 12 11.5625 M) = 6 x 7.833^ + 8 x 11.56252 0.0924 M,)R,= 3 x 18 17x18 L3 X 4-52 4.762662, p-hodnota= Sl = Ef=i n, (L; 3x (0-14) x4.5+6x (3-8) x7.833+8x (9-0) x 11.5625 Var [Sl] = tel) £*=1 ni (L/ _ Mi)2 = ^l^[3x(0-14)2 + 6x(3-8)2 +8 x (9-O)2] = 187.99732 ZL = 408.5/187.9973 = 2.172903, p-hodnota=0.0298 Skw - (Íl)2 = 4.762662 - 2.1729032 = 0.0412, p-hodnota=0.8392 Stanislav Katina 408.5 Testy na porovnanie dvoch kriviek prežívania Cenzurované dáta-prehíad Testy na porovnanie dvoch kriviek prežívania Cenzurované dáta-prehíad Formulácie sčítacím procesom • (X,, Si) nahradíme (A/, (t), Y, {t)), kde N, (ŕ) je počet pozorovaných udalostí v intervale (0,ŕ) v jednotke /', Y, (t) 1 jednotka / je v riziku v čase t 0 inak • V/(ŕ) = / ({T/ > ŕ}) a N, (ŕ) = /({!, <_ŕ, 5 (x) = 1}) • agregovaný_proces_Y (ŕ) =_£, Y, (t), A/_(ŕ) = E/ N, (ŕ), c/A/ (ŕ) = A/V (ŕ) = N(t)- NJt~), kde A/ (ŕ) je suma udalostí do času t vrátane, Y (t) je počet jednotiek v riziku v čase t Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Cenzurované dáta-prehíad Potom asymptoticky T(W t)-E0[T(W,t)] ^Var0[T(W,t)] A/(0,1) Stanislav Katina Testovacia štatistika na porovnanie dvoch kriviek prežívania i {w,ľ) - j0 vv[S)- (s)+- (s) ^ - (s) - (s) j os, kde W(s) je funkcia váh a s > 0 • stredná hodnota E0[T(W,t)] = 0 • rozptyl Var0[T(W,t)] ' rW2(s) ÍJ,(s)Y2_(s) ' f-Jo y,(s) \Y,(s) + Y2(s) í AŇ1(s) + AA72(s)-l\ ^ Y,{s) + Y,{s)-^ J f AÄ/-I (s) +AW2(S)N Yi(s) +Vi(s) OřS Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Cenzurované dáta-prehíad Tarone a Ware (1977) trieda váh O konštantný rozdiel v logitovej škále f (p) = In potom váhy W(t) = 1 O konštantný rozdiel v aritmetickej škále: f (p) = p, potom váhy W(t) = (^z^))-1 = Y\t)/(Y(t) - 1) « 7(ř) O konštantný rozdiel v arcsin škále: f (p) = arcsin ,/p, potom sú váhy rovné W{t) = « y/í7(ŕ) kde př = 1/7(ŕ) a Y(ŕ) počet osôb v riziku v združenom výbere v čase t Vo všeobecnosti môžeme váhy zapísat ako W(t) = g(Y(t)/n) Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Cenzurované dáta-prehíad Testy na porovnanie dvoch kriviek prežívania Cenzurované dáta-prehíad Harrington a Fleming (1982) trieda váh W(t) = Sp(t),p> 0 O p = 0, a teda W{t) = Š°(ŕ) = 1 (Cox-Mantel test alebo log-rank test; Cox, 1972; Mantel, 1966) O p = 1, a teda W{t) = Š1 {t) = S(ŕ) (Gehan-Wilcoxon test alebo Peto-Peto-Wilcoxon test, Gehan, 1965; Peto a Peto, 1972) Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Zovšeobecnený Wilcoxon test Označenia • majme dva náhodné výbery s rozsahmi n-\ a n2 m prvý výber má r-\ cenzúr a n-\ - ^ zlyhaní • druhý výber r2 cenzúr a n2 - r2 zlyhaní • realizácie v(c) xi\ ' i Xj2 'i • • • i Xjr^ *//}■+1 i xjrj+2i • • • ? xjnj J — 1 ? 2 cenzury zlyhania Testovacia štatistika = J2i i Uij> kde í í l +1, akx-i; < x2;- vx-i/ < x±j {o) akx-iž X2/Vxif V X akx-i, > x2j vx)/ > x2; í;' <*2, VX1;. >X^;. (c) (c) kde sumujeme cez všetkých nyn2 porovnávaní Stanislav Katina Formálna formulácia Testovacia štatistika na porovnanie dvoch kriviek prežívania T(WJ) = ^Wj(dv-djfy, potom • stredná hodnota E0[T(W,t)] = 0 • rozptyl Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Zovšeobecnený Wilcoxon test Ak nemáme prítomné cenzúry, sa redukuje na Wilcoxonovu štatistiku Sw, Mann-Whitneho štatistiku SMW a Kendallovo t nasledovne S$ = n2(n^ + n2 + 1) - 2SW = 2SMW - nyn2 = t, kde Sw je suma poradí druhého výberu v združenom výbere, SMW = #(x-i/,X2y), xy, < x2j a posledná rovnost platí aj v prítomnosti zhôd je možné zjednodušit ak všetky cenzúry nastali v čase tn (Halperin, 1960), teda S$ =2S^ + ryr2-nyn2, kde S(z) = sjjj, - ry[n2 - r2), je Mann-Whitney štatistika počítaná len na základe (a?-i - r{) + (n2 - r2) zlyhaní Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Zovšeobecnený Wilcoxon test Testy na porovnanie dvoch kriviek prežívania Zovšeobecnený Wilcoxon test Za platnosti H0 platí • stredná hodnota E0[S[v • rozptyl 1 = 0 Var0[S$] n^n2 ňíňZT) í £7=1 djDj.,(Dj_, + 1) + Ey=i MM + 1) E}=1 (dj(n - Dj - Lj_i)(n - 3D;_i - d, - Ls_, - 1)) }, kde m + n2 = n, D, = £J=1 cř;, D0 = 0, L, = £y=i /y, L0 c/y je počet necenzurovaných pozorovaní (zlyhaní) s poradím j (Rj) v združenom výbere zlyhaní, /, je počet cenzúr z intervalu (Rj, ) 0, Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Zovšeobecnený Wlcoxon test vs. kontingenčné tabufky Majme špeciálny prípad, kde • di = d, h = /, dj = lj = 0 pre / ^ 0 • pozorovania usporiadajme do kontingenčnej tabuľky (KT) 2 x 2 s fixovanými marginálnymi početnosíami ako výsledok náhodného výberu z hypergeometrického rozdelenia • následne sa testovacia štatistika redukuje na rozdiel súčinov diagonálnych a mimodiagonálnych prvkov, čo je v úzkom vzíahu ku pomeru šancí rozptyl Var0[S$] /c/n-! n2 "1+^2-1 Stanislav Katina Ak nemáme žiadne zhody, d\ = ... = dL = '\, ako aj žiadne cenzúry, A| =...=/,= 0 (L = n), potom 1 3" Ak nemáme zhody a všetky cenzúry sa objavili po L = (a?1 - ri) + (n2 - r2) zlyhaniach, teda oř! = ... = cřL = 1,/i = ... = /L_i = 0,/L = n + r2. \Zar0[S£}] = -A7iA72(A7+1) (c), _ n,n2(n-r, - r2) í 1 , m/J- ^T^í) ^Ci+^J+g (n - r-, - r2)2 - 1 Ak máme len zlyhania, ale vyskytujú sa zhody, dostaneme (Hemelhjk, 1952) Var0[S$] = VarQ[SMW|t] = VarQ [SMW] Stanislav Katina nin2Y 1 Stanislav Katina Pre každé ŕ(/), 1 < / < /, môžeme dáta zapísat do KT 2 x 2 výber/status zlyhanie v ŕ(/) nažive v spolu v ŕ(/) 1 d m 3m n m 2 d2i a2i n2i spolu di a, n. • n<\j = # subjektov v prvom NV, ktorí boli v riziku tesne pred časom ŕ(/), n2i = # subjektov v druhom NV, ktorí boli v riziku tesne pred časom ŕ(/), n, = a?i/ + n2, m dy, = # zlyhaní z prvého NV, d2i = # zlyhaní z druhého NV, d i = dy, + d2l m a, = n-, - d, = ay + a2/ = # subjektov, ktorí ostali nažive v čase ŕ(/) • # zlyhaní do času ŕ(/) vrátane d = Eyŕy<ŕ(í) ds Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Prehfad testov • marginálne početnosti v tab., nv, n2, a d, sú náhodné premenné závislé iba na minulosti pred časom ŕ(/) • Mantel a Haenzel (1959): rozdelenie pozorovaní (realizácií) v bunkách KT podmienené pozorovanými marginálnymi početnosťami (d/, a,-, ny,, n2j) za platnosti • to implikuje rozdelenie iba jednej bunky, dy, pretože ostatné početnosti sú íahko odvoditeľné od marginálnych • za platnosti nulovej hypotézy, H0, rozdelenie dy je hypergeometňcké, teda Pr(cŕi,) = • v tejto forme H0 o rovnosti kriviek prežívania implikuje nezávislosť výberu a statusu (nažive alebo zlyhanie) Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Prehfad testov Za platnosti H0 • očakávaná (stredná) hodnota £0 [du] = • rozptyl VarQ [dv] n i - riu nun^aid. J /j?(/í,-1) Informáciu o KT v čase ŕ(/) nám dá nasledovný vzíah Xi 2 _ [du - E0 [dv]Y Var0 [du] x2df, df = 1 Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Príklad Example Majme dáta z klinickej štúdie zhrnuté v nasledovnej tabuíke. ti ni di du E0[du] du-E0[du] Var0[du] 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 X 2 -2.312/1.0302 = 5.179674 p-hodnota=0.0229 z = V2.312/1.0302 = 2.275890 p-hodnota=2 x 0.0114 = 0.0229 Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Prehfad testov Pre všetky tabuíky (/' = 1,2,...,/) píšeme U = ^(dv-E0[du])., /=1 E0 [U] = 0, Var0 [U] = ]T Var0 [du] = ]T nun2iaidi /=1 /=1 ' v ' ' Ak máme fixované d-,, nr,, n2h potom platí (Mantel a Haenzel, 1959) 2 Qmh EL, (cŕií-Eofa,]) í;2 Var0 [cd,] Var0 [Ľ] x2df,df = 1, Q = -r i n- ~ X1 > Z = /,, ril1 ~ A/ (0, 1 ) Var0 [Ľ] Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Prehfad testov Majme váhy w-, asociované s KT v čase ŕ,, potom U = 'ž2wi(dM-E0[du]), /=1 E0 [U] = 0, Var0 [U] = ]T wfVarQ [du] = ]T %: 2 nun2iaidi /=1 Ak máme fixované c/,, n-i,-, n2h potom platí (Mantel a Haenzel, 1959) (U - Ep [L/])2 Var0 [Ľ] 2 z_ U-EoM A/(0,1) Stanislav Katina Testy na porovnanie dvoch kriviek prežívania Prehfad testov Testy na porovnanie dvoch kriviek prežívania Prehfad testov Podía výberu váh w-, rozoznávame nasledovné typy testov: • ak w, = n-,, Q = QG, ide o Gehan-Wilcoxon test (zovšeobecnený Wilcoxonov test; Gehan, 1965), ktorý môžeme zredukoval na Wilcoxonovu štatistiku pri absencii cenzúr • ak Wj> = 1, Q = QMH, ide o Cox-Mantel test (log-rank test; Mantel a Haenzel, 1959) • ak w, = y/rľj, Q = Qtw, ide o Tarone-Ware test (Tarone a Ware, 1977) • ak w, = Sp00ied (t,) = Uj-tj -, Q = Qp, ide o ;£fo'-Eofó]]> /=1 /=1 a V* počítame tak ako V, ale s tým rozdielom, že ide o maticu kxk Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Príklad Example Majme experiment (Thomas et all, 1977), kde sme mali 3 rôzne koncentrácie látky (konc^ = 2.0, konc2 = 1.5 a koncz = 0) a hľadali sme jej účinok na pacientov, u ktorých sme sledovali objavenie sa nádoru (pozri Tab.). Zaujíma nás testovanie A) HQ : A-i (t) = A2 (t) = A3 (t) oproti H| : 3 aspoň jedno / < 7, A, (t) ŕ Ay (t) B) H0 : Ai (t) = A2 (t) = A3 (t) oproti Ai (t) = w*\z (t), A2 (t) = w*\z (t), kde wf (test trendu). koncjj =1,2 Pozn.: časy do zlyhania alebo cenzúry; + znamená cenzúra, n0 = # pozorovaní v t0, koncj je koncentrácia látky v skupine j, konCj "0 2.0 10 41 + 41 + 47 47+ 47+ 58 58 58 100+ 117 1.5 10 43+ 44+ 45+ 67 68+ 136 136 150 150 150 0 9 73+ 74+ 75+ 76 76 76+ 99 166 246+ Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Test trendu Testovaciu štatistiku Qtrend budeme potom písaí nasledovne _(E;=iw;ELi[ďy/-Eo[ďj>]])2_ 'rencf a ak sú váhy rozdelené lineárne, w* = j, potom hovoríme o teste lineárneho trendu Platí Qresidual = Qoverall ~ Qtrend~Xk-2- čo vedie ku záveru, že ak je tento test štatisticky signifikantný, potom váhy použité v Qtrend treba nahradit inými, ktoré by lepšie modelovali trend v dátach Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Príklad 47 58 67 76 99 117 136 150 166 kone. 1 1 3 0 0 0 1 0 0 0 kone. 2 0 0 1 0 0 0 2 3 0 kone. 3 0 0 0 2 1 0 0 0 1 47 58 67 76 99 117 136 150 166 kone. 1 7 2 2 2 2 0 0 0 0 kone. 2 7 7 6 5 5 5 3 0 0 kone. 3 9 9 9 4 2 2 2 2 1 47 58 67 76 99 117 136 150 166 kone. 1 8 5 2 2 2 1 0 0 0 kone. 2 7 7 7 5 5 5 5 3 0 kone. 3 9 9 9 6 3 2 2 2 2 Stanislav Katina Testy na porovnanie viacerých kriviek prežívania Príklad Testy na porovnanie viacerých kriviek prežívania Príklad Prvý čas, kedy sa nádor objavil je ŕ(1) = 47 a KT v ŕ(1) je nasledovná dávka udalosí vŕ(1) bez udalosti vŕ(1) v riziku pred 2.0 1 7 8 1.5 0 7 7 0.0 0 9 9 spolu 1 23 24 Čiastkové výpočty v dávka £o[dyi] Ují 2.0 0.33333 1 - 0.33333 = 0.66667 1.5 0.29167 0-0.29167 = -0.29167 0.0 0.37500 Stanislav Katina r príkazy [r je voíne dostupné na http://cran.r-project.org/] Príklad 1 (poznámky sú uvedené za znakom #): library(survival) attach(ami) names(ami) ## odhady pre obe skupiny zvlast # median, priemerný cas prezivania, rozsahy, pocty zlyhani kmfit.aml <- survfit(Surv(cas, status) ^skupina, data=aml) print(kmfit.ami,show.rmean=TRUE) ## odhady v jednotlivých časoch # fcia prezivania,IS, počet v riziku, počet udalosti v case ## testovanie hypotéz # default rho = 0 [log-rank] test.aml<-survdiff(Surv(cas, status) ^skupina, data=aml) # Peto test rho = 1 test.aml<-survdiff(Surv(cas,status) ^skupina, data=aml,rho = Stanislav Katina d 8x(24-8)x1x(24-1] 242(24-i; V('(D 0.222 -0.097 8x7x1x23 242 x 23 7x(24-7)x1x23 24223 0.097 : 0,207 U V 3.209 -0.803 / 1.319 -0.641 \ i _ / 0.859 0.207 -0.641 2.663 J' "^ 0.207 0.425 Qoveraii = UTV"1U = 8.050, p-hodnota=0.018 Ak w* = (0.0,1.5,2.0)T, potom Ufrend = 36.186,wJV*w* = 5.991, QtrBnd = 6.04, p-hodnota=0.014, a naviac Qresiduai = Qovemii - Qtrend = 2.01, p-hodnota=0.156 Stanislav Katina r príkazy [r je voíne dostupné na http://cran.r-project.org/] Príklad 1 (pokrač.): skup <- as.numeric(aml$skupina) casl <- kmfit.ami.l$time; cas2 <- kmfit.ami.2$time; CAS <- c(casl,cas2) ## kumulativne riziko (KM a NA) # (pocitane rovnako ako v priklade 2, handout 1) LAM.KM <- c(Lambda.KM.1,Lambda.KM.2) par(mfcol=c(1,2)) # nastavenie poctu podokien ## obrázok (limitácia osi-argument xlim,ylim, kde beriem do ## úvahy rozsah združeného vektora kum.rizik) ## KM odhad kum.rizika (podobne aj NA odhad) plot(CAS,LAM.KM,xlab="cas do relapsu (v týždňoch)", ylab="kumulativne riziko", type="n", ylim=range(LAM.KM), xlim=c(0,50)) lines(casl,Lambda.KM.1,lty=l,lwd=2,type="s") lines(cas2,Lambda.KM.2,lty=2,lwd=2,type="s") abline(h=0,col="gray") title(main="KM kumulativne riziko pre AML data") legend("topleft",c("skup M","skup NM"),lty=c(1,2)) Stanislav Katina r príkazy [r je voíne dostupné na http://cran.r-project.org/] r príkazy [r je voíne dostupné na http://cran.r-project.org/] KM kumulativně riziko pre AML data NA kumulativně riziko pre AML data - skup M i i — skup NM i i i i r " 1 i i _ i i—' 1 ~ f 1 o 10 20 30 40 cas do relapsu {v týždňoch) 10 20 30 40 50 cas do relapsu {v týždňoch) Stanislav Katina r príkazy [r je voíne dostupné na http://cran.r-project.org/] Krivky prežívania a IS 1 1 L-l - skup M — skup NM l+ ll cas do relapsu {v týždňoch) plain skala Stanislav Katina Príklad 1 (pokrač.): ## odhady funkcie prežívania pre kazdu skupinu zvlast kmfit.aml.1<-survfit(Surv(cas[skup == 1],status[skup — kmfit.aml.2<-survfit(Surv(cas[skup == 2] ,status [skup — ## obrázok (hrubka čiary lwd, conf.int kresli IS) plot(kmfit.aml.1,xlab="cas do relapsu (v týždňoch)", ylab="pravděpodobnost dožitia",conf.int=FALSE,lwd=2) lines(kmfit.aml.2,lty=2,lwd=2) legend("topright",c("skup M","skup NM"),lty=c(1,2)) title(main="Krivky prezivania a IS",sub="plain skala") 1] )< 2] )< 'd 'd Stanislav Katina r príkazy [r je voíne dostupné na http://cran.r-project.org/] Príklad 2: ## odhady funkcie prezivania pre kazdu skupinu zvlast ## data z přednášky konci <- c(41,41,47,47,47,58,58,58,100,117) statusl <- c (0, 0, 1, 0, 0, 1, 1, 1, 0, 1) konc2 <- c(43,44,45,67,68,136,136,150,150,150) status2 <- c (0, 0, 0, 1, 0, 1, 1, 1, 1, 1) konc3 <- c(73,74,75,76,76,76,99,166,246) status3 <- c(0, 0, 0, 1,1,0,1,1,0) data.all <- c(koncl,konc2,konc3) status.all <- c(statusl,status2,status3) skup.all <- c(rep(1,10),rep(2,10),rep (3,9) ) ## testovanie hypotéz # default rho = 0 [log-rank] survdiff(formula = Surv(data.all, status.all)' # Peto test rho = 1 survdiff(formula = Surv(data.all, status.all)' 'Skup. all) 'Skup.all,rho=l) Stanislav Katina r príkazy [r je voíne dostupné na http://cran.r-project.org/] r príkazy [r je voíne dostupné na http://cran.r-project.org/] Príklad 2 (pokrač.): ## NA kumulatívne riziko # (pocitane rovnako ako v priklade 2, handout 1) LAM <- Lambda.NA.1 plot(cas,Lambda.NA,xlab="cas do relapsu (v týždňoch)",ylab="kumulativne riziko", type="n", ylim=range(LAM) xlim=c(0,166)) lines(cas,Lambda.NA.1,lty=l,lwd=2,type="s") lines(cas,Lambda.NA.2,lty=2,lwd=2,type="s") lines(cas,Lambda.NA.3,lty=3,lwd=2,type="s") abline(h=0,col="gray") title(main="NA kumulativne riziko") legend("bottomleft",c("skup l","skup 2","skup 3"),lty=c(1,2,3)) Stanislav Katina r príkazy [r je voíne dostupné na http://cran.r-project.org/] Príklad 2 (pokrač.): ## odhady funkcie prezivania pre kazdu skupinu zvlast kmfit.l <- survfit(Surv(konci,statusl) ~1) kmfit.2 <- survfit(Surv(konc2,status2) ~1) kmfit.3 <- survfit(Surv(konc3,status3) ~1) ## obrázok (hrubka čiary lwd, conf.int kresli IS) plot(kmfit.1,xlab="cas do zlyhania",ylab="pravděpodobnost dožitia",conf.int=FALSE,lwd=2) lines(kmfit.2,lty=2 , lwd=2) lines(kmfit.3,lty=3,lwd=2) legend("bottomleft",c("skup l","skup 2","skup 3"),lty=c(1,2,3)) title(main="Krivky prezivania a IS",sub="plain skala") Stanislav Katina NA kumulativne riziko skup 1 skup 2 skup 3 ~l- 0 50 100 150 cas do relapsu (v týždňoch) Stanislav Katina Krivky prezivania a IS skup 1 skup 2 skup 3 2: 40 so 60 cas do zlyhania plain skala 100 120 Stanislav Katina