¿Cómo interpretar la salida de predic.coxph?

Después de instalar un modelo cox, es posible hacer predicciones y recuperar el riesgo relativo de nuevos datos. Lo que no entiendo es cómo se calcula el riesgo relativo para un individuo y a qué se relaciona (es decir, el promedio de la población). ¿Alguna recomendación sobre recursos para ayudar a entender (no he avanzado mucho en el análisis de supervivencia, así que cuanto más simple, mejor)

predict.coxph()calcula la razón de riesgo en relación con el promedio de la muestra para todas las variables predictoras . Los factores se convierten en predictores ficticios como de costumbre, cuyo promedio se puede calcular. Recuerde que el modelo Cox PH es un modelo lineal para el log-hazard ln h ( t ) :$p$$\ln h(t)$

Donde es el peligro inicial no especificado. De manera equivalente, el peligro h ( t ) se modela como h ( t ) = h 0 ( t ) ⋅ e β 1 X 1 + ⋯ + β p X p = h 0 ( t ) ⋅ e X β . La razón de riesgo entre dos personas i e i ' con valores predictores$h_{0}(t)$$h(t)$$h(t) = h_{0}(t) \cdot e^{\beta_{1} X_{1} + \dots + \beta_{p} X_{p}} = h_{0}(t) \cdot e^{\bf{X} \bf{\beta}}$$i$$i'$ y X i ' es por lo tanto independiente del riesgo de referencia e independiente del tiempot:$\bf{X}_{i}$$\bf{X}_{i'}$$t$

Para la razón de riesgo estimada entre las personas e i ′ , simplemente conectamos los coeficientes estimados b 1 , ... , b p para el β 1 , ... , β p , dando e X i b y e X i ′ b .$i$$i'$$b_{1}, \ldots, b_{p}$$\beta_{1}, \ldots, \beta_{p}$$e^{\bf{X}_{i} \bf{b}}$$e^{\bf{X}_{i'} \bf{b}}$

Como ejemplo en R, utilizo los datos del apéndice de John Fox en el modelo Cox-PH que proporciona un texto introductorio muy agradable. Primero, buscamos los datos y construimos un modelo simple de Cox-PH para el tiempo de arresto de prisioneros liberados ( fin: factor - recibió ayuda financiera con codificación ficticia "no"-> 0, "yes"-> 1 age,: edad al momento de la liberación, prio: número de condenas anteriores):

> URL   <- "http://socserv.mcmaster.ca/jfox/Books/Companion/data/Rossi.txt"
> Rossi <- read.table(URL, header=TRUE)                  # our data
> Rossi[1:3, c("week", "arrest", "fin", "age", "prio")]  # looks like this
  week arrest fin age prio
1   20      1  no  27    3
2   17      1  no  18    8
3   25      1  no  19   13

> library(survival)                                      # for coxph()    
> fitCPH <- coxph(Surv(week, arrest) ~ fin + age + prio, data=Rossi)    # Cox-PH model
> (coefCPH <- coef(fitCPH))                              # estimated coefficients
     finyes         age        prio 
-0.34695446 -0.06710533  0.09689320 

Ahora conectamos los promedios de muestra para nuestros predictores en la fórmula :$e^{\bf{X} \bf{b}}$

meanFin  <- mean(as.numeric(Rossi$fin) - 1)   # average of financial aid dummy
    meanAge  <- mean(Rossi$age)                   # average age
meanPrio <- mean(Rossi$prio)                  # average number of prior convictions
rMean <- exp(coefCPH["finyes"]*meanFin        # e^Xb
           + coefCPH["age"]   *meanAge
           + coefCPH["prio"]  *meanPrio)

Ahora conectamos los valores predictores de las primeras 4 personas en la fórmula .$e^{\bf{X} \bf{b}}$

r1234 <- exp(coefCPH["finyes"]*(as.numeric(Rossi[1:4, "fin"])-1)
           + coefCPH["age"]   *Rossi[1:4, "age"]
           + coefCPH["prio"]  *Rossi[1:4, "prio"])

Ahora calcule el riesgo relativo para las primeras 4 personas contra el promedio de la muestra y compárelo con el resultado de predict.coxph().

> r1234 / rMean
[1] 1.0139038 3.0108488 4.5703176 0.7722002

> relRisk <- predict(fitCPH, Rossi, type="risk")   # relative risk
> relRisk[1:4]
        1         2         3         4 
1.0139038 3.0108488 4.5703176 0.7722002

Si tiene un modelo estratificado, la comparación predict.coxph()es con los promedios de estratos, esto se puede controlar a través de la referenceopción que se explica en la página de ayuda.