Übung zu Konfidenzintervallen

Zum Datensatz

Für diese Übung brauchen wir den Datensatz einer klinischen Studie. 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 “weightloss” liegt als csv-Datei auf Moodle.

Aufgabe

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

Lösung

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

library(rio)
wloss <- import("WeightLoss.csv")

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

str(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(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 gwü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 = sd / \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 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(wloss$wl1)
sd_wl1 <- sd(wloss$wl1)
n_wl1 <- length(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 dem Rmisc-Package. Ist das mit der CI() Funktion berechnete Konfidenzintervall identisch? 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 aus dem Rmisc-Package.

library(Rmisc)
CI(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 das 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 (nicht sehr sinnvoll aber eine gute Übung…).

Lösung

Berechne das 77% Konfidenzintervall der Variable wl1 (nicht sehr sinnvoll aber eine gute Übung…). 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(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(wloss$wl1[wloss$group == "Control"])
sd_CG <- sd(wloss$wl1[wloss$group == "Control"])
n_CG <- length(wloss$wl1[wloss$group == "Control"])
se_CG <- sd_CG/sqrt(n_CG)

mean_IG <- mean(wloss$wl1[wloss$group == "DietEx"])
sd_IG <- sd(wloss$wl1[wloss$group == "DietEx"])
n_IG <- length(wloss$wl1[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

Übung 5

Aufgabe

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

Lösung

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

library(Rmisc)
group.CI(wl1 ~ group, data = 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 auswenig 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

library(ggplot2)
library(Hmisc)
ggplot(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 8

Aufgabe

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?

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 9

Aufgabe

Du hast im Auftrag einer Stiftung anhand eines Fragebogens Daten zur Arbeitszufriedenheit von Gesundheitspersnal gesammelt. Das 95% CI für den durchschnittliche 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.