Sebelumnya telah dibahas mengenai uji log-rank untuk membandingkan kurva daya tahan antar kelompok. Asumsi-asumsi yang digunakan dalam uji log-rank adalah:

  1. setiap individu yang diamati saling bebas

  2. non-informative censoring: sensoring yang terjadi tidak berkaitan dengan studi yang dilakukan

  3. proporsional hazard (PH): hazard dari satu kelompok seragam/ konstan terhadap kelompok yang lain (jika hazard lebih tinggi maka sepanjang waktu selalu lebih tinggi, dan sebaliknya atau dapat dikatakan bahwa kurva hazard paralel tidak saling memotong)

Salah satu cara untuk memeriksa asumsi proporsional hazard adalah dengan log-log plot, yaitu plot antara \(log(t)\) dan \(-log(-log(S(t)))\). Jika asumsi proporsional hazard tidak terpenuhi, maka biasanya dilakukan transformasi terhadap respon (peubah waktu), biasanya transformasi logaritma, yang mengarah ke accelerated failure time (AFT). Transformasi yang dilakukan mengarahkan kita untuk membuat asumsi bahwa respon mengikuti sebaran tertentu.

Dalam analisis daya tahan ketika waktu daya tahan (respon) diasumsikan memiliki sebaran tertentu maka digunakan model parametrik, misalnya model Eksponensial, Weibull, Logistik, Lognormal, dll. Karena waktu daya tahan diasumsikan memiliki sebaran tertentu, maka fungsi kepekatan peluangnya diketahui, yaitu \(f(t)\), yang dinyatakan sebagai suatu fungsi atas parameter-parameter yang tidak diketahui nilainya. Ingat kembali bahwa fungsi daya tahan dapat dinyatakan sebagai fungsi atas fungsi kepekatan peluang, yaitu: \(S(t)=P(T>t)=\int_t^\infty f(u) du\), dan fungsi hazard dapat dinyatakan sebagai \(h(t)=\frac{-d[S(t)/dt]}{S(t)}\).

Berikut adalah fungsi daya tahan dan fungsi hazard untuk beberapa sebaran:

Distribusi \(S(t)\) \(h(t)\)
Eksponensial \(exp(-\lambda t)\) \(\lambda\)
Weibull \(exp(-\lambda t^p)\) \(\lambda p t^{p-1}\)
Log-logistik \(\frac{1}{1+\lambda t^p}\) \(\frac{\lambda p t^{p-1}}{1+\lambda t^p}\)

Parameter \(\lambda\) merupakan reparameterisasi dari kovariat dan parameter model regresi.

Distribusi Eksponensial

Untuk model eksponensial fungsi hazardnya merupakan hazard konstan, \(h(t)=\lambda\). Hal ini berimplikasi bahwa hazard relatif antara dua sebaran daya tahan eksponensial , yaitu \(\frac{h_i(t)}{h_j(t)}\), konstan untuk semua waktu. Oleh karena itu, model eksponensial merupakan model yang memenuhi proporsinal hazard sekaligus accelerated failure time. Jika model eksponensial didefinisikan sebagai model proporsional hazard maka \(\lambda\) diparameterisasi sebagai \(exp(\beta_0+\sum_{i=1}^p \beta_iX_i)\), sementara jika model eksponensial didefinisikan sebagai model accelerated failure time maka \(\lambda\) diparameterisasi sebagai \(exp(-(\beta_0+\sum_{i=1}^p \beta_iX_i))\).

Ilustrasi

Pada ilustrasi ini akan digunakan data percobaan dari 42 pasien leukimia. Sebagian pasien diberi obat yang mengandung 6-mercaptopurine, sementara lainnya diberi placebo. Waktu yang diamati adalah waktu sampai pasien remisi.

library(MASS)
library(survival)
library(dplyr)
data("gehan")
glimpse(gehan)
## Rows: 42
## Columns: 4
## $ pair  <int> 1, 1, 2, 2, 3, 3, 4, 4, 5, 5, 6, 6, 7, 7, 8, 8, 9, 9, 10, 10, 11…
## $ time  <int> 1, 10, 22, 7, 3, 32, 12, 23, 8, 22, 17, 6, 2, 16, 11, 34, 8, 32,…
## $ cens  <int> 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 1, 0…
## $ treat <fct> control, 6-MP, control, 6-MP, control, 6-MP, control, 6-MP, cont…

Pertama akan diperiksa terlebih dahulu apakah terjadi pelanggaran asumsi proporsional hazard.

plot(survfit(Surv(time,cens)~treat, data=gehan), 
     col=c("black", "red"), fun="cloglog",
     xlab='log(t)', ylab='-log(-log(S(t)))')

Berdasarkan plot log-log tampak pola plot yang paralel, sehingga tidak terjadi pelanggaran asumsi proporsional hazard.

Selanjutnya dilakukan pemodelan regresi eksponensial.

gehan$treat=relevel(gehan$treat,'control')
fit=survreg(Surv(time,cens)~treat, data=gehan, dist='exponential')
summary(fit)
## 
## Call:
## survreg(formula = Surv(time, cens) ~ treat, data = gehan, dist = "exponential")
##             Value Std. Error    z       p
## (Intercept) 2.159      0.218 9.90 < 2e-16
## treat6-MP   1.527      0.398 3.83 0.00013
## 
## Scale fixed at 1 
## 
## Exponential distribution
## Loglik(model)= -108.5   Loglik(intercept only)= -116.8
##  Chisq= 16.49 on 1 degrees of freedom, p= 4.9e-05 
## Number of Newton-Raphson Iterations: 4 
## n= 42

Hasil diatas merupakan hasil model AFT, sehingga \(\lambda\) diparameterisasi sebagai \(exp(-(\beta_0+\sum_{i=1}^p \beta_iX_i))\). Jika dinyatakan dalam persamaan, maka dapat dituliskan sebagai berikut:

\(h(t)=exp(-(2.159+1.527treat:6MP))\)

atau,

\(S(t)=exp(2.159+1.527treat:6MP)\)

Sebelum kita interpretasikan hasil penduga parameter yang diperoleh, biasanya diperiksa terlebih dahulu signifikansi dari penduga parameter tersebut. Untuk semua penduga parameter (intercept dan treatcontrol) nilainya signifikan (p<0.05). Selanjutnya untuk menginterpretasikan masing-masing penduga parameter dilakukan dengan menghitung nilai eksponensialnya.

exp(fit$coefficients[1])
## (Intercept) 
##    8.666667
exp(fit$coefficients[2])
## treat6-MP 
##  4.602564
exp(fit$coefficients[1]+fit$coefficients[2])
## (Intercept) 
##    39.88889

Nilai 8.666667 diinterpretasikan sebagai waktu daya tahan pasien yang memperoleh placebo, nilai 39.88889 adalah waktu daya tahan pasien yang memperoleh obat 6-mercaptopurine, sementara nilai 4.602564 diinterpretasikan sebagai rasio waktu daya tahan antara pasien yang memperoleh obat 6-mercaptopurine dan placebo. Sehingga dapat disimpulkan bahwa waktu daya tahan pasien yang memperoleh obat 6-mercaptopurine 4.602564 lebih lama daripada pasien yang memperoleh obat placebo.

exp(-fit$coefficients[1])
## (Intercept) 
##   0.1153846
exp(-fit$coefficients[2])
## treat6-MP 
## 0.2172702
exp(-(fit$coefficients[1]+fit$coefficients[2]))
## (Intercept) 
##  0.02506964

Jika dihitung nilai eksponensial dari negatif penduga parameter, maka diperoleh fungsi hazard. Nilai 0.1153846 diinterpretasikan sebagai laju hazard pasien yang memperoleh placebo, nilai 0.02506964 adalah laju hazard pasien yang memperoleh obat 6-mercaptopurine, dan nilai 0.2172702 diinterpretasikan sebagai rasio hazard antara pasien yang memperoleh obat 6-mercaptopurine dan placebo.