Multikolinearita

Načítanie dát

econ <- read.csv("economics.csv")
str(econ)

'data.frame':   574 obs. of  7 variables:
 $ X       : int  1 2 3 4 5 6 7 8 9 10 ...
 $ date    : chr  "1967-07-01" "1967-08-01" "1967-09-01" "1967-10-01" ...
 $ pce     : num  507 510 516 512 517 ...
 $ pop     : num  198712 198911 199113 199311 199498 ...
 $ psavert : num  12.6 12.6 11.9 12.9 12.8 11.8 11.7 12.3 11.7 12.3 ...
 $ uempmed : num  4.5 4.7 4.6 4.9 4.7 4.8 5.1 4.5 4.1 4.6 ...
 $ unemploy: int  2944 2945 2958 3143 3066 3018 2878 3001 2877 2709 ...

summary(econ)

       X             date                pce               pop            psavert      
 Min.   :  1.0   Length:574         Min.   :  506.7   Min.   :198712   Min.   : 2.200  
 1st Qu.:144.2   Class :character   1st Qu.: 1578.3   1st Qu.:224896   1st Qu.: 6.400  
 Median :287.5   Mode  :character   Median : 3936.8   Median :253060   Median : 8.400  
 Mean   :287.5                      Mean   : 4820.1   Mean   :257160   Mean   : 8.567  
 3rd Qu.:430.8                      3rd Qu.: 7626.3   3rd Qu.:290291   3rd Qu.:11.100  
 Max.   :574.0                      Max.   :12193.8   Max.   :320402   Max.   :17.300  
    uempmed          unemploy    
 Min.   : 4.000   Min.   : 2685  
 1st Qu.: 6.000   1st Qu.: 6284  
 Median : 7.500   Median : 7494  
 Mean   : 8.609   Mean   : 7771  
 3rd Qu.: 9.100   3rd Qu.: 8686  
 Max.   :25.200   Max.   :15352  

---

udaje <- read.csv("economics.csv")
model <- lm(pce ~ unemploy + psavert + pop, 
            data = udaje)

summary(model)

Call:
lm(formula = pce ~ unemploy + psavert + pop, data = udaje)

Residuals:
     Min       1Q   Median       3Q      Max 
-1045.64  -360.58   -71.61   392.67  1308.28 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -2.547e+04  4.344e+02 -58.617   <2e-16 ***
unemploy    -1.085e-01  1.175e-02  -9.238   <2e-16 ***
psavert      2.019e+02  1.477e+01  13.667   <2e-16 ***
pop          1.143e-01  1.467e-03  77.907   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 490.8 on 570 degrees of freedom
Multiple R-squared:  0.9811,    Adjusted R-squared:  0.981 
F-statistic:  9841 on 3 and 570 DF,  p-value: < 2.2e-16

Všetky premenné sú extrémne významné. Populácia a úspory majú silný pozitívny vzťah so spotrebou. Nezamestnanosť má slabý, ale štatisticky významný negatívny vplyv.

Vizualizácia vybraných premenných

plot(econ$unemploy, type="l", main="Nezamestnanosť")

plot(econ$pce, type="l", main="Osobná spotreba (PCE)")

Krivka nezamestnanosti má „cik-cak“ vzor. Môžme vidieť sezónnosť.

Krivka osobnej spotreby ide neustále hore – od hodnoty približne 1000 až po viac ako 12 000. To znamená, že spotreba domácností v USA dlhodobo rástla. Žiadna sezónnosť - Čiarkovaná línia, nemá žiadny opakujúci sa „cik-cak“ vzor narozdiel od krivky nezamestnanosti.

Jednoduchý regresný model

Vysvetlím miera nezamestnanosti (unemploy) pomocou osobnej spotreby (pce).

model1 <- lm(unemploy ~ pce, data=econ)
summary(model1)

Call:
lm(formula = unemploy ~ pce, data = econ)

Residuals:
    Min      1Q  Median      3Q     Max 
-3215.3 -1546.3  -388.8  1485.8  5493.2 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) 5571.1416   146.7269   37.97   <2e-16 ***
pce            0.4565     0.0245   18.63   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2086 on 572 degrees of freedom
Multiple R-squared:  0.3776,    Adjusted R-squared:  0.3765 
F-statistic: 347.1 on 1 and 572 DF,  p-value: < 2.2e-16

Koeficient pce je významný, ale pravdepodobne zachytáva len spoločný trend, nie reálny vzťah. R² je nízke (≈0.38), model nevysvetľuje nezamestnanosť dobre. Interpretácia je ekonomicky nelogická.

Viacnásobný regresný model

Vysvetlím nezamestnanosť pomocou viacerých premenných.

model2 <- lm(unemploy ~ pce + pop + uempmed, data=econ)
summary(model2)

Call:
lm(formula = unemploy ~ pce + pop + uempmed, data = econ)

Residuals:
    Min      1Q  Median      3Q     Max 
-3182.4  -733.1   -38.1   645.8  3749.1 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -2.659e+04  1.585e+03  -16.78   <2e-16 ***
pce         -1.510e+00  8.132e-02  -18.56   <2e-16 ***
pop          1.406e-01  7.531e-03   18.67   <2e-16 ***
uempmed      6.370e+02  1.560e+01   40.84   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1030 on 570 degrees of freedom
Multiple R-squared:  0.8488,    Adjusted R-squared:  0.848 
F-statistic:  1067 on 3 and 570 DF,  p-value: < 2.2e-16

Model ukazuje, že počet nezamestnaných najlepšie vysvetľuje dĺžka nezamestnanosti a populácia, zatiaľ čo negatívna väzba s PCE je pravdepodobne dôsledkom trendu v časovej rade, nie skutočného ekonomického vzťahu.

Multikolinearita – VIF

library(car)
vif(model2)

      pce       pop   uempmed 
45.187660 41.216070  2.216113 

Výstup ukazuje, že model má extrémne vysokú multikolinearitu medzi premennými pce a pop, čo spôsobuje, že ich koeficienty sú veľmi nespoľahlivé a ich vplyv na nezamestnanosť nemožno z modelu dôveryhodne interpretovať.

Korelačná matica

vars <- econ[, c("pce", "pop", "uempmed", "unemploy")]
round(cor(vars, use="complete.obs"), 3)

           pce   pop uempmed unemploy
pce      1.000 0.987   0.727    0.615
pop      0.987 1.000   0.695    0.634
uempmed  0.727 0.695   1.000    0.869
unemploy 0.615 0.634   0.869    1.000

Pce(spotreba domácností) a pop(populácia) sú veľmi silne korelované (logické – väčšia populácia spotrebuje viac).

Uempmed(doba nezamestnanosti) a unemploy(počet nezamestnaných) sú tiež veľmi silne korelované (dlhšia nezamestnanosť súvisí s väčším počtom nezamestnaných).

Všetky premenné sú medzi sebou pozitívne korelované, žiadna nie je negatívne korelovaná.

pairs(vars,
      main = "Scatterplotová matica – premenné pce, pop, uempmed, unemploy")

Condition Number

library(MASS)
x <- model.matrix(model2)[, -1]
sv <- svd(x)$d
condition_number <- max(sv) / min(sv)
condition_number

[1] 91271.39

Číslo podmienky ~ 91 271 je extrémne vysoké. Prakticky to znamená že matica je silne kolineárna – niektoré premenné sú takmer lineárne závislé.

Alternatívne riešenia multikolinearity

• centrovanie (odčítanie priemeru)- môže znížiť numerické problémy a stabilizovať výpočty. Aj keď neznižuje samotnú koreláciu, pomáha pri výpočte regresných koeficientov.

# Centrovanie
cent <- scale(vars, center = TRUE, scale = FALSE)
round(cor(cent), 3)

           pce   pop uempmed unemploy
pce      1.000 0.987   0.727    0.615
pop      0.987 1.000   0.695    0.634
uempmed  0.727 0.695   1.000    0.869
unemploy 0.615 0.634   0.869    1.000

Centrovanie nezmenilo lineárne vzťahy, len posunulo premenné na stred 0

PCA

vars <- econ[, c("pce", "pop", "uempmed", "unemploy")]

pca_res <- prcomp(vars, center = TRUE, scale. = TRUE)  

summary(pca_res)

Importance of components:
                          PC1    PC2    PC3     PC4
Standard deviation     1.8072 0.7743 0.3572 0.08321
Proportion of Variance 0.8165 0.1499 0.0319 0.00173
Cumulative Proportion  0.8165 0.9664 0.9983 1.00000

pca_res$rotation  

               PC1        PC2         PC3        PC4
pce      0.5120962  0.4817677 -0.08517838  0.7059759
pop      0.5099582  0.4850154  0.22550790 -0.6736831
uempmed  0.5026852 -0.4187833 -0.73740369 -0.1678212
unemploy 0.4743498 -0.5977294  0.63097277  0.1399470

head(pca_res$x)

           PC1        PC2        PC3       PC4
[1,] -2.803222 0.15410408 -0.6711345 0.1294587
[2,] -2.775348 0.13653351 -0.7056593 0.1182991
[3,] -2.781611 0.14724649 -0.6834953 0.1205157
[4,] -2.709410 0.07695554 -0.6918826 0.1137444
[5,] -2.744368 0.11794864 -0.6733345 0.1154366
[6,] -2.737427 0.12175595 -0.7019615 0.1074157

Dáta sú vysoko korelované. PC1 vysvetľuje väčšinu variability (81,6%) – dá sa považovať za hlavnú dimenziu dát. PC2 zachytáva kontrast medzi ekonomickými veličinami (pce, pop) a nezamestnanosťou (uempmed, unemploy). PC3 a PC4 sú zanedbateľné, takže dá sa znížiť dimenzionalita z 4 na 2 bez veľkej straty informácií.