R Anleitungen

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

R Code
# 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.

Oft wollen wir statt den unstandardisierten Regressionskoeffizienten, die standardisierten Koeffizienten haben. In R können wir das Package 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.

R Code
# 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

R Code
# 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.

R Code
# 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.

R Code
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

Anstatt p-Werte zu interpretieren, schauen wir uns die Konfidenzintervalle des Bootstrappings an. Wenn ein 95%-Konfidenzintervall 0 nicht mit einschließt, können wir davon ausgehen, dass der Koeffizient auf einem von α = .05 signifikant ist.

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.

R Code
# 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.

R Code
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.

Aufgrund verschiedener US-Patente ist eine Metrik (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.
R Code
# 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))