Übungsblatt 7

Regressionsanalytisches Modell für zwei Prädiktoren und Regressionsdiagnostik

Sie betrachten ein multiples Regressionsmodell für zwei Prädiktoren $X_1$ und $X_2$.
1. Welche Auswirkungen hat ein großer Standardfehler für B₁ bzw. B₂ auf die Länge der Konfidenzintervalle für $\beta_{1}$ und $\beta_{2}$?
  
  Lösung
  
  Je größer der Standardfehler, desto breiter das entsprechende Konfidenzintervall.
2. Wovon hängt es ab, wie hoch der Standardfehler von B₁ bzw. B₂ ausfällt?
  
  Hinweis:
  
  \[SE\left( B_{1} \right) = \sqrt{Var\left( B_{1} \right)} = \sqrt{\frac{1}{1 - r_{x1x2}^{2}} \cdot \frac{\sigma_{}^{2}}{\sum_{i = 1}^{n}\left( x_{i1} - {\bar{x}}_{1} \right)^{2}}}\]
  
  Lösung
  
  Der Standardfehler von B₁ wird größer,...
  ...je höher die quadrierte Korrelation zwischen x₁ und x₂ ausfällt, je höher also die beiden Prädiktoren miteinander korrelieren.
  ...je größer die Fehlervarianz ist.
  ...je kleiner die Streuung von x₁ ist.
  ...je kleiner die Stichprobe ist.
  
  Analoges gilt für den Standardfehler von B₂
Sie haben für den Sommer 2017 in 200 Städten auf der ganzen Welt die Anzahl der gekauften Eiskugeln in Millionen sowie die Anzahl der Fahrradunfälle erhoben. Laden Sie den SPSS-Datensatz herunter, lesen Sie diesen mit der Funktion read.spss() aus dem foreign package in R ein (Hinweis: Bei diesem Befehl sollten Sie das Argument to.data.frame = TRUE benutzen). Führen Sie dann eine einfache lineare Regression mit der Anzahl an Fahrradunfällen als Kriterium und der Anzahl an verkauften Eiskugeln als Prädiktor durch. (Hinweis: Die Modellvoraussetzungen sind gegeben!)
Lösung
library(foreign) data <- read.spss("Eis.sav", to.data.frame = TRUE)
fit <- lm(Fahrradunfälle ~ Eiskugeln, data = data) summary(fit)

Call: lm(formula = Fahrradunfälle ~ Eiskugeln, data = data) Residuals: Min 1Q Median 3Q Max -63.569 -13.659 0.171 15.084 54.263 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 280.242 4.433 63.216 < 2e-16 *** Eiskugeln 26.010 4.160 6.253 2.43e-09 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 20.08 on 198 degrees of freedom Multiple R-squared: 0.1649, Adjusted R-squared: 0.1607 F-statistic: 39.1 on 1 and 198 DF, p-value: 2.43e-09
1. Berechnen Sie mit R das Konfidenzintervall für β und interpretieren Sie es inhaltlich.
  Lösung
  confint(fit)
  
  2.5 % 97.5 % (Intercept) 271.50029 288.98464 Eiskugeln 17.80712 34.21256
  
  Die plausiblen Werte für die Fahrradunfälle, die pro Million verkaufter Eiskugeln in der Population erwartet werden, liegen zwischen 17.81 und 34.21.
2. Zusätzlich wurde an jedem Strand die Durchschnittstemperatur in Grad Celsius für den Sommer 2017 erhoben. Geben Sie die Modellgleichung für eine multiple lineare Regression mit der Anzahl an Fahrradunfällen als Kriterium und der Anzahl an verkauften Eiskugeln sowie der Durchschnittstemperatur als Prädiktoren an.
  
  Lösung
  
  Y: Anzahl Fahrradunfälle, X₁: Anzahl verkaufter Eiskugeln in Millionen, X₂: Temperatur in Grad Celsius
  
  \[{Y}_{i} = \alpha + {\beta}_{Eis} \cdot {X}_{Eis_i} + {\beta}_{Temp} \cdot{X}_{Temp_i} + \varepsilon_{i}\]
3. Wie würden Sie in diesem Modell ein β_Eis von 0 interpretieren?
  
  Lösung
  
  Bei konstanter Temperatur erwartet man pro Million verkaufter Eiskugeln durchschnittlich einen Anstieg von 0 Fahrradunfällen. Das heißt, bei konstanter Temperatur hängt die Anzahl der erwarteten Fahrradunfälle nicht mehr von der Anzahl verkaufter Eiskugeln ab.
4. Testen Sie mit R, ob β_Eis in der Population 0 beträgt. Dazu müssen Sie die Funktion lm() um einen weiteren Prädiktor erweitern: lm( AV ~ UV1 + UV2, data = Datensatz). Wie interpretieren Sie das Ergebnis des Tests inhaltlich?
  Lösung
  fit2 <- lm(Fahrradunfälle ~ Eiskugeln + Temperatur, data = data) summary(fit2)
  
  Call: lm(formula = Fahrradunfälle ~ Eiskugeln + Temperatur, data = data) Residuals: Min 1Q Median 3Q Max -60.573 -11.940 0.737 12.603 46.418 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 179.393 22.043 8.138 4.46e-14 *** Eiskugeln -3.460 7.459 -0.464 0.643 Temperatur 5.242 1.125 4.661 5.78e-06 *** --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Residual standard error: 19.11 on 197 degrees of freedom Multiple R-squared: 0.2479, Adjusted R-squared: 0.2402 F-statistic: 32.46 on 2 and 197 DF, p-value: 6.534e-13
  
  Bei einem Signifikanzniveau von 0.005 wird die Nullhypothese, dass β_Eis in der Population 0 beträgt, beibehalten (p > 0.005). Das heißt, die Anzahl der verkauften Eiskugeln scheint über die Temperatur hinaus (bzw. bei konstanter Temperatur) keinen Einfluss auf die durchschnittliche Anzahl der erwarteten Fahrradunfälle zu haben
Professorin Freud predigt ihren Studierenden immer, dass Psychologiestudierende mit einem durchschnittlichen IQ von 100, die mindestens 30 Stunden auf die Psychoanalyseklausur gelernt haben, diese mit hoher Sicherheit (sie wählt ein 0.95-Konfidenzniveau) bestehen. Um diese Aussage nächstes Semester mit empirischen Ergebnissen untermauern zu können, hat sie von 200 zufällig gezogenen Psychologiestudierenden die Punktzahl in der Klausur, den IQ-Wert aus einem Intelligenztest und die Anzahl an Stunden, die die Studierenden auf die Klausur gelernt haben, gesammelt. In jeder Klausur von Professorin Freud benötigt man 40 Punkte, um diese zu bestehen. Laden Sie den SPSS-Datensatz herunter und lesen Sie ihn in R ein. Sprechen die Daten dafür, dass die Aussage von Professorin Freud zutreffend ist?
Lösung
daten <- read.spss("Freud.sav", to.data.frame = TRUE)
fit3 <- lm(Klausurpunkte ~ IQ + Lerndauer, data = daten) daten_neu <- data.frame(IQ = 100, Lerndauer = 30) predict(fit3, daten_neu, interval = "prediction", level = 0.95)

fit lwr upr 1 50.21753 40.68859 59.74647
Für einen IQ von 100 und eine Lerndauer von 30 Stunden ergibt sich ein 95%-Vorhersage-Konfidenzintervall von 40.69 bis 59.75 Klausurpunkten. Die Aussage von Professorin Freud ist also insofern zutreffend, als dass man für eine Student*in mit einem IQ von 100 und einer Lerndauer von 30 Stunden mit hoher Sicherheit erwartet, dass diese mehr als 40 Klausurpunkte erreicht, also die Klausur besteht.

Bemerkung: Es sollte jedoch beachtet werden, dass die Ergebnisse aufgrund des nicht- experimentellen Designs der Studie nicht kausal interpretiert werden dürfen! Es ist also nicht sichergestellt, dass eine längere Vorbereitung auf die Klausur wirklich die Wahrscheinlichkeit, die Klausur zu bestehen, erhöht.
Im Übungsblatt 6 haben Sie mit dem Datensatz “Beziehungen.csv” eine Regressionsanalyse zur Vorhersage der Dauer von Beziehungen aus der Altersdifferenz der beiden Partner*innen durchgeführt. Der Output wurde als Objekt gespeichert, das wir im Folgenden “fit” nennen.
1. Führen Sie für den Datensatz eine Regressionsdiagnostik durch.
  
  Folgende neue Befehle benötigen Sie hierfür:
  - Mit der Funktion plot(AV ~ UV, data = Daten) können Sie das Streudiagramm anzeigen lassen. Die Funktion abline(fit) zeichnet in das Streudiagramm die geschätzte Regressionsgerade ein.
  - Mit der Funktion rstandard(fit) können Sie die standardisierten Residuen als Objekt, z.B. „res” speichern und für weitere Überprüfungen verwenden.
  - Mit der Funktion plot(fit$fitted.values, res) können Sie sich die vorhergesagten Werte aus dem Objekt „fit” auf der x-Achse und die Residuen auf der y-Achse in einer Grafik anzeigen lassen.
  - Mit der Funktion which(cooks.distance(fit) > Cut-Off-Wert) können Sie sich die Personen anzeigen lassen, die einen bestimmten Cut-Off-Wert für die Cook’s Distanz überschreiten. Dieser Wert muss für die vorliegenden Daten noch berechnet und bei „Cut-Off-Wert” eingesetzt werden, ein Beispiel wäre: which(cooks.distance(fit) > 0.022).
  - Mit der Funktion plot(cooks.distance(fit)) können Sie sich die Cook’s Distanzen für alle Personen anzeigen lassen. Mit abline(h = Cut-Off-Wert) können Sie in diese Grafik eine Linie auf einer bestimmten Höhe einzeichnen, hier bietet sich der Cut-Off-Wert an, z.B. abline(h = 0.022)
  Lösung
  Daten <- read.csv2("Beziehungen.csv") fit <- lm(Beziehungsdauer ~ Altersdifferenz, data = Daten)
  
  Überprüfung der Linearitätsannahme:
  
  plot(Beziehungsdauer~Altersdifferenz,data=Daten)
  
  Interpretation: Die Linearitätsannahme scheint zu gelten.
  
  Überprüfung der Normalverteilungsannahme:
  
  res <- rstandard(fit) hist(res)
  
  Interpretation: Es kann davon ausgegangen werden, dass die Fehler ε_i einer Normalverteilung folgen.
  
  Überprüfung der Homoskedastizitätsannahme:
  
  plot(fit$fitted.values,res)
  
  Interpretation: Die Homogenitätsannahme scheint zu gelten.
  
  Ausreißeranalyse: Cut-off-Wert Cook’s Distanz: $\frac{4}{n} = \frac{4}{1000} = 0.004$
  
  # Cut-off-Wert Cook’s Distanz: 0.004 plot(cooks.distance(fit)) abline(h=0.004)
  
  which(cooks.distance(fit)>0.004)
  
  24 31 36 58 98 109 114 141 150 190 196 216 225 273 283 300 303 338 339 349 24 31 36 58 98 109 114 141 150 190 196 216 225 273 283 300 303 338 339 349 361 384 422 443 461 472 473 479 538 551 563 573 589 618 637 639 657 659 678 688 361 384 422 443 461 472 473 479 538 551 563 573 589 618 637 639 657 659 678 688 695 703 786 788 812 813 848 875 878 881 886 938 999 695 703 786 788 812 813 848 875 878 881 886 938 999
  
  # Prozentsatz von Cook`s Distance > 0.004` length(which(cooks.distance(fit)>0.004))/nrow(Daten)
  
  [1] 0.053
  
  Wird der Cut-off-Wert als Maßstab zur Bewertung herangezogen, dann weist der Datensatz ca. 5% Einflusswerte auf. Da aber im Modell explizit festgelegt wird, dass die Fehler einer Normalverteilung folgen, sind auch sehr hohe bzw. niedrige Variablenwerte zu erwarten.
2. Warum sind Einflusswerte bei regressionsanalytischen Verfahren problematisch? Wie gehen Sie mit identifizierten Einflusswerten um?
  
  Lösung
  
  Hohe Einflusswerte wirken sich aufgrund ihrer Hebelwirkung auf die Regressions-koeffizienten aus.
  
  Nach der Identifizierung von Einflusswerten wird nach der Ursache für die hohen Ausprägungen gesucht: Liegen Dateneingabefehler vor? Da das statistische Modell außergewöhnlich hohe bzw. niedrige Werte zulässt, sollten Ausreißerwerte nicht entfernt werden (es sei denn es liegt ein Dateneingabefehler vor, der nicht korrigiert werden kann): Es wird in den Modellannahmen festgelegt, dass die Verteilung der Beobachtungswerte einer Normalverteilung folgen, es wird daher erwartet, dass einige wenige Werte außergewöhnlich hoch bzw. niedrig sind.
3. Tragen Sie in die folgende Graphik einen Einflusswert ein und veranschaulichen Sie dessen Wirkung, indem Sie eine zweite Regressionsgerade einzeichnen.
  
  Lösung
Überprüfen Sie für das Regressionsmodell aus Aufgabe 3 dieses Übungsblattes die Modellannahmen der Linearität, Normalverteilung und Homoskedastizität. (Hinweis: Verwenden Sie für die Partiellen Residuen-Plots die Funktion crPlots aus dem car package.)
Lösung
library(car) crPlots(fit3)

res <- rstandard(fit3) hist(res)

plot(fit3$fitted.values, res)
Alle Modellannahmen scheinen erfüllt zu sein. Hinweis: Die pinke Linie im Partiellen Residuen-Plot ist sehr anfällig für Ausreißer und sollte nur mit großer Vorsicht interpretiert werden!