Gráfico de previsão diferente de sobrevivência coxph e rms cph

Criei minha própria versão ligeiramente aprimorada do termplot que utilizo neste exemplo, você pode encontrá-la aqui . Eu publiquei anteriormente no SO, mas quanto mais penso nisso, acredito que isso provavelmente esteja mais relacionado à interpretação do modelo de riscos proporcionais de Cox do que com a codificação real.

O problema

Quando eu olho para um gráfico de Hazard Ratio, espero ter um ponto de referência em que o intervalo de confiança naturalmente seja 0 e este é o caso quando eu uso o cph () do rms packagemas não quando eu uso o coxph () do survival package. O comportamento correto é de coxph () e, em caso afirmativo, qual é o ponto de referência? Além disso, a variável dummy no coxph () tem um intervalo e o valor é diferente de ? $e^0$

Exemplo

Aqui está o meu código de teste:

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

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

dá esse enredo:

cph () termplot2

Os gráficos 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")

dá esse enredo:

coxph () termplot2

Atualizar

Como o @Frank Harrell sugeriu e depois de ajustar a sugestão em seu comentário recente, recebi:

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

Isso deu esse enredo muito legal:

Gráfico de treliça

Eu olhei para o contrast.rms novamente após o comentário e tentei esse código que deu uma plotagem ... embora provavelmente haja muito mais que possa ser feito :-)

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

Deu este lote:

O gráfico de contraste

ATUALIZAÇÃO 2

Thernau teve a gentileza de comentar sobre a falta de confiança na trama:

Os splines de suavização em coxph, como os de gam, são normalizados para que soma (previsão) = 0. Portanto, não tenho um único ponto fixo para o qual a variação é muito pequena.

Embora eu ainda não esteja familiarizado com o GAM, isso parece responder à minha pergunta: isso parece ser uma questão de interpretação.

r survival cox-model Max Gordon
fonte

Vários comentários. Leia primeiro biostat.mc.vanderbilt.edu/Rrms para obter diferenças entre os pacotes rms e Design. Segundo, use plot () em vez de plot.Predict para economizar trabalho. Terceiro, você pode facilmente gerar gráficos para ambos os sexos, por exemplo, usando Predict (ajuste, idade, sexo, diversão = exp) # exp = anti-log; depois plot (resultado) ou plot (resultado, ~ idade | sexo). Você não usa "x = NA" no Predict. O rms usa gráficos de treliça, de modo que os parâmetros par graphics comuns e o mfrow não se aplicam. Veja exemplos no folheto do meu curso rms em biostat.mc.vanderbilt.edu/rms . Para contrast.rms, estude mais a documentação.

precisa

Muito obrigado pela sua contribuição. Atualizei o código com melhores exemplos e adicionei prof. A resposta de Thernau. PS Eu estou realmente animado que o seu planejando uma nova versão do livro, ampliando a seção de viés ponto de corte será muito útil como uma referência

Max Gordon

Você pode usar plote em contrastvez de plot.Predicte contrast.rms. Eu usaria byou lengthdentro, em seqvez de, timese daria contrastduas listas para que você especifique exatamente o que está sendo contrastado. Você também pode usar sombreamento com xYplotfaixas de confiança.

Frank Harrell

Obrigado. Prever porque, então, recebo a ajuda certa no RStudio - algo que, no meu caso, é muito mais vital do que o tempo necessário para escrever o nome completo da função (usando o preenchimento automático (guia), na verdade não perder tanto tempo).

Max