1 Úvod

Tento notebook vysvetľuje tri nadväzujúce oblasti diagnostiky a nápravy lineárneho regresného modelu:

  1. Breusch–Paganov test pre heteroskedasticitu,
  2. Breusch–Godfreyho test pre autokoreláciu,
  3. robustné štandardné chyby a GLS odhad ako nápravné metódy.

Pedagogická logika notebooku je zámerne nasledovná:

  • uvedieme základný teoretický koncept
  • potom odhadneme model metódou OLS,
  • potom overíme predpoklady modelu,
  • následne ukážeme, ako opraviť inferenciu pomocou robustných štandardných chýb,
  • napokon ukážeme, ako zmeniť samotný odhad modelu pomocou GLS.

2 Teoretické východisko

Uvažujme lineárny regresný model v maticovom tvare:

\[ y = X\beta + u \]

Pri klasických predpokladoch platí:

\[ E(u) = 0, \qquad \mathrm{Var}(u) = \sigma^2 I \]

Za týchto predpokladov je OLS odhad daný vzťahom:

\[ \hat{\beta}_{OLS} = (X'X)^{-1}X'y \]

a jeho kovariančná matica je:

\[ \mathrm{Var}(\hat{\beta}_{OLS}) = \sigma^2 (X'X)^{-1} \]

Ak však neplatí

\[ \mathrm{Var}(u) = \sigma^2 I, \]

ale všeobecnejšie

\[ \mathrm{Var}(u) = \Omega, \]

kde \(\Omega\) nie je násobkom jednotkovej matice, potom:

  • OLS odhady môžu zostať nestranné,
  • ale už nemusia byť efektívne,
  • a klasické štandardné chyby nemusia byť správne.

Preto existujú dve základné stratégie:

  1. ponechať OLS koeficienty a opraviť iba štandardné chyby,
  2. zmeniť samotný spôsob odhadu modelu pomocou GLS.

3 Breusch–Paganov test (o heteroskedasticite)

3.1 Myšlienka testu

Uvažujme lineárny regresný model:

\[ y_i = \beta_0 + \beta_1 x_{1i} + \cdots + \beta_k x_{ki} + u_i \]

Za klasického predpokladu homoskedasticity platí:

\[ \mathrm{Var}(u_i) = \sigma^2 \]

To znamená, že rozptyl náhodnej zložky je konštantný pre všetky pozorovania. Breusch–Paganov test skúma, či sa tento rozptyl nemení v závislosti od vysvetľujúcich premenných.

3.2 Hypotézy

\[ H_0: \mathrm{Var}(u_i) = \sigma^2 \]

\[ H_1: \mathrm{Var}(u_i) \neq \sigma^2 \]

3.3 Postup testu

Najprv odhadneme pôvodný model metódou OLS a získame rezíduá:

\[ e_i = y_i - \hat{y}_i \]

Následne odhadneme pomocnú regresiu:

\[ e_i^2 = \alpha_0 + \alpha_1 z_{1i} + \alpha_2 z_{2i} + \cdots + \alpha_m z_{mi} + v_i \]

V praxi sa často volí:

\[ z_{ji} = x_{ji} \]

3.4 Testovacia štatistika

Breusch–Paganov test je založený na LM štatistike:

\[ LM = nR^2 \]

kde \(R^2\) pochádza z pomocnej regresie. Za platnosti nulovej hypotézy platí:

\[ LM \sim \chi^2(m) \]

kde \(m\) je počet regresorov v pomocnej regresii bez interceptu.

3.5 Rozhodovanie

  • veľká p-hodnota \(\Rightarrow\) nulovú hypotézu (o homoskedasticite) nezamietame,
  • malá p-hodnota \(\Rightarrow\) zamietame nulovú hypotézu a usudzujeme na heteroskedasticitu.

4 Breusch–Godfreyho test (o autokorelácii)

4.1 Myšlienka testu

Uvažujme model:

\[ y_t = \beta_0 + \beta_1 x_{1t} + \cdots + \beta_k x_{kt} + u_t \]

Predpokladajme, že chyba môže mať autoregresnú štruktúru rádu \(p\):

\[ u_t = \rho_1 u_{t-1} + \rho_2 u_{t-2} + \cdots + \rho_p u_{t-p} + \varepsilon_t \]

4.2 Hypotézy

\[ H_0: \rho_1 = \rho_2 = \cdots = \rho_p = 0 \]

\[ H_1: \text{aspoň jedno } \rho_j \neq 0 \]

4.3 Postup testu

Najprv odhadneme pôvodný model a získame rezíduá \(e_t\). Potom odhadneme pomocnú regresiu:

\[ e_t = \gamma_0 + \gamma_1 x_{1t} + \cdots + \gamma_k x_{kt} + \rho_1 e_{t-1} + \cdots + \rho_p e_{t-p} + v_t \] kde \(v_t\) je náhodná zložka pomocnej regresie.

4.4 Testovacia štatistika

Aj tu používame LM princíp:

\[ LM = nR^2 \]

kde \(n\) je počet pozorovaní v predchádzajúcej testovanej rovnici (pričom tu dochádza k redukcii o \(p\)), pričom za nulovej hypotézy platí:

\[ LM \sim \chi^2(p) \]

4.5 Rozhodovanie

  • veľká p-hodnota \(\Rightarrow\) nulovú hypotézu o neprítomnosti autokorelácie nezamietame,
  • malá p-hodnota \(\Rightarrow\) zamietame nulovú hypotézu a usudzujeme na autokoreláciu.

5 Robustné štandardné chyby

5.1 Základná idea

Ak sú OLS koeficienty stále konzistentné, ale problém je v nesprávnej kovariančnej matici, môžeme ponechať odhad \(\hat{\beta}_{OLS}\) a korigovať iba štandardné chyby.

Všeobecný robustný tvar kovariančnej matice je:

\[ \widehat{\mathrm{Var}}(\hat{\beta}) = (X'X)^{-1} X' \hat{\Omega} X (X'X)^{-1} \]

kde \(\hat{\Omega}\) zachytáva heteroskedasticitu alebo aj autokoreláciu.

Pripomínam, že pri zanedbaní narušenia vyššie uvedených predpokladov by sme túto maticu odhadovali vzorcom

\[ \widehat{\mathrm{Var}}(\hat{\beta}) = \sigma^2(X'X)^{-1} \]

a dopustili by sme sa chyby v odhade štandardných chýb, čo by mohlo viesť k nesprávnym záverom o štatistickej významnosti koeficientov.

5.2 Whiteove robustné štandardné chyby

Whiteove odhady sú robustné voči heteroskedasticite. Koeficienty regresie sa nemenia, menia sa teda len:

  • štandardné chyby,
  • t-štatistiky,
  • p-hodnoty.

5.3 HAC odhady

Ak chceme korigovať štandardné chyby súčasne na heteroskedasticitu aj autokoreláciu, používame HAC odhady, napríklad Newey–West.

6 GLS odhad

6.1 Generalized Least Squares (GLS) – teoretický základ

6.1.1 1. Motivácia

Uvažujme lineárny regresný model v maticovom tvare:

\[ y = X\beta + u \]

Pri klasických predpokladoch metódy najmenších štvorcov (OLS) platí:

\[ \mathrm{Var}(u) = \sigma^2 I \]

To znamená, že:

  • rozptyl chýb je konštantný (homoskedasticita),
  • chyby nie sú medzi sebou korelované.

Ak je však tento predpoklad porušený, teda:

\[ \mathrm{Var}(u) = \Omega \]

kde \(\Omega \neq \sigma^2 I\), potom:

  • OLS odhad zostáva (za určitých podmienok) nestranný,
  • ale nie je efektívny,
  • a jeho kovariančná matica je nesprávne špecifikovaná.

6.1.2 2. Myšlienka transformácie

Cieľom GLS je transformovať pôvodný model tak, aby náhodná zložka v transformovanom modeli spĺňala klasické predpoklady.

Uvažujme maticu \(P\), pre ktorú platí:

\[ P'P = \Omega^{-1} \]

Následne vynásobíme pôvodný model zľava touto maticou:

\[ Py = PX\beta + Pu \]

Označme:

\[ y^* = Py, \quad X^* = PX, \quad u^* = Pu \]

Dostávame transformovaný model:

\[ y^* = X^*\beta + u^* \]


6.1.3 3. Vlastnosti transformovaného modelu

Pre náhodnú zložku transformovaného modelu platí:

\[ \mathrm{Var}(u^*) = P \Omega P' \]

Keďže platí \(P'P = \Omega^{-1}\), dostávame:

\[ \mathrm{Var}(u^*) = I \]

Teda:

  • transformované chyby sú nekorelované,
  • majú konštantný rozptyl,

a preto spĺňajú klasické predpoklady metódy najmenších štvorcov.


6.1.4 4. GLS odhad

Na transformovaný model môžeme aplikovať OLS. Výsledný GLS odhad má tvar:

\[ \hat{\beta}_{GLS} = (X'\Omega^{-1}X)^{-1}X'\Omega^{-1}y \]


6.1.5 5. Intuícia

  • GLS „váži“ jednotlivé pozorovania podľa variability chýb
  • pozorovania s vyšším rozptylom dostávajú menšiu váhu
  • zohľadňuje aj koreláciu medzi chybami

6.1.6 6. Praktická poznámka

V praxi matica \(\Omega\) zvyčajne nie je známa, preto sa používa:

  • FGLS (Feasible GLS)\(\Omega\) sa najprv odhaduje
  • implementácia napr. pomocou gls() z balíka nlme

Interpretácia :
GLS nemení iba štandardné chyby ako robustné metódy, ale mení samotný odhad modelu.
Využíva informáciu o štruktúre chýb, čím dosahuje efektívnejšie odhady než OLS. ## Intuícia GLS

GLS možno chápať ako transformáciu modelu tak, aby v transformovanom modeli mala chyba „sférický“ tvar. Inými slovami, problémová kovariančná štruktúra sa zohľadní už pri odhade, nie až dodatočne pri výpočte štandardných chýb.

6.2 OLS s robustnými SE verzus GLS

  • OLS + White/HAC SE: koeficienty sa nemenia, opravuje sa iba inferencia.
  • GLS: menia sa aj samotné odhady koeficientov, pretože modelujeme štruktúru chýb.

7 Použité balíky

library(lmtest)
library(sandwich)
library(nlme)

8 Príklad 1: Heteroskedasticita

Vytvoríme dáta s heteroskedasticitou, kde rozptyl chyby rastie s premennou \(x\).

set.seed(123)

n <- 200
x <- runif(n, 1, 10)
u <- rnorm(n, mean = 0, sd = 0.5 * x)
y <- 2 + 3 * x + u

data_bp <- data.frame(y, x)

m_bp <- lm(y ~ x, data = data_bp)
summary(m_bp)
## 
## Call:
## lm(formula = y ~ x, data = data_bp)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -7.3828 -1.5748 -0.1864  1.3994  9.6110 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  2.16353    0.46233    4.68 5.32e-06 ***
## x            2.96535    0.07607   38.98  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.646 on 198 degrees of freedom
## Multiple R-squared:  0.8847, Adjusted R-squared:  0.8841 
## F-statistic:  1520 on 1 and 198 DF,  p-value: < 2.2e-16

8.1 Breusch–Paganov test

bp_res <- bptest(m_bp)
bp_res
## 
##  studentized Breusch-Pagan test
## 
## data:  m_bp
## BP = 24.741, df = 1, p-value = 6.557e-07

8.1.1 Interpretácia

Ak je p-hodnota testu menšia než zvolená hladina významnosti, napríklad \(0.05\), zamietame nulovú hypotézu homoskedasticity. To znamená, že rozptyl rezíduí pravdepodobne závisí od vysvetľujúcej premennej a v modeli je prítomná heteroskedasticita.

8.2 Pomocná regresia

ehat2 <- residuals(m_bp)^2
aux_bp <- lm(ehat2 ~ x, data = data_bp)
summary(aux_bp)
## 
## Call:
## lm(formula = ehat2 ~ x, data = data_bp)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -14.098  -6.205  -2.066   2.126  78.775 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   -2.615      1.975  -1.324    0.187    
## x              1.718      0.325   5.287 3.27e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 11.3 on 198 degrees of freedom
## Multiple R-squared:  0.1237, Adjusted R-squared:  0.1193 
## F-statistic: 27.95 on 1 and 198 DF,  p-value: 3.273e-07
n * summary(aux_bp)$r.squared
## [1] 24.7411

8.2.1 Interpretácia

Ak pomocná regresia vysvetľuje podstatnú časť variability \(e_i^2\), znamená to, že rozptyl chýb nie je konštantný. Vysvetľujúca premenná teda pomáha predpovedať veľkosť chýb.

8.3 Whiteove robustné štandardné chyby

8.3.1 Praktická implementácia v R

white_hc0 <- coeftest(m_bp, vcov = vcovHC(m_bp, type = "HC0"))
white_hc1 <- coeftest(m_bp, vcov = vcovHC(m_bp, type = "HC1"))

white_hc0
## 
## t test of coefficients:
## 
##             Estimate Std. Error t value  Pr(>|t|)    
## (Intercept) 2.163533   0.373375  5.7945 2.664e-08 ***
## x           2.965347   0.080514 36.8301 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
white_hc1
## 
## t test of coefficients:
## 
##             Estimate Std. Error t value  Pr(>|t|)    
## (Intercept)  2.16353    0.37526  5.7655 3.088e-08 ***
## x            2.96535    0.08092 36.6455 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

8.3.2 Interpretácia

Koeficienty regresie zostávajú rovnaké ako pri OLS, ale štandardné chyby sa upravia tak, aby boli robustné voči heteroskedasticite. Ak koeficient zostáva významný aj po tejto úprave, jeho významnosť je robustnejšia.

8.4 GLS pri heteroskedasticite

V balíku nlme môžeme modelovať heteroskedasticitu napríklad pomocou štruktúry varPower.

gls_het <- gls(
  y ~ x,
  data = data_bp,
  weights = varPower(form = ~ x)
)

summary(gls_het)
## Generalized least squares fit by REML
##   Model: y ~ x 
##   Data: data_bp 
##        AIC      BIC   logLik
##   918.1181 931.2711 -455.059
## 
## Variance function:
##  Structure: Power of variance covariate
##  Formula: ~x 
##  Parameter estimates:
##    power 
## 0.750207 
## 
## Coefficients:
##                Value  Std.Error  t-value p-value
## (Intercept) 2.142065 0.26107494  8.20479       0
## x           2.969510 0.06265433 47.39513       0
## 
##  Correlation: 
##   (Intr)
## x -0.85 
## 
## Standardized residuals:
##        Min         Q1        Med         Q3        Max 
## -2.0121760 -0.6537245 -0.1084361  0.6083592  3.1140187 
## 
## Residual standard error: 0.7054642 
## Degrees of freedom: 200 total; 198 residual

8.4.1 Interpretácia

Na rozdiel od Whiteových štandardných chýb teraz nemeníme iba inferenciu, ale aj samotný odhad modelu. GLS priamo zohľadňuje, že rozptyl chýb sa mení s premennou \(x\), a preto môžu byť odhady efektívnejšie.

9 Príklad 2: Autokorelácia

Teraz vytvoríme dáta s autokorelovanými chybami prvého rádu.

set.seed(123)

n <- 200
x <- rnorm(n)
u <- numeric(n)
eps <- rnorm(n, sd = 1)

rho <- 0.7
u[1] <- eps[1]

for (t in 2:n) {
  u[t] <- rho * u[t - 1] + eps[t]
}

y <- 1 + 2 * x + u

data_bg <- data.frame(
  y = y,
  x = x,
  t = 1:n
)

m_bg <- lm(y ~ x, data = data_bg)
summary(m_bg)
## 
## Call:
## lm(formula = y ~ x, data = data_bg)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.7937 -0.9739 -0.0469  0.8434  3.3399 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.15262    0.09073   12.70   <2e-16 ***
## x            2.06205    0.09644   21.38   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.283 on 198 degrees of freedom
## Multiple R-squared:  0.6978, Adjusted R-squared:  0.6963 
## F-statistic: 457.2 on 1 and 198 DF,  p-value: < 2.2e-16

9.1 Breusch–Godfreyho test

bg_res_1 <- bgtest(m_bg, order = 1)
bg_res_1
## 
##  Breusch-Godfrey test for serial correlation of order up to 1
## 
## data:  m_bg
## LM test = 79.796, df = 1, p-value < 2.2e-16
bg_res_2 <- bgtest(m_bg, order = 2)
bg_res_2
## 
##  Breusch-Godfrey test for serial correlation of order up to 2
## 
## data:  m_bg
## LM test = 80.315, df = 2, p-value < 2.2e-16

9.1.1 Interpretácia

Ak je p-hodnota malá, zamietame nulovú hypotézu neprítomnosti autokorelácie. To znamená, že súčasné rezíduum súvisí s minulými rezíduami a klasické OLS štandardné chyby nemusia byť spoľahlivé.

9.2 Pomocná regresia pre BG test

ehat <- residuals(m_bg)

aux_bg_data <- data.frame(
  ehat = ehat[-1],
  x = data_bg$x[-1],
  ehat_lag1 = ehat[-length(ehat)]
)

aux_bg <- lm(ehat ~ x + ehat_lag1, data = aux_bg_data)
summary(aux_bg)
## 
## Call:
## lm(formula = ehat ~ x + ehat_lag1, data = aux_bg_data)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.51584 -0.69443 -0.03918  0.66062  2.39618 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.01436    0.07011  -0.205    0.838    
## x           -0.07870    0.07477  -1.053    0.294    
## ehat_lag1    0.63547    0.05515  11.522   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.989 on 196 degrees of freedom
## Multiple R-squared:  0.4038, Adjusted R-squared:  0.3978 
## F-statistic: 66.39 on 2 and 196 DF,  p-value: < 2.2e-16
nrow(aux_bg_data) * summary(aux_bg)$r.squared
## [1] 80.36491

9.2.1 Interpretácia

Ak je koeficient pri oneskorenom rezíduu významný, znamená to, že minulé chyby pomáhajú vysvetliť súčasné chyby. To je typický znak autokorelácie.

9.3 HAC štandardné chyby

hac_nw <- coeftest(m_bg, vcov = NeweyWest(m_bg))
hac_nw
## 
## t test of coefficients:
## 
##             Estimate Std. Error t value  Pr(>|t|)    
## (Intercept)  1.15262    0.19178  6.0101 8.768e-09 ***
## x            2.06205    0.10274 20.0713 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

9.3.1 Interpretácia

Aj tu zostávajú koeficienty rovnaké ako pri OLS, ale štandardné chyby sú upravené tak, aby boli robustné voči heteroskedasticite aj autokorelácii. Výsledná inferencia je preto dôveryhodnejšia.

9.4 GLS pri autokorelácii

Pri autokorelácii môžeme v gls() špecifikovať napríklad korelačnú štruktúru AR(1).

gls_ar1 <- gls(
  y ~ x,
  data = data_bg,
  correlation = corAR1(form = ~ t)
)

summary(gls_ar1)
## Generalized least squares fit by REML
##   Model: y ~ x 
##   Data: data_bg 
##        AIC      BIC    logLik
##   574.9673 588.1203 -283.4836
## 
## Correlation Structure: AR(1)
##  Formula: ~t 
##  Parameter estimate(s):
##       Phi 
## 0.6476351 
## 
## Coefficients:
##                Value  Std.Error  t-value p-value
## (Intercept) 1.159708 0.19671994  5.89523       0
## x           1.981823 0.06075463 32.62011       0
## 
##  Correlation: 
##   (Intr)
## x 0.007 
## 
## Standardized residuals:
##         Min          Q1         Med          Q3         Max 
## -2.17738142 -0.71702998 -0.03610559  0.66145653  2.58582938 
## 
## Residual standard error: 1.298291 
## Degrees of freedom: 200 total; 198 residual

9.4.1 Interpretácia

GLS teraz priamo modeluje to, že chyby v čase spolu súvisia. Preto sa môžu zmeniť nielen štandardné chyby, ale aj samotné odhadnuté koeficienty.

10 Porovnanie prístupov

10.1 OLS verzus robustné štandardné chyby verzus GLS

coef_ols_bp <- coef(m_bp)
coef_gls_bp <- coef(gls_het)

coef_ols_bg <- coef(m_bg)
coef_gls_bg <- coef(gls_ar1)

coef_ols_bp
## (Intercept)           x 
##    2.163533    2.965347
coef_gls_bp
## (Intercept)           x 
##    2.142065    2.969510
coef_ols_bg
## (Intercept)           x 
##    1.152621    2.062053
coef_gls_bg
## (Intercept)           x 
##    1.159708    1.981823

10.1.1 Interpretácia

Pri OLS a OLS s robustnými štandardnými chybami ostávajú odhadnuté koeficienty rovnaké. Pri GLS sa však môžu zmeniť aj samotné koeficienty, pretože GLS využíva dodatočné informácie o štruktúre chýb.

11 Zhrnutie

11.1 Kedy použiť jednotlivé prístupy?

  • OLS používame ako základný odhad modelu.
  • Breusch–Paganov test používame na test heteroskedasticity.
  • Breusch–Godfreyho test používame na test autokorelácie.
  • Whiteove štandardné chyby použijeme, ak chceme ponechať OLS koeficienty, ale opraviť inferenciu pri heteroskedasticite.
  • HAC štandardné chyby použijeme, ak chceme ponechať OLS koeficienty, ale opraviť inferenciu pri heteroskedasticite aj autokorelácii.
  • GLS použijeme, ak chceme priamo modelovať štruktúru chýb a získať efektívnejší odhad modelu.

12 Záver

Diagnostika regresného modelu nekončí odhadom koeficientov. Po odhade je potrebné overiť, či sú splnené predpoklady o štruktúre chýb, a ak nie sú, zvoliť vhodnú nápravu.

V tomto notebooku sme ukázali tri úrovne práce s problémami rezíduí:

  1. diagnostika pomocou BP a BG testu,
  2. oprava inferencie pomocou Whiteových a HAC štandardných chýb,
  3. zmena samotného odhadu pomocou GLS.

Takto majú študenti k dispozícii ucelenú logickú cestu od diagnostiky až po nápravné metódy.

13 Appendix 1 - HAC

13.1 Whiteove robustné štandardné chyby: HC0, HC1, HC2, HC3

V prítomnosti heteroskedasticity je klasický OLS odhad rozptylu nekonzistentný. Na riešenie tohto problému používame heteroskedasticity-konzistentné (HC) odhady kovariančnej matice.


13.1.1 1. Základ: HC0 (Whiteov odhad)

Uvažujme lineárny regresný model:

\[ y = X\beta + u \]

Whiteov (1980) robustný odhad kovariančnej matice má tvar:

\[ \widehat{\mathrm{Var}}(\hat{\beta}) = (X'X)^{-1} \left( \sum_{i=1}^n x_i x_i' e_i^2 \right) (X'X)^{-1} \]

Tento odhad sa nazýva HC0.

  • Konzistentný pri heteroskedasticite
  • Neobsahuje korekciu pre malé výbery

13.1.2 2. Motivácia pre HC1, HC2, HC3

Aj keď je HC0 asymptoticky konzistentný, v malých výberoch má tendenciu:

  • Podhodnocovať skutočný rozptyl
  • Viesť k častejšiemu zamietaniu nulovej hypotézy

Na nápravu sa používajú alternatívne odhady, ktoré upravujú štvorce rezíduí \(e_i^2\).


13.1.3 3. Rôzne HC odhady

13.1.3.1 HC0

\[ e_i^2 \]

  • Bez úpravy
  • Čisto asymptotický odhad

13.1.3.2 HC1 (korekcia stupňov voľnosti)

\[ \frac{n}{n - k} e_i^2 \]

  • Zohľadňuje stupne voľnosti
  • Často predvolený v softvéroch

13.1.3.3 HC2 (korekcia leverage)

\[ \frac{e_i^2}{1 - h_i} \]

kde \(h_i\) je leverage (diagonálny prvok tzv. hat matice).

  • Koriguje skutočnosť, že pozorovania s vysokým leverage majú tendenciu mať menšie rezíduá

13.1.3.4 HC3 (silnejšia korekcia leverage)

\[ \frac{e_i^2}{(1 - h_i)^2} \]

  • Silnejšia korekcia než HC2
  • Často preferovaný pri malých výberoch
  • Súvisí s jackknife odhadom rozptylu

13.1.3.5 HC4 a vyššie

  • Používajú ešte agresívnejšie úpravy založené na leverage
  • Určené pre extrémne hodnoty leverage
  • Menej často používané v bežnej praxi

13.1.4 4. Zhrnutie

Odhad Úprava Typické použitie
HC0 Žiadna Veľké výbery
HC1 Stupne voľnosti Bežná prax
HC2 Korekcia leverage Stredné výbery
HC3 Silná korekcia leverage Malé výbery (odporúčaný)
HC4+ Extrémne leverage Špeciálne prípady

13.1.5 5. Praktické odporúčania

  • Vo veľkých výberoch: všetky HC odhady dávajú podobné výsledky
  • V malých výberoch: HC3 je často preferovaný
  • V aplikovanej práci:
    • HC1 sa často používa (napr. default v Stata)
    • HC3 je konzervatívnejší a robustnejší

Interpretácia :

Whiteove robustné štandardné chyby nemenia odhady koeficientov, ale upravujú ich štandardné chyby.
Ak sa po tejto úprave zmení štatistická významnosť, pôvodné OLS štandardné chyby boli nespoľahlivé.

14 Príloha 2 - Autokorelačná funkcia rezíduí

14.1 Autokorelácia rezíduí (príloha)

14.1.1 1. Základná myšlienka

Uvažujme lineárny regresný model:

\[ y_t = \beta_0 + \beta_1 x_{1t} + \cdots + \beta_k x_{kt} + u_t \]

Klasický predpoklad OLS je, že chyby sú nekorelované:

\[ \mathrm{Cov}(u_t, u_s) = 0 \quad \text{pre } t \neq s \]

Ak tento predpoklad neplatí, hovoríme o autokorelácii rezíduí.


14.1.2 2. Autokorelačný proces

Najčastejšie sa predpokladá autoregresná štruktúra chýb, napríklad AR(1):

\[ u_t = \rho u_{t-1} + \varepsilon_t, \quad |\rho| < 1 \]

kde:

  • \(\rho\) je parameter autokorelácie
  • \(\varepsilon_t \sim \text{i.i.d.}(0, \sigma^2)\)

14.1.3 3. Vlastnosti autokorelácie

Pre AR(1) proces platí:

\[ \mathrm{Cov}(u_t, u_{t-1}) = \rho \sigma^2 \]

a všeobecnejšie:

\[ \mathrm{Cov}(u_t, u_{t-s}) = \rho^s \sigma^2 \]

To znamená, že:

  • chyby sú medzi sebou korelované v čase
  • korelácia klesá geometricky so vzdialenosťou

14.1.4 4. Kovariančná matica chýb

Pre model s autokoreláciou má kovariančná matica tvar:

\[ \Omega = \begin{pmatrix} \sigma^2 & \rho\sigma^2 & \rho^2\sigma^2 & \cdots \\ \rho\sigma^2 & \sigma^2 & \rho\sigma^2 & \cdots \\ \rho^2\sigma^2 & \rho\sigma^2 & \sigma^2 & \cdots \\ \vdots & \vdots & \vdots & \ddots \end{pmatrix} \]

Táto matica:

  • nie je diagonálna
  • zachytáva koreláciu medzi chybami

14.1.5 5. Dôsledky pre OLS

Pri autokorelácii platí:

  • OLS odhady \(\hat{\beta}\) zostávajú nestranné (pri splnení \(E(u|X)=0\))
  • ale nie sú efektívne
  • štandardné chyby sú nesprávne → testy sú nespoľahlivé

14.1.6 6. Riešenia problému

Existujú dva hlavné prístupy:

14.1.6.1 (a) Oprava inferencie

Použijeme HAC (Newey–West) štandardné chyby:

coeftest(m_bg, vcov = NeweyWest(m_bg))
## 
## t test of coefficients:
## 
##             Estimate Std. Error t value  Pr(>|t|)    
## (Intercept)  1.15262    0.19178  6.0101 8.768e-09 ***
## x            2.06205    0.10274 20.0713 < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

alebo GLS

library(nlme)

gls_ar1 <- gls(
  y ~ x,
  data = data_bg,
  correlation = corAR1(form = ~ t)
)

summary(gls_ar1)
## Generalized least squares fit by REML
##   Model: y ~ x 
##   Data: data_bg 
##        AIC      BIC    logLik
##   574.9673 588.1203 -283.4836
## 
## Correlation Structure: AR(1)
##  Formula: ~t 
##  Parameter estimate(s):
##       Phi 
## 0.6476351 
## 
## Coefficients:
##                Value  Std.Error  t-value p-value
## (Intercept) 1.159708 0.19671994  5.89523       0
## x           1.981823 0.06075463 32.62011       0
## 
##  Correlation: 
##   (Intr)
## x 0.007 
## 
## Standardized residuals:
##         Min          Q1         Med          Q3         Max 
## -2.17738142 -0.71702998 -0.03610559  0.66145653  2.58582938 
## 
## Residual standard error: 1.298291 
## Degrees of freedom: 200 total; 198 residual