Probleme Linearer Regression
Auf dieser Seite findest du Übungen zur Bewältigung von Kollinearitätsproblemen in R.
Der gesamte Code, den du hier findest, kann von dir 1:1 in R umgesetzt werden. Die größte Lernwirkung hat diese Einführung, wenn du möglichst alles einmal selbst in R durchführst.
01 Vorbereitungen
Bevor wir loslegen können, müssen wir die Packages laden, die wir für unser Vorhaben benötigen. Zunächst ist das nur das SfL Package. Falls du dieses noch nicht geladen hast, hole dies nun nach:
library(SfL)Nun laden wir die Datensätze, die wir nutzen werden:
data(data_a)
data(data_c)
data(data_s)
data(data_t)
data(data_v)02 Feststellung Problematischer Korrelation
Je mehr Variablen in einem Regressionsmodell berücksichtigt werden sollen, desto höher ist die Wahrscheinlichkeit, dass einige dieser Variablen stark miteinander korreliert sind. Da Korrelation zwischen Prediktorvariablen zu Kollinearität in linearen Modellen führt, und solche Kollinearität zur Zuverlässigkeit solcher Modelle signifikant verschlechtert, muss dieses Problem unbedingt berücksichtigt werden.
Hierzu muss zunächst überprüft werden, welche Variablen miteinander korreliert sind. Dies funktioniert z.B. mit folgendem Befehl aus dem SfL Package:
library(SfL)
correlation_matrix(data = data_s,
variables = c("baseDurLog", "speakingRate", "googleFreq",
"age", "speaker", "typeOfS", "folSeg", "folType", "item", "transcription", "biphoneProbSum", "biphoneProbSumBin", "monoMultilingual", "pauseBin", "preC"))
Betrachten wir die Korrelationskoeffizienten näher, finden wir folgende problematische Fälle:
folSegvs.folTypebiphoneProbSumBinvs.biphoneProbSumbiphoneProbSumBinvs.itembiphoneProbSumBinvs.transcriptionsbiphoneProbSumvs.itembiphoneProbSumvs.transcriptionitemvs.transcription
03 Vergleich und Ausschluss
Idee: Die stark korrelierten Variablen werden in simple Modelle eingesetzt und diese dann vergliechen. Die Variable mit der jeweils größeren Vorhersagekraft wird behalten.
Potentielles Problem: Der Verlust von Variablen kann einem Verlust an Vorhersagekraft entsprechen.
Vorgehen: Mit der predictor_competition() Funktion des SfL Packages kann herausgefunden werden, welche Variable der jeweils bessere Prediktor ist. Alternativ kann dieses Vorgehen auch manuell durchgeführt werden.
Schritt 1
folSeg vs. folType
predictor_competition(data = data_s,
dependent = "sDurLog",
independent1 = "folSeg",
independent2 = "folType")Analysis of Variance Table
Model 1: sDurLog ~ folSeg
Model 2: sDurLog ~ folType
Res.Df RSS Df Sum of Sq F Pr(>F)
1 135 22.221
2 145 24.453 -10 -2.2321 1.356 0.2074Wie vom Output angedeutet, sind beide Variablen gleich gute Prediktoren. In diesem Fall können wir wählen, welchen wir für die weitere Analyse behalten möchten.
Wir behalten folType. Weitere Vergleiche mit folSeg entfallen.
Schritt 2
biphoneProbSumBin vs. biphoneProbSum
predictor_competition(data = data_s,
dependent = "sDurLog",
independent1 = "biphoneProbSumBin",
independent2 = "biphoneProbSum")Analysis of Variance Table
Model 1: sDurLog ~ biphoneProbSumBin
Model 2: sDurLog ~ biphoneProbSum
Res.Df RSS Df Sum of Sq F Pr(>F)
1 148 24.264
2 148 23.937 0 0.32695 Wie vom Output angedeutet, sind beide Variablen gleich gute Prediktoren. In diesem Fall können wir wählen, welchen wir für die weitere Analyse behalten möchten.
Wir behalten biphoneProbSumBin. Weitere Vergleiche mit biphoneProbSum entfallen.
Schritt 3
item vs. transcription
predictor_competition(data = data_s,
dependent = "sDurLog",
independent1 = "item",
independent2 = "transcription")Analysis of Variance Table
Model 1: sDurLog ~ item
Model 2: sDurLog ~ transcription
Res.Df RSS Df Sum of Sq F Pr(>F)
1 105 16.919
2 114 19.843 -9 -2.9244 2.0165 0.04442 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Wie vom Output angedeutet, ist eine Variable ein besserer Prediktor als die andere Variable. Anhand des RSS-Wertes stellen wir fest, die bessere Prediktorvariable ist item.
Weitere Vergleiche mit transcription entfallen.
Schritt 4
biphoneProbSumBin vs. item
predictor_competition(data = data_s,
dependent = "sDurLog",
independent1 = "biphoneProbSumBin",
independent2 = "item")Analysis of Variance Table
Model 1: sDurLog ~ biphoneProbSumBin
Model 2: sDurLog ~ item
Res.Df RSS Df Sum of Sq F Pr(>F)
1 148 24.264
2 105 16.919 43 7.3446 1.06 0.3961Wie vom Output angedeutet, sind beide Variablen gleich gute Prediktoren. In diesem Fall können wir wählen, welchen wir für die weitere Analyse behalten möchten.
Wir behalten item. Weitere Vergleiche mit biphoneProbSumBin entfallen.
Weiteres Vorgehen
Mit den Gewinnern der vorherigen Schritte sowie den Variablen, die keine problematischen Korrelationskoeffizienten gezeigt haben, können nun regulär lineare Modelle erstellt werden. Zum Beispiel:
mdl <- lm(sDurLog ~
folType + # Gewinner
item + # Gewinner
baseDurLog + # alle weiteren Variablen
speakingRateLog +
googleFreqLog +
age +
speaker +
typeOfS +
monoMultilingual +
pauseBin +
preC,
data = data_s)04 Principal Component Analysis
Möchten wir nach Möglichkeit keine Vielzahl von Variablen verlieren, können wir eine PCA durchführen. Allerdings funktionieren diese nur beschränkt gut bei kategorischen Variablen. Daher wird im folgenden Verlauf auf ein fiktives Beispiel mit rein numerischen Variablen gesetzt.
Zunächst werden die Beispieldaten erstellt:
sigma <- rbind(c(1,-0.8,-0.7), c(-0.8,1, 0.9), c(-0.7,0.9,1))
mu <- c(10, 5, 2)
df1 <- as.data.frame(mvrnorm(n=1000, mu=mu, Sigma=sigma))
sigma <- rbind(c(1,-0.8,-0.7), c(-0.8,1, 0.9), c(-0.7,0.9,1))
mu <- c(9, 4, 1)
df2 <- as.data.frame(mvrnorm(n=1000, mu=mu, Sigma=sigma))
sigma <- rbind(c(1,-0.8,-0.7), c(-0.8,1, 0.9), c(-0.7,0.9,1))
mu <- c(11, 6, 3)
df3 <- as.data.frame(mvrnorm(n=1000, mu=mu, Sigma=sigma))
df <- cbind(df1, df2, df3)
colnames(df) <- c("var1", "var2", "var3", "var4", "var5", "var6", "var7", "var8", "var9")Bitte beachte, dass du beim Umsetzen der PCA minimal andere Ergebnisse erhalten wirst, da die hier erstellten Variablen zufällige Werte enthalten. Die einzelnen Schritte beeinflusst das nicht.
Schritt 1
Zuerst wird ein Data Frame erstellt, der nur die Variablen enthält, die problematische Korrelationen aufweisen. Dafür müssen wir zunächst die Korrelationen überprüfen:
correlation_matrix(df)
Offenbar ist jede Variable mit mindestens einer anderen Variable stark korreliert. Daher müssen wir keinen neuen Data Frame erstellen, sondern können den jetzigen für weitere Schritte nutzen.
Schritt 2
Mit der PCA Funktion von R - pc() - berechnen wir die PCA:
pc <- princomp(df, cor=TRUE, score=TRUE)Schritt 3
Als nächstes können wir uns die Eigenvalues der neu erstellten Variablen anschauen. Hierzu nutzen wir das factoextra Package.
install.packages("factoextra")library("factoextra")
get_eigenvalue(pc) eigenvalue variance.percent cumulative.variance.percent
Dim.1 2.64833838 29.4259820 29.42598
Dim.2 2.62724036 29.1915596 58.61754
Dim.3 2.53887662 28.2097402 86.82728
Dim.4 0.33286865 3.6985406 90.52582
Dim.5 0.30974824 3.4416471 93.96747
Dim.6 0.28595582 3.1772869 97.14476
Dim.7 0.09197875 1.0219862 98.16674
Dim.8 0.08800664 0.9778515 99.14459
Dim.9 0.07698653 0.8554059 100.00000Je höher die Eigenvalue einer neuen Variable - hier Dim.x genannt - desto mehr Varianz kann diese Variable in den Daten erklären.
Schritt 4
Anhand üblicher Kriterien muss nun entschieden werden, welche der neuen Variablen für weitere Analysen behalten werden.
Kriterium 1
Die Eigenvalue einer Variable muss mindestens 1 betragen, da eine Variable mit Eigenvalue 1 mindestens sich selbst erklären kann.
get_eigenvalue(pc) eigenvalue variance.percent cumulative.variance.percent
Dim.1 2.64833838 29.4259820 29.42598
Dim.2 2.62724036 29.1915596 58.61754
Dim.3 2.53887662 28.2097402 86.82728
Dim.4 0.33286865 3.6985406 90.52582
Dim.5 0.30974824 3.4416471 93.96747
Dim.6 0.28595582 3.1772869 97.14476
Dim.7 0.09197875 1.0219862 98.16674
Dim.8 0.08800664 0.9778515 99.14459
Dim.9 0.07698653 0.8554059 100.00000Die folgenden PCA Components erfüllen das Kriterium: Dim.1, Dim.2, Dim.3
Kriterium 2
Die ‚cumulative variance explained‘ der neuen Variablen sollte mindestens 70 % betragen.
get_eigenvalue(pc) eigenvalue variance.percent cumulative.variance.percent
Dim.1 2.64833838 29.4259820 29.42598
Dim.2 2.62724036 29.1915596 58.61754
Dim.3 2.53887662 28.2097402 86.82728
Dim.4 0.33286865 3.6985406 90.52582
Dim.5 0.30974824 3.4416471 93.96747
Dim.6 0.28595582 3.1772869 97.14476
Dim.7 0.09197875 1.0219862 98.16674
Dim.8 0.08800664 0.9778515 99.14459
Dim.9 0.07698653 0.8554059 100.00000Die folgenden PCA Components erfüllen das Kriterium: Dim.1, Dim.2, Dim.3
Kriterium 3
Es sollte keine großen Unterschiede zwischen den Eigenvalues der neuen Variablen geben.
get_eigenvalue(pc) eigenvalue variance.percent cumulative.variance.percent
Dim.1 2.64833838 29.4259820 29.42598
Dim.2 2.62724036 29.1915596 58.61754
Dim.3 2.53887662 28.2097402 86.82728
Dim.4 0.33286865 3.6985406 90.52582
Dim.5 0.30974824 3.4416471 93.96747
Dim.6 0.28595582 3.1772869 97.14476
Dim.7 0.09197875 1.0219862 98.16674
Dim.8 0.08800664 0.9778515 99.14459
Dim.9 0.07698653 0.8554059 100.00000Die folgenden PCA Components erfüllen das Kriterium: Dim.1, Dim.2, Dim.3
Kriterium 4
Neue Variablen sind nur hilfreich, wenn sie auch interpretierbar sind. Obwohl es sich hier um fiktive Variablen handelt, können wir dieses Kriterium heranziehen. Components, die “ein bisschen von allem” in ihren Loadings enthalten, sind nicht interpretierbar. Enthalten Components allerdings große Anteile einiger bestimmter Variablen, sind die i.d.R. interpretierbar.
pc$loadings[,1:3] Comp.1 Comp.2 Comp.3
var1 0.2383359 0.40961666 0.2877117
var2 -0.2595344 -0.45275418 -0.2968547
var3 -0.2635848 -0.42026070 -0.2935514
var4 0.1430101 -0.36227606 0.3882110
var5 -0.1440861 0.38539583 -0.4389340
var6 -0.1408366 0.36496130 -0.4273685
var7 0.4807232 -0.09797024 -0.2606992
var8 -0.5118988 0.12084916 0.2837557
var9 -0.5021733 0.11591147 0.2615927Die folgenden PCA Components erfüllen das Kriterium: Dim.1, Dim.2, Dim.3
Entscheidung
Nach allen Kriterien sind folgende Components als neue Variablen zu behalten: Dim.1, Dim.2, Dim.3
Schritt 5
Nun müssen die Variablen, die behalten werden, extrahiert werden:
pcframe <- as.data.frame(pc$scores)[1:3]Bei einem echten Anwendungsfall würde man diese als zusätzliche Spalten zum eigentlichen Data Frame hinzufügen, damit sie im weiteren Analyseprozess zur Verfügung stehen.
Schritt 6
Nun steht die Interpretation der neuen Variablen an. In Schritt 4 haben wir die Loadings der einzelnen Components schon kurz gesehen. Hier schauen wir uns diese nun im Detail an.
pc$loadings[,1:3] Comp.1 Comp.2 Comp.3
var1 0.2383359 0.40961666 0.2877117
var2 -0.2595344 -0.45275418 -0.2968547
var3 -0.2635848 -0.42026070 -0.2935514
var4 0.1430101 -0.36227606 0.3882110
var5 -0.1440861 0.38539583 -0.4389340
var6 -0.1408366 0.36496130 -0.4273685
var7 0.4807232 -0.09797024 -0.2606992
var8 -0.5118988 0.12084916 0.2837557
var9 -0.5021733 0.11591147 0.2615927Die Loadings einer PCA Component stellen ihre Korrelation mit den ursprünglichen Variablen dar. Je höher die Korrelation, desto mehr der ursprünglichen Variable ist in einer Component enthalten.
Comp.1 umfasst die Effekte von var7, var8 und var9 am stärksten, da hier die höchsten Korrelationen festzustellen sind. Dabei geht der Effekt von var7 in die entgegengesetzte Richtung zu den Effekten von var8 und var9.
Comp.2 umfasst die Effekte von var1, var2 und var3 am stärksten, da hier die höchsten Korrelationen festzustellen sind. Dabei geht der Effekt von var1 in die entgegengesetzte Richtung zu den Effekten von var2 und var3.
Comp.3 umfasst die Effekte von var4, var5 und var6 am stärksten, da hier die höchsten Korrelationen festzustellen sind. Dabei geht der Effekt von var4 in die entgegengesetzte Richtung zu den Effekten von var5 und var6.
Weiteres Vorgehen
Mit den neuen Variablen sowie den Variablen, die nicht Teil der PCA waren (falls es solche gab), können nun regulär lineare Modelle erstellt werden. Zum Beispiel:
# Erstellen einer abhängigen Variable zur Veranschaulichung
pcframe$numbers2 <- c(rep(1:10, 50)*-1)
mdl <- lm(numbers2 ~ Comp.1 + Comp.2 + Comp.3, pcframe)
anova(mdl)Analysis of Variance Table
Response: numbers2
Df Sum Sq Mean Sq F value Pr(>F)
Comp.1 1 11.5 11.490 1.4009 0.236857
Comp.2 1 69.3 69.268 8.4454 0.003741 **
Comp.3 1 0.1 0.141 0.0172 0.895672
Residuals 996 8169.1 8.202
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 105 Aufgaben
Checke bei den unten angebenen Modellen, ob es hohe Korrelationskoeffizienten zwischen Prediktoren gibt. Falls ja, nutze eine der beiden Strategien zur Vermeidung von Kollinearitätsproblemen.
Aus data_s:
sDurLog~bsex+age+locationsDurLog~typeOfS+pauseBinsDurLog~baseDurLog+speakingRatesDurLog~biphoneProbSumBin+folType+preCsDurLog~bsex+age+location+typeOfS+pauseBin+baseDurLog+speakingRate+biphoneProbSumBin+folType+preC
Aus data_a:
height~bsex+hairheight~hair+eyesheight~bsex+hair+hair+eyes
Aus data_v:
duration ~ vowel + structure + rate