Prediciendo el riesgo absoluto usando la regresión de Cox

Question

Estoy tratando de usar R para predecir el riesgo absoluto de desarrollar eventos adversos en una cohorte y compararlo con el resultado observado. ¿Debo usar survreg o coxph para hacer esto? ¿Alguien tan amable de explicar cómo hacer esto con código R?

El período medio de seguimiento de mi cohorte es solo de hasta 6 años, ¿así que puedo predecir el riesgo absoluto para cada individuo al final del período de seguimiento (incluyendo datos censurados y no censurados)?

Answer 1

Entonces, por "riesgo absoluto" voy a asumir que te refieres a la probabilidad acumulativa de un evento en el tiempo t, o la amenaza en el tiempo t.

En resumen, no, un modelo de riesgos proporcionales de Cox realmente no te devuelve esa información: el modelo en sí no calcula la amenaza subyacente, solo la diferencia relativa en la amenaza entre los valores de las covariables. Esto te libera bastante de tener que especificar la amenaza subyacente de tu resultado, que en muchos casos es desconocida, no particularmente de interés en primer lugar, o difícil de especificar usando un modelo paramétrico.

Si quieres estimar la función de riesgo subyacente, necesitas usar modelos de supervivencia paramétricos, como los utilizados en survreg. Hay una gran cantidad de tutoriales en línea, que incluyen código, para survreg y modelos paramétricos.

Este es uno de mis favoritos, ya que incluye un tratamiento teórico y bastante código. La documentación del paquete 'Survival' también es un buen lugar para comenzar.

Answer 2

Estoy de acuerdo con las respuestas anteriores de que Cox no fue principalmente diseñado para cálculos de riesgo absoluto, pero creo que fue más por razones históricas (las personas se centraban en riesgos relativos y no absolutos), y no veo una razón legítima para no usarlo para predicciones. Ya hay trabajos en modelado de predicción que utilizan con éxito Cox, pero creo que prefieren el paquete glmnet para la familia Cox, con flexibilidad en los términos de regularización y el paquete hdnom para las probabilidades de supervivencia (hdnom:::glmnet_survcurve).

Volviendo a la pregunta, yendo en ciclos para el modelo Cox básico de survfit, llegué a este código para obtener las probabilidades individuales estimadas de un evento:

bh=basehaz(cox_model_fit, newdata = data, centered = FALSE) #función de riesgo acumulado base para todas las covariables establecida en 0 (centered = false)
lps <- predict(cox_model_fit, data, type = "lp") #predictores lineales, betas x covariables para cada observación

tiempo_de_interes = 6 #si te interesa la supervivencia en t=6

i_tiempo_de_interes = match(1, round(bh$time,1) == tiempo_de_interes, nomatch = -100) #molestamente, no hay argumento de tiempo en la función basehaz, así que encuentra ese tiempo en la lista (se construye en todos los tiempos censurados y de eventos, así que elige un punto así o necesitas interpolar entre los puntos disponibles)

prob_evento_t = 1-exp(-bh[i_tiempo_de_interes,1])^exp(lps) #calcula el riesgo de evento a partir del riesgo base y los predictores lineales. devuelve un vector de probabilidades para cada persona en los "datos"

Luego puedes verificar prob_evento_t por separado para casos y controles/censurados y probar la diferencia, etc.

O bien, usa esta función para calcular la estadística de concordancia para datos de supervivencia, esto da la probabilidad de que dos observaciones aleatorias de casos y controles tengan las probabilidades en el orden correcto (menores posibilidades de evento para controles que para los casos). Para resultado binario, el índice de concordancia es igual al área bajo la curva ROC, por lo que es una medida muy común de qué tan bueno es un modelo para clasificación/modelos binarios.

survConcordance(Surv(data$survtimed, data$outcome) ~ lps)$concordance