# NENAČÍTAVAM žiadne balíky cez library(), aby to nespadlo.
# Len si overím, či náhodou nie sú dostupné.
have_lmtest   <- requireNamespace("lmtest",   quietly = TRUE)
have_sandwich <- requireNamespace("sandwich", quietly = TRUE)

have_lmtest
## [1] FALSE
have_sandwich
## [1] FALSE

1 Údaje a ich príprava

V tomto cvičení používam údaje z Eurostatu
súbor: une_rt_a__custom_18708117_linear.csv.
Zameriavam sa na ročnú mieru nezamestnanosti na Slovensku v rokoch 2014–2024.

# načítanie údajov
udaje <- read.csv("une_rt_a__custom_18708117_linear.csv",
                  stringsAsFactors = FALSE)

str(udaje)
## 'data.frame':    33 obs. of  11 variables:
##  $ DATAFLOW   : chr  "ESTAT:UNE_RT_A(1.0)" "ESTAT:UNE_RT_A(1.0)" "ESTAT:UNE_RT_A(1.0)" "ESTAT:UNE_RT_A(1.0)" ...
##  $ LAST.UPDATE: chr  "11/09/25 23:00:00" "11/09/25 23:00:00" "11/09/25 23:00:00" "11/09/25 23:00:00" ...
##  $ freq       : chr  "Annual" "Annual" "Annual" "Annual" ...
##  $ age        : chr  "From 15 to 74 years" "From 15 to 74 years" "From 15 to 74 years" "From 15 to 74 years" ...
##  $ unit       : chr  "Percentage of population in the labour force" "Percentage of population in the labour force" "Percentage of population in the labour force" "Percentage of population in the labour force" ...
##  $ sex        : chr  "Total" "Total" "Total" "Total" ...
##  $ geo        : chr  "Czechia" "Czechia" "Czechia" "Czechia" ...
##  $ TIME_PERIOD: int  2014 2015 2016 2017 2018 2019 2020 2021 2022 2023 ...
##  $ OBS_VALUE  : num  6.1 5.1 4 2.9 2.2 2 2.6 2.8 2.2 2.6 ...
##  $ OBS_FLAG   : chr  "" "" "" "" ...
##  $ CONF_STATUS: logi  NA NA NA NA NA NA ...
head(udaje)
##              DATAFLOW       LAST.UPDATE   freq                 age
## 1 ESTAT:UNE_RT_A(1.0) 11/09/25 23:00:00 Annual From 15 to 74 years
## 2 ESTAT:UNE_RT_A(1.0) 11/09/25 23:00:00 Annual From 15 to 74 years
## 3 ESTAT:UNE_RT_A(1.0) 11/09/25 23:00:00 Annual From 15 to 74 years
## 4 ESTAT:UNE_RT_A(1.0) 11/09/25 23:00:00 Annual From 15 to 74 years
## 5 ESTAT:UNE_RT_A(1.0) 11/09/25 23:00:00 Annual From 15 to 74 years
## 6 ESTAT:UNE_RT_A(1.0) 11/09/25 23:00:00 Annual From 15 to 74 years
##                                           unit   sex     geo TIME_PERIOD
## 1 Percentage of population in the labour force Total Czechia        2014
## 2 Percentage of population in the labour force Total Czechia        2015
## 3 Percentage of population in the labour force Total Czechia        2016
## 4 Percentage of population in the labour force Total Czechia        2017
## 5 Percentage of population in the labour force Total Czechia        2018
## 6 Percentage of population in the labour force Total Czechia        2019
##   OBS_VALUE OBS_FLAG CONF_STATUS
## 1       6.1                   NA
## 2       5.1                   NA
## 3       4.0                   NA
## 4       2.9                   NA
## 5       2.2                   NA
## 6       2.0                   NA

Vyberiem len Slovensko a vytvorím časový rad nezamestnanosti:

# filtrovanie na Slovensko a výber stĺpcov
udaje_svk <- udaje[udaje$geo == "Slovakia",
                   c("TIME_PERIOD", "OBS_VALUE")]

# premenujem stĺpce
names(udaje_svk) <- c("Year", "Unemployment")

# zoradím podľa rokov
udaje_svk <- udaje_svk[order(udaje_svk$Year), ]

udaje_svk
##    Year Unemployment
## 23 2014         13.1
## 24 2015         11.5
## 25 2016          9.6
## 26 2017          8.1
## 27 2018          6.5
## 28 2019          5.7
## 29 2020          6.7
## 30 2021          6.8
## 31 2022          6.1
## 32 2023          5.8
## 33 2024          5.3

Základné štatistiky:

summary(udaje_svk$Unemployment)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   5.300   5.950   6.700   7.745   8.850  13.100

A jednoduchý graf vývoja nezamestnanosti v čase (base R):

plot(udaje_svk$Year, udaje_svk$Unemployment,
     type = "b",
     xlab = "Rok",
     ylab = "Nezamestnanosť (%)",
     main = "Miera nezamestnanosti na Slovensku")

2 Základný lineárny regresný model

Ako jednoduchý základný model používam regresiu miery nezamestnanosti na časový trend:

\[ Unemployment_t = \beta_0 + \beta_1 \, Year_t + e_t. \]

model <- lm(Unemployment ~ Year, data = udaje_svk)
summary(model)
## 
## Call:
## lm(formula = Unemployment ~ Year, data = udaje_svk)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -2.0455 -0.6850  0.3918  0.8605  1.9591 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1378.8300   255.0566   5.406 0.000430 ***
## Year          -0.6791     0.1263  -5.376 0.000447 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.325 on 9 degrees of freedom
## Multiple R-squared:  0.7625, Adjusted R-squared:  0.7361 
## F-statistic:  28.9 on 1 and 9 DF,  p-value: 0.0004471

Do grafu doplním regresnú priamku:

plot(udaje_svk$Year, udaje_svk$Unemployment,
     xlab = "Rok",
     ylab = "Nezamestnanosť (%)",
     main = "Nezamestnanosť na Slovensku a lineárny trend")
points(udaje_svk$Year, udaje_svk$Unemployment, pch = 16)
abline(model, lty = 2, col = "blue")

V texte viem stručne opísať:

3 Reziduá a autokorelácia

Pri časových radoch ma zaujíma, či sú reziduá nezávislé v čase.

udaje_svk$residuals <- residuals(model)
udaje_svk$fitted    <- fitted(model)

head(udaje_svk)
##    Year Unemployment  residuals    fitted
## 23 2014         13.1  1.9590909 11.140909
## 24 2015         11.5  1.0381818 10.461818
## 25 2016          9.6 -0.1827273  9.782727
## 26 2017          8.1 -1.0036364  9.103636
## 27 2018          6.5 -1.9245455  8.424545
## 28 2019          5.7 -2.0454545  7.745455

3.1 Graf rezíduí v čase a ACF

# reziduá v čase
plot(udaje_svk$Year, udaje_svk$residuals,
     type = "b",
     xlab = "Rok",
     ylab = "Reziduá",
     main = "Reziduá základného modelu v čase")
abline(h = 0, col = "red")

# ACF rezíduí (funkcia acf je v base/stats)
acf(udaje_svk$residuals,
    main = "ACF rezíduí – základný model")

Z grafov viem vizuálne posúdiť, či reziduá naznačujú autokoreláciu.

3.2 Durbin–Watsonov štatistický ukazovateľ (bez balíka)

Keďže balík lmtest nemusí byť dostupný, vypočítam si Durbin–Watsonov štatistický ukazovateľ sama z definície:

\[ DW = \frac{\sum_{t=2}^{T} (e_t - e_{t-1})^2}{\sum_{t=1}^{T} e_t^2}. \]

e <- udaje_svk$residuals
num <- sum((e[2:length(e)] - e[1:(length(e)-1)])^2)
den <- sum(e^2)
DW_stat <- num / den
DW_stat
## [1] 0.4733047

Hodnotu DW môžem slovne okomentovať (≈ 2 bez autokorelácie, < 2 pozitívna, > 2 negatívna autokorelácia).

Ak by bol balík lmtest náhodou k dispozícii, viem použiť aj klasický test:

lmtest::dwtest(model)

4 Dynamický model s oneskorenou nezamestnanosťou

Ak autokorelácia vyzerá ako problém, jednou z možností je pridať do modelu
oneskorenú hodnotu nezamestnanosti. Vytvorím lag o 1 rok:

\[ Unemployment_t = \alpha_0 + \alpha_1 \, Unemployment_{t-1} + \alpha_2 \, Year_t + u_t. \]

# oneskorená nezamestnanosť – jednoduchý lag o 1 rok
Unemp_lag1 <- c(NA, head(udaje_svk$Unemployment, -1))

udaje_dyn <- data.frame(
  Year         = udaje_svk$Year,
  Unemployment = udaje_svk$Unemployment,
  Unemp_lag1   = Unemp_lag1
)

# odstránim prvý riadok s NA
udaje_dyn <- udaje_dyn[!is.na(udaje_dyn$Unemp_lag1), ]

udaje_dyn
##    Year Unemployment Unemp_lag1
## 2  2015         11.5       13.1
## 3  2016          9.6       11.5
## 4  2017          8.1        9.6
## 5  2018          6.5        8.1
## 6  2019          5.7        6.5
## 7  2020          6.7        5.7
## 8  2021          6.8        6.7
## 9  2022          6.1        6.8
## 10 2023          5.8        6.1
## 11 2024          5.3        5.8

Odhadnem dynamický model:

model_dyn <- lm(Unemployment ~ Unemp_lag1 + Year,
                data = udaje_dyn)

summary(model_dyn)
## 
## Call:
## lm(formula = Unemployment ~ Unemp_lag1 + Year, data = udaje_dyn)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.8265 -0.2940 -0.1893  0.4054  1.1145 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)   
## (Intercept)  54.18923  297.51125   0.182  0.86063   
## Unemp_lag1    0.70369    0.17212   4.088  0.00464 **
## Year         -0.02605    0.14673  -0.178  0.86413   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6549 on 7 degrees of freedom
## Multiple R-squared:  0.9142, Adjusted R-squared:  0.8897 
## F-statistic: 37.29 on 2 and 7 DF,  p-value: 0.000185

Skontrolujem ACF rezíduí aj pri dynamickom modeli:

acf(residuals(model_dyn),
    main = "ACF rezíduí – dynamický model")

Ak je lmtest dostupný, môžem znovu urobiť testy:

lmtest::dwtest(model_dyn)
lmtest::bgtest(model_dyn, order = 1)

5 Newey–West robustné štandardné chyby (ak sú balíky dostupné)

Ak sú nainštalované balíky lmtest a sandwich, viem získať
Newey–West robustné štandardné chyby:

# robustné štandardné chyby pre základný model
lmtest::coeftest(model, vcov. = sandwich::NeweyWest(model))

# robustné štandardné chyby pre dynamický model
lmtest::coeftest(model_dyn, vcov. = sandwich::NeweyWest(model_dyn))

Ak balíky nie sú dostupné, tento chunk sa automaticky preskočí a knit nespadne.

6 Stručný záver

V cvičení som:

Na záver viem stručne zhodnotiť: