Comment faire une analyse ROC en R avec un modèle de Cox

10

J'ai créé quelques modèles de régression de Cox et j'aimerais voir à quel point ces modèles fonctionnent et j'ai pensé qu'une courbe ROC ou une statistique c pourrait être utile de la même manière que ces articles:

JN Armitage et JH van der Meulen, «Identification de la comorbidité chez les patients chirurgicaux à l'aide de données administratives avec le Royal College of Surgeons Charlson Score», British Journal of Surgery, vol. 97, num. 5, art. 772-781, Maj 2010.

Armitage a utilisé la régression logistique, mais je me demande s'il est possible d'utiliser un modèle du package de survie, le survivalROC donne un indice que cela est possible, mais je ne sais pas comment faire fonctionner cela avec une régression de Cox régulière.

Je serais reconnaissant si quelqu'un me montrait comment faire une analyse ROC sur cet exemple:

library(survival)
data(veteran)

attach(veteran)
surv <- Surv(time, status)
fit <- coxph(surv ~ trt + age + prior, data=veteran)
summary(fit)

Si possible, j'apprécierais à la fois la sortie c-statique brute et un joli graphique

Merci!

Mettre à jour

Merci beaucoup pour les réponses. @Dwin: Je voudrais juste être sûr d'avoir bien compris avant de sélectionner votre réponse.

Le calcul tel que je le comprends selon la suggestion de DWin:

library(survival)
library(rms)
data(veteran)

fit.cph <- cph(surv ~ trt + age + prior, data=veteran, x=TRUE, y=TRUE, surv=TRUE)

# Summary fails!?
#summary(fit.cph)

# Get the Dxy
v <- validate(fit.cph, dxy=TRUE, B=100)
# Is this the correct value?
Dxy = v[rownames(v)=="Dxy", colnames(v)=="index.corrected"]

# The c-statistic according to the Dxy=2(c-0.5)
Dxy/2+0.5

Je ne suis pas familier avec la fonction de validation et le bootstrap mais après avoir regardé le prof. La réponse de Frank Harrel ici sur R-help J'ai pensé que c'était probablement le moyen d'obtenir le Dxy. L'aide pour valider les états:

... La corrélation de rang Dxy de Somers à calculer à chaque rééchantillonnage (cela prend un peu plus de temps que les statistiques basées sur la vraisemblance). Les valeurs correspondant à la ligne Dxy sont égales à 2 * (C - 0,5) où C est l'indice C ou la probabilité de concordance.

Je suppose que je suis surtout confus par les colonnes. J'ai pensé que la valeur corrigée était celle que je devais utiliser mais je n'ai pas vraiment compris la sortie de validation:

      index.orig training    test optimism index.corrected   n
Dxy      -0.0137  -0.0715 -0.0071  -0.0644          0.0507 100
R2        0.0079   0.0278  0.0037   0.0242         -0.0162 100
Slope     1.0000   1.0000  0.2939   0.7061          0.2939 100
...

Dans la question R-help, j'ai compris que je devrais avoir "surv = TRUE" dans le cph si j'ai des strates mais je ne sais pas quel est le but du paramètre "u = 60" dans la fonction de validation. Je vous serais reconnaissant de bien vouloir m'aider à les comprendre et à vérifier que je n'ai commis aucune erreur.

Max Gordon
la source
2
Je voudrais probablement jeter un œil au paquet rms et à sa cph()commande.
chl
2
index.correctedc'est ce qu'il faut souligner. Ce sont des estimations des performances futures probables. u=60n'est pas nécessaire validatecar vous n'avez pas de strates. Si vous aviez des strates, les courbes de survie peuvent se croiser et vous devez spécifier un moment particulier pour obtenir la zone ROC généralisée.
Frank Harrell

Réponses:

2

@chl a indiqué une réponse spécifique à votre question. La cphfonction du package «rms» produira un Somers-D qui peut être transformé trivialement en un c-index. Cependant, Harrell (qui a introduit le c-index dans la pratique biostatistique) pense que ce n'est pas judicieux comme stratégie générale pour évaluer les mesures pronostiques, car il a un faible pouvoir de discrimination entre les alternatives. Au lieu de compter sur la littérature chirurgicale pour vos conseils méthodologiques, il serait plus sage de rechercher la sagesse accumulée dans le texte de Harrell, "Regression Modeling Strategies" ou Steyerberg's "Clinical Prediction Models".

DWin
la source
4
Merci pour la note. Je pense que et ne sont pas mauvais pour décrire la discrimination prédictive d'un seul modèle pré-spécifié. Mais comme vous l'avez dit, ils n'ont pas le pouvoir de faire plus que cela. CDxyC
Frank Harrell
Merci pour vos réponses, ma situation est que j'ai trois scores différents que je veux comparer et voir comment ils fonctionnent. Je n'ai pas eu le temps d'examiner la partie Somers-D et je reviendrai une fois que j'aurai eu le temps (j'ai jeté un coup d'œil rapide et je n'ai rien trouvé d'utile). J'ai également commandé le livre @FrankHarrell, "Regression Modeling Strategies", ISBN 13: 978-0387952321, et j'espère qu'il me guidera dans mes choix.
Max Gordon
2
Puisque Dxy = 2 * (c- 0,5), le calcul de c Dxy donné devrait être trivial.
DWin
3

Selon vos besoins, l'intégration d'un modèle dans un modèle plus grand et un test de rapport de vraisemblance "chunk" pour la valeur ajoutée des variables supplémentaires vous donneront un test puissant. Mon livre parle d'un indice issu de cette approche (l '"indice d'adéquation").χ2

Frank Harrell
la source
+1 pour m'avoir guidé dans la bonne direction. Je viens de terminer la statistique C et le score plus détaillé que je regarde avait une statistique C de 0,4365081 tandis que l'autre avait 0,4414625 (je suppose que je devrais compter 0,5-Dxy / 2 dans mon cas). J'ai pris un bon moment à faire le calcul sur mon échantillon de 140 000; J'ai dû abaisser les bootstraps à 10 et je ne sais pas quel est l'impact de cela. J'ai hâte de lire votre livre (il est par la poste) et j'espère qu'il m'aidera à mieux comprendre la méthodologie et à comparer la statistique C avec l'indice d'adéquation.
Max Gordon
Bien. Il n'est pas facile de dire si 0,44 contre 0,43 signifie beaucoup sans regarder les distributions des valeurs prédites.
Frank Harrell
Je comprends qu'il est difficile de commenter des chiffres comme celui-là. Je vais essayer de regarder la distribution. Ma principale interprétation du résultat est que très peu est expliqué par mon modèle et même s'il y a une petite différence, ce n'est probablement pas très important. Il serait intéressant de savoir à quoi s'attendre dans un contexte de survie - atteindre une valeur de 0,8 comme ils l'ont fait dans l'analyse que j'ai référencée dans ma question semble assez loin ... mais là encore ma survie est la survie d'une prothèse implantée et pas la survie des patients. Ils ont également utilisé une régression logistique qui modifie peut-être l'estimation.
Max Gordon
La régression logistique ne fonctionnerait pas si le temps est important ou si le temps de suivi varie selon les sujets. De retour à la question d'origine, les risques prévus auront une distribution étroite si très peu de variations sont expliquées par le modèle.
Frank Harrell
Je viens de recevoir votre livre ... J'ai eu un verrouillage rapide sur la partie survie mais quand j'essaie votre étude de cas au chapitre 20, mais je reçois une erreur sur la partie impute (w, sz): 'variable sz n'a pas un attribut names () '. J'ai suivi chapt. 8: a chargé la base de données avec getHdata (prostate) (n'a pas pu trouver le site Web dans le livre), a fait le w <- transcan (~ sz + sg + ap + sbp + dbp + age + wt + hg + ekg + pf + bm + hx, imputé = T, transformé = T, imcat = "arbre", données = prostate) mais je n'ai rien trouvé sur la dénomination ...
Max Gordon