Survival Analysis (Modul: Lebensdaueranalyse) R OLAND R AU Universität Rostock, Sommersemester 2015 16. Juni 2015 © Roland Rau Survival Analysis 1 / 18 Teams, Themen, & Termine Zwei Termine mit jeweils 3 “slots” stehen zur Auswahl: 07. Juli & 14. Juli # 1 2 3 4 5 6 Thema Versicherungsstatus & Sterblichkeit SES & Sterblichkeit BMI & Sterblichkeit Rauchersterblichkeit Kosten für med. Hilfe & Sterblichkeit Self rated health & Sterblichkeit © Roland Rau Team Grützmacher, Smilde-Becker Lepzien, Schulz Köpke, Liewert, Sauerberg Dronia, Eckardt Epperlein Wenau Survival Analysis Termin / “Slot” 1/2 2/1 1/1 2/3 1/3 2/2 2 / 18 Proportional Hazards Modelle Grundsätzliche Form: hi (t|X) = h0 (t)eβ1 x1,i +β2 x2,i +... Trennung in Baseline Hazard (h0 (t)) und Effekt/e (β1 , β2 , . . .) der Kovariaten (x1 , x2 , . . .). Falls h0 (t) parametrisch spezifiziert ⇒ parametrische Survival-Regression Falls h0 (t) ⇒ Cox-Modell © Roland Rau Survival Analysis 3 / 18 Zur Erinnerung I © Roland Rau Survival Analysis 4 / 18 Zur Erinnerung II wolpertinger <- read.table("wolpertinger-geschlecht-2015.txt", header=TRUE) geschlecht <- ifelse(wolpertinger$sex==1, 1, 0) loglike <- function(theta, t, covariate,delta) { lambda <- theta[1] beta <- theta[2] return(-(sum( (log(lambda) + beta*covariate)*delta) + sum(-lambda * exp(beta * covariate) * t))) } ergebnis <- optim(par=c(0.5, 1.2), fn=loglike, t=wolpertinger$zeit, delta=wolpertinger$cens, covariate=geschlecht) ## ## ## ## Warning: Warning: Warning: Warning: NaNs NaNs NaNs NaNs produced produced produced produced ergebnis ## ## ## ## ## ## ## ## ## ## ## ## ## $par [1] 0.1584 0.2470 $value [1] 929.5 $counts function gradient 71 NA $convergence [1] 0 © Roland Rau Survival Analysis 5 / 18 Einfacher: I library(survival) s.obj <- Surv(wolpertinger$zeit, wolpertinger$cens) survreg(s.obj ~ geschlecht, dist="exponential") ## ## ## ## ## ## ## ## ## ## ## ## Call: survreg(formula = s.obj ~ geschlecht, dist = "exponential") Coefficients: (Intercept) geschlecht 1.8423 -0.2469 Scale fixed at 1 Loglik(model)= -929.5 Loglik(intercept only)= -932.1 Chisq= 5.22 on 1 degrees of freedom, p= 0.022 n= 1000 © Roland Rau Survival Analysis 6 / 18 Zur Erinnerung: Cox-Modell I coxmodel <- coxph(s.obj ~ geschlecht) coxmodel ## ## ## ## ## ## ## ## Call: coxph(formula = s.obj ~ geschlecht) coef exp(coef) se(coef) z p geschlecht 0.249 1.28 0.109 2.3 0.022 Likelihood ratio test=5.29 on 1 df, p=0.0214 © Roland Rau n= 1000, number of events= 343 Survival Analysis 7 / 18 Numerisches Beispiel: Wolpertinger (reloaded) Wie interpretieren wir nun diese Ergebnisse? round(coef(coxmodel),4) ## geschlecht ## 0.2492 Da wir die Variable x so gewählt haben, dass für das weibliche Geschlecht x = 0 gilt, so erkennen wir: βx h(t) = h0 e β0 = h0 e 0 = h0 e = h0 dass der Baseline-Hazard λ =0.247 dem Hazard der Frauen entspricht. Dieser wird bei Männern um den Faktor βx h(t) = h0 e β1 = h0 e 0.2492 = h0 e = h0 × 1.282999 erhöht. Das Hazard Ratio beträgt damit: HR = β x h0 (t)e 1 ♂ h0 β x (t)e 1 ♀ = β x h0 (t)e 1 ♂ h0 (t) β x = e 1 ♂ = 1.282999 Das bedeutet, dass in unserem Modell Männer in jeder Altersstufe ein in etwa 30% höheres Sterberisiko als Frauen besitzen (28.2999%). © Roland Rau Survival Analysis 8 / 18 Das Cox-Modell mit unseren Daten I Zuerst müssen wir natürlich wieder unsere Daten einlesen und leicht überarbeiten: newdata <- read.table("smalldataset2015.txt", header=TRUE, sep=",") newdata2 <- subset(newdata, status %in% c(0,1)) newdata3 <- subset(newdata2, age %in% 18:84) newdata4 <- subset(newdata3, yob %in% 1914:1980) eintrittsalter <- newdata4$age delta <- newdata4$status austrittsalter <- ifelse(newdata4$status == 0, 2007 - newdata4$yob, # also falls rechtszensiert newdata4$year - newdata4$yob+0.5) # falls Ereignis surv.obj <- Surv(eintrittsalter, austrittsalter, delta) © Roland Rau Survival Analysis 9 / 18 Ein erstes Cox-Modell mit “unseren” Daten Wir berechnen nun mit unseren Daten ein Modell analog zu den Wolpertingern: coxmodel1 <- coxph(surv.obj ~ as.factor(newdata4$sex)) coxmodel1 ## ## ## ## ## ## ## ## Call: coxph(formula = surv.obj ~ as.factor(newdata4$sex)) coef exp(coef) se(coef) z p as.factor(newdata4$sex)2 -0.436 0.646 0.0356 -12.3 0 Likelihood ratio test=149 on 1 df, p=0 n= 33364, number of events= 3234 Wir sehen also, dass das Geschlecht 2 (=Frauen) ein niedrigeres Sterberisiko hat als die Referenzkategorie Männer (β = −0.4362). Dies können wir als relatives Risiko interpretieren in der Form von β = e−0.4362 = 0.6465 Wenn die Annahmen unseres Modells also stimmen, so haben Frauen in jeder Altersstufe ein 1 − 0.6465=0.3535 =35.35% niedrigeres Sterberisiko. © Roland Rau Survival Analysis 10 / 18 Ein zweites Cox-Modell mit “unseren” Daten Welchen Einfluss hat das Rauchen auf die Sterblichkeit? Wir operationalisieren diese Frage mittels der Variable eversmoke. table(newdata4$eversmoke) ## ## 1 2 ## 15898 17325 7 86 8 8 9 47 Wir können uns nun entscheiden, ob wir die Kategorien 7–9 ausschliessen oder aber als Residualkategorie zusammenfassen. newdata4$raucher <- "Residual" newdata4$raucher[newdata4$eversmoke==1] <- "Raucher" newdata4$raucher[newdata4$eversmoke==2] <- "Nichtraucher" table(newdata4$raucher) ## ## Nichtraucher ## 17325 Raucher 15898 Residual 141 Es bietet sich zudem an, die Variable explizit als kategoriale Variable zu definieren mit einer Referenzgruppe. newdata4$raucher <- as.factor(newdata4$raucher) newdata4$raucher <- relevel(newdata4$raucher, ref="Nichtraucher") © Roland Rau Survival Analysis 11 / 18 Ein zweites Cox-Modell mit “unseren” Daten coxmodel2 <- coxph(surv.obj ~ as.factor(newdata4$sex) + newdata4$raucher) coxmodel2 ## ## ## ## ## ## ## ## ## ## Call: coxph(formula = surv.obj ~ as.factor(newdata4$sex) + newdata4$raucher) coef exp(coef) se(coef) z p as.factor(newdata4$sex)2 -0.293 0.746 0.0366 -8.012 1.1e-15 newdata4$raucherRaucher 0.571 1.769 0.0382 14.932 0.0e+00 newdata4$raucherResidual -0.023 0.977 0.3030 -0.076 9.4e-01 Likelihood ratio test=380 on 3 df, p=0 © Roland Rau n= 33364, number of events= 3234 Survival Analysis 12 / 18 Kleine Aufgabe Wie hoch ist mit den gegebenen Schätzwerten der vorherigen Folie das relative Sterberisiko von rauchenden Frauen im Vergleich zu nichtrauchenden Männern? von rauchenden Männern im Vergleich zu nichtrauchenden Frauen? HR = h(t)Frau, raucht ho (t)eβ1 x1,Frau +β2 x2,raucht 1.3205 e−0.293·1+0.571·1 = ≈ ≈ 1.3205 = β x +β x 1 1,Mann 2 2,rauchtnicht h(t)Mann, raucht nicht e0·0+0·0 1 ho (t)e HR = h(t)Mann, raucht ho (t)eβ1 x1,Mann +β2 x2,raucht e0·0+0.571·1 1.77004 = −0.293·1+0·0 ≈ = ≈ 2.3726 β x +β x h(t)Frau, raucht nicht e 0.74602 ho (t)e 1 1,Frau 2 2,rauchtnicht © Roland Rau Survival Analysis 13 / 18 Wie können wir graphisch diese Ergebnisse darstellen? © Roland Rau Survival Analysis 14 / 18 Nächste Veranstaltung: Das Cox-Modell mit unseren Daten Wie überprüfe ich die Proportionalitätsannahme? Was tun, wenn die Proportionalitätsannahme verletzt ist? Wie plotte ich Ergebnisse in R? © Roland Rau Survival Analysis 15 / 18 Vielen Dank für Ihre Aufmerksamkeit! © Roland Rau Survival Analysis 16 / 18 Lizenz This open-access work is published under the terms of the Creative Commons Attribution NonCommercial License 2.0 Germany, which permits use, reproduction & distribution in any medium for non-commercial purposes, provided the original author(s) and source are given credit. Für ausführlichere Informationen: http://creativecommons.org/licenses/by-nc/2.0/de/ (Deutsch) http://creativecommons.org/licenses/by-nc/2.0/de/deed.en (English) © Roland Rau Survival Analysis 17 / 18 Kontakt Universität Rostock Institut für Soziologie und Demographie Lehrstuhl für Demographie Ulmenstr. 69 18057 Rostock Germany Tel.: +49-381-498 4044 Fax.: +49-381-498 4395 Email: [email protected] Sprechstunde im Sommersemester 2015: Mittwochs, 09:00–10:00 (und nach Vereinbarung) © Roland Rau Survival Analysis 18 / 18
© Copyright 2024 ExpyDoc