Anggota Kelompok:

1. Moch. Wahyu Septiadi (5003211019)

2. Abdillah Ilham (5003211069)

3. Erlin Wahyuni (5003211097)

4. Rafi Aditya Aryaputra (5003211129)

5. Calista Deva Salfatah (5003211137)

Import Library dan Data

library(survival)
## Warning: package 'survival' was built under R version 4.3.3
leu = read.csv("leu2.csv")
leu
##    Time Censor Treatment logWBC
## 1     6      1         1   2.31
## 2     6      1         1   4.06
## 3     6      1         1   3.28
## 4     7      1         1   4.43
## 5    10      1         1   2.96
## 6    13      1         1   2.88
## 7    16      1         1   3.60
## 8    22      1         1   2.32
## 9    23      1         1   2.57
## 10    6      0         1   3.20
## 11    9      0         1   2.80
## 12   10      0         1   2.70
## 13   11      0         1   2.60
## 14   17      0         1   2.16
## 15   19      0         1   2.05
## 16   20      0         1   2.01
## 17   25      0         1   1.78
## 18   32      0         1   2.20
## 19   32      0         1   2.53
## 20   34      0         1   1.47
## 21   35      0         1   1.45
## 22    1      1         0   2.80
## 23    1      1         0   5.00
## 24    2      1         0   4.91
## 25    2      1         0   4.48
## 26    3      1         0   4.01
## 27    4      1         0   4.36
## 28    4      1         0   2.42
## 29    5      1         0   3.49
## 30    5      1         0   3.97
## 31    8      1         0   3.52
## 32    8      1         0   3.05
## 33    8      1         0   2.32
## 34    8      1         0   3.26
## 35   11      1         0   3.49
## 36   11      1         0   2.12
## 37   12      1         0   1.50
## 38   12      1         0   3.06
## 39   15      1         0   2.30
## 40   17      1         0   2.95
## 41   22      1         0   2.73
## 42   23      1         0   1.97

Tranforming into Survival Object

y = Surv(leu$Time, leu$Censor)

Modelling Cox Proportional-Hazard

Model 1: Treatment

model1 <- coxph(y ~ leu$Treatment)
model1
## Call:
## coxph(formula = y ~ leu$Treatment)
## 
##                  coef exp(coef) se(coef)      z        p
## leu$Treatment -1.5721    0.2076   0.4124 -3.812 0.000138
## 
## Likelihood ratio test=16.35  on 1 df, p=5.261e-05
## n= 42, number of events= 30

Model 2: With White Blood Cell

model2 <- coxph(y ~ leu$Treatment + leu$logWBC)
model2
## Call:
## coxph(formula = y ~ leu$Treatment + leu$logWBC)
## 
##                  coef exp(coef) se(coef)      z       p
## leu$Treatment -1.3861    0.2501   0.4248 -3.263  0.0011
## leu$logWBC     1.6909    5.4243   0.3359  5.034 4.8e-07
## 
## Likelihood ratio test=46.71  on 2 df, p=7.187e-11
## n= 42, number of events= 30

Model 3: With Interactions

model3 <- coxph(y ~ leu$Treatment + leu$logWBC + leu$Treatment*leu$logWBC)
model3
## Call:
## coxph(formula = y ~ leu$Treatment + leu$logWBC + leu$Treatment * 
##     leu$logWBC)
## 
##                              coef exp(coef) se(coef)      z        p
## leu$Treatment            -2.37491   0.09302  1.70547 -1.393    0.164
## leu$logWBC                1.55489   4.73459  0.39866  3.900 9.61e-05
## leu$Treatment:leu$logWBC  0.31752   1.37372  0.52579  0.604    0.546
## 
## Likelihood ratio test=47.07  on 3 df, p=3.356e-10
## n= 42, number of events= 30

Assumption Checking

One-at-a-time Strategy (Graphic Method)

One-at-a-time Strategy: Treatment

y = Surv(leu$Time,leu$Censor)
minusloglog = function(x){-log(-log(x))}

kmfit_leu = survfit(y ~ leu$Treatment)
plot(kmfit_leu, fun = minusloglog, xlab ="Time",
       ylab = "-log-log S", lty = c("solid","dashed"))
legend("topright", c("Placebo", "6-MP treatment"), lty = c("solid", "dashed"))

Grafik \(-log(-log(\hat{S(t)}))\) survival function menunjukkan bahwa kelompok placebo (garis penuh) memiliki tingkat kelangsungan hidup yang lebih rendah dibandingkan kelompok 6-Mercaptopurine (6-MP) treatment (garis putus-putus), terlihat dari kurva kelompok placebo yang menurun lebih cepat. Namun yang sesungguhnya perlu diperhatikan untuk pengujian asumsi proportional pada Cox Proportional Hazard Model adalah proporsionalitas dari \(-log(-log(\hat{S(t)}))\) untuk setiap treatment. Pada plot di atas, antara treatment placebo dan MP-6 kurang lebih sejajar. Sehingga, Treatment tidak melanggar asumsi proporsionalitas.

One-at-a-time strategy: logWBC

leu$logWBC_cat = cut(leu$logWBC, breaks = c(0,2.3, 3, 100), 
                      labels = c("Low", "Medium","High"))
kmfit_leu2 = survfit(y ~ leu$logWBC_cat)
names(kmfit_leu2)
##  [1] "n"         "time"      "n.risk"    "n.event"   "n.censor"  "surv"     
##  [7] "std.err"   "cumhaz"    "std.chaz"  "strata"    "type"      "logse"    
## [13] "conf.int"  "conf.type" "lower"     "upper"     "t0"        "call"
plot(kmfit_leu2, fun = minusloglog, xlab ="Time",
     ylab = "-log-log S", col = c("black", "red", "blue"))
legend("topright", c("Low", "Medium", "High"),
       lty=c("solid"), col = c("black", "red", "blue"))

Pada prediktor LogWBC, hasil yang didapakan untuk setiap kategorisasi kondisi White Blood Cell cenderung parallel, sehingga bisa dikatakan bahwa prediktor Kategorisasi LogWBC memenuhi asumsi proporsional.

Goodness-of-fit Method

gof_test = cox.zph(model2)
gof_test
##                  chisq df    p
## leu$Treatment 8.27e-05  1 0.99
## leu$logWBC    4.00e-02  1 0.84
## GLOBAL        4.02e-02  2 0.98

Pada Goodness of fit, kita melakukan uji apakah suatu prediktor memiliki kesesuaian dengan asumsi proporsionalitas dalam model Cox Proporsional Hazard. Dengan hipotesis uji untuk goodness of fit berikut. \[H_0: \rho_i = 0; i=1,2, ..., p\] \[H_1: \rho \ne 0; i=1,2, ..., p\] Dimana gagal tolak \(H_0\) artinya prediktor memnuhi asumsi Proportional Hazard dan tolak \(H_1\) sebaliknya untuk setiap prediktor. Dari hasil tersebut, semua prediktor memenuhi asumsi proportional hazard karena \(P-Value>\alpha (0,05)\).