Cvičenie 6 – Diagnostika regresného modelu

Cieľ cvičenia

Cieľom cvičenia je odhadnúť lineárny regresný model a overiť, či sú splnené základné predpoklady klasického lineárneho regresného modelu. Skúmame, ako ziskovosť firmy, veľkosť firmy a zadlženosť vplývajú na hodnotu ESG indexu.

Uvažovaný model má tvar:

\[ ESG_i = \beta_0 + \beta_1 ROA_i + \beta_2 SIZE_i + \beta_3 DEBT_i + \varepsilon_i \]

kde:

\(ESG_i\) je ESG index firmy,
\(ROA_i\) je rentabilita aktív,
\(SIZE_i\) je veľkosť firmy,
\(DEBT_i\) je zadlženosť firmy,
\(\varepsilon_i\) je náhodná zložka.

Ekonomicky očakávame:

\(\beta_1 > 0\): ziskovejšie firmy majú viac zdrojov na ESG aktivity,
\(\beta_2 > 0\): väčšie firmy sú pod väčším tlakom investorov a regulácie,
\(\beta_3 < 0\): viac zadlžené firmy môžu ESG aktivity odkladať z dôvodu finančných obmedzení.

Načítanie balíkov

library(lmtest)
library(car)
library(tseries)
library(dplyr)
library(janitor)

Načítanie a úprava dát

udaje <- read.csv2("data.csv", header = TRUE, sep = ";", dec = ".")

udaje <- udaje %>%
  clean_names() %>%
  rename(
    market_cap = market_capitalization,
    ROA = return_on_assets,
    DEBT = debt_to_asset,
    SIZE = firm_size,
    SDI = social_disclosure_index,
    EDI = environmental_disclosure_index,
    GDI = governance_disclosure_index,
    ESG = esg_index
  ) %>%
  mutate(
    DEBT = as.numeric(gsub(",", "", DEBT)),
    ROA = as.numeric(gsub(",", "", ROA)),
    market_cap = as.numeric(gsub(",", "", market_cap)),
    SIZE = as.numeric(gsub(",", "", SIZE)),
    SDI = as.numeric(gsub(",", "", SDI)),
    EDI = as.numeric(gsub(",", "", EDI)),
    GDI = as.numeric(gsub(",", "", GDI)),
    ESG = as.numeric(gsub(",", "", ESG))
  )

# ak zadanie vyžaduje iba jeden rok, stačilo by odkomentovať nasledujúci riadok
# udaje <- udaje %>% filter(years == 2013)

udaje %>% summarise_all(~ sum(is.na(.)))

##   years companies exchange_sector primary_business tobin_q market_cap ROA DEBT
## 1     1         0               0                0       0          2   0    2
##   SIZE SDI EDI GDI ESG
## 1    2   1   1   2   1

Vo výsledku sa objavilo niekoľko chýbajúcich hodnôt, preto bolo pri odhade modelu automaticky vyradených 5 pozorovaní. Napriek tomu zostáva rozsah vzorky dostatočný na odhad regresného modelu.

Odhad regresného modelu

model <- lm(ESG ~ ROA + SIZE + DEBT, data = udaje)
summary(model)

## 
## Call:
## lm(formula = ESG ~ ROA + SIZE + DEBT, data = udaje)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.29734 -0.08720 -0.01507  0.06931  0.55865 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -3.297e-01  3.680e-02  -8.961  < 2e-16 ***
## ROA         -3.649e-07  1.090e-05  -0.033  0.97331    
## SIZE         8.731e-02  4.893e-03  17.842  < 2e-16 ***
## DEBT        -2.835e-04  1.030e-04  -2.753  0.00606 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.122 on 751 degrees of freedom
##   (5 observations deleted due to missingness)
## Multiple R-squared:  0.3128, Adjusted R-squared:   0.31 
## F-statistic: 113.9 on 3 and 751 DF,  p-value: < 2.2e-16

Na základe výstupu modelu dostávame tieto hlavné výsledky:

premenná SIZE má štatisticky významný pozitívny vplyv na ESG index,
premenná DEBT má štatisticky významný negatívny vplyv na ESG index,
premenná ROA nie je štatisticky významná.

Konkrétne platí, že pri zvýšení premennej SIZE o jednu jednotku sa ESG index v priemere zvýši približne o 0.0907 bodu, pri nezmenených ostatných premenných. Naopak, vyššia zadlženosť je spojená s miernym poklesom ESG indexu. Premenná ROA sa v tomto modeli ako významný determinant ESG nepotvrdila.

Koeficient determinácie je približne \(R^2 = 0.3209\), čo znamená, že model vysvetľuje asi 32.1 % variability ESG indexu. Model ako celok je štatisticky významný, keďže p-hodnota F-testu je menšia ako 0.001.

Diagnostické grafy

Residuals vs Fitted

plot(model, which = 1)

Residuals vs Fitted

Reziduály sú rozložené približne okolo nulovej hodnoty, čo naznačuje, že model ako celok nemá výrazné systematické skreslenie. Väčšina bodov sa nachádza v pomerne úzkom pásme okolo nuly. Červená LOESS krivka však vykazuje mierne zakrivenie, čo môže naznačovať, že lineárny funkčný tvar nemusí úplne zachytiť vzťah medzi vysvetľujúcimi premennými a ESG indexom. Zároveň možno pozorovať niekoľko bodov s vyššími rezíduami, ktoré môžu predstavovať odľahlé alebo vplyvné pozorovania.

Q-Q plot

plot(model, which = 2)

Normal Q-Q plot

V centrálnej časti Q-Q grafu ležia body relatívne blízko referenčnej priamky, čo naznačuje, že stredná časť rozdelenia rezíduí je približne normálna. Na koncoch, najmä v pravom chvoste, sa však body od priamky výraznejšie odchyľujú. To znamená, že reziduály nie sú dokonale normálne rozdelené a v dátach sa vyskytujú extrémnejšie kladné reziduály.

Scale-Location plot

plot(model, which = 3)

Scale-Location plot

Scale-Location graf ukazuje, že rozptyl rezíduí je približne stabilný, no červená krivka mierne rastie. To naznačuje slabú heteroskedasticitu, teda mierny nárast variability rezíduí pri vyšších vyrovnaných hodnotách. Nejde však o veľmi silné porušenie predpokladu konštantného rozptylu.

Residuals vs Leverage

plot(model, which = 5)

Residuals vs Leverage

Väčšina pozorovaní má nízky leverage, čo je priaznivé. V grafe sa však objavuje jedno výraznejšie pozorovanie s vysokou hodnotou leverage, ktoré môže byť potenciálne vplyvné. Hoci žiadny bod zjavne neprekračuje najvyššie hranice Cookovej vzdialenosti, toto pozorovanie si zaslúži dodatočnú pozornosť.

Test normality rezíduí

Shapiro-Wilkov test

shapiro.test(residuals(model))

## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(model)
## W = 0.95033, p-value = 2.972e-15

Shapiro-Wilkov test testuje nulovú hypotézu, že reziduály pochádzajú z normálneho rozdelenia. V našom prípade je p-hodnota veľmi malá (menšia ako 0.05), preto nulovú hypotézu zamietame. Reziduály teda nie sú normálne rozdelené.

Jarque-Bera test

jarque.bera.test(residuals(model))

## 
##  Jarque Bera Test
## 
## data:  residuals(model)
## X-squared = 217.2, df = 2, p-value < 2.2e-16

Aj Jarque-Bera test vedie k rovnakému záveru, keďže p-hodnota je taktiež veľmi malá. Predpoklad normality rezíduí teda nie je splnený. Výsledok je v súlade s Q-Q grafom, ktorý ukázal odchýlky najmä v chvostoch rozdelenia.

Test heteroskedasticity

Breusch-Pagan test

bptest(model)

## 
##  studentized Breusch-Pagan test
## 
## data:  model
## BP = 13.976, df = 3, p-value = 0.002937

Breusch-Pagan test skúma nulovú hypotézu homoskedasticity. Keďže p-hodnota testu je 0.002937, nulovú hypotézu zamietame a usudzujeme, že v modeli je prítomná heteroskedasticita.

Test autokorelácie

Durbin-Watson test

dwtest(model)

## 
##  Durbin-Watson test
## 
## data:  model
## DW = 0.57013, p-value < 2.2e-16
## alternative hypothesis: true autocorrelation is greater than 0

Durbin-Watson test dáva hodnotu približne 0.57013, čo je výrazne menej ako 2. To by štandardne naznačovalo kladnú autokoreláciu rezíduí.

Breusch-Godfrey test

bgtest(model, order = 1)

## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  model
## LM test = 385.91, df = 1, p-value < 2.2e-16

Aj Breusch-Godfrey test potvrdzuje autokoreláciu rezíduí, keďže p-hodnota je menšia ako 0.001.

Poznámka: tieto testy sa interpretujú najmä pri časových alebo panelových dátach s prirodzeným usporiadaním pozorovaní. Ak ide o čisto prierezové dáta firiem, interpretácia autokorelácie je menej prirodzená a tento výsledok treba vnímať opatrne.

Odľahlé a vplyvné pozorovania

Test na odľahlé hodnoty

outlierTest(model)

##     rstudent unadjusted p-value Bonferroni p
## 620 4.646116         3.9946e-06    0.0030159
## 40  4.366480         1.4401e-05    0.0108720

Test identifikoval najmenej dve pozorovania s výrazne veľkými študentizovanými rezíduami. Tieto pozorovania možno považovať za potenciálne odľahlé hodnoty a bolo by vhodné ich vecne skontrolovať v pôvodných dátach.

Cookova vzdialenosť

cd <- cooks.distance(model)
head(sort(cd, decreasing = TRUE), 10)

##        163        744        745        570        569        568        743 
## 2.96745709 0.03707414 0.03189280 0.02051884 0.01890510 0.01834270 0.01333937 
##         20        100         70 
## 0.01251919 0.01092030 0.01045708

Najvyššiu Cookovu vzdialenosť má pozorovanie 163 s hodnotou približne 2.96745709, čo je veľmi vysoká hodnota. To naznačuje, že toto pozorovanie môže mať výrazný vplyv na odhad regresných koeficientov. Ostatné hodnoty Cookovej vzdialenosti sú už podstatne menšie. Pri finálnej interpretácii modelu je preto vhodné zvážiť aj kontrolný odhad bez tohto pozorovania.

Krátke teoretické poznámky

Štandardizované rezíduá sú rezíduá vydelené odhadovanou smerodajnou odchýlkou, vďaka čomu sú navzájom porovnateľné.

Leverage vyjadruje, ako veľmi sa pozorovanie odlišuje v priestore vysvetľujúcich premenných od ostatných pozorovaní. Súvisí s diagonálnymi prvkami tzv. hat matice:

\[ H = X(X'X)^{-1}X' \]

pričom platí:

\[ \hat{y} = Hy \]

Pre i-te pozorovanie je leverage definovaný ako:

\[ h_i = x_i'(X'X)^{-1}x_i \]

LOESS krivka v diagnostických grafoch predstavuje lokálne vyhladený trend. Ak je približne horizontálna, model je pravdepodobne špecifikovaný vhodne. Jej zakrivenie môže naznačovať nelinearitu alebo inú systematickú chybu modelu.

Záver

Cieľom cvičenia bolo odhadnúť regresný model pre ESG index a overiť základné predpoklady lineárnej regresie. Z výsledkov vyplynulo, že štatisticky významný pozitívny vplyv na ESG má veľkosť firmy, zatiaľ čo zadlženosť pôsobí negatívne. Premenná ROA sa v modeli ako významná nepotvrdila.

Pri diagnostike modelu sa ukázalo, že nie všetky predpoklady sú úplne splnené. Reziduály nie sú úplne normálne rozdelené, objavila sa mierna heteroskedasticita a v dátach sa nachádzajú aj niektoré odľahlé alebo vplyvné pozorovania.

Celkovo však model zachytáva základné vzťahy medzi ESG indexom, veľkosťou firmy a zadlženosťou, aj keď výsledky treba interpretovať s určitou opatrnosťou.