Zu den Daten

Die folgende 2x2 Tabelle zeigt die Daten einer Studie, welche durchgeführt wurde, um die Wirksamkeit von Aspirin zur Prävention eines frühzeitigen Todes nach einem Herzinfarkt zu evaluieren. Hier findest du weitere Informationen zum den Daten.

##          outcome
## treatment survived death
##   placebo      354    52
##   aspirin      725    85

Du kannst dir diese Tabelle in R selber “nachbauen”:

# Matrix anhand der gegebenen Werte erstellen
m <- matrix(c(354, 52, 725, 85), byrow = TRUE, ncol = 2)
# Reihen und Spalten benenen
rownames(m) <- c("placebo", "aspirin")
colnames(m) <- c("survived", "died")
# aus einer Matrix eine Tabelle machen
tab <- as.table(m)
tab
##         survived died
## placebo      354   52
## aspirin      725   85

Manchmal möchte man auf Grundlage einer 2x2 Tabelle die nicht aggregierten Daten zurückgewinnen (z.B. um Grafiken zu erstellen). Dazu gibt es eine praktische Funktion aus dem Paket rstatix.

library(rstatix)
d.aspirin <- counts_to_cases(tab)
# Variablen des erstellen DF's benennen
names(d.aspirin) <- c("gruppe", "outcome")
head(d.aspirin)
##    gruppe  outcome
## 1 placebo survived
## 2 placebo survived
## 3 placebo survived
## 4 placebo survived
## 5 placebo survived
## 6 placebo survived

Proportionen pro Gruppe

Aufgaben

Verwende für die Aufgaben die oben gezeigte Tabelle.

  1. Berechne anhand der Normalapproximation für beide Gruppen die Proportion und ein dazugehörigens 95% CI der Personen, welche gestorben sind.

  2. Vergleiche die manuellen Berechnungen mit dem Resultat der in R implementierten Funktion prop.test(). Interpretiere!

Lösungen

Wir basieren die Berechnungen auf den folgenden Angaben:

tab
##         survived died
## placebo      354   52
## aspirin      725   85

Anhand der Tabelle berechnen wir zuerst die Proportion, den Standardfehler und schliesslich das 95% CI. Beachten Sie dass

\[ SE_{\hat{\pi}} = \sqrt{\frac{\hat{\pi}(1-\hat{\pi})}{n}}. \]

cases <- tab[, 2] # Todesfälle pro Gruppe
n <- tab[, 1] + tab[, 2] # n pro Gruppe
prop <- cases/n # Proportion pro Gruppe
se.prop <- sqrt((prop * (1-prop))/n) # Standardfehler pro Gruppe
CI.lower <- prop-1.96*se.prop
CI.upper <- prop+1.96*se.prop
cbind(CI.lower, CI.upper)
##           CI.lower  CI.upper
## placebo 0.09557231 0.1605853
## aspirin 0.08383222 0.1260443

Mit prop.test() erhalten wir die folgenden Resultate.

Placebo:

prop.test(x = cases["placebo"], n = n["placebo"])
## 
##  1-sample proportions test with continuity correction
## 
## data:  cases["placebo"] out of n["placebo"], null probability 0.5
## X-squared = 223.16, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
##  0.09793556 0.16545553
## sample estimates:
##         p 
## 0.1280788

Aspirin:

prop.test(x = cases["aspirin"], n = n["aspirin"])
## 
##  1-sample proportions test with continuity correction
## 
## data:  cases["aspirin"] out of n["aspirin"], null probability 0.5
## X-squared = 504.1, df = 1, p-value < 2.2e-16
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
##  0.08510476 0.12860969
## sample estimates:
##         p 
## 0.1049383

Die Proportionen betragen 12.81% (Placebo), bzw. 10.49% (Aspirin).

Weil das n relativ gross ist, führt die Berechnung unter der Normalapproximation fast zum gleichen Ergebnis wie die prop.test() Funktion. Wie sehen, dass sich die 95% CIs der beiden Gruppen überlappen. Aus diesem Grund wissen wir bereits jetzt, dass auf dem 5%-Niveau keine Evidenz vorliegt, dass sich die beiden Proportionen unterscheiden.

Standardmässig testet prop.test(), ob sich die Proportion von 0.5 unterscheidet. In diesem Fall interessiert uns dieser Test (und folglich der P-Wert nicht). Mit dem Argument p = … könnte man \(H_0\) individuell anpassen. Es kommt bei Proportionen jedoch häufig vor, dass man primär an einem 95% CI interessiert ist.

Proportionsdifferenzen

Aufgaben

Es wird weiterhin mit der 2x2 Tabelle von oben gerechnet.

  1. Berechnen Sie die Proportionsdifferenz zwischen der Placebo- und der Aspiringruppe und berechnen Sie mit Hilfe der Normalapproximation ein 95% CI und einen P-Wert für diese Differenz.

  2. Wiederholen Sie die Berechnungen mit Hilfe der prop.test() Funktion. Vergleichen und Interpretieren Sie.

Lösungen

Wir berechnen die Differenz und den Standardfehler für diese Differenz:\[ SE_{\hat{\pi_1}-\hat{\pi_2}} = \sqrt{SE_{\hat{\pi_1}}^2 + SE_{\hat{\pi_2}}^2} \]

Ob Sie Placebo-Aspirin oder Aspirin-Placebo rechnen, spielt keine Rolle, solange die Interpretation sinngemäss erfolgt.

prop.diff <- prop["placebo"] - prop["aspirin"]
se.prop.diff <- sqrt(se.prop["placebo"]^2 + se.prop["aspirin"]^2)
CI.diff.lower <- prop.diff-1.96*se.prop.diff
CI.diff.upper <- prop.diff+1.96*se.prop.diff
cbind(CI.diff.lower, CI.diff.upper)
##         CI.diff.lower CI.diff.upper
## placebo   -0.01561689    0.06189798

Den P-Wert berechnen wir via z-Transformation:

z <- prop.diff/se.prop.diff
p <- 2*pnorm(-abs(z))
p
##   placebo 
## 0.2419047

Mit der prop.test() Funktion sind die Ergebniss sehr ähnlich:

prop.test(x = cases, n = n)
## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  cases out of n
## X-squared = 1.2264, df = 1, p-value = 0.2681
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.01746498  0.06374608
## sample estimates:
##    prop 1    prop 2 
## 0.1280788 0.1049383

Die Proportion der Personen, welche gestorben sind, ist in der Aspiringruppe 2.3% tiefer. Dies Entspricht der Risk Difference (RD). Wir vertrauen zu 95% darauf, dass die wahre RD im Intervall [-0.017; 0.064] liegt. Das 95% CI sagt uns, dass die RD sowohl negativ wie positiv sein kann. Die Evidenz gegen \(H_0\) ist also schwach. Dies bestätigt der hohe P-Wert: Die Wahrscheinlichkeit für einen Unterschied von 2.3% (oder einen extremeren) unter der Annahme, dass \(H_0\) stimmt, ist 26.81%.