Analyse de survie dans R

Paquets

Certains packages que nous utiliserons aujourd’hui incluent:

• lubridate
• survival
• survminer

library(survival)
library(survminer)
library(lubridate)

Qu’est-ce que les données de survie?

Données de temps sur événement qui consistent en une heure de début et une heure de fin distinctes.

Exemples de cancer

• Temps de la chirurgie à la mort
• Délai entre le début du traitement et la progression
• Délai entre la réponse et la récidive

Exemples d’autres domaines

Les données temporelles sont courantes dans de nombreux domaines, y compris, mais sans s’y limiter

• Temps entre l’infection à VIH et le développement du SIDA
• Temps de crise cardiaque
• Délai d’apparition de la toxicomanie
• Temps pour commencer l’activité sexuelle
• Temps de dysfonctionnement de la machine

Alias ​​pour l’analyse de survie

Parce que l’analyse de survie est courante dans de nombreux autres domaines, elle porte également d’autres noms

• Analyse de fiabilité
• Analyse de durée
• Analyse de l’historique des événements
• Analyse des événements

L’ensemble de données pulmonaires

le lung jeu de données est disponible à partir du survival paquet dans R. Les données contiennent des sujets atteints d’un cancer du poumon avancé du North Central Cancer Treatment Group. Certaines variables que nous utiliserons pour démontrer les méthodes aujourd’hui incluent

• temps: Temps de survie en jours
• statut: censure statut 1 = censuré, 2 = mort
• sexe: Homme = 1 Femme = 2

Qu’est-ce que la censure?

RICH JT, NEELY JG, PANIELLO RC, VOELKER CCJ, NUSSENBAUM B, WANG EW. UN GUIDE PRATIQUE POUR COMPRENDRE LES COURBES KAPLAN-MEIER. Chirurgie de la tête et du cou en oto-rhino-laryngologie: journal officiel de l’American Academy of Otolaryngology Head and Neck Surgery. 2010; 143 (3): 331-336. doi: 10.1016 / j.otohns.2010.05.007.

Types de censure

Un sujet peut être censuré en raison de:

• Perte de suivi
• Retrait de l’étude
• Aucun événement à la fin de la période d’étude fixe

Plus précisément, ce sont des exemples de droite censure.

La censure gauche et la censure d’intervalle sont également possibles, et des méthodes existent pour analyser ce type de données, mais cette formation sera limitée à la censure droite.

Données de survie censurées

Dans cet exemple, comment calculer la proportion de personnes sans événement à 10 ans?

Les sujets 6 et 7 étaient sans événement à 10 ans. Les sujets 2, 9 et 10 avaient événement avant 10 ans. Les sujets 1, 3, 4, 5 et 8 ont été censuré avant 10 ans, donc nous ne savons pas s’ils ont organisé l’événement ou non d’ici 10 ans – comment incorporer ces sujets dans notre estimation?

Répartition du temps de suivi

• Les sujets censurés fournissent toujours des informations et doivent donc être correctement inclus dans l’analyse
• La distribution des temps de suivi est biaisée et peut différer entre les patients censurés et ceux qui ont des événements
• Les délais de suivi sont toujours positifs

Composantes des données de survie

Pour sujet (je):

1. Heure de l’évènement (T_i )

2. Temps de censure (C_i )

3. Indicateur d’événement ( delta_i ):

   • 1 si événement observé (c.-à-d. (T_i leq C_i ))
   • 0 si censuré (c.-à-d. (T_i> C_i ))

4. Temps observé (Y_i = min (T_i, C_i) )

Les heures observées et un indicateur d’événement sont fournis dans le lung Les données

• temps: Temps de survie en jours ( (Y_i ))
• statut: censure statut 1 = censuré, 2 = mort ( ( delta_i ))

Gérer les dates en R

Les données seront souvent accompagnées de dates de début et de fin plutôt que de durées de survie pré-calculées. La première étape consiste à s’assurer que celles-ci sont formatées en tant que dates dans R.

Créons un petit exemple de jeu de données avec des variables sx_date pour la date de la chirurgie et last_fup_date pour la dernière date de suivi.

date_ex <- 
  tibble(
    sx_date = c("2007-06-22", "2004-02-13", "2010-10-27"), 
    last_fup_date = c("2017-04-15", "2018-07-04", "2016-10-31")
    )

date_ex

## # A tibble: 3 x 2
##   sx_date    last_fup_date
##                 
## 1 2007-06-22 2017-04-15   
## 2 2004-02-13 2018-07-04   
## 3 2010-10-27 2016-10-31

Nous voyons que ce sont deux variables de caractère, ce qui sera souvent le cas, mais nous avons besoin qu'elles soient formatées en dates.

Formatage des dates - base R

date_ex %>% 
  mutate(
    sx_date = as.Date(sx_date, format = "%Y-%m-%d"), 
    last_fup_date = as.Date(last_fup_date, format = "%Y-%m-%d") 
    )

## # A tibble: 3 x 2
##   sx_date    last_fup_date
##               
## 1 2007-06-22 2017-04-15   
## 2 2004-02-13 2018-07-04   
## 3 2010-10-27 2016-10-31

• Notez qu'en base R le format doit inclure le séparateur ainsi que le symbole. par exemple. si votre date est au format m / d / Y alors vous auriez besoin format = "%m/%d/%Y"
• Consultez la liste complète des symboles de format de date sur https://www.statmethods.net/input/dates.html

Formatage des dates - paquet de lubrification

Nous pouvons également utiliser le lubridate package pour formater les dates. Dans ce cas, utilisez le ymd fonction

date_ex %>% 
  mutate(
    sx_date = ymd(sx_date), 
    last_fup_date = ymd(last_fup_date)
    )

## # A tibble: 3 x 2
##   sx_date    last_fup_date
##               
## 1 2007-06-22 2017-04-15   
## 2 2004-02-13 2018-07-04   
## 3 2010-10-27 2016-10-31

• Notez que contrairement à la base R option, les séparateurs n'ont pas besoin d'être spécifiés
• La page d'aide de ?dmy affichera toutes les options de format.

Calcul des temps de survie - base R

Maintenant que les dates sont formatées, nous devons calculer la différence entre l'heure de début et l'heure de fin dans certaines unités, généralement des mois ou des années. À la base R, utilisation difftime pour calculer le nombre de jours entre nos deux dates et le convertir en une valeur numérique en utilisant as.numeric. Convertissez ensuite en années en divisant par 365.25, le nombre moyen de jours dans une année.

date_ex %>% 
  mutate(
    os_yrs = 
      as.numeric(
        difftime(last_fup_date, 
                 sx_date, 
                 units = "days")) / 365.25
    )

## # A tibble: 3 x 3
##   sx_date    last_fup_date os_yrs
##                 
## 1 2007-06-22 2017-04-15      9.82
## 2 2004-02-13 2018-07-04     14.4 
## 3 2010-10-27 2016-10-31      6.01

Calcul des temps de survie - lubrification

En utilisant le lubridate paquet, l'opérateur %--% désigne un intervalle de temps, qui est ensuite converti en nombre de secondes écoulées à l'aide as.duration et finalement converti en années en divisant par dyears(1), qui donne le nombre de secondes dans une année.

date_ex %>% 
  mutate(
    os_yrs = 
      as.duration(sx_date %--% last_fup_date) / dyears(1)
    )

## # A tibble: 3 x 3
##   sx_date    last_fup_date os_yrs
##                 
## 1 2007-06-22 2017-04-15      9.82
## 2 2004-02-13 2018-07-04     14.4 
## 3 2010-10-27 2016-10-31      6.02

• Remarque: nous devons charger le lubridate package à l'aide d'un appel à library afin de pouvoir accéder aux opérateurs spéciaux (similaire à la situation avec les tuyaux)

Indicateur d'événement

Pour les composantes des données de survie, j'ai mentionné l'indicateur d'événement:

Indicateur d'événement ( delta_i ):

• 1 si événement observé (c.-à-d. (T_i leq C_i ))
• 0 si censuré (c.-à-d. (T_i> C_i ))

Cependant, R les Surv La fonction acceptera également TRUE / FALSE (TRUE = événement) ou 1/2 (2 = événement).

dans le lung données, nous avons:

• statut: censure statut 1 = censuré, 2 = mort

Fonction de survie

La probabilité qu'un sujet survive au-delà d'une période donnée

[S

(S

En théorie, la fonction de survie est fluide; en pratique, nous observons les événements sur une échelle de temps discrète.

Probabilité de survie

• Probabilité de survie à un certain moment, (S
• Peut être estimé comme le nombre de patients en vie sans perte de suivi à ce moment, divisé par le nombre de patients en vie juste avant ce moment
• le Kaplan-Meier L'estimation de la probabilité de survie est le produit de ces probabilités conditionnelles jusqu'à cette date
• Au temps 0, la probabilité de survie est 1, c'est-à-dire (S (t_0) = 1 )

Création d'objets de survie

La méthode de Kaplan-Meier est le moyen le plus courant d'estimer les temps de survie et les probabilités. Il s'agit d'une approche non paramétrique qui se traduit par une fonction de pas, où il y a un pas vers le bas chaque fois qu'un événement se produit.

• le Surv fonction du survival package crée un objet de survie à utiliser comme réponse dans une formule de modèle. Il y aura une entrée pour chaque sujet qui est le temps de survie, qui est suivie d'un + si le sujet était censuré. Regardons les 10 premières observations:

Surv(lung$time, lung$status)[1:10]

##  [1]  306   455  1010+  210   883  1022+  310   361   218   166

Estimation des courbes de survie avec la méthode de Kaplan-Meier

• le survfit La fonction crée des courbes de survie basées sur une formule. Générons la courbe de survie globale pour toute la cohorte, affectons-la à l'objet f1et regardez le names de cet objet:

f1 <- survfit(Surv(time, status) ~ 1, data = lung)
names(f1)

##  [1] "n"         "time"      "n.risk"    "n.event"   "n.censor"  "surv"     
##  [7] "std.err"   "cumhaz"    "std.chaz"  "type"      "logse"     "conf.int" 
## [13] "conf.type" "lower"     "upper"     "call"

Quelques éléments clés de cette survfit Les objets qui seront utilisés pour créer des courbes de survie comprennent:

• time, qui contient les points de début et de fin de chaque intervalle de temps
• surv, qui contient la probabilité de survie correspondant à chaque time

Parcelle de Kaplan-Meier - base R

Maintenant, nous traçons le survfit objet dans la base R pour obtenir le tracé de Kaplan-Meier.

plot(survfit(Surv(time, status) ~ 1, data = lung), 
     xlab = "Days", 
     ylab = "Overall survival probability")

• Le tracé par défaut dans la base R montre la fonction de pas (ligne continue) avec les intervalles de confiance associés (lignes pointillées)
• Les lignes horizontales représentent la durée de survie pour l'intervalle
• Un intervalle se termine par un événement
• La hauteur des lignes verticales montre le changement de probabilité cumulative
• Les observations censurées, indiquées par des graduations, réduisent la survie cumulée entre les intervalles. (Remarque les graduations pour les patients censurés ne sont pas affichées par défaut, mais pourraient être ajoutées à l'aide de l'option mark.time = TRUE)

Parcelle de Kaplan-Meier - ggsurvplot

Alternativement, le ggsurvplot fonction du survminer le package est construit sur ggplot2, et peut être utilisé pour créer des tracés Kaplan-Meier. Vérifiez cheatsheet pour le survminer paquet.

ggsurvplot(
    fit = survfit(Surv(time, status) ~ 1, data = lung), 
    xlab = "Days", 
    ylab = "Overall survival probability")

• Le tracé par défaut utilisant ggsurvplot montre la fonction de pas (ligne continue) avec les bandes de confiance associées (zone ombrée).
• Les graduations des patients censurés sont affichées par défaut, obscurcissant quelque peu la ligne elle-même dans cet exemple, et pourraient être supprimées à l'aide de l'option censor = FALSE

Estimation (X)-survie d'années

Une quantité souvent intéressante dans une analyse de survie est la probabilité de survivre au-delà d'un certain nombre ((X)) d'années.

Par exemple, pour estimer la probabilité de survie à (1) année, utilisation summary avec le times argument (Remarque les time variable dans le lung les données sont en fait en jours, nous devons donc utiliser times = 365.25)

summary(survfit(Surv(time, status) ~ 1, data = lung), times = 365.25)

## Call: survfit(formula = Surv(time, status) ~ 1, data = lung)
## 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##   365     65     121    0.409  0.0358        0.345        0.486

Nous constatons que le (1)de probabilité de survie à 5 ans dans cette étude est de 41%.

Les limites inférieures et supérieures associées de l'intervalle de confiance à 95% sont également affichées.

(X)la survie sur une année et la courbe de survie

le (1)-la probabilité de survie annuelle est le point sur l'axe des y qui correspond à (1) année sur l'axe des x pour la courbe de survie.

(X)ans de survie est souvent estimée incorrectement

Que se passe-t-il si vous utilisez une estimation «naïve»?

121 des 228 patients sont décédés (1) année donc:

[Big(1 - frac{121}{228}Big) times 100 = 47%] - Vous obtenez un Incorrect estimation de la (1)probabilité de survie sur une année lorsque vous ignorez le fait que 42 patients ont été censurés avant (1) an.

• Rappelez-vous correct estimation de la (1)-la probabilité de survie à un an était de 41%.

Impact sur (X)ans de survie en ignorant la censure

• Imaginez deux études, chacune avec 228 sujets. Il y a 165 décès dans chaque étude. Pas de censure dans l'un (ligne orange), 63 patients censurés dans l'autre (ligne bleue)
• Ignorer la censure conduit à une surestimer de la probabilité de survie globale, car les sujets censurés ne fournissent que des informations pour partie du temps de suivi, puis tomber hors de l'ensemble de risque, réduisant ainsi la probabilité cumulée de survie

Comparaison des temps de survie entre les groupes

• Nous pouvons effectuer des tests de signification entre groupes en utilisant un test de log-rank
• Le test du log-rank pondère également les observations sur toute la durée du suivi et est le moyen le plus courant de comparer les temps de survie entre les groupes
• Il existe des versions qui pèsent plus lourdement le suivi précoce ou tardif qui pourraient être plus appropriées en fonction de la question de recherche (voir ?survdiff pour différentes options de test)

Nous obtenons la valeur de p de log-rank en utilisant le survdiff fonction. Par exemple, nous pouvons tester s'il y avait une différence de temps de survie selon le sexe dans le lung Les données

survdiff(Surv(time, status) ~ sex, data = lung)

## Call:
## survdiff(formula = Surv(time, status) ~ sex, data = lung)
## 
##         N Observed Expected (O-E)^2/E (O-E)^2/V
## sex=1 138      112     91.6      4.55      10.3
## sex=2  90       53     73.4      5.68      10.3
## 
##  Chisq= 10.3  on 1 degrees of freedom, p= 0.001

Le modèle de régression de Cox

Nous pouvons vouloir quantifier la taille d'un effet pour une seule variable, ou inclure plus d'une variable dans un modèle de régression pour tenir compte des effets de plusieurs variables.

Le modèle de régression de Cox est un modèle semi-paramétrique qui peut être utilisé pour ajuster des modèles de régression univariable et multivariable qui ont des résultats de survie.

[h(t|X_i) = h_0

(h

Quelques hypothèses clés du modèle:

• censure non informative
• risques proportionnels

Remarque: des modèles de régression paramétrique pour les résultats de survie sont également disponibles, mais ils ne seront pas abordés dans cette formation

Nous pouvons ajuster les modèles de régression pour les données de survie en utilisant le coxph fonction, qui prend une Surv objet sur le côté gauche et a une syntaxe standard pour les formules de régression dans R sur le côté droit.

coxph(Surv(time, status) ~ sex, data = lung)

## Call:
## coxph(formula = Surv(time, status) ~ sex, data = lung)
## 
##        coef exp(coef) se(coef)      z       p
## sex -0.5310    0.5880   0.1672 -3.176 0.00149
## 
## Likelihood ratio test=10.63  on 1 df, p=0.001111
## n= 228, number of events= 165

Formatage des résultats de régression de Cox

Nous pouvons voir une version ordonnée de la sortie en utilisant le tidy fonction du broom paquet:

broom::tidy(
  coxph(Surv(time, status) ~ sex, data = lung), 
  exp = TRUE
  ) %>% 
  kable()

Ou utiliser tbl_regression du gtsummary paquet

coxph(Surv(time, status) ~ sex, data = lung) %>% 
  gtsummary::tbl_regression(exp = TRUE) 

Rapports de risques

• La quantité d'intérêt d'un modèle de régression de Cox est un rapport de risque (HR). Le HR représente le rapport des dangers entre deux groupes à un moment donné.

• Le HR est interprété comme le taux instantané de survenue de l'événement d'intérêt chez ceux qui sont toujours à risque pour l'événement. C'est ne pas un risque, bien qu'il soit communément interprété comme tel.

• Si vous avez un paramètre de régression (bêta) (de la colonne estimate dans notre coxph) puis HR = ( exp ( beta) ).

• Une RH < 1 indicates reduced hazard of death whereas a HR > 1 indique un risque accru de décès.

• Ainsi, notre HR = 0,59 implique qu'environ 0,6 fois plus de femmes meurent que d'hommes, à un moment donné.

Exemple: réponse tumorale

Exemple: La survie globale est mesurée depuis le début du traitement et l'intérêt porte sur l'association entre la réponse complète au traitement et la survie.

• Anderson et al (JCO, 1983) ont expliqué pourquoi les méthodes traditionnelles telles que les tests de log-rank ou la régression de Cox sont biaisés en faveur des répondants dans ce scénario et ont proposé l'approche historique.
• le hypothèse nulle dans l'approche du point de repère, la survie à partir du point de repère ne dépend pas du statut de la réponse au point de repère.

Anderson, J., Cain, K. et Gelber, R. (1983). Analyse de la survie par réponse tumorale. Journal of Clinical Oncology: Journal officiel de l'American Society of Clinical Oncology, 1 (11), 710-9.

Autres exemples

D'autres covariables d'intérêt possibles pour la recherche sur le cancer qui peuvent ne pas être mesurées au départ comprennent:

• échec de la transplantation
• maladie du greffon contre l'hôte
• deuxième résection
• Thérapie adjuvante
• conformité
• événements indésirables

Exemple de données - l'ensemble de données BMT du package SemiCompRisks

Données sur 137 patients transplantés de moelle osseuse. Les variables d'intérêt comprennent:

• T1 temps (en jours) avant le décès ou le dernier suivi
• delta1 indicateur de décès; 1 mort, 0 vivant
• TA temps (en jours) de la maladie aiguë du greffon contre l'hôte
• deltaA indicateur de maladie aiguë du greffon contre l'hôte; 1-Maladie aiguë du greffon contre l'hôte développée, 0-Maladie aiguë du greffon contre l'hôte jamais développée

Chargeons les données à utiliser dans des exemples tout au long

data(BMT, package = "SemiCompRisks")

Méthode historique

1. Sélectionnez une heure fixe après la référence comme heure de référence. Remarque: cela doit être fait sur la base des informations cliniques, avant l'inspection des données
2. Sous-ensemble de la population pour ceux suivis au moins jusqu'à l'heure du point de repère. Remarque: signalez toujours le nombre exclu en raison d'un événement d'intérêt ou de censure avant l'heure du point de repère.
3. Calculez le suivi à partir du point de repère et appliquez les tests de log-rank traditionnels ou la régression de Cox

dans le BMT l'intérêt des données réside dans l'association entre la maladie aiguë du greffon contre l'hôte (aGVHD) et la survie. Mais aGVHD est évalué après la greffe, qui est notre temps de référence, ou le début du suivi.

Étape 1 Sélectionnez l'heure du point de repère

En règle générale, aGVHD survient dans les 90 premiers jours suivant la transplantation, nous utilisons donc un point de repère de 90 jours.

L'intérêt porte sur l'association entre la maladie aiguë du greffon contre l'hôte (aGVHD) et la survie. Mais aGVHD est évalué après la greffe, qui est notre temps de référence, ou le début du suivi.

Étape 2 Sous-ensemble de la population pour ceux suivis au moins jusqu'à l'heure du point de repère

lm_dat <- 
  BMT %>% 
  filter(T1 >= 90) 

Cela réduit notre taille d'échantillon de 137 à 122.

• Les 15 patients exclus sont décédés avant le repère de 90 jours

L'intérêt porte sur l'association entre la maladie aiguë du greffon contre l'hôte (aGVHD) et la survie. Mais aGVHD est évalué après la greffe, qui est notre temps de référence, ou le début du suivi.

Étape 3 Calculez le temps de suivi à partir du point de repère et appliquez les méthodes traditionnelles.

lm_dat <- 
  lm_dat %>% 
  mutate(
    lm_T1 = T1 - 90
    )

lm_fit <- survfit(Surv(lm_T1, delta1) ~ deltaA, data = lm_dat)

Exemple de repère de régression de Cox utilisant des données BMT

Dans la régression de Cox, vous pouvez utiliser le subset option dans coxph d'exclure les patients qui n'ont pas été suivis pendant la période historique

coxph(
  Surv(T1, delta1) ~ deltaA, 
  subset = T1 >= 90, 
  data = BMT
  ) %>% 
  gtsummary::tbl_regression(exp = TRUE)

Covariable dépendante du temps

Une alternative à une analyse historique est l'incorporation d'une covariable dépendante du temps. Cela peut être plus approprié lorsque

• la valeur d'une covariable change avec le temps
• il n'y a pas de moment historique évident
• l'utilisation d'un point de repère entraînerait de nombreuses exclusions

Configuration des données covariables en fonction du temps

Analyse des covariables dépendantes du temps dans R nécessite la configuration d'un ensemble de données spécial. Voir le document détaillé à ce sujet de l'auteur du survival paquet Utilisation de covariables dépendantes du temps et de coefficients dépendants du temps dans le modèle de Cox.

Il n'y avait pas de variable ID dans le BMT données, qui sont nécessaires pour créer le jeu de données spécial, alors créez-en un appelé my_id.

bmt <- rowid_to_column(BMT, "my_id")

Utilisez le tmerge fonctionner avec le event et tdc options de fonction pour créer l'ensemble de données spécial.

• tmerge crée un long ensemble de données avec plusieurs intervalles de temps pour les différentes valeurs de covariable pour chaque patient
• event crée le nouvel indicateur d'événement pour aller avec les intervalles de temps nouvellement créés
• tdc crée l'indicateur de covariable dépendant du temps pour aller avec les intervalles de temps nouvellement créés

td_dat <- 
  tmerge(
    data1 = bmt %>% select(my_id, T1, delta1), 
    data2 = bmt %>% select(my_id, T1, delta1, TA, deltaA), 
    id = my_id, 
    death = event(T1, delta1),
    agvhd = tdc(TA)
    )

Covariée dépendante du temps - exemple d'un seul patient

Pour voir ce que cela fait, regardons les données des 5 premiers patients individuels.

Les variables d'intérêt dans les données originales ressemblaient à

##   my_id   T1 delta1   TA deltaA
## 1     1 2081      0   67      1
## 2     2 1602      0 1602      0
## 3     3 1496      0 1496      0
## 4     4 1462      0   70      1
## 5     5 1433      0 1433      0

Le nouvel ensemble de données pour ces mêmes patients ressemble à

##   my_id   T1 delta1 id tstart tstop death agvhd
## 1     1 2081      0  1      0    67     0     0
## 2     1 2081      0  1     67  2081     0     1
## 3     2 1602      0  2      0  1602     0     0
## 4     3 1496      0  3      0  1496     0     0
## 5     4 1462      0  4      0    70     0     0
## 6     4 1462      0  4     70  1462     0     1
## 7     5 1433      0  5      0  1433     0     0

Covariée dépendante du temps - régression de Cox

Maintenant, nous pouvons analyser cette covariable dépendante du temps comme d'habitude en utilisant la régression de Cox avec coxph et une modification de notre utilisation de Surv d'inclure des arguments à la fois time et time2

coxph(
  Surv(time = tstart, time2 = tstop, event = death) ~ agvhd, 
  data = td_dat
  ) %>% 
  gtsummary::tbl_regression(exp = TRUE)

Résumé

Nous constatons que la maladie aiguë du greffon contre l'hôte n'est pas significativement associée au décès en utilisant une analyse historique ou une covariable dépendante du temps.

Souvent, on voudra utiliser l'analyse des points de repère pour la visualisation d'une seule covariable et la régression de Cox avec une covariable dépendante du temps pour la modélisation univariable et multivariable.

Paquets

Le paquet principal à utiliser dans les analyses des risques concurrents est

library(cmprsk)

Quels sont les risques concurrents?

Lorsque les sujets ont plusieurs événements possibles dans un cadre temporel

Exemples:

• récurrence
• décès par maladie
• décès d'autres causes
• réponse au traitement

Tout ou partie de ceux-ci (entre autres) peuvent être des événements possibles dans une étude donnée.

Donc quel est le problème?

La dépendance non observée entre les événements est le problème fondamental qui conduit à la nécessité d'une attention particulière.

Par exemple, on peut imaginer que les patients qui récidivent sont plus susceptibles de mourir, et donc les temps de récidive et les temps de décès ne seraient pas des événements indépendants.

Contexte des risques concurrents

Deux approches de l'analyse en présence de multiples résultats potentiels:

1. Risque lié à la cause d'un événement donné: il représente le taux par unité de temps de l'événement parmi ceux qui n'ont pas échoué à d'autres événements
2. Incidence cumulée d'un événement donné: cela représente le taux par unité de temps de l'événement ainsi que l'influence des événements concurrents

Chacune de ces approches ne peut éclairer qu'un aspect important des données tout en en masquant éventuellement d'autres, et l'approche choisie devrait dépendre de la question d'intérêt.

Un tas de notes et de références supplémentaires

• Lorsque les événements sont indépendants (presque jamais vrais), les dangers spécifiques à la cause ne sont pas biaisés
• Lorsque les événements sont dépendants, une variété de résultats peut être obtenue selon le réglage
• L'incidence cumulée en utilisant Kaplan-Meier est toujours> = incidence cumulée en utilisant des méthodes de risques concurrents, donc ne peut conduire qu'à une surestimation de l'incidence cumulée, le degré de surestimation dépend des taux d'événements et de la dépendance entre les événements
• Pour établir qu'une covariable agit effectivement sur l'événement d'intérêt, les risques spécifiques à la cause peuvent être préférés pour le traitement ou les tests d'effet marqueur pronostique
• Pour établir le bénéfice global, les risques de sous-distribution peuvent être préférés pour la construction de nomogrammes pronostiques ou la prise en compte des effets économiques sur la santé pour mieux comprendre l'influence du traitement et d'autres covariables à l'échelle absolue

Dignam JJ, Zhang Q, Kocherginsky M. L'utilisation et l'interprétation des modèles de régression des risques concurrents. Clin Cancer Res. 2012; 18 (8): 2301-8.

Kim HT. Incidence cumulée dans les données sur les risques concurrents et l'analyse de régression des risques concurrents. Clin Cancer Res. 15 janv.2007; 13 (2 pt 1): 559-65.

Satagopan JM, Ben-Porat L, Berwick M, Robson M, Kutler D, Auerbach AD. Une note sur les risques concurrents dans l'analyse des données de survie. Br J Cancer. 2004; 91 (7): 1229-35.

Austin, P. et Fine, J. (2017). Recommandations pratiques pour la communication des analyses de modèle Fine-Gray pour les données de risque concurrentes. Statistiques en médecine, 36 (27), 4391-4400.

Incidence cumulée des risques concurrents

• Estimation non paramétrique de l'incidence cumulée
• Estime l'incidence cumulative de l'événement d'intérêt
• À tout moment, la somme de l'incidence cumulée de chaque événement est égale à l'incidence cumulée totale de tout événement (ce n'est pas vrai dans le contexte spécifique à la cause)
• Le test de Gray est un test du chi carré modifié utilisé pour comparer 2 groupes ou plus

Exemple de données sur le mélanome du package MASS

Nous utilisons le Melanoma les données du MASS pour illustrer ces concepts. Il contient des variables:

• time temps de survie en jours, éventuellement censuré.
• status 1 est mort d'un mélanome, 2 vivants, 3 morts d'autres causes.
• sex 1 = mâle, 0 = femelle.
• age Age en années.
• year De fonctionnement.
• thickness épaisseur de la tumeur en mm.
• ulcer 1 = présence, 0 = absence.

data(Melanoma, package = "MASS")

Incidence cumulée dans les données sur le mélanome

Estimer l'incidence cumulée dans le contexte de risques concurrents à l'aide du cuminc fonction.

Remarque: dans le Melanoma données, les patients censurés sont codés comme (2 ) pour status, nous ne pouvons donc pas utiliser le cencode option par défaut de (0 )

cuminc(Melanoma$time, Melanoma$status, cencode = 2)

## Estimates and Variances:
## $est
##           1000       2000       3000      4000      5000
## 1 1 0.12745714 0.23013963 0.30962017 0.3387175 0.3387175
## 1 3 0.03426709 0.05045644 0.05811143 0.1059471 0.1059471
## 
## $var
##             1000         2000         3000        4000        5000
## 1 1 0.0005481186 0.0009001172 0.0013789328 0.001690760 0.001690760
## 1 3 0.0001628354 0.0002451319 0.0002998642 0.001040155 0.001040155

Tracer l'incidence cumulée - base R

Générer une base R tracer avec tous les paramètres par défaut.

ci_fit <- 
  cuminc(
    ftime = Melanoma$time, 
    fstatus = Melanoma$status, 
    cencode = 2
    )

plot(ci_fit)

Dans la légende:

• Le premier chiffre indique le groupe, dans ce cas il n'y a qu'une estimation globale donc il est (1) pour les deux
• Le deuxième chiffre indique le type d'événement, dans ce cas la ligne continue est (1) pour la mort d'un mélanome et la ligne pointillée est (3 ) pour la mort d'autres causes

Tracer l'incidence cumulée - ggcompetingrisks

Nous pouvons également tracer l’incidence cumulée en utilisant ggscompetingrisks fonction du survminer paquet.

Dans ce cas, nous obtenons un panneau étiqueté en fonction du groupe et une légende étiquetée événement, indiquant le type d'événement pour chaque ligne.

Remarques

• vous pouvez utiliser l'option multiple_panels = FALSE d'avoir tous les groupes tracés sur un seul panneau
• contrairement à la base R l'axe y ne passe pas à 1 par défaut, vous devez donc le modifier manuellement

ggcompetingrisks(ci_fit)

Comparer l'incidence cumulative entre les groupes

Dans cuminc Le test de Gray est utilisé pour les tests entre groupes.

Par exemple, comparez les Melanoma résultats selon ulcer, la présence ou l'absence d'ulcération. Les résultats des tests se trouvent dans Tests.

ci_ulcer <- 
  cuminc(
    ftime = Melanoma$time, 
    fstatus = Melanoma$status, 
    group = Melanoma$ulcer,
    cencode = 2
    )

ci_ulcer[["Tests"]]

##        stat           pv df
## 1 26.120719 3.207240e-07  1
## 3  0.158662 6.903913e-01  1

Tracer l'incidence cumulée selon le groupe - ggcompetingrisks

ggcompetingrisks(
  fit = ci_ulcer, 
  multiple_panels = FALSE,
  xlab = "Days",
  ylab = "Cumulative incidence of event",
  title = "Death by ulceration",
  ylim = c(0, 1)
)

Tracer l'incidence cumulée par groupe - manuellement

Remarque Je trouve personnellement le ggcompetingrisks fonction manque de personnalisation, en particulier par rapport à ggsurvplot. Je fais généralement mon propre tracé, en créant d'abord un ensemble de données bien rangé du cuminc ajuster les résultats, puis tracer les résultats. Voir le code source de cette présentation pour plus de détails sur le code sous-jacent.

Tracer un type d'événement unique - manuellement

Souvent, un seul des types d'événements sera intéressant, même si nous voulons toujours tenir compte de l'événement en compétition. Dans ce cas, l'événement d'intérêt peut être tracé seul. Encore une fois, je le fais manuellement en créant d'abord un ensemble de données bien rangé du cuminc ajuster les résultats, puis tracer les résultats. Voir le code source de cette présentation pour plus de détails sur le code sous-jacent.

Ajouter les nombres à risque table

Vous voudrez peut-être ajouter les nombres de table de risque à un graphique d'incidence cumulée, et il n'y a pas de moyen facile de le faire à ma connaissance. Voir le code source de cette présentation pour un exemple

1. Obtenez un complot de la base R, ggcompetingrisks, ou ggplot
2. Obtenez le nombre à risque table à partir d'un ggsurvplot en utilisant le survfit où tous les événements comptent comme un seul point de terminaison composite
   • Forcer les axes à avoir les mêmes limites et ruptures et titres
   • Assurez-vous que les couleurs / types de ligne correspondent aux étiquettes de groupe
   • Essayez d'obtenir la même taille de police
3. Combinez ensuite l'intrigue et la table des risques. Je utilise le plot_grid fonction du cowplot paquet pour cela
   • Je ne sais PAS comment modifier la taille du texte «Numéro à risque»…

Régression des risques concurrents

Deux approches:

1. Dangers spécifiques à la cause
   • taux d'occurrence instantané du type d'événement donné chez des sujets actuellement sans événement
   • estimée à l'aide de la régression de Cox (coxph fonction)
2. Dangers de sous-distribution
   • taux instantané d'occurrence du type d'événement donné chez des sujets qui n'ont pas encore vécu un événement de ce type
   • estimée par régression Fine-Gray (crr fonction)

Régression des risques concurrents dans les données sur le mélanome - approche du risque de sous-distribution

Imaginons que nous souhaitons étudier l'effet de l'âge et du sexe sur la mort due au mélanome, la mort due à d'autres causes étant un événement concurrent.

Remarques:

• crr nécessite la spécification de covariables comme matrice
• Si plusieurs événements présentent un intérêt, vous pouvez demander des résultats pour un autre événement en utilisant le failcode par défaut, les résultats sont renvoyés pour failcode = 1

shr_fit <- 
  crr(
    ftime = Melanoma$time,
    fstatus = Melanoma$status,
    cov1 = Melanoma[, c("sex", "age")],
    cencode = 2
    )

shr_fit

## convergence:  TRUE 
## coefficients:
##     sex     age 
## 0.58840 0.01259 
## standard errors:
## [1] 0.271800 0.009301
## two-sided p-values:
##  sex  age 
## 0.03 0.18

Dans l'exemple précédent, les deux sex et age ont été codés comme des variables numériques. le crr La fonction ne peut naturellement pas gérer les variables de caractères, et vous obtiendrez une erreur, donc si des variables de caractères sont présentes, nous devons créer des variables factices en utilisant model.matrix

# Create an example character variable
chardat <- 
  Melanoma %>% 
  mutate(
    sex_char = ifelse(sex == 0, "Male", "Female")
  )

# Create dummy variables with model.matrix
# The [, -1] removes the intercept
covs1 <- model.matrix(~ sex_char + age, data = chardat)[, -1]

# Now we can pass that to the cov1 argument, and it will work
crr(
  ftime = chardat$time,
  fstatus = chardat$status,
  cov1 = covs1,
  cencode = 2
  )

Formatage des résultats de crr

Sortie de crr n'est pas pris en charge par broom::tidy() ou gtsummary::tbl_regression() en ce moment. Comme alternative, essayez le (pas flexible, mais mieux que rien?) mvcrrres de mon ezfun paquet

ezfun::mvcrrres(shr_fit) %>% 
  kable()

Régression des risques concurrents dans les données sur le mélanome - approche du risque par cause

Censurer tous les sujets qui n'ont pas eu d'intérêt, dans ce cas la mort d'un mélanome, et utiliser coxph comme avant. Ainsi, les patients décédés d'autres causes sont désormais censurés pour l'approche du risque spécifique à la cause des risques concurrents.

Les résultats peuvent être formatés avec broom::tidy() ou gtsummary::tbl_regression()

chr_fit <- 
  coxph(
    Surv(time, ifelse(status == 1, 1, 0)) ~ sex + age, 
    data = Melanoma
    )

broom::tidy(chr_fit, exp = TRUE) %>% 
  kable()

gtsummary::tbl_regression(chr_fit, exp = TRUE)

Ce que nous avons couvert

• Les bases de l'analyse de survie, y compris la fonction de survie de Kaplan-Meier et la régression de Cox
• Analyse des points de repère et covariables dépendantes du temps
• Incidence cumulée et régression pour les analyses des risques concurrents

Ce qui reste?

Une variété de morceaux qui peuvent arriver et être utiles à savoir:

1. Évaluation de l'hypothèse des risques proportionnels
2. Faire un complot de survie en douceur basé sur (X)survie à 5 ans selon une covariable continue
3. Survie conditionnelle

Évaluation des risques proportionnels

Une hypothèse du modèle de régression des risques proportionnels de Cox est que les risques sont proportionnels à chaque instant dans le suivi. Comment pouvons-nous vérifier si nos données répondent à cette hypothèse?

Utilisez le cox.zph fonction du package de survie. Il en résulte deux choses principales:

1. Un test d'hypothèse pour savoir si l'effet de chaque covariable diffère selon le temps, et un test global de toutes les covariables à la fois.
   • Ceci est fait par testiung pour un effet d'interaction entre la covariable et le log (temps)
   • Une valeur p significative indique que l'hypothèse des risques proportionnels est violée
2. Tracés des résidus de Schoenfeld
   • La déviation d'une ligne à pente nulle prouve que l'hypothèse des risques proportionnels est violée

mv_fit <- coxph(Surv(time, status) ~ sex + age, data = lung)
cz <- cox.zph(mv_fit)
print(cz)

##        chisq df    p
## sex    2.608  1 0.11
## age    0.209  1 0.65
## GLOBAL 2.771  2 0.25

plot(cz)

Tracé de survie en douceur - quantile de survie

Parfois, vous souhaiterez visualiser une estimation de survie selon une variable continue. le sm.survival fonction du sm package vous permet de le faire pour un quantile de la distribution des données de survie. Le quantile par défaut est p = 0.5 pour la survie médiane.

library(sm)

sm.options(
  list(
    xlab = "Age (years)",
    ylab = "Time to death (years)")
  )

sm.survival(
  x = lung$age,
  y = lung$time,
  status = lung$status,
  h = sd(lung$age) / nrow(lung)^(-1/4)
  )

• Les x représentent des événements
• Les o représentent la censure
• La ligne est une estimation lissée de la survie médiane selon l'âge
  • Dans ce cas, trop lisse!

L'option h est le paramètre de lissage. Cela devrait être lié à l'écart type de la covariable continue, (X). Suggéré pour commencer ( frac {sd (x)} {n ^ {- 1/4}} ) puis réduisez de (1/2 ), (1/4 ), etc. pour obtenir une bonne quantité de lissage. Le tracé précédent était trop lisse, alors réduisons-le de (1/4 )

sm.survival(
  x = lung$age,
  y = lung$time,
  status = lung$status,
  h = (1/4) * sd(lung$age) / nrow(lung)^(-1/4)
  )

Survie conditionnelle

Parfois, il est intéressant de générer des estimations de survie parmi un groupe de patients qui ont déjà survécu pendant un certain temps.

[S(y|x) = frac{S(x + y)}{S(x)}]

• (y ): nombre d'années de survie supplémentaires d'intérêt
• (X): nombre d'années de survie d'un patient

Zabor, E., Gonen, M., Chapman, P., et Panageas, K. (2013). Pronostic dynamique utilisant des estimations conditionnelles de survie. Cancer, 119 (20), 3589-3592.

Estimations conditionnelles de survie

Les estimations sont faciles à générer avec les mathématiques de base par vous-même.

Alternativement, j'ai un package simple en développement appelé condsurv pour générer des estimations et des graphiques liés à la survie conditionnelle. Nous pouvons utiliser le conditional_surv_est pour obtenir des estimations et des intervalles de confiance à 95%. Condition de survie à 6 mois

remotes::install_github("zabore/condsurv")

library(condsurv)

fit1 <- survfit(Surv(time, status) ~ 1, data = lung)

prob_times <- seq(365.25, 182.625 * 5, 182.625)

purrr::map_df(
  prob_times, 
  ~conditional_surv_est(
    basekm = fit1, 
    t1 = 182.625, 
    t2 = .x) 
  ) %>% 
  mutate(months = round(prob_times / 30.4)) %>% 
  select(months, everything()) %>% 
  kable()

Tracés de survie conditionnelle

Nous pouvons également visualiser des données de survie conditionnelles basées sur différentes durées de survie. le condsurv::condKMggplot la fonction peut aider à cela.

cond_times <- seq(0, 182.625 * 4, 182.625)

gg_conditional_surv(
  basekm = fit1, 
  at = cond_times,
  main = "Conditional survival in lung data",
  xlab = "Days"
  ) +
  labs(color = "Conditional time")

Le tracé résultant a une courbe de survie pour chaque fois que nous nous conditionnons. Dans ce cas, la première ligne est la courbe de survie globale car elle conditionne au temps 0.

knit_exit()