Unterschiedliches Vorhersageplot von Survival Coxph und RMS Cph

Ich habe meine eigene, leicht verbesserte Version des Termplots erstellt, die ich in diesem Beispiel verwende. Sie finden sie hier . Ich habe zuvor auf SO gepostet, aber je mehr ich darüber nachdenke, desto mehr hängt dies wahrscheinlich mit der Interpretation des Cox Proportional Hazards-Modells zusammen als mit der tatsächlichen Codierung.

Das Problem

Wenn ich mir ein Hazard Ratio-Diagramm anschaue, erwarte ich einen Referenzpunkt, an dem das Konfidenzintervall natürlich 0 ist, und dies ist der Fall, wenn ich cph () aus dem verwende, rms packageaber nicht, wenn ich coxph () aus dem verwende survival package. Ist korrektes Verhalten von coxph () und wenn ja, was ist der Bezugspunkt? Außerdem hat die Dummy-Variable in coxph () ein Intervall und der Wert ist anders als ?$e^0$

Beispiel

Hier ist mein Testcode:

# 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)

Die cph-Diagramme

Dieser Code:

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")

gibt diese Handlung:

Die Coxph-Diagramme

Dieser Code:

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")

gibt diese Handlung:

Aktualisieren

Wie @Frank Harrell vorschlug und nachdem ich den Vorschlag in seinem letzten Kommentar angepasst hatte, bekam ich:

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)))

Dies gab diese sehr schöne Handlung:

Ich habe mir nach dem Kommentar noch einmal die kontrast.rms angesehen und diesen Code ausprobiert, der eine Handlung ergab ... obwohl wahrscheinlich noch viel mehr getan werden kann :-)

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")

Gab diese Handlung:

UPDATE 2

Prof. Thernau war so freundlich, die Handlungen zu kommentieren, denen es an Selbstvertrauen mangelt:

Die Glättungssplines in Coxph werden wie die in Gam so normalisiert, dass die Summe (Vorhersage) = 0 ist. Ich habe also keinen festen Einzelpunkt, für den die Varianz besonders klein ist.

Obwohl ich mit GAM noch nicht vertraut bin, scheint dies meine Frage zu beantworten: Dies scheint eine Frage der Interpretation zu sein.

Ich denke, es sollte definitiv einen Punkt geben, an dem das Konfidenzintervall die Breite Null hat. Sie können auch einen dritten Weg ausprobieren, bei dem ausschließlich Effektivfunktionen verwendet werden. Unter der Hilfedatei befindet sich ein Beispiel für kontrast.rms, um ein Hazard-Ratio-Diagramm zu erhalten. Es beginnt mit dem Kommentar # zeigt separate Schätzungen nach Behandlung und Geschlecht. Sie müssen Anti-Log, um das Verhältnis zu erhalten.