In acest curs ne propune sa interpretam rezultatele obtinute prin modelarea Cox, pentru setul de date veteran. Vom studia problema 2, pag 120, cartea A self learning book, Kleimaum.
Vom studia urmatoarele obiective:
sa se construiasca un model Cox avand variabilele tratament, celltype, karno, diagtime, age si prior;
folosindu-se modelul Cox anterior, sa se compare riscul de deces la persoanele cu celule de tip adeno, fata de persoanele cu celule de tip large;
De asemenea, sa se compare riscul de deces la persoanele cu celule de tip adeno, fata de persoanele cu celule de tip squamous;
pe baza modelului Cox de la punctul a), se poate vedea un efect al tratamentului asupra timpului de supravietuire?
sa se deduca o expresie a functiei de supravietuire pentru o persoana careia i s-a administrat tratament si care are celule de tip squamous (celelalte variabile vor fi ajustate la o valoare medie);
sa se constuiasca un model Cox care are in plus fata de modelul Cox de la a) teremenii: treatment × performance status; treatment × disease duration; treatment × age; si treatment × prior therapy. Pentru acest model revazut, sa se dea o expresie pentru riscul de deces determinat de administrarea tratamentului (celelate variabile fiind ajustate la o valoare medie).
Vom incarca setul de date veteran si vom construi 3 variabilele dummy asociate variabilei celltype (aceasta variabila este o variabila de tip categorie, cu 4 nivele). Trebuie mentionat ca cea dea a patra categorie, cu celltype=large, va fi considerara categorie de referinta. Iata mai jos comenzile care vor defini variabilele \(z_3,z_4,z_5\) (reamintim ca variabila dummy \(z_3\), are valoarea 1, daca celltype=adeno si in rest are valoarea 0; asemanator pentru restul variabilelor).
library(survival)
attach(veteran)
head(veteran)## trt celltype time status karno diagtime age prior
## 1 1 squamous 72 1 60 7 69 0
## 2 1 squamous 411 1 70 5 64 10
## 3 1 squamous 228 1 60 3 38 0
## 4 1 squamous 126 1 60 9 63 10
## 5 1 squamous 118 1 70 11 65 10
## 6 1 squamous 10 1 20 5 49 0S=Surv(veteran$time,veteran$status)
z3=(celltype=="adeno")*1
z3## [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [38] 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [75] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1
## [112] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0z4=(celltype=="smallcell")*1
z4## [1] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
## [38] 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [75] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0
## [112] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 z5=(celltype=="squamous")*1
z5## [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [38] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1
## [75] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## [112] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0Vom construi modelul cerut in problema la punctul a). Iata mai jos comanda care ne va da coeficientii \(\beta_i\) ce corespund variabilelor din model.
m=coxph(S~trt+z3+z4+z5+karno+diagtime+age+prior)
m## Call:
## coxph(formula = S ~ trt + z3 + z4 + z5 + karno + diagtime + age +
## prior)
##
## coef exp(coef) se(coef) z p
## trt 2.946e-01 1.343e+00 2.075e-01 1.419 0.15577
## z3 7.948e-01 2.214e+00 3.029e-01 2.624 0.00869
## z4 4.603e-01 1.584e+00 2.662e-01 1.729 0.08383
## z5 -4.013e-01 6.695e-01 2.827e-01 -1.420 0.15574
## karno -3.282e-02 9.677e-01 5.508e-03 -5.958 2.55e-09
## diagtime 8.132e-05 1.000e+00 9.136e-03 0.009 0.99290
## age -8.706e-03 9.913e-01 9.300e-03 -0.936 0.34920
## prior 7.159e-03 1.007e+00 2.323e-02 0.308 0.75794
##
## Likelihood ratio test=62.1 on 8 df, p=1.799e-10
## n= 137, number of events= 128Modelul Cox asociat este
\(h(T,X)=h_0(t) \exp(0.29*trt+0.79 z_3+0.46*z_4-0.4*z_5-0.32*karno+0.00*diagtime-0.087*age+0.007*prior)\)
\(HR=\frac{h(t, celltype=ADENO)}{h(t, celltype=large)}=\exp(0.79)=2.2\)
Deci, riscul de deces in categoria persoanelor cu celule de tip adeno este de 2.2 ori mai mare decat al celor din categoria large.
exp(0.79)## [1] 2.203396In mod asemanator vom riscul de deces intre categoriile de persoane cu celule de tip adeno fata de cele cu celule de tip squamous, vom calcula
\(HR=\frac{h(t, celltype=adeno)}{h(t, celltype=squamous)}=\exp(0.79-(-0.4))=\exp(1.19)\approx 3.3\)
exp(1.19)## [1] 3.287081Deci, riscul de deces in categoria persoanelor cu celule de tip adeno este de 3.3 ori mai mare decat al celor din categoria squamous.
In output-ul de la modelul a) se observa ca variabila trt are p-value=0.155>0.05, deci variabila este nesemnificativa. Deci, putem concluziona ca variabila tratament nu are efect asupra timpului de supravietuire.
Folosind relatia dintre functia de supravietuire si functia de hazard, avem ca
\(\hat{S}(t,X)=[\hat{S}_0(t,X)]^{2*\beta_1+\beta_5-0.32*\bar{karno}+0.00*\bar{diagtime}-0.087*\bar{age}+0.007*\bar{prior}}\)
mf=coxph(S~trt+z3+z4+z5+karno+diagtime+age+prior+trt*karno+trt*diagtime+trt*age+trt*prior)
mf## Call:
## coxph(formula = S ~ trt + z3 + z4 + z5 + karno + diagtime + age +
## prior + trt * karno + trt * diagtime + trt * age + trt *
## prior)
##
## coef exp(coef) se(coef) z p
## trt 2.065596 7.890003 1.450990 1.424 0.15457
## z3 0.795478 2.215499 0.303171 2.624 0.00869
## z4 0.520981 1.683679 0.271207 1.921 0.05474
## z5 -0.374699 0.687496 0.289041 -1.296 0.19485
## karno -0.014985 0.985127 0.018879 -0.794 0.42734
## diagtime -0.006551 0.993470 0.042786 -0.153 0.87830
## age 0.008530 1.008567 0.029257 0.292 0.77062
## prior 0.135326 1.144910 0.073353 1.845 0.06506
## trt:karno -0.011222 0.988841 0.011090 -1.012 0.31158
## trt:diagtime 0.004920 1.004932 0.023127 0.213 0.83153
## trt:age -0.014903 0.985207 0.018781 -0.794 0.42748
## trt:prior -0.090822 0.913181 0.048532 -1.871 0.06130
##
## Likelihood ratio test=67.81 on 12 df, p=8.219e-10
## n= 137, number of events= 128Deci, modelul Cox este de forma:
\(h(t,X)=h_0(t)\exp[\sum\limits_{i}\beta_i X_i]\), unde \(\beta_i\) sunt coeficientii variabilelor \(X_i\) din output R.