Datensatz

Der Datensatz “melanom.csv” besteht aus 205 Zeilen und 7 Spalten. Es geht um das Überleben von Patient*innen nach einer Operation auf Grund eins malignen Melanoms.

Codebook:

  • V1: a numeric vector, row number (not meaningful)
  • no: a numeric vector, patient code
  • status: a numeric vector code, survival status; 1: dead from melanoma, 2: alive, 3: dead from other cause
  • days: a numeric vector, observation time
  • ulc: a numeric vector code, ulceration; 1: present, 2: absent
  • thick: a numeric vector, tumor thickness (mm)
  • sex: a numeric vector code; 1: female, 2: male

Source:
P.K. Andersen, Ø. Borgan, R.D. Gill, and N. Keiding (1991), Statistical Models Based on Counting Processes, Appendix 1, Springer-Verlag


Datenimport und Datenvorbereitung

Aufgaben

  1. Lade den Datensatz herunter, lade ihn in R und schau dir die Struktur des Datensatzes an.
  2. Ändere die Variablen sex und ulc in Faktoren um.
  3. Die Variable status ist als numerische Variable hinterlegt. Erstelle einen logischen Vektor (TRUE/FALSE) mit dem Namen died. Die Werte 1 und 3 sollen TRUE annehmen und der Wert 2 FALSE. Hier ein Beispiel, was mit einem logischen Vektor gemeint ist:
x <- c(1,2,1,2)
x
## [1] 1 2 1 2
x1 <- x == 1 # Alle 1er werden zu TRUE
x1
## [1]  TRUE FALSE  TRUE FALSE
x2 <- x != 1 # Alle 2er werden zu TRUE
x2
## [1] FALSE  TRUE FALSE  TRUE


Die Variablen x1 und x2 sind nun logische Vektoren:

class(c(x1, x2))
## [1] "logical"

Lösungen

  1. Lade den Datensatz herunter, lade ihn in R und schau dir die Struktur des Datensatzes an.
library(rio)
meldata  <- import("../../Data/melanom.csv")
str(meldata)
## 'data.frame':    205 obs. of  7 variables:
##  $ V1    : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ no    : int  789 13 97 16 21 469 685 7 932 944 ...
##  $ status: int  3 3 2 3 1 1 1 1 3 1 ...
##  $ days  : int  10 30 35 99 185 204 210 232 232 279 ...
##  $ ulc   : int  1 2 2 2 1 1 1 1 1 1 ...
##  $ thick : int  676 65 134 290 1208 484 516 1288 322 741 ...
##  $ sex   : int  2 2 2 1 2 2 2 2 1 1 ...

  1. Ändere die Variablen sex und ulc in Faktoren um.
meldata$sex <- factor(meldata$sex, levels = c(1,2), labels = c("female", "male"))
meldata$ulc <- factor(meldata$ulc, levels = c(1,2), labels = c("present", "absent"))
str(meldata)
## 'data.frame':    205 obs. of  7 variables:
##  $ V1    : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ no    : int  789 13 97 16 21 469 685 7 932 944 ...
##  $ status: int  3 3 2 3 1 1 1 1 3 1 ...
##  $ days  : int  10 30 35 99 185 204 210 232 232 279 ...
##  $ ulc   : Factor w/ 2 levels "present","absent": 1 2 2 2 1 1 1 1 1 1 ...
##  $ thick : int  676 65 134 290 1208 484 516 1288 322 741 ...
##  $ sex   : Factor w/ 2 levels "female","male": 2 2 2 1 2 2 2 2 1 1 ...

  1. Die Variable status ist als numerische Variable hinterlegt. Erstelle einen logischen Vektor (TRUE/FALSE) mit dem Namen died. Die Werte 1 und 3 sollen TRUE annehmen und der Wert 2 FALSE.
meldata$died <- meldata$status != 2 # somit wird aus 1 und 3 TRUE und aus 2 FALSE
meldata$died <- meldata$status == 1 | meldata$status == 3 # äquivalent, nur länger
meldata$died
##   [1]  TRUE  TRUE FALSE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE
##  [13]  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE
##  [25]  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE
##  [37]  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE  TRUE FALSE  TRUE FALSE
##  [49] FALSE FALSE  TRUE  TRUE FALSE  TRUE FALSE  TRUE FALSE  TRUE FALSE  TRUE
##  [61] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE  TRUE FALSE FALSE  TRUE
##  [73] FALSE FALSE  TRUE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
##  [85] FALSE FALSE FALSE FALSE  TRUE FALSE FALSE FALSE FALSE FALSE FALSE  TRUE
##  [97] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [109] FALSE FALSE  TRUE  TRUE FALSE  TRUE FALSE  TRUE FALSE  TRUE FALSE FALSE
## [121] FALSE FALSE FALSE FALSE FALSE  TRUE FALSE FALSE FALSE FALSE  TRUE FALSE
## [133] FALSE FALSE FALSE FALSE  TRUE FALSE FALSE FALSE FALSE FALSE  TRUE FALSE
## [145] FALSE FALSE FALSE FALSE  TRUE FALSE FALSE FALSE FALSE  TRUE FALSE FALSE
## [157] FALSE FALSE FALSE  TRUE FALSE  TRUE FALSE FALSE FALSE FALSE FALSE FALSE
## [169] FALSE FALSE  TRUE FALSE FALSE FALSE FALSE FALSE FALSE  TRUE FALSE FALSE
## [181] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [193] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [205] FALSE

Kontrolle, ob ein logischer Vektor erstellt wurde:

class(meldata$died)
## [1] "logical"

Manuelle Berechnung von Raten

Aufgaben

  1. Berechne die Anzahl follow-up Tage und die Anzahl Todesfälle für die Personen mit und ohne Ulzeration. Speichere die Daten in zwei neuen Objekten ab (fupdays und events)
  2. Berechne die Rate pro 1000 Personentage für die Personen mit und ohne Ulzeration.
  3. Berechne für beide Raten das 95% CI. Benütze dazu die cipoisson() Funktion aus dem survival Package.
  4. Berechne das Rate Ratio und interpretiere dieses.

Lösungen

  1. Berechne die Anzahl follow-up Tage und die Anzahl Todesfälle für die Personen mit und ohne Ulzeration. Speichere die Daten in zwei neuen Objekten ab (fupdays und events)

fupdays:

fupdays <- tapply(meldata$days, meldata$ulc, sum)
fupdays
## present  absent 
##  163603  277721


rates:

events <- tapply(meldata$died, meldata$ulc, sum)
events
## present  absent 
##      48      23

  1. Berechne die Rate pro 1000 Personentage für die Personen mit und ohne Ulzeration.
rates <- 1000*events/fupdays
rates
##    present     absent 
## 0.29339315 0.08281693

  1. Berechne für beide Raten das 95% CI. Benütze dazu die cipoisson() Funktion aus dem survival Package.
library(survival)
1000*cipoisson(events, fupdays)
##              lower     upper
## present 0.21632495 0.3889968
## absent  0.05249883 0.1242661

  1. Berechne das Rate Ratio und interpretiere dieses.
rates[1]/rates[2]
##  present 
## 3.542671

Es bedeutet, dass die Todesrate bei Personen mit Ulzeration um den Faktor 3.54 höher ist, als bei Personen ohne Ulzeration.


Raten via Poisson Regression

Aufgaben

  1. Berechne die gleichen Raten und das gleiche Rate Ratio wie oben, aber mit der glm() Funktion. Lass dir die summary deines Modells ausgeben. Definiere bei der Variable ulc die Ausprägung “absent” als Referenz.
  2. Berechne anhand des Outputs der Poisson Regression die Todesrate pro 1000 Personentage für die Gruppe mit Ulzeration. Dieser Wert sollte gleich sein, wie denjenigen welchen du oben schon manuell berechnet hast.
  3. Bestimme anhand des Outputs der Poisson Regression das Rate Ratio. Auch dieses sollte mit dem oben manuell berechneten Wert übereinstimmen.
  4. Rechne das Modell mit dem Geschlecht und der Dicke des Tumors (thick) als weitere Prädiktorvariablen und lass dir die summary ausgeben. Interpretiere den Output.

Lösungen

  1. Berechne die gleichen Raten und das gleiche Rate Ratio wie oben, aber mit der glm() Funktion. Lass dir die summary deines Modells ausgeben. Definiere bei der Variable ulc die Ausprägung “absent” als Referenz.
meldata$ulc<-relevel(meldata$ulc,ref="absent") # absent wird als Referenz gebraucht.
glm <- glm(died ~ ulc + offset(log(days)), data = meldata, family = poisson(link = "log"))
summary(glm)
## 
## Call:
## glm(formula = died ~ ulc + offset(log(days)), family = poisson(link = "log"), 
##     data = meldata)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.6235  -0.8107  -0.5672   0.9774   3.1623  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -9.3989     0.2085 -45.076  < 2e-16 ***
## ulcpresent    1.2649     0.2536   4.988 6.11e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 274.48  on 204  degrees of freedom
## Residual deviance: 247.34  on 203  degrees of freedom
## AIC: 393.34
## 
## Number of Fisher Scoring iterations: 6

  1. Berechne anhand des Outputs der Poisson Regression die Todesrate pro 1000 Personentage für die Gruppe mit Ulzeration. Dieser Wert sollte gleich sein, wie denjenigen welchen du oben schon manuell berechnet hast.
exp(glm$coefficients[1] + glm$coefficients[2])*1000
## (Intercept) 
##   0.2933932

  1. Bestimme anhand des Outputs der Poisson Regression das Rate Ratio. Auch dieses sollte mit dem oben manuell berechneten Wert übereinstimmen.
exp(glm$coefficients[2])
## ulcpresent 
##   3.542671

  1. Rechne das Modell mit dem Geschlecht und der Dicke des Tumors (thick) als weitere Prädiktorvariablen und lass dir die Summary ausgeben. Interprtiere den Output.
glm3 <- glm(died ~ ulc + sex + thick + offset(log(days)), 
            data = meldata, family = poisson(link = "log"))
summary(glm3)
## 
## Call:
## glm(formula = died ~ ulc + sex + thick + offset(log(days)), family = poisson(link = "log"), 
##     data = meldata)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.0901  -0.8270  -0.5218   0.8690   3.1175  
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -9.7746901  0.2443700 -40.000  < 2e-16 ***
## ulcpresent   0.9644074  0.2710134   3.559 0.000373 ***
## sexmale      0.4496266  0.2396143   1.876 0.060592 .  
## thick        0.0010366  0.0003484   2.976 0.002924 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 274.48  on 204  degrees of freedom
## Residual deviance: 235.44  on 201  degrees of freedom
## AIC: 385.44
## 
## Number of Fisher Scoring iterations: 6

Damit es etwas schöner aussieht, kann man den Output auch direkt runden. Hier am Beispiel der CIs.

round(cbind(RR = exp(glm3$coefficients), exp(confint(glm3))),3)
##                RR 2.5 % 97.5 %
## (Intercept) 0.000 0.000  0.000
## ulcpresent  2.623 1.558  4.529
## sexmale     1.568 0.979  2.514
## thick       1.001 1.000  1.002

Interpretation

Das Rate Ratio bzgl. Ulzeration verändert sich von 3.54 auf 2.62 mit Adjustieren nach Geschlecht und Dicke des Tumors. Ulzeration correliert demnach mit mindestens einer der anderen beiden Variablen. Fügt man die Variablen schrittweise dem Modell hinzu, sieht man, dass Ulzeration insbesondere mit der Dicke des Tumors korreliert.

glm2 <- glm(died ~ ulc + sex + offset(log(days)), data = meldata, family = poisson(link = "log"))
summary(glm2)
## 
## Call:
## glm(formula = died ~ ulc + sex + offset(log(days)), family = poisson(link = "log"), 
##     data = meldata)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.8542  -0.8437  -0.5376   0.8957   3.0577  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -9.5860     0.2319 -41.332  < 2e-16 ***
## ulcpresent    1.2044     0.2552   4.718 2.38e-06 ***
## sexmale       0.5133     0.2389   2.148   0.0317 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 274.48  on 204  degrees of freedom
## Residual deviance: 242.76  on 202  degrees of freedom
## AIC: 390.76
## 
## Number of Fisher Scoring iterations: 6
glm3 <- glm(died ~ ulc + sex + thick + offset(log(days)), data = meldata, family = poisson(link = "log"))
summary(glm3)
## 
## Call:
## glm(formula = died ~ ulc + sex + thick + offset(log(days)), family = poisson(link = "log"), 
##     data = meldata)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -2.0901  -0.8270  -0.5218   0.8690   3.1175  
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -9.7746901  0.2443700 -40.000  < 2e-16 ***
## ulcpresent   0.9644074  0.2710134   3.559 0.000373 ***
## sexmale      0.4496266  0.2396143   1.876 0.060592 .  
## thick        0.0010366  0.0003484   2.976 0.002924 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 274.48  on 204  degrees of freedom
## Residual deviance: 235.44  on 201  degrees of freedom
## AIC: 385.44
## 
## Number of Fisher Scoring iterations: 6
library(lmtest)
lrtest(glm, glm2, glm3)
## Likelihood ratio test
## 
## Model 1: died ~ ulc + offset(log(days))
## Model 2: died ~ ulc + sex + offset(log(days))
## Model 3: died ~ ulc + sex + thick + offset(log(days))
##   #Df  LogLik Df  Chisq Pr(>Chisq)   
## 1   2 -194.67                        
## 2   3 -192.38  1 4.5792   0.032363 * 
## 3   4 -188.72  1 7.3201   0.006819 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Wenn wir die Variablen einzeln hinzunehmen sehen wir, dass sich der Effekt von ulc hauptsächlich durch die Hinzunahme mit der Dicke des Tumors verändert. Demnach korrelieren diese beiden Variablen zu einem gewissen Grad. Der lrtest() zeigt uns, dass die Hinzuname beider Variablen das Modell signifikant verbessert.

Was sehen wir sonst noch?

Männer scheinen ein höheres Risiko zu haben (RR = 1.57, bei gleichem Status bezüglich Ulzeration und bei gleicher Tumordicke), jedoch schliesst das 95% CI die 1 knapp mit ein. Der Output zeigt ebenfalls, dass das Risiko, um zu sterben mit der Dicke des Tumors steigt: pro mm steigt das Risiko um den Faktor 1.0010 (bei gleichem Geschlecht und gleichem Status bzgl. Ulzeration).



Survival Curves

Um Überlebenskurven zu rechnen, brauchen wir mindestens zwei Informationen:

  1. Wie lange eine Person beobachtet wurde
  2. Ob das Outcome bei der Person auftrat

Mit der Funktion Surv() aus dem Package survival kann man ein Objekt erstellen, welches beide diese Informationen beinhaltet.

Aufgaben

  1. Erstelle ein Surv-Objekt melanom.surv für den Datensatz “melanom”. Lass dir die ersten zehn Beobachtungen anzeigen. Wie viele der ersten 10 Beobachtungen hatten das Outcome? Tipp: Surv() aus dem Paket survival.
  2. Berechne mit der survfit() Funktion die Überlebenskurve für die ganze Gruppe. Speichere die berechneten Werte im Objekt surv.all. Plotte die Überlebenskurve. Tipp: plot()
  3. Lass dir die Überlebenstabelle anzeigen. Wie viel Prozent der Personen leben noch nach 667 Tagen?
  4. Berechne die Überlebenskurve für die Personen mit und ohne Ulzeration separat. Stelle beide Kurven in einer Graphik dar. Lass dir auch die separaten Überlebenstabellen anzeigen.
  5. Unterscheiden sich die beiden Kurven statistisch signifikant? Führe den log-rank test durch, um die Differenz der beiden Überlebenskurven zu testen. Interpretiere das Resultat. Tipp: survdiff()

Lösungen

  1. Erstelle ein Surv-Objekt melanom.surv für den Datensatz “melanom”. Lass dir die ersten zehn Beobachtungen anzeigen. Wie viele der ersten 10 Beobachtungen hatten das Outcome?
library(survival)
melanom.surv <- Surv(meldata$days, meldata$died)
melanom.surv[1:10] # Alternative: head(melanom.surv, n = 10)
##  [1]  10   30   35+  99  185  204  210  232  232  279

Es sind 9, alle ausser die dritte Person, diese Beobachtung nennt man “zensiert”.


  1. Berechne mit der survfit() Funktion die Überlebenskurve für die ganze Gruppe. Speichere die berechneten Werte im Objekt surv.all. Plotte die Überlebenskurve.
surv.all <- survfit(melanom.surv ~ 1 )
plot(surv.all, main = "Survival Curve for the whole group", xlab = "Days", ylab = "Relative survival")

  1. Lass dir die Überlebenstabelle anzeigen. Wie viel Prozent der Personen leben noch nach 667 Tagen?
summary(surv.all)
## Call: survfit(formula = melanom.surv ~ 1)
## 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##    10    205       1    0.995 0.00487        0.986        1.000
##    30    204       1    0.990 0.00686        0.977        1.000
##    99    202       1    0.985 0.00840        0.969        1.000
##   185    201       1    0.980 0.00968        0.962        1.000
##   204    200       1    0.976 0.01081        0.955        0.997
##   210    199       1    0.971 0.01181        0.948        0.994
##   232    198       2    0.961 0.01357        0.935        0.988
##   279    196       1    0.956 0.01436        0.928        0.984
##   295    195       1    0.951 0.01510        0.922        0.981
##   355    194       1    0.946 0.01580        0.916        0.978
##   386    193       1    0.941 0.01646        0.910        0.974
##   426    192       1    0.936 0.01709        0.903        0.970
##   469    191       1    0.931 0.01769        0.897        0.967
##   493    190       1    0.927 0.01826        0.891        0.963
##   529    189       1    0.922 0.01881        0.885        0.959
##   621    188       1    0.917 0.01934        0.880        0.955
##   629    187       1    0.912 0.01985        0.874        0.952
##   659    186       1    0.907 0.02034        0.868        0.948
##   667    185       1    0.902 0.02081        0.862        0.944
##   718    184       1    0.897 0.02127        0.856        0.940
##   752    183       1    0.892 0.02171        0.851        0.936
##   779    182       1    0.887 0.02214        0.845        0.932
##   793    181       1    0.882 0.02255        0.839        0.928
##   817    180       1    0.877 0.02295        0.834        0.924
##   826    179       1    0.873 0.02334        0.828        0.920
##   833    178       1    0.868 0.02372        0.822        0.915
##   858    177       1    0.863 0.02409        0.817        0.911
##   869    176       1    0.858 0.02444        0.811        0.907
##   872    175       1    0.853 0.02479        0.806        0.903
##   967    174       1    0.848 0.02513        0.800        0.899
##   977    173       1    0.843 0.02546        0.795        0.895
##   982    172       1    0.838 0.02578        0.789        0.890
##  1041    171       1    0.833 0.02609        0.784        0.886
##  1055    170       1    0.828 0.02639        0.778        0.882
##  1062    169       1    0.824 0.02669        0.773        0.878
##  1075    168       1    0.819 0.02697        0.767        0.873
##  1156    167       1    0.814 0.02725        0.762        0.869
##  1228    166       1    0.809 0.02753        0.757        0.865
##  1252    165       1    0.804 0.02779        0.751        0.860
##  1271    164       1    0.799 0.02805        0.746        0.856
##  1312    163       1    0.794 0.02831        0.741        0.852
##  1427    162       1    0.789 0.02855        0.735        0.847
##  1435    161       1    0.784 0.02879        0.730        0.843
##  1506    159       1    0.779 0.02903        0.725        0.838
##  1516    155       1    0.774 0.02928        0.719        0.834
##  1525    154       1    0.769 0.02951        0.714        0.829
##  1548    152       1    0.764 0.02975        0.708        0.825
##  1560    150       1    0.759 0.02999        0.703        0.820
##  1584    148       1    0.754 0.03022        0.697        0.816
##  1621    146       1    0.749 0.03045        0.692        0.811
##  1667    137       1    0.743 0.03072        0.686        0.806
##  1690    134       1    0.738 0.03098        0.680        0.801
##  1726    131       1    0.732 0.03125        0.673        0.796
##  1860    117       1    0.726 0.03161        0.667        0.791
##  1933    110       1    0.719 0.03200        0.659        0.785
##  2061     95       1    0.712 0.03255        0.651        0.779
##  2062     94       1    0.704 0.03307        0.642        0.772
##  2085     92       1    0.697 0.03359        0.634        0.766
##  2103     90       1    0.689 0.03409        0.625        0.759
##  2108     88       1    0.681 0.03459        0.617        0.752
##  2256     80       1    0.673 0.03519        0.607        0.745
##  2388     75       1    0.664 0.03585        0.597        0.738
##  2467     69       1    0.654 0.03659        0.586        0.730
##  2565     63       1    0.644 0.03746        0.574        0.721
##  2782     57       1    0.632 0.03846        0.561        0.712
##  3042     52       1    0.620 0.03960        0.547        0.703
##  3154     46       1    0.607 0.04097        0.531        0.692
##  3182     44       1    0.593 0.04229        0.515        0.682
##  3338     35       1    0.576 0.04435        0.495        0.670
##  3458     28       1    0.555 0.04729        0.470        0.656

Es sind 90.2 %.


  1. Berechne die Überlebenskurve für die Personen mit und ohne Ulzeration separat. Stelle beide Kurven in einer Graphik dar.
surv.ulc <- survfit(melanom.surv ~ meldata$ulc)
plot(surv.ulc, main = "Survival Curve for the whole group", xlab = "Days", ylab = "Relative survival", conf.int = TRUE, col = c(2,3))
legend("bottomleft", lty = c(1,1), col = c(2,3), legend = c("Ulceration", "No ulceration"))

Hier noch die Life-tables pro Gruppe:

summary(surv.ulc)
## Call: survfit(formula = melanom.surv ~ meldata$ulc)
## 
##                 meldata$ulc=absent 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##    30    115       1    0.991 0.00866        0.974        1.000
##    99    113       1    0.983 0.01224        0.959        1.000
##   817    112       1    0.974 0.01495        0.945        1.000
##   858    111       1    0.965 0.01720        0.932        0.999
##   869    110       1    0.956 0.01915        0.919        0.994
##   872    109       1    0.947 0.02089        0.907        0.989
##  1041    108       1    0.939 0.02246        0.896        0.984
##  1427    107       1    0.930 0.02390        0.884        0.978
##  1435    106       1    0.921 0.02523        0.873        0.972
##  1548    101       1    0.912 0.02658        0.861        0.966
##  1560     99       1    0.903 0.02786        0.850        0.959
##  1584     97       1    0.893 0.02909        0.838        0.952
##  1690     85       1    0.883 0.03059        0.825        0.945
##  1933     71       1    0.871 0.03259        0.809        0.937
##  2061     62       1    0.856 0.03496        0.791        0.928
##  2062     61       1    0.842 0.03710        0.773        0.918
##  2085     60       1    0.828 0.03904        0.755        0.909
##  2103     58       1    0.814 0.04090        0.738        0.898
##  2388     49       1    0.798 0.04331        0.717        0.887
##  2782     41       1    0.778 0.04642        0.692        0.875
##  3154     31       1    0.753 0.05126        0.659        0.860
##  3182     29       1    0.727 0.05568        0.626        0.845
##  3458     19       1    0.689 0.06457        0.573        0.828
## 
##                 meldata$ulc=present 
##  time n.risk n.event survival std.err lower 95% CI upper 95% CI
##    10     90       1    0.989  0.0110        0.967        1.000
##   185     89       1    0.978  0.0155        0.948        1.000
##   204     88       1    0.967  0.0189        0.930        1.000
##   210     87       1    0.956  0.0217        0.914        0.999
##   232     86       2    0.933  0.0263        0.883        0.986
##   279     84       1    0.922  0.0282        0.869        0.979
##   295     83       1    0.911  0.0300        0.854        0.972
##   355     82       1    0.900  0.0316        0.840        0.964
##   386     81       1    0.889  0.0331        0.826        0.956
##   426     80       1    0.878  0.0345        0.813        0.948
##   469     79       1    0.867  0.0358        0.799        0.940
##   493     78       1    0.856  0.0371        0.786        0.931
##   529     77       1    0.844  0.0382        0.773        0.923
##   621     76       1    0.833  0.0393        0.760        0.914
##   629     75       1    0.822  0.0403        0.747        0.905
##   659     74       1    0.811  0.0413        0.734        0.896
##   667     73       1    0.800  0.0422        0.721        0.887
##   718     72       1    0.789  0.0430        0.709        0.878
##   752     71       1    0.778  0.0438        0.696        0.869
##   779     70       1    0.767  0.0446        0.684        0.859
##   793     69       1    0.756  0.0453        0.672        0.850
##   826     68       1    0.744  0.0460        0.660        0.840
##   833     67       1    0.733  0.0466        0.647        0.831
##   967     66       1    0.722  0.0472        0.635        0.821
##   977     65       1    0.711  0.0478        0.623        0.811
##   982     64       1    0.700  0.0483        0.611        0.801
##  1055     63       1    0.689  0.0488        0.600        0.791
##  1062     62       1    0.678  0.0493        0.588        0.782
##  1075     61       1    0.667  0.0497        0.576        0.772
##  1156     60       1    0.656  0.0501        0.564        0.761
##  1228     59       1    0.644  0.0505        0.553        0.751
##  1252     58       1    0.633  0.0508        0.541        0.741
##  1271     57       1    0.622  0.0511        0.530        0.731
##  1312     56       1    0.611  0.0514        0.518        0.721
##  1506     55       1    0.600  0.0516        0.507        0.710
##  1516     53       1    0.589  0.0519        0.495        0.700
##  1525     52       1    0.577  0.0521        0.484        0.689
##  1621     51       1    0.566  0.0523        0.472        0.678
##  1667     50       1    0.555  0.0525        0.461        0.668
##  1726     49       1    0.543  0.0526        0.449        0.657
##  1860     42       1    0.530  0.0529        0.436        0.645
##  2108     32       1    0.514  0.0538        0.419        0.631
##  2256     29       1    0.496  0.0548        0.400        0.616
##  2467     24       1    0.475  0.0563        0.377        0.600
##  2565     19       1    0.450  0.0586        0.349        0.581
##  3042     16       1    0.422  0.0613        0.318        0.561
##  3338     12       1    0.387  0.0655        0.278        0.539

  1. Unterscheiden sich die beiden Kurven statistisch signifikant? Führe den log-rank test durch, um die Differenz der beiden Überlebenskurven zu testen. Interpretiere das Resultat.
survdiff(melanom.surv ~ meldata$ulc)
## Call:
## survdiff(formula = melanom.surv ~ meldata$ulc)
## 
##                       N Observed Expected (O-E)^2/E (O-E)^2/V
## meldata$ulc=absent  115       23     44.5      10.4      27.9
## meldata$ulc=present  90       48     26.5      17.3      27.9
## 
##  Chisq= 27.9  on 1 degrees of freedom, p= 1e-07

Der p-Wert ist sehr klein. Sollte \(H_0\) wahr sein, dass sich die beiden Kurven nicht unterscheiden, dann ist die Wahrscheinlichkeit, um zufällig den beobachteten Unterschied oder einen extremeren zu erhalten, 0.00001%. Es liegt also starke Evidenz gegen \(H_0\) vor.