Breviar teoretic*

Metodele actuariale se folosesc mai ales in situatia in care numarul inregistrarilor este mare, iar momentele de timp in care are loc evenimentul sunt masurate precis, astfel setul de date continand multe valori distincte.

In acest caz, metoda Kalpan-Meier produce tabele greu de citit, iar curba de supravietuire este foarte “fragmentata”.

In acest caz, se recomanda metodele actuariale si se va calcula estimatorul actuarial.

Timpii de supravietuire urmariti sunt grupati in intervale, de obicei de lungimi egale. Pentru fiecare interval \([t_i,t_{i+1})\) avem urmatoarele notatii:

-\(d_i\) numarul de persoane care au decedat in acest interval;

-\(c_i\) numarul de persoane cenzurate in acest interval;

-\(n_i\) numarul de persoane in viata la inceputul acestui interval;

-\(\tilde{n_i}=n_i-c_i/2\) numarul de persoane aflate la risc in timpul acestui interval; (aici se face presupunerea ca observatiile cenzurate sunt la risc doar pe jumatate din interval)

-rata de supravietuire \(1-d_i/\tilde{n_i}\).

Atunci, formula de limita-produs pentru estimarea probabilitatilor de supravietuire este:

\[\hat{S_i}=\prod\limits_{j=1}^i(1-\frac{d_j}{\tilde{n}_j})\]

Aplicatie (Exemplu 3.5/pag 44, din cartea *)

Un medicament nou este testat pe pacientii cu leucemie. Variabila care inregistreaztimpul (in sap.) pana boala revine contine urmatoarele observatii:

1,1,2,2,3,4,6,8,8,9,9+,10,10,11,11,11+,12,13+,14,14+,15,15,15+,16,17,17,19,20+,23,25+,26,27+,29.

Ne propunem estimarea probabilitatii de supravietuire.

Metoda actuariala analitica

Pentru calculul probabilitatilor de supravietuire pe intervale va rugam sa consultati materialul Prezentare Curs 4 de pe CV-UPT.

Metoda actuariala in R

Datele inregistrate se vor impartii in 6 intervale egale: [0,5),[5,10),[10,15),[15,20),[20,25),[25,30).

Pentru a aplica metodele actuarile este necesar sa incarcam urmatoarele pachete:

library(survival)
library(KMsurv)

Vom incarca in R momentele de timp, astfel:

-variabila tis inregistreaza intervalele de timp predefinite

-variabila time1 inregistreaza momentele de timp cu eveniment

-variabila time2 inregistreaza momentele de timp cenzurate

tis<-c(0,5,10,15,20,25,30)
time1<- c(1,1,2,2,3,4,6,8,8,9,10,10,11,11,12,14,15,15,16,17,17,19,23,26,29)
time2<-c(9,11,13,14,15,20,25,27)

Ne propunem sa numaram cate evenimente sunt inregistrate pe fiecare interval. Pentru a face aceasta numarare, se va defini functia nevent:

nsubs<-length(time1)+length(time2)
 index<-seq(1,length(tis)-1)
 nevent<-sapply(index,function(index){length(time1[time1<tis[index+1]&time1>=tis[index]])})
 nevent
## [1] 6 4 6 6 1 2

Ne propunem sa numaram cate evenimente sunt cenzurate pe fiecare interval. Pentru a face aceasta numarare, se va defini functia nlost astfel:

 nlost<-sapply(index, function(index){length(time2[time2<tis[index+1]&time2>=tis[index]])})
nlost
## [1] 0 1 3 1 1 2

Se va apela functia lifetable ce ne va da tabelul probabilitatilor de supravietuire pe intervale.

mytable<-lifetab(tis,nsubs,nlost,nevent)
 mytable
##       nsubs nlost nrisk nevent      surv         pdf     hazard    se.surv
## 0-5      33     0  33.0      6 1.0000000 0.036363636 0.04000000 0.00000000
## 5-10     27     1  26.5      4 0.8181818 0.024699828 0.03265306 0.06714081
## 10-15    22     3  20.5      6 0.6946827 0.040664352 0.06857143 0.08054306
## 15-20    13     1  12.5      6 0.4913609 0.047170648 0.12631579 0.09010497
## 20-25     6     1   5.5      1 0.2555077 0.009291188 0.04000000 0.08376358
## 25-30     4     2   3.0      2 0.2090517          NA         NA 0.08039061
##            se.pdf  se.hazard
## 0-5   0.013428163 0.01624808
## 5-10  0.011558838 0.01627204
## 10-15 0.014736482 0.02757976
## 15-20 0.016360417 0.04892942
## 20-25 0.008939145 0.03979950
## 25-30          NA         NA

Se va extrage din tabelul anterior coloana 5 care contine vectorul \(r\) al probabilitatilor de supravietuire, pe fiecare interval.

r<-c(mytable[,5],mytable[,5][length(mytable[,5])])

Se va reprezenta grafic acest vector.

r
## [1] 1.0000000 0.8181818 0.6946827 0.4913609 0.2555077 0.2090517 0.2090517
plot(tis,r,axes=TRUE, type="s",xlab="Timpul de suprav pe intervale",ylab="Valoarea functiei de supravietuire")

Metoda Kaplan-Meier clasica

time<-c(time1,time2)
time
##  [1]  1  1  2  2  3  4  6  8  8  9 10 10 11 11 12 14 15 15 16 17 17 19 23 26 29
## [26]  9 11 13 14 15 20 25 27
status1<-rep(1, times=length(time1),)
status2<-rep(0, times=length(time2))
status<-c(status1,status2)

Vom crea un tabel cu inregistrarile problemei

status
##  [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0
df<-data.frame(cbind(indiv=1:(length(time1)+length(time2)),time,status))
head(df)
##   indiv time status
## 1     1    1      1
## 2     2    1      1
## 3     3    2      1
## 4     4    2      1
## 5     5    3      1
## 6     6    4      1

Vom crea obiectul tip din Analiza Supravietuirii, ce contine inregistrarile din studiul analizat, iar observatiile cenzurate sunt marcate cu +

my.surv<-Surv(time,status) 
my.surv
##  [1]  1   1   2   2   3   4   6   8   8   9  10  10  11  11  12  14  15  15  16 
## [20] 17  17  19  23  26  29   9+ 11+ 13+ 14+ 15+ 20+ 25+ 27+

Vom aplica metoda Kaplan-Meier pentru a estima probabilitatile de supravietuire pentru fiecare moment de timp

my.fit<-survfit(my.surv~1) 
my.fit
## Call: survfit(formula = my.surv ~ 1)
## 
##       n events median 0.95LCL 0.95UCL
## [1,] 33     25     15      10      19
summary(my.fit)
## Call: survfit(formula = my.surv ~ 1)
## 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##     1     33       2    0.939  0.0415       0.8614        1.000
##     2     31       2    0.879  0.0568       0.7742        0.998
##     3     29       1    0.848  0.0624       0.7346        0.980
##     4     28       1    0.818  0.0671       0.6966        0.961
##     6     27       1    0.788  0.0712       0.6600        0.940
##     8     26       2    0.727  0.0775       0.5901        0.896
##     9     24       1    0.697  0.0800       0.5566        0.873
##    10     22       2    0.634  0.0843       0.4881        0.822
##    11     20       2    0.570  0.0870       0.4229        0.769
##    12     17       1    0.537  0.0881       0.3890        0.740
##    14     15       1    0.501  0.0892       0.3533        0.710
##    15     13       2    0.424  0.0906       0.2788        0.644
##    16     10       1    0.381  0.0909       0.2391        0.609
##    17      9       2    0.297  0.0883       0.1656        0.532
##    19      7       1    0.254  0.0852       0.1318        0.491
##    23      5       1    0.203  0.0820       0.0924        0.448
##    26      3       1    0.136  0.0778       0.0441        0.417
##    29      1       1    0.000     NaN           NA           NA

Vom reprezenta grafic curba probabilitatilor de supravietuire

plot(survfit(my.surv~1),xlim = c(0, 33), ylim = c(0, 1),mark.time=TRUE, xlab="timp",ylab="probabilitati de supravietuire", conf.int=FALSE) 

title('Curba Kaplan-Maier')

Prin comparearea celor doua grafice se observa ca o analiza statistica riguroasa se poate face mai simplu prin analiza pe intervale; de altfel, exista multe studii clinice care au datele deja organizate pe intervale( in special in cazul in care inregistrarile se fac in anumite saptamani ale studiului clinic, conform unei proceduri stabilite anterior, prin planul statistic).

*Referinta: Clinical Statistics, Olga Korosteleva