📊 Eksplorasi & Deskriptif
#Eksplorasi & Deskriptif
library(survival)
library(survminer)
## Loading required package: ggplot2
## Loading required package: ggpubr
##
## Attaching package: 'survminer'
## The following object is masked from 'package:survival':
##
## myeloma
library(ggpubr)
library(ggplot2)
#1.Data set yang digunakan
data(ovarian)
## Warning in data(ovarian): data set 'ovarian' not found
head(ovarian)
## futime fustat age resid.ds rx ecog.ps
## 1 59 1 72.3315 2 1 1
## 2 115 1 74.4932 2 1 1
## 3 156 1 66.4658 2 1 2
## 4 421 0 53.3644 2 2 1
## 5 431 1 50.3397 2 1 1
## 6 448 0 56.4301 1 1 2
#2. Deskripsi Data
summary(ovarian)
## futime fustat age resid.ds
## Min. : 59.0 Min. :0.0000 Min. :38.89 Min. :1.000
## 1st Qu.: 368.0 1st Qu.:0.0000 1st Qu.:50.17 1st Qu.:1.000
## Median : 476.0 Median :0.0000 Median :56.85 Median :2.000
## Mean : 599.5 Mean :0.4615 Mean :56.17 Mean :1.577
## 3rd Qu.: 794.8 3rd Qu.:1.0000 3rd Qu.:62.38 3rd Qu.:2.000
## Max. :1227.0 Max. :1.0000 Max. :74.50 Max. :2.000
## rx ecog.ps
## Min. :1.0 Min. :1.000
## 1st Qu.:1.0 1st Qu.:1.000
## Median :1.5 Median :1.000
## Mean :1.5 Mean :1.462
## 3rd Qu.:2.0 3rd Qu.:2.000
## Max. :2.0 Max. :2.000
#Distribusi Status (fustat)
table(ovarian$fustat)
##
## 0 1
## 14 12
📈 Analisis Kaplan–Meier
#Kurva Kaplan-Meier
km_fit <- survfit(Surv(futime, fustat) ~ rx, data = ovarian)
ggsurvplot(
km_fit,
data = ovarian,
pval = TRUE,
risk.table = TRUE,
surv.median.line = "hv",
title = "Kaplan–Meier Curve by Treatment (rx)"
)
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## ℹ The deprecated feature was likely used in the ggpubr package.
## Please report the issue at <https://github.com/kassambara/ggpubr/issues>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
## Ignoring unknown labels:
## • colour : "Strata"
#Uji Log-Rank
survdiff(Surv(futime, fustat) ~ rx, data = ovarian)
## Call:
## survdiff(formula = Surv(futime, fustat) ~ rx, data = ovarian)
##
## N Observed Expected (O-E)^2/E (O-E)^2/V
## rx=1 13 7 5.23 0.596 1.06
## rx=2 13 5 6.77 0.461 1.06
##
## Chisq= 1.1 on 1 degrees of freedom, p= 0.3
⚙️ Pemodelan Cox Proportional Hazards
#5. Estimasi Model
fit <- coxph(
Surv(futime, fustat) ~ age + resid.ds + rx + ecog.ps,
data = ovarian
)
summary(fit)
## Call:
## coxph(formula = Surv(futime, fustat) ~ age + resid.ds + rx +
## ecog.ps, data = ovarian)
##
## n= 26, number of events= 12
##
## coef exp(coef) se(coef) z Pr(>|z|)
## age 0.12481 1.13294 0.04689 2.662 0.00777 **
## resid.ds 0.82619 2.28459 0.78961 1.046 0.29541
## rx -0.91450 0.40072 0.65332 -1.400 0.16158
## ecog.ps 0.33621 1.39964 0.64392 0.522 0.60158
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## age 1.1329 0.8827 1.0335 1.242
## resid.ds 2.2846 0.4377 0.4861 10.738
## rx 0.4007 2.4955 0.1114 1.442
## ecog.ps 1.3996 0.7145 0.3962 4.945
##
## Concordance= 0.807 (se = 0.068 )
## Likelihood ratio test= 17.04 on 4 df, p=0.002
## Wald test = 14.25 on 4 df, p=0.007
## Score (logrank) test = 20.81 on 4 df, p=3e-04
#6. Uji Asumsi Proportional Hazards
ph_test <- cox.zph(fit)
ph_test
## chisq df p
## age 0.170 1 0.680
## resid.ds 1.155 1 0.282
## rx 0.595 1 0.440
## ecog.ps 2.928 1 0.087
## GLOBAL 4.455 4 0.348
ggcoxzph(ph_test)
#7. Kurva Survival Terprediksi
ggadjustedcurves(
fit,
data = ovarian,
variable = "rx",
legend.title = "Treatment",
title = "Adjusted Survival Curves by Treatment"
)
## Ignoring unknown labels:
## • fill : "Treatment"
## • linetype : "Treatment"
## • shape : "Treatment"
