Cvičenie 10 – Multikolinearita
Vanessa Vargová
knitr::opts_chunk$set(echo = TRUE, message = FALSE, warning = FALSE)1. Úvod
V tomto zadaní znova pracujem s údajmi World Happiness Report za rok 2019 z kaggle (https://www.kaggle.com/datasets/unsdsn/world-happiness), sú uložené v zložke “udaje” pod názvom 2019.csv
2. Údaje
udaje_raw <- read.csv("udaje/2019.csv", header = TRUE)
# zjednotenie názvov stĺpcov
names(udaje_raw) <- make.names(names(udaje_raw))
# výber relevantných premenných
udaje <- udaje_raw[, c(
"Score",
"GDP.per.capita",
"Social.support",
"Healthy.life.expectancy"
)]
# imputácia chýbajúcich hodnôt mediánom
column_medians <- sapply(udaje, median, na.rm = TRUE)
for (col in names(udaje)) {
udaje[[col]][is.na(udaje[[col]])] <- column_medians[col]
}3. Základný regresný model
Vysvetľujem hodnotu Score pomocou premenných GDP per capita, Social support a Healthy life expectancy. Všetky tri regresné koeficienty sú *štatisticky významné na hladine významnosti 1 % a majú kladný vplyv na skóre šťastia. Upravený koeficient determinácie je približne 0,72, čo znamená, že model vysvetľuje približne 72 % variability závislej premennej.
model <- lm(
Score ~ GDP.per.capita + Social.support + Healthy.life.expectancy,
data = udaje
)
summary(model)##
## Call:
## lm(formula = Score ~ GDP.per.capita + Social.support + Healthy.life.expectancy,
## data = udaje)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.7018 -0.4155 -0.0520 0.4535 1.3369
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.1350 0.2116 10.088 < 2e-16 ***
## GDP.per.capita 0.8098 0.2358 3.434 0.000766 ***
## Social.support 1.3219 0.2483 5.324 3.58e-07 ***
## Healthy.life.expectancy 1.2977 0.3661 3.544 0.000523 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.588 on 152 degrees of freedom
## Multiple R-squared: 0.7263, Adjusted R-squared: 0.7209
## F-statistic: 134.5 on 3 and 152 DF, p-value: < 2.2e-164. Korelačná matica a scatterploty
Z korelačnej matice vidíme pomerne silné kladné korelácie medzi vysvetľujúcimi premennými, najmä medzi GDP per capita a Healthy life expectancy (≈ 0,84). Tieto hodnoty naznačujú možnú prítomnosť multikolinearity, aj keď korelácie ešte nepresahujú kritickú hranicu 0,9.
xvars <- udaje[, c("GDP.per.capita", "Social.support", "Healthy.life.expectancy")]
round(cor(xvars), 3)## GDP.per.capita Social.support Healthy.life.expectancy
## GDP.per.capita 1.000 0.755 0.835
## Social.support 0.755 1.000 0.719
## Healthy.life.expectancy 0.835 0.719 1.000pairs(
xvars,
main = "Scatterplotová matica – vysvetľujúce premenné (2019)"
)
5. VIF
Žiadna hodnota nepresahuje hranicu 5, podľa VIF nejde o závažnú multikolinearitu.
library(car)
vif(model)## GDP.per.capita Social.support Healthy.life.expectancy
## 3.956068 2.473462 3.5227426. Condition Number
Výsledok 14,6 signalizuje miernu multikolinearitu. Tento ukazovateľ potvrdzuje výsledky korelačnej matice a VIF.
X <- model.matrix(model)[, -1]
XtX <- t(X) %*% X
eig <- eigen(XtX)
condition_number <- sqrt(max(eig$values) / min(eig$values))
condition_number## [1] 14.637417. Riešenie 1: vynechanie premennej
Po vynechaní premennej GDP per capita a Healthy life expectancy zostávajú oba zjednodušené modely štatisticky významné. Upravený koeficient determinácie sa však mierne znižuje v porovnaní so základným modelom, čo naznačuje stratu vysvetľovacej schopnosti. Vynechanie premennej teda síce znižuje multikolinearitu, ale za cenu horšej kvality modelu.
model_no_gdp <- lm(
Score ~ Social.support + Healthy.life.expectancy,
data = udaje
)
summary(model_no_gdp)##
## Call:
## lm(formula = Score ~ Social.support + Healthy.life.expectancy,
## data = udaje)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.41862 -0.41895 -0.05979 0.46693 1.45281
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.8664 0.2035 9.173 2.95e-16 ***
## Social.support 1.6661 0.2350 7.089 4.65e-11 ***
## Healthy.life.expectancy 2.1051 0.2904 7.248 1.95e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6084 on 153 degrees of freedom
## Multiple R-squared: 0.7051, Adjusted R-squared: 0.7012
## F-statistic: 182.9 on 2 and 153 DF, p-value: < 2.2e-16model_no_health <- lm(
Score ~ GDP.per.capita + Social.support,
data = udaje
)
summary(model_no_health)##
## Call:
## lm(formula = Score ~ GDP.per.capita + Social.support, data = udaje)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.00390 -0.40373 -0.03284 0.44986 1.35217
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.3298 0.2120 10.991 < 2e-16 ***
## GDP.per.capita 1.3465 0.1875 7.182 2.80e-11 ***
## Social.support 1.5375 0.2496 6.159 6.18e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6098 on 153 degrees of freedom
## Multiple R-squared: 0.7037, Adjusted R-squared: 0.6998
## F-statistic: 181.7 on 2 and 153 DF, p-value: < 2.2e-168. Riešenie 2: škálovanie premenných
Po centrovaní a škálovaní vysvetľujúcich premenných sa hodnota Condition Number výrazne znížila na približne 4, čo signalizuje nízku multikolinearitu. Hodnoty VIF zostali nezmenené, čo je očakávané. Model si zároveň zachoval rovnakú vysvetľovaciu schopnosť ako pôvodný model, avšak regresné koeficienty stratili priamu interpretáciu v pôvodných jednotkách.
udaje$gdp_c <- scale(udaje$GDP.per.capita, center = TRUE, scale = TRUE)
udaje$soc_c <- scale(udaje$Social.support, center = TRUE, scale = TRUE)
udaje$hl_c <- scale(udaje$Healthy.life.expectancy, center = TRUE, scale = TRUE)
model_centered <- lm(
Score ~ gdp_c + soc_c + hl_c,
data = udaje
)
summary(model_centered)##
## Call:
## lm(formula = Score ~ gdp_c + soc_c + hl_c, data = udaje)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.7018 -0.4155 -0.0520 0.4535 1.3369
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.40710 0.04708 114.848 < 2e-16 ***
## gdp_c 0.32262 0.09394 3.434 0.000766 ***
## soc_c 0.39550 0.07428 5.324 3.58e-07 ***
## hl_c 0.31420 0.08865 3.544 0.000523 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.588 on 152 degrees of freedom
## Multiple R-squared: 0.7263, Adjusted R-squared: 0.7209
## F-statistic: 134.5 on 3 and 152 DF, p-value: < 2.2e-16vif(model_centered)## gdp_c soc_c hl_c
## 3.956068 2.473462 3.522742Xc <- model.matrix(model_centered)[, -1]
XtXc <- t(Xc) %*% Xc
eig_c <- eigen(XtXc)
condition_number_c <- sqrt(max(eig_c$values) / min(eig_c$values))
condition_number_c## [1] 3.970619. Riešenie 3: úprava jednotiek (GDP v tisícoch)
Po prepočítaní premennej GDP per capita na tisíce sa regresné koeficienty dostali do porovnateľných rádov a zachovali si interpretovateľnosť. Hodnoty VIF ostali rovnaké ako v základnom modeli, avšak Condition Number výrazne vzrástol, čo naznačuje, že samotná zmena jednotiek neodstraňuje numerickú nestabilitu matice \(\mathbf{X}^T\mathbf{X}\). Tento prístup je však vhodný z hľadiska interpretácie koeficientov.
udaje$GDPpc_1000 <- udaje$GDP.per.capita / 1000
model_gdp_1000 <- lm(
Score ~ GDPpc_1000 + Social.support + Healthy.life.expectancy,
data = udaje
)
summary(model_gdp_1000)##
## Call:
## lm(formula = Score ~ GDPpc_1000 + Social.support + Healthy.life.expectancy,
## data = udaje)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.7018 -0.4155 -0.0520 0.4535 1.3369
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.1350 0.2116 10.088 < 2e-16 ***
## GDPpc_1000 809.8204 235.8097 3.434 0.000766 ***
## Social.support 1.3219 0.2483 5.324 3.58e-07 ***
## Healthy.life.expectancy 1.2977 0.3661 3.544 0.000523 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.588 on 152 degrees of freedom
## Multiple R-squared: 0.7263, Adjusted R-squared: 0.7209
## F-statistic: 134.5 on 3 and 152 DF, p-value: < 2.2e-16vif(model_gdp_1000)## GDPpc_1000 Social.support Healthy.life.expectancy
## 3.956068 2.473462 3.522742Xg <- model.matrix(model_gdp_1000)[, -1]
XtXg <- t(Xg) %*% Xg
eig_g <- eigen(XtXg)
condition_number_g <- sqrt(max(eig_g$values) / min(eig_g$values))
condition_number_g## [1] 6766.17210. Zhrnutie
Na základe výsledkov môžeme konštatovať, že v dátach World Happiness Report 2019 je prítomná mierna multikolinearita medzi vysvetľujúcimi premennými. Ako najefektívnejšie riešenie z hľadiska stability modelu sa ukázalo škálovanie premenných, zatiaľ čo úprava jednotiek je vhodná najmä vtedy, ak je prioritou zachovanie interpretovateľnosti regresných koeficientov.
Ďakujem za dnešné cvičenie Vanessa Vargová