Übung zu Konfidenzintervallen

Zum Datensatz

Für diese Übung brauchen wir den Datensatz einer klinischen Studie (weight-loss.csv). Die 22 Personen wurden randomisiert in eine der beiden Gruppen zugeteilt:

Kontrollgruppe (keine Intervention)
Diät + Exercise

Das Outcome ist der Gewichtsverlust in Kilogramm nach einem (wl1) und nach 3 Monaten (wl3). Der Datensatz liegt als csv-Datei auf Moodle.

Aufgabe

Lade den Datensatz herunter, speichere ihn ab und importiere ihn in R. Kontrolliere, ob der Datensatz korrekt eingelesen wurde (verschiedene Varianten möglich).

Lösung

Working directory definieren
library(rio) laden
Datensatz importieren (auf richtig Benennung achten)

library(rio)
df.wloss <- import("../data/weight-loss.csv")

Kontrolliere, ob der Datensatz korrekt eingelesen wurde (verschiedene Varianten möglich).

str(df.wloss)

## 'data.frame':    22 obs. of  4 variables:
##  $ ID   : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ group: chr  "Control" "Control" "Control" "Control" ...
##  $ wl1  : int  4 4 4 3 5 6 6 5 5 3 ...
##  $ wl3  : int  3 3 1 1 2 4 4 1 1 2 ...

summary(df.wloss)

##        ID           group                wl1             wl3       
##  Min.   : 1.00   Length:22          Min.   :3.000   Min.   :1.000  
##  1st Qu.: 6.25   Class :character   1st Qu.:4.000   1st Qu.:1.000  
##  Median :11.50   Mode  :character   Median :5.000   Median :2.000  
##  Mean   :16.95                      Mean   :5.273   Mean   :2.136  
##  3rd Qu.:28.75                      3rd Qu.:6.000   3rd Qu.:3.000  
##  Max.   :34.00                      Max.   :9.000   Max.   :4.000

Berechnung des Konfidenzintervalls

Berechnung der unteren (LL) und oberen (UL) Schranke eines Konfidenzintervalls:

\(LL = \bar{x} - z_{1-(p/2)} SE\)

\(UL = \bar{x} + z_{1-(p/2)} SE\)

Wobei \(z_{1-(p/2)}\) den gewünschten Prozentsatz des Intervalls definiert. Bei einem 95% CI wäre \(p\) also 0.05, was auf beide Seiten aufgeteilt wird (darum \(p/2\)).

Der Standardfehler berechnet sich wie folgt:

\(SE = \frac{s}{\sqrt{n}}\)

Übung 1

Aufgabe

Berechne das 95% Konfidenzintervall der Variable wl1, unabhängig von der Gruppe. Für diese Übung soll das Konfidenzintervall manuell unter Verwendung der z-Verteilung berechnet werden.

Lösung

Berechne das 95% Konfidenzintervall der Variable wl1, unabhängig von der Gruppe. Für diese Übung soll das Konfidenzintervall manuell berechnet werden.

mean.wl1 <- mean(df.wloss$wl1)
sd.wl1 <- sd(df.wloss$wl1)
n.wl1 <- length(df.wloss$wl1)
se.wl1 <- sd.wl1/sqrt(n.wl1)
z <- qnorm(1-(0.05/2))

LL <- mean.wl1 - z * se.wl1
UL <- mean.wl1 + z * se.wl1

Untere Schranke:

LL

## [1] 4.497039

Obere Schranke:

UL

## [1] 6.048416

Übung 2

Aufgabe

Berechne das gleiche Konfidenzintervall wie in Übung 1 erneut, dieses Mal aber mit der CI() Funktion aus library(Rmisc). Ist das mit der CI() Funktion berechnete Konfidenzintervall identisch mit dem manuell berechneten? Wann ja, warum? Wenn nicht, warum nicht?

Lösung

Berechne das gleiche Konfidenzintervall wie in Übung 1 erneut, dieses Mal aber mit der CI() Funktion library(Rmisc).

library(Rmisc)
CI(df.wloss$wl1)

##    upper     mean    lower 
## 6.095769 5.272727 4.449685

Das mit der CI() Funktion berechnete Konfidenzintervall ist etwas breiter als das manuell berechnete, weil diese Funktion die t-Verteilung und nicht die z-Verteilung braucht. Bei kleinen Stichproben ist die t-Verteilung genauer. Je grösser \(n\), desto näher ist die t-Verteilung an der z-Verteilung.

Hier ein Vergleich der z-Verteilung mit der t-Verteilungen mit verschiedenen Freiheitsgeraden (df). Man sieht, dass je kleiner das df, desto breiter die Verteilung:

range <- seq(-4, 4, length = 100)
             
my.df <- data.frame(x <- rep(range, 5),
                    df = rep(c("df_03", "df_10", "df_30", "df_50", "normal"), each = 100),
                    value = c(dt(range, 3),
                              dt(range, 10),
                              dt(range, 30),
                              dt(range, 50),
                              dnorm(range)))

library(ggplot2)
ggplot(my.df, aes(x = x, y = value, colour = df)) + 
  geom_line() +
  xlab("") +
  ylab("Density") +
  ggtitle("z and t distributions")

Übung 3

Aufgabe

Berechne das 77% Konfidenzintervall der Variable wl1 unter Verwendung der z-Verteilung.

Lösung

Berechne das 77% Konfidenzintervall der Variable wl1.

Manuell berechnet:

z <- qnorm(1-(0.23/2)) # p = 23%
LL <- mean.wl1 - z * se.wl1
UL <- mean.wl1 + z * se.wl1

Untere Schranke:

LL

## [1] 4.797665

Obere Schranke:

UL

## [1] 5.747789

Mit der CI-Funktion:

CI(df.wloss$wl1, ci = 0.77)

##    upper     mean    lower 
## 5.762010 5.272727 4.783444

Übung 4

Aufgabe

Berechne nun das 95% CI für die Mittelwertsdifferenz der beiden Gruppen zum ersten Zeitpunkt (manuelle Berechnung).

Lösung

Berechne nun das 95% CI für die Mittelwertsdifferenz der beiden Gruppen zum ersten Zeitpunkt (manuelle Berechnung).

mean.CG <- mean(df.wloss$wl1[df.wloss$group == "Control"])
sd.CG <- sd(df.wloss$wl1[df.wloss$group == "Control"])
n.CG <- length(df.wloss$wl1[df.wloss$group == "Control"])
se.CG <- sd.CG/sqrt(n.CG)

mean.IG <- mean(df.wloss$wl1[df.wloss$group == "DietEx"])
sd.IG <- sd(df.wloss$wl1[df.wloss$group == "DietEx"])
n.IG <- length(df.wloss$wl1[df.wloss$group == "DietEx"])
se.IG <- sd.IG/sqrt(n.IG)

mean.diff <- mean.CG - mean.IG
se.mean.diff <- sqrt(se.CG^2+se.IG^2)

z <- qnorm(1-0.05/2)

LL <- mean.diff - z * se.mean.diff
UL <- mean.diff + z * se.mean.diff

Untere Schranke:

LL

## [1] -3.205477

Obere Schranke

UL

## [1] -0.1945231

Anmerkung: Weil die Varianzen nicht homogen sind (und die Gruppen nicht genau gleich gross sind) würde die Berechnung mit der gepoolten Standardabweichung zu einem noch genaueren Resultat führen:

\[ s_{pooled}^2 = \frac{1}{df} * ((n - 1) * s_X^2 + (m-1) * s_Y^2) \] und

\[ SE_{X - Y} = \sqrt{\frac{s_{pooled}^2}{n} + \frac{s_{pooled}^2}{m}} = s_{pooled} * \sqrt{\frac{1}{n} + \frac{1}{m}}. \]

Übung 5

Aufgabe

Berechne die 95% CI’s der Variable wl1 für jede Gruppe sepparat. Brauche dafür die group.CI() Funktion aus library(Rmisc).

Lösung

Berechne die 95% CI’s der Variable wl1 für jede Gruppe sepparat. Brauche dafür die group.CI() Funktion aus library(Rmisc).

library(Rmisc)
group.CI(wl1 ~ group, data = df.wloss)

##     group wl1.upper wl1.mean wl1.lower
## 1 Control  5.135370      4.5  3.864630
## 2  DietEx  7.810215      6.2  4.589785

Übung 6

Aufgabe

Stelle die beiden oben berechneten 95% CI’s mit Hilfe des ggplot2-Packages graphisch dar. Es ist nicht die Meinung, dass du diesen Code auswendig kennst. Du kannst ihn einfach aus den Workshops kopieren und die relevanten Variablen anpassen. Nicht vergessen, das ggplot2 und das Hmisc Package zuerst zu laden.

Warum sind die Intervalle relativ breit?

Lösung

Stelle die beiden oben berechneten 95% CI’s mit Hilfe von library(ggplot2) graphisch dar. Es ist nicht die Meinung, dass du diesen Code auswendig kennst. Du kannst ihn einfach aus den Workshops kopieren und die relevanten Variablen anpassen. Nicht vergessen, das ggplot2 und das Hmisc Package zuerst zu laden.

library(ggplot2)
library(Hmisc)
ggplot(df.wloss, aes(x = group, y = wl1)) +
  stat_summary(fun = mean, geom = "point") +
  stat_summary(fun.data = mean_cl_normal, geom = "errorbar", width = 0.1) +
    xlab("")

Warum sind die Intervalle relativ breit? Weil das \(n\) relativ “tief” und demnach der SE “gross” ist (Wurzel-n-Gesetz)

Übung 7

Aufgabe

Nehmen wir auf Grund einer Umfrage an, die beobachtete mittlere Schuhgrösse von NBA Spielern (n = 50) ist 50.2 mit [95% CI = 49.7 bis 50.7]. Wie gross ist die Standardabweichung in dieser Stichprobe?

Lösung

Nehmen wir auf Grund einer Umfrage an, die mittlere Schuhgrösse von NBA Spielern (n = 50) ist 50.2 [95% CI = 49.7 bis 50.7]. Wie gross ist die Standardabweichung in dieser Stichprobe?

z <- qnorm(1-(0.05/2))
sd <- ((50.7 - 50.2) / z) * sqrt(50) #Umformung der Formel für SE
sd

## [1] 1.803877

Übung 8

Aufgabe

Du hast im Auftrag einer Stiftung anhand eines Fragebogens Daten zur Arbeitszufriedenheit von Gesundheitspersonal gesammelt. Das 95% CI für den durchschnittlichen Score (Skala von 0 bis 100) geht von 39 bis 75.

Wie ist ein CI überhaupt zu interpretieren?
Die Stiftung möchte eine genauere Schätzung des Mittelwerts. Was könntest du tun?

Lösung

Bei wiederholt gezogenen Stichproben würden 95 % der auf die gleiche Art berechneten Vertrauensintervalle dem wahren Mittelwert enthalten.
Das einfachste wäre, wenn man mehr Leute befragt. Das \(n\) hat den grössten Einfluss auf den Standardfehler (Wurzel-n-Gesetz), welcher wiederum die Breite des CI’s bestimmt. In diesem Fall müsste man wohl aber auch der Frage nachgehen, warum die Streuung der Befragten so gross ist.