Pengecekan asumsi Proporsional Hazard dengan R


Pada uji log-rank, diperlukan asumsi Proportional Hazard (PH) untuk mendapatkan hasil pengujian yang valid.
Apakah yang dimaksud dengan asumsi PH?
Misalkan ingin dilihat survival experience sejumlah pasien. Pada pasien tersebut diberikan 2 jenis treatment, sebut treatment A dan B. Sebelumnya, untuk mengetahui apakah ada perbedaan survival experience antara kedua kelompok pasien tersebut, dapat dilakukan dengan uji log-rank.
Apakah hasil dari uji log-rank ini valid? Ya, jika asumsi PH terpenuhi.
Misal:
  • Pasien 1 diberi treatment A, dengan hazard hA(t)
  • Pasien 2 diberi treatment B, dengan hazard hB(t)
Maka, rasio hazard dari kedua pasien tersebut adalah:
HR=hA(t)hB(t)=c
untuk setiap titik waktu t, dengan c adalah suatu konstanta.
Sehingga, hA(t)=c×hB(t); yang menyatakan bahwa, hazard dari suatu individu adalah proportional terhadap hazard individu yang lain, berdasarkan atribut mereka (dalam hal ini treatment).
Bagaimana melakukan pengecekan asumsi PH?
Secara visual, dapat dilakukan dengan plot log(H(t)) terhadap waktu t dengan H(t) adalah fungsi kumulatif hazard.
Sebagai contoh, akan digunakan data leukemia.
library(survival) 
library(MASS) 

data(gehan) 
summary(gehan)
##       pair         time            cens            treat   
##  Min.   : 1   Min.   : 1.00   Min.   :0.0000   6-MP   :21  
##  1st Qu.: 6   1st Qu.: 6.00   1st Qu.:0.0000   control:21  
##  Median :11   Median :10.50   Median :1.0000               
##  Mean   :11   Mean   :12.88   Mean   :0.7143               
##  3rd Qu.:16   3rd Qu.:18.50   3rd Qu.:1.0000               
##  Max.   :21   Max.   :35.00   Max.   :1.0000
 
Pada data ini, waktu survival menyatakan lamanya durasi waktu dari dimulainya pengobatan (saat pasien dinyatakan sudah sembuh total atau sebagian) hingga penyakit tersebut kambuh kembali (dalam bulan).
Apakah ada perbedaan survival experience antara kelompok pasien yang diberi treatment 6-MP dengan kelompok pasien yang diberi placebo?
Bisa dilihat dari grafik fungsi survival kedua kelompok tersebut.
fit.gehan <- survfit(Surv(time, cens) ~ treat, data = gehan,
     conf.type = "log-log")

plot(fit.gehan, lty = 1:2, col=1:2,
     main=expression(paste("Estimasi Kaplan-Meier ", hat(S)(t))),
     xlab="Waktu survival (T)", ylab="Fungsi survival, S(t)")
legend(20, 1, c("treat: 6-MP", "treat: placebo"), lty = 1:2, col=1:2, bty="n")
 
l (T)", ylab="Fungsi survival, S(t)") legend(20, 1, c("treat: 6-MP", "treat: placebo"), lty = 1:2, col=1:2, bty="n")
Terlihat jelas bahwa pemberian 6-MP lebih membantu para pasien tersebut dibandingkan dengan placebo. Hal ini karena fungsi survival kelompok 6-MP selalu lebih tinggi daripada fungsi survival kelompok placebo; yang berarti pada tiap waktu t, peluang pasien dengan 6-MP untuk melewati waktu tersebut dengan tidak kambuh adalah lebih besar daripada peluang pada pasien dengan placebo.
Uji formalnya, dapat dilakukan dengan uji log-rank, sebagai berikut.
survdiff(Surv(time, cens) ~ treat, data = gehan)
## Call:
## survdiff(formula = Surv(time, cens) ~ treat, data = gehan)
## 
##                N Observed Expected (O-E)^2/E (O-E)^2/V
## treat=6-MP    21        9     19.3      5.46      16.8
## treat=control 21       21     10.7      9.77      16.8
## 
##  Chisq= 16.8  on 1 degrees of freedom, p= 4e-05
Hasil di atas mengonfirmasi visualisasi dari grafik taksiran survival yang kita peroleh sebelumnya, yaitu ada perbedaan survival experience dari kedua kelompok pasien tersebut (p-value=4e-5<0.05).
Apakah hasil pengujian ini valid?
fit.gehan <- survfit(Surv(time, cens) ~ treat, data = gehan,
     conf.type = "log-log")

plot(fit.gehan, lty = 1:2, col=1:2,
     main=expression(paste("Estimasi Kaplan-Meier ", hat(S)(t))),
     xlab="Waktu survival (T)", ylab="Fungsi survival, S(t)")
legend(20, 1, c("treat: 6-MP", "treat: placebo"), lty = 1:2, col=1:2, bty="n")
Terlihat jelas bahwa pemberian 6-MP lebih membantu para pasien tersebut dibandingkan dengan placebo. Hal ini karena fungsi survival kelompok 6-MP selalu lebih tinggi daripada fungsi survival kelompok placebo; yang berarti pada tiap waktu t, peluang pasien dengan 6-MP untuk melewati waktu tersebut dengan tidak kambuh adalah lebih besar daripada peluang pada pasien dengan placebo.
Uji formalnya, dapat dilakukan dengan uji log-rank, sebagai berikut.
survdiff(Surv(time, cens) ~ treat, data = gehan)
## Call:
## survdiff(formula = Surv(time, cens) ~ treat, data = gehan)
## 
##                N Observed Expected (O-E)^2/E (O-E)^2/V
## treat=6-MP    21        9     19.3      5.46      16.8
## treat=control 21       21     10.7      9.77      16.8
## 
##  Chisq= 16.8  on 1 degrees of freedom, p= 4e-05
Hasil di atas mengonfirmasi visualisasi dari grafik taksiran survival yang kita peroleh sebelumnya, yaitu ada perbedaan survival experience dari kedua kelompok pasien tersebut (p-value=4e-5<0.05).
Apakah hasil pengujian ini valid?
plot( survfit( Surv(time, cens ) ~ treat, data=gehan),
fun="cloglog", lty=1:2, mark.time=FALSE,
xlab="Waktu survival, T", ylab="log(H(t))" )
legend( 1,1, lty=1:2, legend=c("Control (N=21)", "6-MP (N=21)"),
bty="n" )
Terlihat bahwa plot logaritma dari fungsi kumulatif hazard untuk kedua kelompok tersebut paralel. Hal ini menyatakan bahwa asumsi PH yang diperlukan dalam uji log-rank terpenuhi, sehingga hasil dari uji tersebut valid.

Tidak ada komentar:

Posting Komentar

Relevant & Respectful Comments Only.