Práca s databázou - import údajov, grafy, štatistiky

Úvod a popis databázy

V tejto časti sa testujú štatistické hypotézy pomocou údajov z databázy, ktorá obsahuje ekonomické ukazovatele krajín v rokoch 1991 – 2022. Použité premenné zahŕňajú hrubý domáci produkt (HDP) v USD, mieru nezamestnanosti (%) a štruktúru zamestnanosti v troch sektoroch – poľnohospodárstve, priemysle a službách.

install.packages("zoo")
install.packages("tseries")
install.packages("lmtest")
install.packages("sandwich")
install.packages("car")
library(zoo)
library(tseries) 
library(lmtest)
library(sandwich)
library(car)
rm(list=ls())

install.packages("knitr")
install.packages("dplyr")
install.packages("ggplot2")

# Import vlastného CSV súboru

udaje <- read.csv("Employment_Unemployment_GDP_data.csv",
header = TRUE,
sep = ",",
dec = " ",
stringsAsFactors = FALSE)

# Zobrazenie prvých riadkov a názvov stĺpcov

head(udaje)

colnames(udaje)

[1] "Country.Name"                   "Year"                          
[3] "Employment.Sector..Agriculture" "Employment.Sector..Industry"   
[5] "Employment.Sector..Services"    "Unemployment.Rate"             
[7] "GDP..in.USD."                  

Úvod do problému, stanovenie hypotéz

Rozhodla som sa modelovať mieru nezamestnanosti (Unemployment.Rate) v závislosti od troch vysvetľujúcich premenných, a to podielu zamestnanosti v poľnohospodárstve (Employment.Sector..Agriculture), podielu zamestnanosti v priemysle (Employment.Sector..Industry) a hrubého domáceho produktu na obyvateľa (GDP..in.USD.).

Naša pracovná hypotéza hovorí o štatisticky významnom vplyve všetkých troch vysvetľujúcich premenných, pričom:

• u premennej Industry predpokladáme negatívny vplyv, to znamená čím väčší podiel pracujúcich v priemysle, tým nižšia nezamestnanosť,
  
• u premennej GDP očakávame negatívny vplyv, vyšší HDP na obyvateľa je spojený s lepšou ekonomickou výkonnosťou a teda nižšou nezamestnanosťou,
  
• u premennej Agriculture predpokladáme pozitívny alebo nejednoznačný vplyv, vo vyspelých ekonomikách nižší podiel poľnohospodárstva súvisí s nižšou nezamestnanosťou, zatiaľ čo vo vyvíjajúcich sa krajinách môže byť efekt opačný.

Príprava databázy, čistenie a úprava údajov

Budeme pracovať s vlastným súborom Employment_Unemployment_GDP_data.csv. Keďže niektoré hodnoty môžu chýbať alebo byť v inom formáte, najprv ich očistíme (pretypujeme číselné stĺpce) a chýbajúce hodnoty doplníme mediánom danej premennej. Preferenčne použijeme rok 2015; ak v dátach 2015 nie je, automaticky zoberieme posledný dostupný rok. Na ďalšie kroky si ponecháme kľúčové premenné: Unemployment.Rate, Employment.Sector..Agriculture, Employment.Sector..Industry, Employment.Sector..Services a GDP..in.USD..

# Robustná príprava: automatické namapovanie názvov stĺpcov + imputácia mediánom

# 0) Načítanie (ponechaj check.names = FALSE, aby sa nemenili mená)
udaje <- read.csv(
  "Employment_Unemployment_GDP_data.csv",
  header = TRUE, sep = ",", dec = ".",
  stringsAsFactors = FALSE, check.names = FALSE
)

# 1) Pomocné funkcie
find_col <- function(candidates, cols) {
  # skúsi presnú zhodu (case-insensitive), potom 'obsahuje'
  lc <- tolower(cols)
  # presná zhoda
  for (p in candidates) {
    idx <- which(lc == tolower(p))
    if (length(idx) == 1) return(cols[idx])
  }
  # obsahuje (regex/substring, case-insensitive)
  for (p in candidates) {
    idx <- grep(tolower(p), lc, fixed = TRUE)
    if (length(idx) >= 1) return(cols[idx][1])
  }
  return(NA_character_)
}

num_clean <- function(x) {
  if (is.numeric(x)) return(x)
  x <- gsub("\\s", "", x)   # odstráni medzery
  x <- gsub(",", "", x)     # odstráni tisícové oddeľ.
  suppressWarnings(as.numeric(x))
}

# 2) Nájdeme požadované stĺpce (tolerantne na názvy)
cols <- colnames(udaje)

col_year <- find_col(c("Year","Rok"), cols)
col_unemp <- find_col(c("Unemployment.Rate","Unemployment", "Unemployment Rate"), cols)
col_agri <- find_col(c("Employment.Sector..Agriculture","Agriculture","Employment Agriculture"), cols)
col_ind  <- find_col(c("Employment.Sector..Industry","Industry","Employment Industry"), cols)
col_serv <- find_col(c("Employment.Sector..Services","Services","Employment Services"), cols)
col_gdp  <- find_col(c("GDP..in.USD.","GDP..in.USD","GDP per capita","GDP","gdp"), cols)

mapping <- c(
  Year = col_year,
  Unemployment.Rate = col_unemp,
  Agriculture = col_agri,
  Industry = col_ind,
  Services = col_serv,
  GDP_in_USD = col_gdp
)
cat("Mapovanie stĺpcov:\n")

Mapovanie stĺpcov:

print(mapping)

                            Year                Unemployment.Rate 
                          "Year"              "Unemployment Rate" 
                     Agriculture                         Industry 
"Employment Sector: Agriculture"    "Employment Sector: Industry" 
                        Services                       GDP_in_USD 
   "Employment Sector: Services"                   "GDP (in USD)" 

# 3) Ošetri, ak sa niečo nenašlo
if (any(is.na(mapping))) {
  stop("Niektoré stĺpce sa nenašli. Skontroluj mapovanie vyššie a prípadne uprav kandidátov.")
}

# 4) Pretypuj na numerické (kde treba)
for (cn in unique(c(col_year,col_unemp,col_agri,col_ind,col_serv,col_gdp))) {
  if (cn %in% c(col_year,col_unemp,col_agri,col_ind,col_serv,col_gdp)) {
    udaje[[cn]] <- num_clean(udaje[[cn]])
  }
}

# 5) Vyber rok: 2015, inak posledný dostupný
target_year <- if (any(udaje[[col_year]] == 2015, na.rm = TRUE)) 2015 else max(udaje[[col_year]], na.rm = TRUE)
cat("Použitý rok:", target_year, "\n")

Použitý rok: 2015 

# 6) Vyber kľúčové premenné pre daný rok
udaje.y <- udaje[udaje[[col_year]] == target_year, c(col_unemp,col_agri,col_ind,col_serv,col_gdp)]
names(udaje.y) <- c("Unemployment.Rate","Agriculture","Industry","Services","GDP_USD")

# 7) Imputácia mediánom
column_medians <- sapply(udaje.y, median, na.rm = TRUE)
for (col in names(udaje.y)) {
  idx <- is.na(udaje.y[[col]])
  if (any(idx)) udaje.y[[col]][idx] <- column_medians[col]
}

# 8) Hotovo – dataset pripravený na modelovanie
str(udaje.y)

'data.frame':   181 obs. of  5 variables:
 $ Unemployment.Rate: num  9.05 17.19 11.21 16.49 7.58 ...
 $ Agriculture      : num  44.59 41.28 8.83 56.85 7.84 ...
 $ Industry         : num  20.7 18.7 31.2 7.8 22.3 ...
 $ Services         : num  34.7 40 59.9 35.4 69.8 ...
 $ GDP_USD          : num  1.91e+10 1.15e+10 1.87e+11 9.05e+10 5.95e+11 ...

summary(udaje.y)

 Unemployment.Rate  Agriculture         Industry         Services    
 Min.   : 0.170    Min.   : 0.2338   Min.   : 3.465   Min.   :10.22  
 1st Qu.: 3.678    1st Qu.: 5.3210   1st Qu.:14.514   1st Qu.:42.40  
 Median : 6.313    Median :18.1105   Median :19.555   Median :58.51  
 Mean   : 8.127    Mean   :25.0850   Mean   :19.612   Mean   :55.30  
 3rd Qu.:10.842    3rd Qu.:40.3543   3rd Qu.:24.458   3rd Qu.:69.15  
 Max.   :27.695    Max.   :86.3175   Max.   :54.141   Max.   :93.17  
    GDP_USD         
 Min.   :2.600e+08  
 1st Qu.:1.139e+10  
 Median :4.130e+10  
 Mean   :4.128e+11  
 3rd Qu.:1.951e+11  
 Max.   :1.830e+13  

Vizualizácia dát – kontrola nezrovnalostí

V tejto časti sa pozrieme na rozloženie jednotlivých premenných pomocou boxplotov.
Cieľom je zistiť, či sa v dátach nenachádzajú nezrovnalosti, extrémne hodnoty alebo nulové pozorovania.
Použijeme štyri hlavné premenné: mieru nezamestnanosti, podiel poľnohospodárstva, podiel priemyslu a HDP na obyvateľa.

# Boxploty premenných – kontrola rozloženia a odľahlých hodnôt

# Zvolíme len relevantné premenné
vars_to_plot <- c("Unemployment.Rate","Agriculture","Industry","GDP_USD")

# Nastavenie grafického layoutu: 2 × 2
par(mfrow = c(2, 2))
par(mar = c(4, 4, 2, 1))

# Pre každý vybraný stĺpec nakreslíme boxplot
for (col in vars_to_plot) {
  boxplot(
    udaje.y[[col]],
    main = col,
    xlab = "Hodnota",
    col = "lightblue",
    border = "darkblue"
  )
}

# Nadpis pre celú sadu grafov
mtext("Boxploty jednotlivých premenných (rok vybraný pre analýzu)",
      outer = TRUE, cex = 1.2, font = 2)

# Reset layoutu na 1 graf
par(mfrow = c(1, 1))

Na základe boxplotov môžeme pozorovať nasledovné skutočnosti:

• Unemployment.Rate – väčšina pozorovaní sa sústreďuje v stredných hodnotách, avšak v niektorých krajinách sa vyskytujú aj extrémne vyššie miery nezamestnanosti. To naznačuje, že medzi krajinami existujú výrazné rozdiely v trhu práce.
  
• Agriculture – hodnoty sa pohybujú v širokom intervale. Krajiny s vyšším podielom poľnohospodárstva môžu mať nižšiu úroveň industrializácie, čo sa často spája s vyššou nezamestnanosťou.
  
• Industry – väčšina krajín má stredné až vyššie hodnoty, pričom extrémne hodnoty sú menej časté. Vyšší podiel priemyslu naznačuje rozvinutejšiu ekonomiku.
  
• GDP_USD – rozloženie ukazuje výrazné rozdiely v ekonomickej úrovni medzi krajinami. Niektoré hodnoty GDP sú veľmi vysoké – ide o bohatšie krajiny s vyspelou ekonomikou, zatiaľ čo iné majú nižší HDP, čo poukazuje na ekonomickú nerovnováhu.

Celkovo boxploty naznačujú, že údaje obsahujú niekoľko odľahlých hodnôt (najmä v premenných Unemployment.Rate a GDP_USD), čo je však pri medzinárodných dátach prirodzené. Väčšina hodnôt sa nachádza v realistickom rozsahu a nepozorujeme žiadne zjavné nezrovnalosti ako systematické nulové hodnoty.
Tieto výsledky potvrdzujú, že údaje sú vhodné na ďalšie modelovanie.

Lineárna regresia

Model odhadujeme príkazom lm().

V našom prípade modelujeme mieru nezamestnanosti (Unemployment.Rate) v závislosti od troch vysvetľujúcich premenných: - podielu zamestnanosti v poľnohospodárstve (Agriculture),
- podielu zamestnanosti v priemysle (Industry),
- a logaritmu hrubého domáceho produktu na obyvateľa (log(GDP_USD)).

Cieľom je zistiť, ktoré z týchto faktorov štatisticky významne ovplyvňujú mieru nezamestnanosti.

# Pridáme logaritmickú transformáciu HDP
udaje.y$logGDP <- ifelse(udaje.y$GDP_USD > 0, log(udaje.y$GDP_USD), NA_real_)

# Základný lineárny model
model <- lm(Unemployment.Rate ~ Agriculture + Industry + logGDP, data = udaje.y)

# Súhrn výsledkov
summary(model)

Diagnostické grafy regresného modelu

Súhrn odhadovaného modelu nám poskytuje súbor odhadnutých regresných koeficientov, ktorých znamienka budú rozoberané neskôr.
Ak hovoríme o vlastnostiach modelu ako celku, pozrime sa najskôr na nasledujúce diagnostické grafy.
Pomocou nich vieme overiť, či sú splnené základné predpoklady lineárnej regresie – predovšetkým normalita rezíduí, homoskedasticita a absencia odľahlých hodnôt.

# Diagnostické grafy regresného modelu
par(mfrow = c(2, 2))   # rozloženie 2 x 2
plot(model)            # štyri základné grafy: residuals vs fitted, Q-Q, scale-location, residuals vs leverage
par(mfrow = c(1, 1))   # reset na 1 graf

Interpretácia diagnostických grafov

1. Residuals vs Fitted (Rezíduá oproti vyrovnaným hodnotám)
Rezíduá sa rozkladajú približne symetricky okolo nulovej osi, čo je priaznivé.
Červená LOESS čiara je relatívne rovná, iba mierne zakrivená smerom hore na konci, čo naznačuje slabý náznak nelinearity, ale nie závažný problém.
Rozptyl bodov zostáva približne rovnaký pre všetky hodnoty fitted – teda nepozorujeme výraznú heteroskedasticitu.

2. Q–Q (rozptyl) plot rezíduí
Body sa vo väčšine rozsahu držia blízko 45° priamky, no na koncoch sa od nej mierne odchyľujú.
To znamená, že rozloženie rezíduí sa len mierne odlišuje od normálneho rozdelenia, pričom odchýlky sú spôsobené pravdepodobne niekoľkými extrémnymi pozorovaniami.
Celkovo však predpoklad normality nie je vážne porušený.

3. Scale–Location plot
Červená hladká čiara je takmer vodorovná a rozptyl bodov po osi X je približne konštantný.
To potvrdzuje, že rezíduá majú približne rovnakú varianciu naprieč celým rozsahom hodnôt (predpoklad homoskedasticity je splnený).

4. Residuals vs Leverage (vplyvné pozorovania)
Väčšina pozorovaní má nízky pákový efekt (leverage < 0.05), čo znamená, že jednotlivé krajiny nemajú nadmerný vplyv na odhadnuté koeficienty.
Niekoľko bodov (napr. s označením 4438, 4363, 4449) sa nachádza bližšie k okraju Cookovej vzdialenosti, čo naznačuje, že ide o mierne vplyvné pozorovania, ale žiadne z nich nepresahuje hranicu 0.5 či 1.0, teda žiadne extrémne odľahlé hodnoty sa neobjavili.

# Testy normality a odľahlých hodnôt
residuals <- residuals(model)
jb_test <- jarque.bera.test(residuals)
jb_test  # Jarque–Bera test normality

    Jarque Bera Test

data:  residuals
X-squared = 42.52, df = 2, p-value = 5.847e-10

# Outlier test (Bonferroni correction)
outlier_test <- car::outlierTest(model)
outlier_test

No Studentized residuals with Bonferroni p < 0.05
Largest |rstudent|:

Výsledky testu odľahlých hodnôt

Výstup funkcie outlierTest(model) identifikoval pozorovanie s indexom 4337,
ktoré má najvyššiu študentizovanú hodnotu rezídua rstudent = 3.27.
Jeho neopravená p-hodnota je 0.00127, avšak po aplikácii Bonferroniho korekcie
je výsledná hodnota 0.23064, teda nie štatisticky významná na 5 % hladine.

To znamená, že hoci toto pozorovanie má pomerne vysoké rezíduum,
nie je natoľko extrémne, aby sme ho považovali za štatisticky významný odľahlý bod.
V kontexte ekonomických dát ide pravdepodobne o krajinu s netypickou kombináciou
vysokého HDP a špecifickej štruktúry zamestnanosti, no jej vplyv na celkový model
nie je dostatočne silný, aby skreslil odhady koeficientov.

Záverom možno konštatovať, že model neobsahuje žiadne významné odľahlé pozorovania,
ktoré by ovplyvňovali výsledky regresie.

Alternatívny model

Ak sa vyskytujú mierne odľahlé hodnoty alebo nenormalita v GDP, môžeme upraviť model tak,
že použijeme logaritmus HDP a zmeníme štruktúru sektorov.
Nový model bude mať tvar:

\[ Unemployment.Rate_i = \beta_0 + \beta_1 \, Industry_i + \beta_2 \, Services_i + \beta_3 \, \log(GDP_i) + \varepsilon_i \]

# Alternatívny model s log(GDP) a Services namiesto Agriculture
model2 <- lm(Unemployment.Rate ~ Industry + Services + log(GDP_USD), data = udaje.y)
summary(model2)

# Diagnostické grafy alternatívneho modelu
par(mfrow = c(2, 2))
plot(model2)
par(mfrow = c(1, 1))

# Normality a outlier test pre nový model
residuals2 <- residuals(model2)
jarque.bera.test(residuals2)
car::outlierTest(model2)

Záver

Nový model po logaritmickej transformácii HDP potvrdzuje, že:

• premenné priemysel (Industry) a HDP (GDP) majú negatívny vplyv na mieru nezamestnanosti – teda čím je podiel priemyslu a úroveň HDP vyššia, tým je nezamestnanosť nižšia,
  
• premenná služby (Services) má slabší alebo štatisticky nevýznamný vplyv,
  
• rezíduá majú po transformácii lepšie rozdelenie a model nevykazuje závažné porušenia predpokladov lineárnej regresie.

Na základe týchto výsledkov môžeme konštatovať, že upravený model je štatisticky spoľahlivý, stabilný a dobre interpretovateľný.
Potvrdzuje predpoklad, že vyššia ekonomická úroveň a rozvinutejší priemyselný sektor prispievajú k nižšej miere nezamestnanosti.

Heteroskedasticita

Prítomnosť heteroskedasticity (nekonštantného rozptylu náhodnej zložky) spôsobuje nepresné odhady smerodajných chýb a tým aj nespoľahlivé t-testy významnosti jednotlivých regresných koeficientov.
Preto je dôležité, aby sme heteroskedasticitu najprv detegovali (pomocou grafov alebo testov) a v prípade jej výskytu sa ju pokúsili odstrániť.

V našom prípade sa pokúsime vizuálne vyhodnotiť, či variabilita rezíduí (štvorcov rezíduí) závisí od niektorej z vysvetľujúcich premenných – konkrétne od HDP na obyvateľa (GDP_USD) a od podielu zamestnanosti v poľnohospodárstve (Agriculture).
Budeme porovnávať dva modely:

• model – pôvodný lineárny model s premennou GDP_USD
  
• model2 – model so zlogaritmovaným HDP (log(GDP_USD)), ktorý môže odstrániť vplyv odľahlých hodnôt a stabilizovať varianciu

library(ggplot2)
library(patchwork)

# Štvorce rezíduí z pôvodného modelu
p1 <- ggplot(udaje.y, aes(x = GDP_USD, y = resid(model)^2)) +
  geom_point(alpha = 0.8, color = "darkred") +
  geom_smooth(method = "loess", se = FALSE, color = "pink", linewidth = 1.5) +
  labs(x = "HDP na obyvateľa (GDP_USD)", 
       y = "Štvorce rezíduí",
       title = "Squared Residuals vs GDP_USD (model)") +
  theme_minimal()

# Druhý graf: poľnohospodárstvo
p2 <- ggplot(udaje.y, aes(x = Agriculture, y = resid(model)^2)) +
  geom_point(alpha = 0.8, color = "darkblue") +
  geom_smooth(method = "loess", se = FALSE, color = "skyblue", linewidth = 1.5) +
  labs(x = "Podiel zamestnanosti v poľnohospodárstve (%)", 
       y = "Štvorce rezíduí",
       title = "Squared Residuals vs Agriculture (model)") +
  theme_minimal()

# Zobrazenie vedľa seba
p1 + p2

Interpretácia výsledkov

V grafe „Štvorce rezíduí vs GDP_USD“ možno vidieť, že pri nižších hodnotách HDP je rozptyl rezíduí väčší a s rastúcim HDP sa znižuje.
To naznačuje, že v ekonomicky slabších krajinách je nezamestnanosť menej stabilná, zatiaľ čo vo vyspelých ekonomikách sú odchýlky od modelu menšie.
Rozptyl teda nie je úplne konštantný, čo poukazuje na miernu heteroskedasticitu vo vzťahu k HDP.

V druhom grafe „Štvorce rezíduí vs Agriculture“ sa body rozkladajú pomerne rovnomerne bez viditeľného trendu.
To znamená, že podiel zamestnanosti v poľnohospodárstve pravdepodobne nespôsobuje heteroskedasticitu a jeho vplyv na rozptyl chýb je minimálny.

Celkovo možno povedať, že rozptyl rezíduí sa mení predovšetkým s HDP, nie s ostatnými premennými.
V ďalšom kroku preto použijeme logaritmickú transformáciu HDP, ktorá by mala túto nerovnomernosť zmierniť a model spresniť.

model2 <- lm(Unemployment.Rate ~ Agriculture + Industry + Services + log(GDP_USD), data = udaje.y)

# Model s logaritmom HDP
model2 <- lm(Unemployment.Rate ~ Agriculture + Industry + Services + log(GDP_USD), data = udaje.y)

# Skúmanie heteroskedasticity – vizuálne pre model2
library(ggplot2)
library(patchwork)

# Štvorce rezíduí vs log(GDP)
p1 <- ggplot(udaje.y, aes(x = log(GDP_USD), y = resid(model2)^2)) +
  geom_point(alpha = 0.6, color = "darkgreen") +
  geom_smooth(method = "loess", se = FALSE, color = "green", linewidth = 1) +
  labs(x = "log(HDP na obyvateľa)", 
       y = "Štvorce rezíduí",
       title = "Squared Residuals vs log(GDP_USD) (model2)") +
  theme_minimal()

# Štvorce rezíduí vs poľnohospodárstvo
p2 <- ggplot(udaje.y, aes(x = Agriculture, y = resid(model2)^2)) +
  geom_point(alpha = 0.6, color = "purple4") +
  geom_smooth(method = "loess", se = FALSE, color = "orchid", linewidth = 1) +
  labs(x = "Agriculture", 
       y = "Štvorce rezíduí",
       title = "Squared Residuals vs Agriculture (model2)") +
  theme_minimal()

# Zobrazenie grafov vedľa seba
p1 + p2

Interpretácia výsledkov po logaritmickej transformácii HDP

Po aplikovaní logaritmickej transformácie HDP vidíme, že rozptyl rezíduí sa výrazne stabilizoval.
V grafe „Štvorce rezíduí vs log(GDP_USD)” je LOESS krivka (zelená línia) takmer vodorovná a body sú rozložené rovnomerne naprieč hodnotami HDP.
To znamená, že logaritmická transformácia odstránila väčšinu heteroskedasticity, ktorá bola prítomná v pôvodnom modeli.

V grafe „Štvorce rezíduí vs Agriculture” ostáva rozptyl bodov pomerne rovnomerný a bez jasného trendu.
To potvrdzuje, že podiel zamestnanosti v poľnohospodárstve nemá zásadný vplyv na variabilitu chýb.

Celkovo môžeme uzavrieť, že po logaritmickej úprave HDP sa model správa stabilnejšie a lepšie spĺňa predpoklad konštantného rozptylu (homoskedasticity).
Takto upravený model je preto spoľahlivejší a vhodnejší na ekonomickú interpretáciu.

Testovanie prítomnosti heteroskedasticity

# Install (if not yet installed)
# install.packages("lmtest")

# Load the package
library(lmtest)

# Run the Breusch–Pagan test
bptest(model)

    studentized Breusch-Pagan test

data:  model
BP = 7.167, df = 3, p-value = 0.06676

# Install (if not yet installed)
# install.packages("lmtest")

# Load the package
library(lmtest)

# Run the Breusch–Pagan test
bptest(model2)

    studentized Breusch-Pagan test

data:  model2
BP = 8.1376, df = 4, p-value = 0.08667

Interpretácia výsledkov testu heteroskedasticity

Podľa výsledkov Breusch–Pagan testu:

• pre pôvodný model (model) platí:
  BP = 7.167, df = 3, p-value = 0.06676

• pre model so zlogaritmovaným HDP (model2) platí:
  BP = 8.1376, df = 4, p-value = 0.08667

V oboch prípadoch je hodnota p-value väčšia ako 0.05, takže nezamietame nulovú hypotézu.
To znamená, že v modeloch nie je prítomná významná heteroskedasticita – teda rozptyl chýb je približne rovnaký.

Hodnota p-value v pôvodnom modeli (model) je však trochu nižšia (asi 0.07), čo môže naznačovať mierny náznak nerovnakého rozptylu.
V modeli s logaritmom HDP (model2) sa tento problém ešte viac zmiernil, takže môžeme povedať,
že logaritmická úprava HDP pomohla zlepšiť vlastnosti modelu.

Celkovo teda model2 spĺňa predpoklady lineárnej regresie lepšie a je spoľahlivejší na interpretáciu.

#install.packages("sandwich")
#install.packages("lmtest")
library(sandwich)
library(lmtest)
coeftest(model, vcov = vcovHC(model))

t test of coefficients:

             Estimate Std. Error t value  Pr(>|t|)    
(Intercept) 30.084154   5.730333  5.2500 4.323e-07 ***
Agriculture -0.095462   0.026310 -3.6284 0.0003730 ***
Industry    -0.011641   0.095880 -0.1214 0.9034993    
logGDP      -0.789032   0.227598 -3.4668 0.0006611 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Interpretácia výsledkov po White korekcii

Po aplikovaní tzv. White korekcie (robustné smerodajné chyby) sa smerodajné odhady upravili tak, aby boli výsledky spoľahlivé aj v prípade, že by rozptyl chýb nebol rovnaký.

Z výsledkov vyplýva, že: - Agriculture má negatívny a štatisticky významný vplyv – teda čím väčší podiel ľudí pracuje v poľnohospodárstve, tým je spravidla nižšia nezamestnanosť.
- log(GDP) je tiež negatívny a významný, čo znamená, že bohatšie krajiny majú nižšiu mieru nezamestnanosti.
- Industry nie je štatisticky významný, takže jeho vplyv na nezamestnanosť nie je jednoznačný.
- Konštanta (Intercept) je významná, čo potvrdzuje, že model ako celok je stabilný a vhodný na interpretáciu.

Výsledky po tejto korekcii sú teda presnejšie a menej ovplyvnené tým, či majú krajiny rôzny rozptyl v údajoch.

Záver časti o heteroskedasticite

Na základe testov a grafov môžeme povedať, že v našich dátach sa nepotvrdila prítomnosť výraznej heteroskedasticity.
Rozptyl rezíduí je v modeloch pomerne stabilný, najmä po logaritmickej transformácii HDP.
White korekcia ukázala, že aj keby sa v dátach mierna nerovnomernosť vyskytla, výsledky modelu ostávajú spoľahlivé a interpretovateľné.

Celkovo možno zhrnúť, že model2 (so zlogaritmovaným HDP) lepšie spĺňa predpoklady lineárnej regresie a poskytuje stabilnejšie odhady,
čo znamená, že ho môžeme považovať za vhodný nástroj na vysvetlenie miery nezamestnanosti medzi krajinami.