Análisis de Supervivencia

Omar Rodríguez Torres

2024-09-20

Introducción

En este reporte, realizaremos un análisis de supervivencia utilizando el modelo Kaplan-Meier en un conjunto de datos relacionado con el carcinoma hepatocelular (HCC). Exploraremos cómo diferentes variables afectan la supervivencia de los pacientes y visualizaremos las curvas de supervivencia correspondientes.

¿Qué es el estimador de Kaplan-Meier?

El estimador de Kaplan-Meier es una herramienta no paramétrica utilizada para estimar la función de supervivencia a partir de datos de tiempo hasta evento. Es especialmente útil en estudios donde algunos sujetos pueden no haber experimentado el evento de interés al final del período de estudio (datos censurados).

La fórmula básica del estimador de Kaplan-Meier es:

\[ \hat{S}(t) = \prod_{t_i \leq t} \left(1 - \frac{d_i}{n_i}\right) \]

Donde:

\(\hat{S}(t)\) es la estimación de la probabilidad de supervivencia más allá del tiempo \(t\).
\(t_i\) son los tiempos en los que ocurre al menos un evento.
\(d_i\) es el número de eventos que ocurren en el tiempo \(t_i\).
\(n_i\) es el número de individuos en riesgo justo antes del tiempo \(t_i\).

Este estimador es fundamental para entender cómo cambia la probabilidad de supervivencia a lo largo del tiempo en estudios clínicos y epidemiológicos.

Configuración del Entorno

Primero, cargamos las librerías necesarias para el análisis.

# Cargar las librerías necesarias
library(survival)
library(survminer)
library(dplyr)
library(DT) # Librería añadida para tablas dinámicas

Carga y Preparación de Datos

Descargamos los datos del archivo hcc-data.txt desde el siguiente enlace. Este conjunto de datos, que corresponde a pacientes con carcinoma hepatocelular (HCC), fue recolectado en un hospital universitario en Portugal e incluye información clínica real de 165 pacientes diagnosticados con HCC. A continuación, cargamos los datos y asignamos nombres descriptivos a las columnas para facilitar su análisis.

# Leer los datos
data <- read.table("hcc+survival/hcc-survival/hcc-data.txt", header = FALSE, sep = ",")

# Cambiar los nombres de las columnas
colnames_data <- c(
  "gender", "symptoms", "alcohol", "hepatitis_b_surface_antigen", "hepatitis_b_e_antigen",
  "hepatitis_b_core_antibody", "hepatitis_c_virus_antibody", "cirrhosis", "endemic_countries",
  "smoking", "diabetes", "obesity", "hemochromatosis", "arterial_hypertension",
  "chronic_renal_insufficiency", "human_immunodeficiency_virus", "nonalcoholic_steatohepatitis",
  "esophageal_varices", "splenomegaly", "portal_hypertension", "portal_vein_thrombosis",
  "liver_metastasis", "radiological_hallmark", "age_at_diagnosis", "grams_of_alcohol_per_day",
  "packs_of_cigarettes_per_year", "performance_status", "encephalopathy_degree", "ascites_degree",
  "international_normalised_ratio", "alpha_fetoprotein", "haemoglobin", "mean_corpuscular_volume",
  "leukocytes", "platelets", "albumin", "total_bilirubin", "alanine_transaminase",
  "aspartate_transaminase", "gamma_glutamyl_transferase", "alkaline_phosphatase",
  "total_proteins", "creatinine", "number_of_nodules", "major_dimension_of_nodule",
  "direct_bilirubin", "iron", "oxygen_saturation", "ferritin", "class_attribute"
)

colnames(data) <- colnames_data

# Convertir variables numéricas y categóricas adecuadamente
data$class_attribute <- as.numeric(data$class_attribute)
data$gender <- factor(data$gender, levels = c(1, 0), labels = c("Hombre", "Mujer"))
data$cirrhosis <- factor(data$cirrhosis, levels = c(1, 0), labels = c("Sí", "No"))
data = data %>% mutate(
  class_attribute = ifelse(class_attribute == 0, 1, 0)
)

Verificamos los primeros registros de los datos:

# Mostrar una tabla dinámica de los primeros registros
datatable(head(data), options = list(pageLength = 5))

Análisis de Supervivencia Global

Creamos un objeto de supervivencia utilizando la edad al diagnóstico como tiempo y el atributo de clase como evento (1 = muerte, 0 = censura).

# Crear el objeto de supervivencia
surv_object <- Surv(time = data$age_at_diagnosis, event = data$class_attribute)

Ajustamos el modelo Kaplan-Meier y presentamos un resumen.

# Ajustar el modelo Kaplan-Meier
fit <- survfit(surv_object ~ 1)

# Resumen del modelo
summary(fit)

## Call: survfit(formula = surv_object ~ 1)
## 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##    27    162       1   0.9938 0.00615      0.98184        1.000
##    40    159       1   0.9876 0.00873      0.97061        1.000
##    41    157       1   0.9813 0.01070      0.96053        1.000
##    43    155       1   0.9750 0.01237      0.95102        0.999
##    44    153       1   0.9686 0.01383      0.94185        0.996
##    45    152       1   0.9622 0.01514      0.93300        0.992
##    49    148       2   0.9492 0.01750      0.91552        0.984
##    51    143       2   0.9359 0.01961      0.89827        0.975
##    52    139       2   0.9225 0.02152      0.88124        0.966
##    56    131       1   0.9154 0.02248      0.87241        0.961
##    57    126       1   0.9082 0.02344      0.86335        0.955
##    59    120       1   0.9006 0.02444      0.85394        0.950
##    60    117       1   0.8929 0.02541      0.84445        0.944
##    61    113       3   0.8692 0.02818      0.81567        0.926
##    62    108       1   0.8611 0.02905      0.80604        0.920
##    63    101       2   0.8441 0.03088      0.78569        0.907
##    64     95       3   0.8174 0.03352      0.75431        0.886
##    66     85       2   0.7982 0.03538      0.73178        0.871
##    67     80       2   0.7782 0.03720      0.70864        0.855
##    68     73       1   0.7676 0.03819      0.69626        0.846
##    69     70       2   0.7457 0.04012      0.67101        0.829
##    71     66       1   0.7344 0.04108      0.65810        0.819
##    72     58       5   0.6710 0.04627      0.58621        0.768
##    73     52       2   0.6452 0.04796      0.55777        0.746
##    74     45       2   0.6166 0.04993      0.52607        0.723
##    75     39       2   0.5849 0.05214      0.49119        0.697
##    76     35       1   0.5682 0.05326      0.47287        0.683
##    77     29       2   0.5290 0.05633      0.42939        0.652
##    78     26       2   0.4883 0.05889      0.38554        0.619
##    79     21       3   0.4186 0.06276      0.31200        0.562
##    80     17       1   0.3940 0.06372      0.28694        0.541
##    81     14       1   0.3658 0.06508      0.25813        0.518
##    82     12       2   0.3048 0.06701      0.19814        0.469
##    83      9       1   0.2710 0.06759      0.16620        0.442
##    84      8       2   0.2032 0.06550      0.10806        0.382
##    85      6       1   0.1694 0.06273      0.08194        0.350
##    86      5       1   0.1355 0.05862      0.05802        0.316
##    87      4       1   0.1016 0.05285      0.03666        0.282
##    88      2       1   0.0508 0.04460      0.00909        0.284

Como se puede observar la salida de R muestra la tabla de supervivencia, la tabla de riesgos y la tabla de eventos en cada tiempo. En dicha tabla se puede identificar:

time: Tiempo en el que ocurre el evento.
n.risk: Número de individuos en riesgo en ese tiempo.
n.event: Número de eventos en ese tiempo.
survival: Probabilidad de supervivencia.
std.err: Error estándar de la probabilidad de supervivencia.
lower 95% CI: Límite inferior del intervalo de confianza del 95%.
upper 95% CI: Límite superior del intervalo de confianza del 95%.

A continución, graficamos la curva Kaplan-Meier para visualizar estos resultados.

# Graficar la curva Kaplan-Meier
ggsurvplot(fit, data = data, conf.int = TRUE, ggtheme = theme_minimal(),
           title = "Curva de Supervivencia Kaplan-Meier para HCC",
           xlab = "Edad al Diagnóstico (años)", ylab = "Probabilidad de Supervivencia")

Se observa que la probabilidad de supervivencia disminuye a lo largo del tiempo, lo cual es esperado en un análisis de supervivencia. Sin embargo, es importante tener en cuenta que la curva de supervivencia es una estimación y está sujeta a incertidumbre.

Análisis por Género

Ahora podríamos estar interesados en cómo el género afecta la supervivencia de los pacientes con HCC. Ajustamos el modelo Kaplan-Meier por género y graficamos las curvas correspondientes. La manera más directa de hacerlo es filtrar los datos por género y ajustar el modelo por separado.

Hombres

De la tabla data filtramos los datos para hombres y ajustamos el modelo Kaplan-Meier.

# Filtrar datos para hombres
male_data <- data %>% filter(gender == "Hombre")

# Crear objeto de supervivencia
male_surv <- Surv(time = male_data$age_at_diagnosis, event = male_data$class_attribute)

# Ajustar el modelo Kaplan-Meier
male_fit <- survfit(male_surv ~ 1)

# Graficar la curva Kaplan-Meier para hombres
ggsurvplot(male_fit, data = male_data, conf.int = TRUE, ggtheme = theme_minimal(),
           title = "Curva Kaplan-Meier para Hombres",
           xlab = "Edad al Diagnóstico (años)", ylab = "Probabilidad de Supervivencia")

Mujeres

Para las mujeres, realizamos un procedimiento similar.

# Filtrar datos para mujeres
female_data <- data %>% filter(gender == "Mujer")

# Crear objeto de supervivencia
female_surv <- Surv(time = female_data$age_at_diagnosis, event = female_data$class_attribute)

# Ajustar el modelo Kaplan-Meier
female_fit <- survfit(female_surv ~ 1)

# Graficar la curva Kaplan-Meier para mujeres
ggsurvplot(female_fit, data = female_data, conf.int = TRUE, ggtheme = theme_minimal(),
           title = "Curva Kaplan-Meier para Mujeres",
           xlab = "Edad al Diagnóstico (años)", ylab = "Probabilidad de Supervivencia")

Comparación por Género

Sin embargo, para poder realizar esto de manera más eficiente, podemos ajustar el modelo Kaplan-Meier por género directamente. Trabajar de manera conjunta nos permite comparar las curvas de supervivencia de hombres y mujeres en un solo gráfico. Para ello utilizamos la función survfit con la fórmula surv_object ~ gender. Cabe destacar que para este análisis se considera como variable clasificadora el género, por lo que se debe tener en cuenta que la variable género es categórica.

# Ajustar el modelo Kaplan-Meier por género
fit_gender <- survfit(surv_object ~ gender, data = data)

# Resumen del ajuste
summary(fit_gender)

## Call: survfit(formula = surv_object ~ gender, data = data)
## 
##                 gender=Hombre 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##    41    131       1   0.9924  0.0076       0.9776        1.000
##    44    129       1   0.9847  0.0108       0.9638        1.000
##    45    128       1   0.9770  0.0131       0.9516        1.000
##    49    126       2   0.9615  0.0169       0.9289        0.995
##    51    121       2   0.9456  0.0200       0.9072        0.986
##    52    117       2   0.9294  0.0227       0.8860        0.975
##    56    112       1   0.9211  0.0240       0.8753        0.969
##    57    108       1   0.9126  0.0252       0.8645        0.963
##    59    102       1   0.9036  0.0265       0.8532        0.957
##    60     99       1   0.8945  0.0278       0.8417        0.951
##    61     95       3   0.8663  0.0313       0.8070        0.930
##    62     90       1   0.8566  0.0324       0.7954        0.923
##    63     85       2   0.8365  0.0346       0.7713        0.907
##    64     79       3   0.8047  0.0379       0.7338        0.882
##    66     69       2   0.7814  0.0402       0.7064        0.864
##    67     64       2   0.7570  0.0425       0.6781        0.845
##    69     56       2   0.7299  0.0451       0.6467        0.824
##    71     53       1   0.7162  0.0463       0.6310        0.813
##    72     47       4   0.6552  0.0514       0.5618        0.764
##    73     42       1   0.6396  0.0525       0.5446        0.751
##    74     36       2   0.6041  0.0553       0.5049        0.723
##    75     30       2   0.5638  0.0585       0.4601        0.691
##    76     27       1   0.5429  0.0599       0.4373        0.674
##    77     22       1   0.5183  0.0621       0.4098        0.655
##    78     20       1   0.4923  0.0641       0.3814        0.636
##    79     17       3   0.4055  0.0697       0.2894        0.568
##    80     13       1   0.3743  0.0710       0.2581        0.543
##    81     10       1   0.3368  0.0731       0.2201        0.515
##    82      9       1   0.2994  0.0739       0.1845        0.486
##    83      7       1   0.2566  0.0747       0.1450        0.454
##    84      6       2   0.1711  0.0702       0.0766        0.382
##    85      4       1   0.1283  0.0643       0.0480        0.343
##    87      3       1   0.0855  0.0553       0.0241        0.304
## 
##                 gender=Mujer 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##    27     30       1    0.967  0.0328       0.9045        1.000
##    40     28       1    0.932  0.0463       0.8456        1.000
##    43     25       1    0.895  0.0576       0.7888        1.000
##    68     16       1    0.839  0.0765       0.7017        1.000
##    72     11       1    0.763  0.1006       0.5889        0.988
##    73     10       1    0.686  0.1159       0.4930        0.956
##    77      7       1    0.588  0.1346       0.3758        0.921
##    78      6       1    0.490  0.1435       0.2763        0.870
##    82      3       1    0.327  0.1642       0.1221        0.875
##    86      2       1    0.163  0.1418       0.0299        0.895
##    88      1       1    0.000     NaN           NA           NA

# Graficar las curvas Kaplan-Meier por género
ggsurvplot(fit_gender, data = data, conf.int = TRUE, ggtheme = theme_minimal(),
           title = "Curva Kaplan-Meier por Género", pval = TRUE, risk.table = TRUE,
           xlab = "Edad al Diagnóstico (años)", ylab = "Probabilidad de Supervivencia",
           legend.title = "Género")

Se observa que las curvas de supervivencia para hombres y mujeres son similares, lo cual sugiere que el género no tiene un efecto significativo en la supervivencia de los pacientes con HCC en este conjunto de datos.

Análisis por Cirrosis

Podemos realizar un análisis similar para la cirrosis, una condición clínica común en pacientes con HCC. Ajustamos el modelo Kaplan-Meier por cirrosis y graficamos las curvas correspondientes.

# Ajustar el modelo Kaplan-Meier por cirrosis
fit_cirrhosis <- survfit(surv_object ~ cirrhosis, data = data)

# Resumen del ajuste
summary(fit_cirrhosis)

## Call: survfit(formula = surv_object ~ cirrhosis, data = data)
## 
##                 cirrhosis=Sí 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##    40    145       1   0.9931 0.00687       0.9797        1.000
##    41    143       1   0.9862 0.00972       0.9673        1.000
##    43    141       1   0.9792 0.01190       0.9561        1.000
##    44    139       1   0.9721 0.01375       0.9455        0.999
##    45    138       1   0.9651 0.01535       0.9355        0.996
##    49    136       2   0.9509 0.01811       0.9161        0.987
##    51    131       2   0.9364 0.02053       0.8970        0.977
##    52    127       2   0.9216 0.02270       0.8782        0.967
##    56    119       1   0.9139 0.02380       0.8684        0.962
##    57    115       1   0.9059 0.02488       0.8584        0.956
##    59    109       1   0.8976 0.02601       0.8481        0.950
##    60    106       1   0.8891 0.02710       0.8376        0.944
##    61    102       3   0.8630 0.03022       0.8058        0.924
##    62     97       1   0.8541 0.03119       0.7951        0.917
##    63     90       2   0.8351 0.03326       0.7724        0.903
##    64     84       3   0.8053 0.03626       0.7373        0.880
##    66     75       2   0.7838 0.03834       0.7122        0.863
##    67     71       2   0.7617 0.04031       0.6867        0.845
##    68     64       1   0.7498 0.04140       0.6729        0.836
##    69     61       2   0.7253 0.04354       0.6447        0.816
##    71     57       1   0.7125 0.04460       0.6303        0.806
##    72     50       5   0.6413 0.05025       0.5500        0.748
##    73     44       2   0.6121 0.05202       0.5182        0.723
##    74     37       2   0.5790 0.05422       0.4820        0.696
##    75     32       2   0.5429 0.05654       0.4426        0.666
##    76     28       1   0.5235 0.05775       0.4217        0.650
##    77     22       2   0.4759 0.06153       0.3693        0.613
##    79     17       3   0.3919 0.06711       0.2802        0.548
##    81     11       1   0.3563 0.06983       0.2426        0.523
##    83      8       1   0.3117 0.07395       0.1958        0.496
##    84      7       1   0.2672 0.07561       0.1534        0.465
##    85      6       1   0.2227 0.07499       0.1151        0.431
##    86      5       1   0.1781 0.07201       0.0807        0.393
##    87      4       1   0.1336 0.06636       0.0505        0.354
##    88      2       1   0.0668 0.05773       0.0123        0.363
## 
##                 cirrhosis=No 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##    27     15       1    0.933  0.0644       0.8153        1.000
##    78      7       2    0.667  0.1659       0.4094        1.000
##    80      4       1    0.500  0.1906       0.2369        1.000
##    82      3       2    0.167  0.1502       0.0285        0.975
##    84      1       1    0.000     NaN           NA           NA

# Graficar las curvas Kaplan-Meier por cirrosis
ggsurvplot(fit_cirrhosis, data = data, conf.int = TRUE, ggtheme = theme_minimal(),
           title = "Curva Kaplan-Meier por Cirrosis", pval = TRUE, risk.table = TRUE,
           xlab = "Edad al Diagnóstico (años)", ylab = "Probabilidad de Supervivencia",
           legend.title = "Cirrosis")

Se observa que la curva de supervivencia para pacientes sin cirrosis es significativamente más baja que para aquellos sin cirrosis en edades tempranas, esto se debe a que la cirrosis es una condición médica que sucede por un deterioro prologonado del hígado por lo que se presenta en edades avanzadas. Y esto se confirma ya que después de los 50 años la curva se invierte y la supervivencia de los pacientes con cirrosis es menor, factor que se mantiene constante a lo largo del tiempo.

Análisis por Consumo de Alcohol

Sin embargo, anteriormente hemos agrupado a los pacientes por variables categóricas, ahora vamos a realizar un análisis por consumo de alcohol basado en el consumo diario de alcohol en gramos. El consumo de alcohol es una medida continúa por lo que no se puede realizar un análisis directo como en los casos anteriores. Para ello, primero preparamos los datos y clasificamos el consumo de alcohol en tres niveles: bajo, alto y muy alto, basado en los valores de consumo de alcohol

# Reemplazar "?" por NA y convertir a numérico
data_alcohol <- data %>%
  mutate(grams_of_alcohol_per_day = na_if(grams_of_alcohol_per_day, "?"),
         grams_of_alcohol_per_day = as.numeric(grams_of_alcohol_per_day))

Con ello realizamos la transfromación de los datos de un formato string a un formato numérico y reemplazamos los valores faltantes por NA. Sin embargo, los datos siguen siendo continuos por lo que se debe clasificar en niveles para poder realizar un análisis de supervivencia. Para ello, primero observamos un resumen de los datos de consumo de alcohol y obtenemos los cuartiles para clasificar el consumo de alcohol.

summary(data_alcohol$grams_of_alcohol_per_day)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
##    0.00    0.00   75.00   71.01  100.00  500.00      48

Como se observa, los datos de consumo de alcohol varían desde 0 hasta 500 gramos por día, con un promedio de 71.01 gramos por día. A continuación, clasificamos el consumo de alcohol en tres niveles: bajo (menos de 75 gramos por día), alto (entre 75 y 100 gramos por día) y muy alto (más de 100 gramos por día).

Clasificamos el consumo de alcohol según los valores de CAPD.

# Clasificar el consumo de alcohol basado en MCV
data_alcohol <- data_alcohol %>%
  mutate(alcohol = case_when(
    grams_of_alcohol_per_day < 75 | is.na(grams_of_alcohol_per_day)  ~ "Bajo",
    grams_of_alcohol_per_day >= 75 & grams_of_alcohol_per_day < 100 ~ "Alto",
    grams_of_alcohol_per_day >= 100 ~ "Muy alto"
  ))

datatable(head(data_alcohol), options = list(pageLength = 5))

Con ello, se ha clasificado el consumo de alcohol en tres niveles: bajo, alto y muy alto. A continuación, ajustamos el modelo Kaplan-Meier por nivel de alcohol y graficamos las curvas correspondientes.

# Crear el objeto de supervivencia
surv_object_alcohol <- Surv(time = data_alcohol$age_at_diagnosis, event = data_alcohol$class_attribute)

# Ajustar el modelo Kaplan-Meier por nivel de alcohol
fit_alcohol <- survfit(surv_object_alcohol ~ alcohol, data = data_alcohol)

# Graficar las curvas Kaplan-Meier por nivel de alcohol
ggsurvplot(
  fit_alcohol, data = data_alcohol, conf.int = TRUE, 
  ggtheme = theme_minimal(), 
  pval = TRUE,
  title = "Curva Kaplan-Meier - Clasificación de Alcohol por MCV",
  legend.title = "Nivel de Alcohol",
  xlab = "Edad al Diagnóstico (años)", ylab = "Probabilidad de Supervivencia",
  risk.table = TRUE
)

Con ello se observa que la curva de supervivencia para los pacientes con un consumo de alcohol bajo es significativamente más alta que para aquellos con un consumo alto o muy alto. Sin embargo, no se observan diferencias significativas entre los pacientes con un consumo alto y muy alto de alcohol. Esto sugiere que el consumo moderado de alcohol puede estar asociado con una mayor supervivencia en pacientes con HCC. Sin embargo, es importante tener en cuenta que este análisis se basa en los datos de MCV y puede no reflejar completamente el consumo real de alcohol de los pacientes. Además que para el grupo de alto consumo se cuenta con un número de pacientes muy bajo. Y por último, por convenencia, se ha clasificado el consumo de alcohol con valores NA en el grupo de bajo consumo para tener datos con los suficientes pacientes para realizar el análisis, por lo que se debe tener en cuenta que este análisis puede no ser del todo preciso.

Conclusiones

En este análisis, hemos explorado cómo diferentes variables clínicas afectan la supervivencia de pacientes con HCC. Observamos diferencias significativas en las curvas de supervivencia al comparar por género y cirrosis. Sin embargo, el análisis por consumo de alcohol basado en MCV no mostró diferencias estadísticamente significativas, lo cual podría deberse a limitaciones en los datos o en la clasificación utilizada.

Estos hallazgos pueden ayudar a identificar factores de riesgo y guiar estrategias de tratamiento. Es importante considerar que el análisis de supervivencia es una herramienta poderosa en la investigación médica y debe complementarse con estudios adicionales para validar los resultados.

Referencias

Librerías utilizadas: survival, survminer, dplyr, DT
Datos: Archivo hcc-data.txt
Documentación de survival: CRAN Survival Package
Documentación de survminer: CRAN Survminer Package
Documentación de DT: CRAN DT Package
Conjunto de Datos: UCI Machine Learning Repository - HCC Survival Data

Nota: Este análisis es de carácter ilustrativo y se recomienda una revisión exhaustiva de los datos y métodos antes de extraer conclusiones clínicas.