Gráfico de predicción diferente de la supervivencia coxph y rms cph

Question

He creado mi propia versión ligeramente mejorada del termplot que utilizo en este ejemplo, puedes encontrarlo aquí . Ya he publicado anteriormente en SO pero cuanto más lo pienso, creo que esto probablemente esté más relacionado con la interpretación del modelo de riesgos proporcionales de Cox que con la codificación real.

El problema

Cuando miro un gráfico de Hazard Ratio espero tener un punto de referencia donde el intervalo de confianza sea naturalmente 0 y este es el caso cuando uso el cph() del rms package pero no cuando uso el coxph() del survival package . ¿Es correcto el comportamiento de coxph() y si es así cuál es el punto de referencia? Además, la variable ficticia en coxph() tiene un intervalo y el valor es distinto de $e^0$ ?

Ejemplo

Aquí está mi código de prueba:

# Load libs
library(survival)
library(rms)

# Regular survival
survobj <- with(lung, Surv(time,status))

# Prepare the variables
lung$sex <- factor(lung$sex, levels=1:2, labels=c("Male", "Female"))
labels(lung$sex) <- "Sex"
labels(lung$age) <- "Age"

# The rms survival
ddist <- datadist(lung)
options(datadist="ddist")
rms_surv_fit <- cph(survobj~rcs(age, 4)+sex, data=lung, x=T, y=T)

Las parcelas cph

Este código:

termplot2(rms_surv_fit, se=T, rug.type="density", rug=T, density.proportion=.05,
          se.type="polygon", yscale="exponential", log="y", 
          xlab=c("Age", "Sex"), 
          ylab=rep("Hazard Ratio", times=2),
          main=rep("cph() plot", times=2),
          col.se=rgb(.2,.2,1,.4), col.term="black")

da esta trama:

Las parcelas coxph

Este código:

termplot2(surv_fit, se=T, rug.type="density", rug=T, density.proportion=.05, 
          se.type="polygon", yscale="exponential", log="y", 
          xlab=c("Age", "Sex"), 
          ylab=rep("Hazard Ratio", times=2),
          main=rep("coxph() plot", times=2),
          col.se=rgb(.2,.2,1,.4), col.term="black")

da esta trama:

Actualización

Como @Frank Harrell sugirió y después de ajustar a lo largo de la sugerencia en su reciente comentario que tengo:

p <- Predict(rms_surv_fit, age=seq(50, 70, times=20), 
             sex=c("Male", "Female"), fun=exp)
plot.Predict(p, ~ age | sex,
             col="black",
             col.fill=gray(seq(.8, .75, length=5)))

Esto le dio a esta trama muy agradable:

He vuelto a mirar el contrast.rms después del comentario y he probado este código que ha dado una gráfica... aunque probablemente se pueda hacer mucho más :-)

w <- contrast.rms(rms_surv_fit, 
                  list(sex=c("Male", "Female"), 
                       age=seq(50, 70, times=20)))

xYplot(Cbind(Contrast, Lower, Upper) ~ age | sex, 
       data=w, method="bands")

Dio esta parcela:

ACTUALIZACIÓN 2

El profesor Thernau tuvo la amabilidad de comentar la falta de confianza en las parcelas:

Los splines de suavizado en coxph, al igual que los de gam, están normalizados de modo que sum(prediction) =0. Así que no tengo un punto fijo para para el cual la varianza es extra pequeña.

Aunque todavía no estoy familiarizado con la GAM, esto parece responder a mi pregunta: esto parece ser una cuestión de interpretación.

Answer 1

Creo que debería haber un punto en el que el intervalo de confianza sea de ancho cero. También puede intentar una tercera forma que es utilizar únicamente las funciones rms. Hay un ejemplo en el archivo de ayuda para contrast.rms para obtener un gráfico de hazard ratio. Comienza con el comentario # mostrar estimaciones separadas por tratamiento y sexo. Tendrá que antilogarizar para obtener la proporción.