Lernziele und Literatur

Fast alle im Bachelor-Modul “Empirische Gesundheitsökonomik” behandelten Methoden sind Varianten der linearen Regression. Die Studierenden verstehen, dass lineare Regression ein beschreibendes Verfahren zur Schätzung bedingter Mittelwerte ist, und sie können Regressionskoeffizienten in einfachen und multiplen linearen Regressionen, einschließlich solcher mit Interaktionstermen, korrekt interpretieren.

Literatur
  • Gelman & Hill, Kapitel 3.1-3.3
  • James et al. An Introduction to Statistical Learning, Kapitel 3.

Einleitung

Die meisten Einführungen in die lineare Regressionsanalyse beschreiben diese als den Versuch, eine Gerade “optimal” durch eine Punktwolke zu legen. Dies erfolgt mit dem Ziel der bestmöglichen Beschreibung des (linearen) Zusammenhangs zwischen zwei kontinuierlichen Größen. “Optimal” bedeutet, dass die Summe der quadrierten Residuen, d.h. der Abweichungen der Punkte einer Punktwolke von der geschätzten Gerade minimiert wird. Das mathematische Verfahren heißt Methode der kleinsten Quadrate (KQ-Methode), bzw. auf Englisch Ordinary Least Squares (OLS) und ist im Grundsatz aus einführenden Statistikvorlesungen bekannt.

Das folgende Beispiel zeigt eine typische Anwendung, die dieser Motivation folgt. Hier wurden Ausgaben für Lebensmittel und Einkommen (jeweils in logarithmierter Form) in einem Koordinatensystem abgetragen und eine Regressiongerade nach der Methode der kleinsten Quadrate durch die Punktwolke gelegt. Ziel ist die empirische Schätzung einer sogenannten Engel-Kurve, wie sie aus der mikroökonomischen Theorie bekannt ist, also des Zusammenhangs von Einkommen und Konsumausgaben.

Jeder Punkt der Abbildung stellt einen Haushalt und seine Kombination von Einkommen und Ausgaben für Lebensmittel innerhalb eines bestimmten Zeitraums dar (Querschnittsdaten). Der Haushalt ganz rechts hat beispielsweise ein Einkommen von \(\exp(8.5) \approx 4\,900\) Geldeinheiten und verausgabt \(\exp(7.5) \approx 1\,800\) Geldeinheiten für Lebensmittel. Die Gerade durch die Punktwolke ist die geschätzte Engel-Kurve. Sie hat die Form:

\[\mbox{Log (Ausgaben)} = 0.54 + 0.86 \times \mbox{Log (Einkommen)}\]

Der Wert 0.86 kann, da sowohl die abhängige Variable “Ausgaben” als auch der Prädiktor “Einkommen” logarithmiert wurde, direkt als Einkommenselastizität der Nachfrage nach Lebensmitteln interpretiert werden (eine Erklärung findet sich in Kapitel 4).

In der “Empirischen Gesundheitsökonomik” werden wir Regression nicht als Verfahren zur Schätzung funktionaler Zusammenhänge, sondern vielmehr als ein leistungsfähiges Verfahren zur Schätzung von (bedingten) Mittelwerten, sowie deren Vergleich zwischen verschiedenen Gruppen, interpretieren. Dabei müssen – wie wir sehen werden – weder die abhängige Variable noch die Prädiktoren kontinuierlich sein. Ein typisches Beispiel für die Anwendung der Regressionsanalyse in diesem Sinne ist der Vergleich von Versuchs- und Kontrollgruppe einer randomisierten Studie in Bezug auf den Mittelwert der primären Ergebnisvariable. Der Prädiktor “Gruppenzugehörigkeit” ist eine kategoriale Variable mit zwei Kategorien (auch binäre oder dichotome Variable genannt).

Lineare Regression mit einem binären Prädiktor

Zur Illustration der Regressionsanalyse verwenden wir zunächst wieder den HBSC-Datensatz. Das einfachste Beispiel ist die Regression mit einer kontinuierlichen abhängigen Variablen (hier: bodyheight) und einem binären Prädiktor (hier: sex). Zunächst ist nur das Verständnis der Regressionskoeffizienten, d.h. des Achsenabschnitts und der Steigungskoeffizienten, wichtig. Ihre Interpretation bereit erfahrungsgemäß Probleme. Bevor die Interpretation der Koeffizienten bei verschiedenen Modellspezifikationen verstanden wurde, ist eine Beschäftigung mit weiterführenden Themen wie Schätzfehler, Inferenzstatistik, Signifikanz, Identifikation kausaler Effekt, usw., nicht sinnvoll. Wir werden dies später ausführlich behandeln.

In R werden lineare Regressionen mit der Funktion lm() geschätzt. Da es sich um eine “alte” Base-R-Funktion handelt, sind die Daten nicht das erste genannte Argument der Funktion. Stattdessen ist das erste Argument eine Formel bodyheight ~ 1 + sex, d.h. wir schätzen die mittlere Körpergröße als (lineare) Funktion einer Konstanten (1) und des Prädiktors Geschlecht. Die Regressionskonstante für den Achsenabschnitt muss nicht eigens angegeben werden. Die Formel bodyheight ~ sex würde zum gleichen Resultat führen. Wir werden im Folgenden auf Angabe der Konstante in der Formel verzichten.

lm(bodyheight ~ 1 + sex, data = hbsc)
## 
## Call:
## lm(formula = bodyheight ~ 1 + sex, data = hbsc)
## 
## Coefficients:
## (Intercept)      sexGirl  
##     164.130       -4.082

Das Gleichung der berechneten Regressionsgerade lautet:

\[ \widehat{\mbox{bodyheight}} = 164.130 - 4.082 \cdot \mbox{sexGirl} \]

wobei \(\widehat{\mbox{bodyheight}}\) die vorhersagte oder erwartete Körpergröße bei gegebenem Geschlecht eines Kindes ist. Das Dach über bodyheight zeigt der Konvention folgend an, dass es sich um einen Schätzwert handelt. Interpretation:

Alternativ könnten Mittelwerte nach Gruppen auch mit group_by() und summarise() berechnet werden:

hbsc %>% 
  group_by(sex) %>% 
  summarise(bodyheight = mean(bodyheight, na.rm = TRUE))
## `summarise()` ungrouping output (override with `.groups` argument)
## # A tibble: 2 x 2
##   sex   bodyheight
##   <fct>      <dbl>
## 1 Boy         164.
## 2 Girl        160.

Bei der Ausgabe von Tibbles werden Zahlen kompakt dargestellt, hier auf ganze Zahlen gerundet. Berechnet werden die Mittelwerte natürlich genauer. Um in der Ausgabe die gleiche Anzahl an Nachkommastellen zu bekommen wie bei lm(), können wir mit der Funktion options(pillar.sigfig = 6) die Anzahl der bei Ausgabe dargestellten signifikanten Ziffern auf 6 erhöhen (die Grundeinstellung (default) ist 3).

options(pillar.sigfig = 6)
hbsc %>% 
  group_by(sex) %>% 
  summarise(bodyheight = mean(bodyheight, na.rm = TRUE))
## `summarise()` ungrouping output (override with `.groups` argument)
## # A tibble: 2 x 2
##   sex   bodyheight
##   <fct>      <dbl>
## 1 Boy      164.130
## 2 Girl     160.048

Die Ergebnisse entsprechen nun exakt denen von lm dargestellten. Nochmal: das zeigt, dass Regression ein Verfahren zur Schätzung gruppenspezifischer Mittelwerte ist.

Übungsaufgaben mit R
  1. Laden Sie den Datensatz “NHANES” aus dem Paket “NHANES”. NHANES steht für National Health and Nutrition Examination Survey, eine große, regelmäßig durchgeführte Studie zum Gesundheitsstatus und Ernährungszustand von Erwachsenen und Kindern in den USA.
  2. Behalten Sie nur Beobachtungen mit einem Alter von 50 Jahren oder älter und speichern Sie den Dataframe unter einem anderen Namen ab (nur im Arbeitsspeicher, nicht auf der Festplatte). Wieviele Beobachtungen hat der Datensatz jetzt? Wieviele Männer, wieviele Frauen?
  3. Stellen Sie den Zusammenhang zwischen Geschlecht und systolischem Blutdruck in einer geeigneten Abbildung dar, z.B. in zwei Boxplots. Welches Geschlecht hat im Schnitt einen höheren Blutdruck?
  4. Berechnen Sie den mittleren diastolischem und systolischem Blutdruck nach Geschlecht mit group_by() und summarise().
  5. Quantifizieren Sie den Zusammenhang von Geschlecht und systolischem Blutdruck mit Hilfe einer linearen Regression. Vergleichen Sie die Ergebnisse mit dem aus Aufgabe 4. Interpretieren Sie den Achsenabschnitt und den Steigungskoeffizienten inhaltlich und quantitativ.

Faktoren als Prädiktoren

Nun mag man sich fragen, was R eigentlich berechnet, wenn ein binärer Prädiktor mit Ausprägungen wie “Boy” und “Girl” verwendet wird. Damit man rechnen kann, müssen diese Ausprägungen eine numerische Repräsentation haben. Tatsächlich werden – ohne dass man das sieht – die Werte Null und Eins verwendet. Um dies beispielhaft nachzuvollziehen, ändern wir unsere Regression leicht ab. Wir erzeugen zuerst eine neue numerische Variable girl, die den Wert Eins für alle Mädchen und den Wert Null für alle Jungen hat. Dazu bedienen wir uns der Funktion ifelse(). Durch Voranstellen von hbsc$ stellen wir sicher, dass die neue Variable Teil unseres Dataframe hbsc wird. (Eine Alternative stellt die mit # herauskommentierte Verwendung der Pipe %>% und mutate() dar. Diese Alternative wäre hier umständlicher. Sollen aber mehrere Variable neu erzeugt werden, dann kann mutate() eine deutlich übersichtlichere Variante sein.)

hbsc$girl <- ifelse(hbsc$sex=="Girl", 1, 0)
# hbsc <- hbsc %>% 
#   mutate(girl = ifelse(sex=="Girl", 1, 0))
names(hbsc)
## [1] "sex"        "age"        "bodyheight" "bodyweight" "thinkbody" 
## [6] "girl"

names() nach Erzeugen der Variable girl zeigt, dass girl am Ende des Dataframe eingefügt wurde. Durch Kreuztabellieren der alten Variable sex und der neuen Variable girl können wir den Erfolg unserer ifelse()-Transformation überprüfen:

flat_table(hbsc, sex, girl)
##      girl   0   1
## sex              
## Boy       260   0
## Girl        0 240

Offenbar haben wie gewünscht alle Jungen den Wert girl = 0 und alle Mädchen den Wert girl = 1. Nun verwenden wir diese Variable als Prädiktor in einer Regression:

lm(bodyheight ~ girl, data = hbsc)
## 
## Call:
## lm(formula = bodyheight ~ girl, data = hbsc)
## 
## Coefficients:
## (Intercept)         girl  
##     164.130       -4.082

Wir erhalten die gleiche geschätzte Regressiongleichung wie zuvor. Das heißt: R führt bei binären Faktoren in Regressionsmodellen intern die gleiche Transformation durch, die wir hier explizit mit girl durchgeführt haben: für eine Ausprägung des Merkmals wird der Wert numerisch auf Null, für die andere Ausprägung des Merkmals numerisch auf Eins gesetzt.

Exkurs*

Wie R intern mit Faktoren und anderen Ausdrücken in Formeln umgeht, kann man mit Hilfe der Funktion model.matrix() überprüfen. Diese zeigt, wie die Prädiktoren zur Schätzung eines durch Formel angegebenen Regressionsmodells intern aufbereitet werden: Für den Achsenabschnitt gibt es überall eine Eins, beim Geschlecht für Jungen eine Null, für Mädchen eine Eins. Wir geben nur die ersten sechs Beobachtungen aus:

head(model.matrix(bodyheight ~ sex, data = hbsc ))
##   (Intercept) sexGirl
## 1           1       1
## 2           1       1
## 3           1       0
## 4           1       1
## 5           1       0
## 6           1       0

Vorsicht vor Scheingenauigkeit

Bei Interpretation von Ergebnissen wird oft die Genauigkeit, die uns der Output von R vorgaukelt, missverstanden. Selbst erfahrene Forscher scheinen manchmal zu glauben, dass die per Voreinstellung gewählte Anzahl an Nachkommastellen auch in Veröffentlichungen berichtet werden sollte. R hat bei der Darstellung des Regressionsergebnisses ohne unser Zutun drei Nachkommastellen ausgewählt, aber niemand misst die Größe eines Menschen in cm auf drei Nachkommastellen genau. Daher ist es sinnvoll, auf ganze cm zu runden. Man kann sich die Arbeit des Rundens von Hand u.a. durch Anwendung der Funktion print() mit der Option digits abnehmen lassen. Dazu speichern wir das Regressionsergebnis zunächst in einem Objekt namens fit geben dieses Objekt anschließend mit print() aus. Hinweis: das Argument digits der Funktion print() ist nicht als Anzahl an Nachkommastellen zu interpretieren, sondern als Anzahl an signifikanten Stellen, d.h. die Anzahl an Ziffern ungleich Null hinter den führenden Nullen.

fit <- lm(bodyheight ~ sex, data = hbsc)
print(fit, digits = 1)
## 
## Call:
## lm(formula = bodyheight ~ sex, data = hbsc)
## 
## Coefficients:
## (Intercept)      sexGirl  
##         164           -4

Alternativ kann man die beiden Zeilen auch als Pipe spezifizieren und erspart sich damit das Abspeichern des Regressionergebnisses:

lm(bodyheight ~ sex, data = hbsc) %>% print(digits = 1)
## 
## Call:
## lm(formula = bodyheight ~ sex, data = hbsc)
## 
## Coefficients:
## (Intercept)      sexGirl  
##         164           -4

Das Regressionsergebnis lautet nun gerundet:

\[ \widehat{\mbox{bodyheight}} = 164 - 4 \cdot \mbox{sexGirl} \]

Interpretation:

  • Für Jungen (mit dem Wert sexGirl=0) wird eine mittlere Größe von \(164 - 4\cdot 0=164\) cm geschätzt
  • Für Mädchen (mit dem Wert sexGirl=1) wird eine mittlere Größe von \(164 - 4\cdot 1=160\) cm geschätzt
  • Der Unterschied beträgt 4. Die Mädchen sind im Schnitt 4 cm kleiner als die Jungen.

Graphische Darstellung

Zur Illustration verdeutlichen wir das Ergebnis der Regressionanalyse anhand eines Streudiagramms. Zunächst stellen wird die Ursprungsdaten dar:

ggplot(data = hbsc, aes(x = sex, y = bodyheight)) +
  geom_jitter(width = .05, height = .5)

Das Geom geom_jitter() erzeugt eine Variante des Streudiagramms (geom_point()), bei der auf alle Datenpunkte noch eine kleine Zufallszahl mit Mittelwert Null addiert wird, um das sogenannte Overplotting zu vermeiden. Die Abbildung zeigt die Körpergröße aller Jungen und Mädchen im Datensatz. Durch die lineare Regression wird das arithmetische Mittel in jeder Gruppe berechnet und das Ergebnis ist eine Gerade, die diese beiden arithmetischen Mittel verbindet. Wir zeichnen nun zuerst die arithmethischen Mittel (nach Geschlecht) ein. Diese werden durch stat_summary() berechnet und durch große blaue Punkte dargestellt.

ggplot(data = hbsc, aes(x = sex, y = bodyheight)) +
  geom_jitter(width = .05, height = .5) +
  stat_summary(fun.y = "mean", geom = "point",
    size = 5, color = "blue")

Zuletzt verbinden wir die beiden großen blauen Punkte. Die Verbindungsstrecke erzeugen wir durch das Geom geom_smooth(), das verbunden mit dem Argument method = "lm" nach der Methode der kleinsten Quadrate berechnete Regressiongeraden berechnet und in den Graphen einzeichnet.

ggplot(data = hbsc, aes(x = sex, y = bodyheight)) +
  geom_jitter(width = .05, height = 0.5) +
  stat_summary(fun.y = "mean", geom = "point",
    size = 5, color = "blue") +
  geom_smooth(aes(x = as.numeric(sex)), method = "lm", se = FALSE)

Nochmal: die Details der Erstellung dieser Abbildungen mit R sind nicht wichtig. Sie sind hier nur für den interessierten Leser angegeben. Die Darstellung soll vielmehr verdeutlichen: Lineare Regression ist ein Verfahren zur Schätzung gruppenspezifischer Mittelwerte. Im Fall eines binären Prädiktors gibt der Achsenabschnitt den Mittelwerte der Referenzgruppe (hier die Jungen), und der Steigungskoeffizient den Unterschied in den Mittelwerten der jeweils anderen Gruppe im Vergleich zur Referenzgruppe an.

Lineare Regression mit einem kontinuierlichen Prädiktor

Erhöht man die Anzahl möglicher unterschiedlicher Ausprägungen des Prädiktors bis hin zu potenziell unendlich, dann kann man nicht mehr für jede dieser Ausprägungen einen Mittelwert der abhängigen Variable berechnen und sinnvoll. Vielmehr muss man vereinfachende Annahmen darüber treffen, wie sich die bedingten Mittelwerte verhalten, wenn sich der Prädiktor verändert. Eine häufig getroffene Annahme ist, dass jede Veränderung des Prädiktors um eine jeweils Einheit mit einer gleich starken Veränderung des bedingten Mittelwerts der abhängigen Variable verbunden ist. Mit anderen Worten: der bedingte Mittelwert wird als eine lineare Funktion des Prädiktors geschätzt.

Im folgenden wollen wir das mittlere Körpergewicht als lineare Funktion der Körpergröße der Kinder im Datensatz vorhersagen. Zunächst jedoch stellen wir die Daten als Streudiagramm dar, dazu mit stat_summary() das mittlere Körpergewicht bei jeder Ausprägung der Körpergröße:

hbsc %>% 
  ggplot(aes(x = bodyheight, y = bodyweight)) +
  geom_point(size = 1, color = "grey40") +
  stat_summary(color = "blue")

Die Ausgangsdaten werden als kleine, graue Punkte dargestellt, die bedingten Mittelwerte für jeden cm Körpergröße in blau. Man erkennt, dass stat_summary() nur dann Mittelwerte und Fehlermargen einzeichnet, wenn mindestens zwei Datenpunkte vorhanden sind. (Das hängt damit zusammen, dass mit einem Datenpunkt keine Fehlermarge berechnet werden kann.) Jedenfalls sieht man auch, dass eine Verbindung der bedingten Mittelwerte durch eine Gerade keine schlechte Anpassung darstellen würde. Im Gegenteil: betrachtet man etwa die beiden etwas höher gelegenen Punkte bei 171 und 172 cm in der Mitte der Abbildung, dann drängt sich der Eindruck auf, dass es sich um zufällige Abweichungen nach oben handelt. Niemand würde behaupten, dass bei exakt diesen beiden Größen Menschen systematisch schwerer sind, als wenn sie einen cm kleiner oder einen cm größer sind. Die Verbindung der bedingten Mittelwerte durch eine Gerade gleicht solche kleinen Abweichungen in sinnvoller Weise aus.

hbsc %>% 
  ggplot(aes(x = bodyheight, y = bodyweight)) +
  geom_point(size = 1, color = "grey40") +
  stat_summary(color = "blue") +
  geom_smooth(method = "lm")

Hier wurde das Element geom_smooth() nach der Methode lm hinzugefügt. Das ist nichts anderes als die geschätzte Regressionsgerade. Sie stellt die Verbindung aller (mit der Anzahl der Beobachtungen gewichteten) geschätzen Mittelwerte des Körpergewichts bei verschiedenen Körpergrößen dar. Der graue Bereich um die blaue Regressiongerade herum stellt die Unsicherheit dieser auf einer Stichprobe von Kindern beruhenden Schätzung dar. Dazu mehr im Kapitel “Statistische Signifkanz”. Wird dieses “Konfidenzband” als störend empfunden, kann es mit se = FALSE abgeschaltet werden.

Die Parameter der eingezeichneten Regressiongerade (mit der Funktion lm() berechnet) lauten:

lm(bodyweight ~ bodyheight, data = hbsc)
## 
## Call:
## lm(formula = bodyweight ~ bodyheight, data = hbsc)
## 
## Coefficients:
## (Intercept)   bodyheight  
##    -95.5521       0.9094

Die geschätzte Regressiongleichung, wieder unter Verzicht auf unnötige Nachkommastellen, lautet:

\[ \widehat{\mbox{bodyweight}} = -95.6 + 0.9 \cdot \mbox{bodyheight} \]

Interpretation:

Übungsaufgaben mit R
  1. Laden Sie den Datensatz “NHANES” aus dem Paket “NHANES” und behalten Sie nur Beobachtungen mit einem Alter von 50 Jahren oder älter.
  2. Stellen Sie den Zusammenhang von Alter und systolischem Blutdruck in einem Streudiagramm dar. Verwenden Sie erst geom_point und dann geom_jitter() und vergleichen Sie: welche Darstellung gefällt Ihnen besser? Was ist bei 80 Jahren los?
  3. Ergänzen Sie das Streudiagramm um eine Regressiongerade (mit oder ohne Konfidenzband).
  4. Quantifizieren Sie den Zusammenhang von Alter und systolischem Blutdruck mit Hilfe einer linearen Regressionanalyse. Interpretieren Sie Achsenabschnitt und Steigungskoeffizient inhaltlich und quantitativ. Falls eine Interpretation nicht sinnvoll ist, geben Sie an warum.

Lineare Regression mit mehreren Prädiktoren

Die Stärke der Regressionsanalyse besteht darin, dass Mittelwertunterschiede zwischen Gruppen in Abhängigkeit von mehreren Prädiktoren gleichzeitig geschätzt werden können. Allerdings wird zugleich die Interpretation einzelner Koeffizienten schwieriger. Im Regelfall gilt die ceteris paribus-Interpretation: ein Steigungskoeffizient gibt den Mittelwertunterschied in Bezug auf einen Prädiktor an, wenn zugleich alle anderen Prädiktoren konstant gehalten werden. Wir illustrieren dies nun am Beispiel des Zusammenhangs von Körpergröße, Alter und Geschlecht.

lm(bodyheight ~ age + sex, data = hbsc) %>% print(digits = 1)
## 
## Call:
## lm(formula = bodyheight ~ age + sex, data = hbsc)
## 
## Coefficients:
## (Intercept)          age      sexGirl  
##         101            5           -4

Die geschätzte Gleichung für das Körpergewicht mit Körpergröße und Geschlecht als Prädiktoren lautet: \[ \widehat{\mbox{bodyweight}}= 101 + 5 \cdot \mbox{age} - 4 \cdot \mbox{sexGirl} \]

Interpretation:

Graphische Abbildung (Streudiagramm und Regressionsgeraden) des geschätzten Modells:

Auf den ersten Blick fällt auf, dass man hier drei und nicht eine Punktwolke hat: bei 11, 13 und 15 Jahren. Das hat mit einer Besonderheit der HBSC-Studie zu tun, in deren Stichprobe sich – von wenigen Ausnahmen abgesehen – nur 11-, 13- und 15- jährige befinden. Ansonsten zeigt sich aber der – in diesem Alter zu erwartende – positive Zusammenhang zwischen Alter und Körpergröße.

Die Abbildung zeigt nun als Ergebnis der Regression zwei parallele Geraden der mittleren Größe von Jungen und Mädchen in Abhängigkeit vom Alter im vertikalen Abstand von 4 cm. Dies verdeutlicht die Interpretation des Koeffizienten von sexGirl als ceteris paribus: Mädchen sind dieser Schätzung zufolge in jedem Alter 4 cm kleiner als Jungen. Bei genauerer Betrachtung der Punktewolken fällt jedoch auf, dass gerade im Alter 15 die Jungen häufig oberhalb “ihrer” Regressionsgerade liegen und im Alter 13 eher unterhalb. Möglicherweise ist die dem geschätzten Modell innewohnende Annahme, dass Jungen und Mädchen mit jedem Lebensjahr gleich wachsen, falsch. Will man Geraden mit unterschiedlicher Steigung für unterschiedliche Gruppen schätzen, dann greift man oft auf Interaktionsterme zurück.

Interaktionsterme

Interaktionsterme machen die Interpretation einzelner Koeffizienten noch schwieriger, da die quantitative Interpretation dann von anderen Größen im Modell abhängt. Gilt ohne Interaktionsterm, dass ein Steigungskoeffizient den ceteris-paribus-Mittelwertunterschied in Bezug auf einen Prädiktor angibt, so können bei Einschluss von Interaktionstermen die Mittelwertunterschiede in Bezug auf einen Prädiktor von der Ausprägung anderer Prädiktoren abhängen. Wir illustrieren dies wieder am Beispiel des Zusammenhangs von Körpergröße, Alter und Geschlecht, und fügen einen Interaktionsterm als das mathematische Produkt von Alter und Geschlecht hinzu. In R-Formeln wird dies unmittelbar durch den Doppelpunkt erzeugt. Das Produkt kann, muss aber nicht als neue Variable des Dataframe erzeugt werden.

lm(bodyheight ~ age + sex + age:sex, data = hbsc) %>% print(digits = 1)
## 
## Call:
## lm(formula = bodyheight ~ age + sex + age:sex, data = hbsc)
## 
## Coefficients:
## (Intercept)          age      sexGirl  age:sexGirl  
##          89            6           20           -2

Die geschätzte Gleichung mit Interaktion lautet:

\[ \widehat{\mbox{bodyheight}}= 89 + 6 \cdot \mbox{age} + 20 \cdot \mbox{sexGirl} - 2 \cdot \mbox{age} \cdot \mbox{sexGirl} \]

Interpretation:

In der Abbildung sieht man deutlich, dass das Größenwachstum der Mädchen zwischen Alter 13 und 15 kleiner ist als das Größenwachstum der Jungen.

ggplot(hbsc, aes(x = age, y = bodyheight, color = sex)) +
  geom_point() +
  geom_smooth(method = "lm", se = FALSE)
## `geom_smooth()` using formula 'y ~ x'

Bemerkung: dass die Schätzung teilweise Koeffizienten geliefert hat, deren Interpretation nicht sinnvoll ist, ist je nach Anwendung wenig befriedigend. Wir werden später lernen, wie man durch lineare Transformation der Ausgangsvariablen in Modellen mit Interaktionstermen zu sinnvolleren Darstellungen kommen kann.

Getrennte Schätzung statt Interaktionsterm

Wenn einer der Prädiktoren binär ist, dann bietet sich als Alternative zum Einfügen von Interaktionstermen als äquivalente Alternative das getrennte Schätzen zweier Regressionmodelle an, hier für Jungen und Mädchen getrennt. Die Funktion filter() dient dabei der Auswahl von Datenreihen (also Beobachtungen) nach dem Geschlecht.

lm(bodyheight ~ age, data = filter(hbsc, sex == "Boy")) %>% print(digits = 1)
## 
## Call:
## lm(formula = bodyheight ~ age, data = filter(hbsc, sex == "Boy"))
## 
## Coefficients:
## (Intercept)          age  
##          89            6
lm(bodyheight ~ age, data = filter(hbsc, sex == "Girl")) %>% print(digits = 1) 
## 
## Call:
## lm(formula = bodyheight ~ age, data = filter(hbsc, sex == "Girl"))
## 
## Coefficients:
## (Intercept)          age  
##         109            4

Vergleich man die Schätzwerte mit denen oben, dann stellt man folgendes fest: Die Regression für Jungen liefert die gleichen Werte wie oben, wenn man sexGirl auf Null setzt. Die Regression für Mädchen liefert die gleichen Werte wie oben, wenn man sexGirl auf Eins setzt, d.h. wenn man (Intercept) und sexGirl, sowie age und age:sexGirl addiert.

Übungsaufgaben mit R
  1. Laden Sie den Datensatz “NHANES” aus dem Paket “NHANES” und behalten Sie nur Beobachtungen mit einem Alter von 50 Jahren oder älter.
  2. Verwenden Sie Alter und Geschlecht zur Vorhersage des mittleren systolischen Blutdrucks in einer Regressionsanalyse. Interpretieren Sie Achsenabschnitt und Steigungskoeffizienten der Regression inhaltlich und quantitativ. Falls eine Interpretation nicht sinnvoll ist, geben Sie an warum.
  3. Verwenden Sie Alter und Geschlecht zur Vorhersage des mittleren systolischen Blutdrucks in einer Regressionsanalyse mit Interaktionterm. Interpretieren Sie Achsenabschnitt und Steigungskoeffizienten inhaltlich und quantitativ. Falls eine Interpretation nicht sinnvoll ist, geben Sie an warum. Für welches Geschlecht steigt der mittlere systolische Blutdruck stärker mit dem Alter?
  4. Erstellen Sie ein Streudiagramm des Zusammenhangs von Alter und systolischem Blutdruck. Stellen Sie Männer und Frauen in unterschiedlichen Farben dar (Hinweis: aes(color = ...). Zeichnen Sie jeweils Regressiongeraden für Männer und Frauen ein.

Wiederholungsfragen

Geben Sie an, ob die Aussage wahr oder falsch ist, und begründen Sie Ihre Antwort. Beziehen Sie sich dabei auf Inhalte des Kapitels, und nennen Sie zusätzliche Annahmen, die Sie möglicherweise treffen müssen.

  1. Bei der Darstellung von Mittelwerten und Regressionsergebnissen sollte man immer soviele Nachkommastellen angeben, wie von R ausgegeben werden.
  2. Der Achsenabschnitt einer linearen Regression kann bei kontinuierlichen Prädiktoren stets sinnvoll interpretiert werden.
  3. Lineare Regression ist als statistisches Verfahren nicht nur zur Schätzung des linearen Zusammenhangs zwischen kontinuierlichen Variablen geeignet.
  4. Der Koeffizient eines binären Prädiktors gibt den Unterschied im Mittelwert zwischen zwei Gruppen an.
  5. Nach Gruppen getrennte Regressionen sind äquivalent zu einer einzigen Regression mit Indikatorvariablen für die Gruppenzugehörigkeit.