R: Multiple Regression
R bietet eine Vielzahl von Funktionen und Packages mit denen eine ganze Reihe von Regressionverfahren samt Regressionsdiagnostiken durchgeführt werden können.
Ein Regressionsmodell aufstellen
# Beispielmodell mit y als Kriterium und x1, x2 und x3 als Prädiktoren
Modell <- lm(y ~ x1 + x2 + x3, data = Regressionsdaten)
# Zusammenfassung der Regression anzeigen
summary(Modell)
# Call:
# lm(formula = y ~ x1 + x2 + x3, data = Regressionsdaten, data = data)
#
# Residuals:
# Min 1Q Median 3Q Max
# -1.37482 -0.34992 0.00423 0.33156 1.76924
#
# Coefficients:
# Estimate Std. Error t value Pr(>|t|)
# (Intercept) 0.5628660 0.1865514 3.017 0.00274 **
# x1 -0.0008655 0.0012821 -0.675 0.50005
# x2 0.0033282 0.0634318 0.052 0.95818
# x3 0.9031152 0.0314215 28.742 < 2e-16 ***
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
# Residual standard error: 0.5246 on 353 degrees of freedom
# (25 observations deleted due to missingness)
# Multiple R-squared: 0.7038, Adjusted R-squared: 0.7013
# F-statistic: 279.6 on 3 and 353 DF, p-value: < 2.2e-16
In der Ausgabe sehen wir zuerst noch einmal unseren Aufruf (Call
), dann die Residuen (Residuals
) und danach die Koeffizienten und die verbundenen Statistiken (Coefficients
).
In der Tabelle der Koeffizienten ist die erste Spalte der unstandardisierte Regressionskoeffizient, danach der Standardfehler und am Ende der t-Wert mit dem korrespondierenden p-Wert, der uns sagt, ob der Regressionskoeffizient signifikant ist.
lm.beta
und darin die gleichnamige Funktion lm.beta()
verwenden, um die standardisierten Regressionskoeffizienten zu erhalten, in unserem Beispiel über lm.beta(Modell)
.Interaktionen und Funktionen
Wir können beim Aufstellen unseres Regressionsmodells Interaktionen, also das Produkt zwischen zwei Variablen, angeben. Dies geschieht in R allerdings nicht über den üblichen Multiplikationsoperator (*
), sondern über den Doppelpunkt (:
), wie unten in Modell2
.
Wir können auch Funktionen auf Variablen anwenden, ohne Dafür eine neue Variable definieren zu müssen. In Modell3
haben wir so unser Kriterium y
logarithmiert in das Regressionsmodell aufgenommen.
# Modell mit y als Kriterium und x3 sowie die Interaktion von x1 und x2 als Prädiktoren
Modell2 <- lm(y ~ x1:x2 + x3, data = Regressionsdaten)
# Modell mit dem Logarithmus von y als Kriterium und x1, x2 und x3 als Prädiktoren
Modell3 <- lm(log(y) ~ x1 + x2 + x3, data = Regressionsdaten)
Andere nützliche Funktionen
# Regressionskoeffizienten anzeigen
coefficients(Modell)
# 95 % Konfidenzintervall der Regressionskoeffizienten berechnen
confint(Modell, level=0.95)
# Vorhergesagte Werte anzeigen
fitted(Modell)
# Residuen des Regressionsmodells
residuals(Modell)
Voraussetzungen überprüfen
Zum Überprüfung der Multikollinearität (VIF) und der Autokorrelation benötigen wir noch das Paket car
, dass uns zwei Funktionen hierfür zur Verfügung stellt.
# Residuen auf Normalverteilung überprüfen
shapiro.test(residuals(Modell))
# Autokorrelation überprüfen (Durbin-Watson-Test)
car::durbinWatsonTest(Modell)
# Multikollinearität mit VIF (variance influence factor) überprüfen
car::vif(Modell)
Bootstrapping
Bootstrapping ist eine robuste Methode, statistische Modelle zu berechnen. Es berechnet ein statistisches Modell wiederholt (meist einige tausend Mal) mit verschiedenen zufälligen Zusammenstellungen aus dem ursprünglichen Datensatz. Es hat selbst kaum Voraussetzungen und produziert stets normalverteilte Variablen (je höher die Anzahl der Wiederholungen, desto eher entsprechen die Variablen einer Normalverteilung).
Wir können Bootstrapping mit dem R-Paket car
auch für die Berechnung von Regressionsmodellen verwenden.
library(car)
# Bootstrapping mit 5000 Wiederholungen berechnen
bootstrapRegression <- Boot(Modell, R=5000)
# Zusammenfassung
summary(bootstrapRegression)
# 95 % Konfidenzintervalle für das Bootstrapping ausgeben
confint(bootstrapRegression, level = 0.95)
# Bootstrap bca confidence intervals
#
# 2.5 % 97.5 %
# (Intercept) -1.878251610 0.07500605
# x1 0.003333067 0.12184393
# x2 0.080634292 0.49185284
# x3 0.106946437 0.21257247
Schrittweise Regression und Variablenauswahl
Wie der Name „schrittweise Regression“ andeutet, werden bei diesem Verfahren Variablen schrittweise ausgewählt. Unabhängige Variablen werden anhand ihrer statistischen Signifikanz der Reihe nach entweder hinzugefügt oder entfernt. Bei der schrittweisen Auswahl wird entweder die signifikanteste Variable hinzugefügt oder die am wenigsten signifikante Variable entfernt. Es werden nicht alle möglichen Modelle in Betracht gezogen und es wird am Ende des Algorithmus ein einziges Regressionsmodell generiert.
Schrittweise Regression
In R benutzen wir aus dem Package MASS
die Funktion stepAIC
, die die schrittweise Regression für uns berechnet.
# Signifikante Variablen werden dem Modell schrittweise hinzugefügt
MASS::stepAIC(Modell, direction = "forward")
# Nicht-signifikante Variablen werden aus dem Modell schrittweise entfernt
MASS::stepAIC(Modell, direction = "backward")
# Variablen werden dem Modell schrittweise hinzugefügt und entfernt
MASS::stepAIC(Modell, direction = "both")
Alle Variablenkombinationen
Stattdessen können wir aber auch alle Variablenkombinationen im Regressionsmodell berechnen und das Modell auswählen, was beispielsweise die höchste Varianz aufklärt (d.h. den höchsten R²-Wert hat). Genau das macht die Funktion regsubsets
aus dem R-Paket leaps
.
leaps <- leaps::regsubsets(y ~ x1 + x2 + x3, data = Regressionsdaten, nbest=5)
# Zusammenfassung der 5 besten Modelle anzeigen
summary(leaps)
# Die 5 besten Modelle grafisch darstellen (mit R²)
plot(leaps)
Relative Wichtigkeit
Einen ähnlichen Ansatz verfolgt das R-Paket relaimpo. Es beinhaltet sechs verschiedene Metriken, welche die Wichtigkeit der einzelnen Prädiktoren im Regressionsmodell bestimmen sollen. Zusätzlich dazu beinhaltet es Funktionen zur grafischen Darstellung und Bootstrapping der Koeffizienten.
pmvd
) nicht im Paket enthalten, dass über install.packages()
installiert wird. Will man diese Metrik berechnen, muss man das Paket direkt von der Seite der Beuth Hochschule für Technik Berlin herunterladen und manuell installieren.# Relative Wichtigkeit eines Modells berechnen
relaimpo::calc.relimp(Modell, type = c("lmg", "last", "first", "betasq", "pratt"), rela = TRUE )
# Regression mit Bootstrapping und relativer Wichtigkeit berechnen
bootstrapRegression <- relaimpo::boot.relimp(Modell, type = c("lmg", "last", "first", "betasq", "pratt"), rela = TRUE, b = 5000, diff = T, rank = T )
relaimpo::booteval.relimp(bootstrapRegression)
# Grafische Ausgabe der Berechnungen
plot(relaimpo::booteval.relimp(bootstrapRegression, sort=TRUE))