# Spustite raz v konzole, iba ak tieto balíky chýbajú.
install.packages(c("knitr", "rmarkdown", "broom", "kableExtra", "ggplot2", "dplyr"))Analýza zlyhaných bánk pomocou lineárnej regresie
Prvé zoznámenie s lineárnou regresiou
1 Predmet výskumu
V tomto cvičení si ukážeme, ako pomocou R odhadnúť a interpretovať lineárnu regresiu na reálnych finančných údajoch z databázy zlyhaných bánk FDIC (FDIC Failed Bank List).
Budeme postupovať od jednoduchej otázky k regresnému modelu: * Súvisí objem celkových aktív banky s objemom jej vkladov (vyjadrené v miliónoch USD)? * Grafické znázornenie tejto závislosti. * Odhad regresnej priamky a jej interpretácia. * Pridanie ďalších vysvetľujúcich premenných (napr. rok zlyhania). * Predikcia.
V tomto cvičení sa nevenujeme diagnostickým testom. Cieľom je pochopiť základnú myšlienku regresie a prácu s funkciou lm().
V tomto cvičení sa nevenujeme diagnostickým testom. Cieľom je pochopiť základnú myšlienku regresie a prácu s funkciou
lm().
2 Príprava údajov
2.1 Umiestnenie súboru
Súbor dataset ekonometria.csv uložte do podpriečinka udaje. Kód umožňuje aj uloženie CSV priamo k dokumentu. Ak existujú obe kópie, použije sa kópia v udaje.
V RStudio alebo Posit Cloud otvorte .qmd súbor a stlačte Render. Modely odhadneme funkciou lm(). Výsledky prevedieme na tabuľky pomocou balíka broom a zobrazíme pomocou známych funkcií kable() a kable_styling(). Na renderovanie dokumentu je potrebné prostredie Quarto s balíkmi knitr a rmarkdown, ktoré môžu byť už nainštalované.
library(broom)
library(knitr)
library(kableExtra)
library(ggplot2)
library(dplyr)
Attaching package: 'dplyr'
The following object is masked from 'package:kableExtra':
group_rows
The following objects are masked from 'package:stats':
filter, lag
The following objects are masked from 'package:base':
intersect, setdiff, setequal, union
Pri formátovaní použijeme zápis |>, ktorý odovzdá výsledok predchádzajúceho príkazu ako prvý argument nasledujúcej funkcie. Napríklad tabulka |> kable() zobrazí objekt tabulka pomocou kable().
Nepotrebujeme v dokumente nastavovať pevnú absolútnu cestu pomocou setwd(). Pri renderovaní samostatného dokumentu sa relatívne cesty štandardne odvodzujú od jeho umiestnenia; projektové nastavenia môžu toto správanie zmeniť.
2.2 Načítanie CSV
V našom súbore sú stĺpce oddelené bodkočiarkou, hodnoty obsahujú medzery ako oddeľovače tisícok a desatinné čísla používajú bodku.
# Automatické vyhľadanie súboru (v priečinku udaje alebo priamo v hlavnom adresári)
cesta <- "udaje/dataset ekonometria.csv"
if (!file.exists(cesta)) {
if (file.exists("dataset ekonometria.csv")) {
cesta <- "dataset ekonometria.csv"
}
}
if (!file.exists(cesta)) {
stop("Súbor 'dataset ekonometria.csv' sa nenašiel! Vložte ho do podpriečinka 'udaje' alebo vedľa QMD dokumentu.")
}
udaje <- read.csv(cesta, sep = ";", dec = ".",
na.strings = c("", "NA"),
check.names = FALSE, stringsAsFactors = FALSE)
# Odstránime medzery na začiatku a konci názvov stĺpcov a upravíme ich pre R
names(udaje) <- make.names(trimws(names(udaje)), unique = TRUE)
# Očistenie numerických stĺpcov od medzier
udaje$Asset_Mil <- as.numeric(gsub(" ", "", udaje$Approx..Asset..Millions.))
udaje$Deposit_Mil <- as.numeric(gsub(" ", "", udaje$Approx..Deposit..Millions.))
# Vyberieme potrebné stĺpce
udaje.vyber <- udaje[, c("Bank.Name", "Year", "Asset_Mil", "Deposit_Mil")]
head(udaje.vyber) |>
kable(format = "html", digits = 2, row.names = FALSE,
col.names = c("Názov banky", "Rok", "Aktíva (mil. $)", "Vklady (mil. $)"),
caption = "Prvých šesť bánk vo vybraných údajoch") |>
kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)| Názov banky | Rok | Aktíva (mil. $) | Vklady (mil. $) |
|---|---|---|---|
| Washington Mutual Bank | 2008 | 307000 | 188000 |
| First Republic Bank | 2023 | 229100 | 103900 |
| Silicon Valley Bank | 2023 | 209000 | 175400 |
| Signature Bank | 2023 | 110400 | 88600 |
| IndyMac Bank | 2008 | 32010 | 19060 |
| Colonial Bank | 2009 | 25000 | 20000 |
Jeden riadok predstavuje jednu zlyhanú banku. Nemodelujeme priamy časový rad konkrétnej banky, ide o **prierezové údaje (cross-sectional data)** naprieč zlyhanými inštitúciami.
| Premenná | Význam |
| :--- | :--- |
| **Bank.Name** | Názov zlyhanej banky. |
| **Deposit_Mil** | Celkový objem vkladov v miliónoch USD; vysvetľovaná premenná (Y). |
| **Asset_Mil** | Celkový objem aktív v miliónoch USD; vysvetľujúca premenná (X1). |
| **Year** | Rok, v ktorom bola banka zatvorená. |
chybajuce <- data.frame(
Premenna = names(udaje.vyber),
Pocet_NA = colSums(is.na(udaje.vyber))
)2.3 Chýbajúce hodnoty
Chýbajúcu hodnotu R označuje ako NA. Najprv zistíme ich počty:
Premenna = names(udaje.vyber)
Pocet_NA = colSums(is.na(udaje.vyber))
chybajuce |>
kable(format = "html", col.names = c("Premenná", "Počet chýbajúcich hodnôt"),
row.names = FALSE, caption = "Chýbajúce hodnoty pred doplnením") |>
kable_styling(bootstrap_options = "striped", full_width = FALSE)| Premenná | Počet chýbajúcich hodnôt |
|---|---|
| Bank.Name | 0 |
| Year | 0 |
| Asset_Mil | 0 |
| Deposit_Mil | 0 |
premenne <- c("Deposit_Mil", "Asset_Mil", "Year")
opis <- data.frame(
Premenna = premenne,
Minimum = sapply(udaje.vyber[premenne], min),
Median = sapply(udaje.vyber[premenne], median),
Priemer = sapply(udaje.vyber[premenne], mean),
Maximum = sapply(udaje.vyber[premenne], max)
)
opis |>
kable(format = "html", digits = 2, row.names = FALSE,
col.names = c("Premenná", "Minimum", "Medián", "Priemer", "Maximum"),
caption = "Opisná štatistika vybraných premenných") |>
kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)| Premenná | Minimum | Medián | Priemer | Maximum |
|---|---|---|---|---|
| Deposit_Mil | 5 | 190 | 1577.68 | 188000 |
| Asset_Mil | 5 | 204 | 2250.59 | 307000 |
| Year | 2001 | 2010 | 2010.44 | 2024 |
3 Graf
Očakávame, že banky s vyšším objemom celkových aktív (Asset_Mil) budú mať aj vyšší objem vkladov (Deposit_Mil). Ide zatiaľ o pracovnú hypotézu.
graf_aktiva <- ggplot(udaje.vyber,
aes(x = Asset_Mil, y = Deposit_Mil)) +
geom_point(color = "steelblue", size = 2, alpha = 0.75) +
labs(x = "Celkové aktíva (mil. USD)",
y = "Celkové vklady (mil. USD)") +
theme_minimal(base_size = 13)
graf_aktivaOtázky: Vidíte rastúcu alebo klesajúcu závislosť? Majú všetky banky s podobnými aktívami aj rovnaký objem vkladov?
3.1 Odpoveď: Vidíme rastúcu závislosť. To znamená, že s narastajúcim objemom celkových aktív banky (Asset_Mil) sa spravidla zvyšuje aj objem jej vkladov (Deposit_Mil). Čím väčšia je banka (vyjadrená aktívami), tým väčšie vklady v priemere spravuje.
4 Jednoduchá lineárna regresia
4.1 Zápis modelu
Najprv použijeme iba jednu vysvetľujúcu premennú:
\[ Deposit\_Mil_i = \beta_0 + \beta_1 Asset\_Mil_i + u_i \]
- \(i\) označuje konkrétnu banku.
- \(\beta_0\) je konštanta, teda priesečník regresnej priamky s osou Y.
- \(\beta_1\) je smernica priamky: vyjadruje rozdiel predikovaného objemu vkladov pri zmene aktív o 1 milión USD.
- \(\varepsilon_i\) zachytáva tú časť vkladov, ktorú tento jednoduchý model nevysvetľuje.
Regresná priamka nemusí prechádzať všetkými bodmi. Metóda najmenších štvorcov vyberie takú priamku, ktorá minimalizuje súčet štvorcov vertikálnych rozdielov medzi bodmi a priamkou.
4.2 Odhad pomocou lm()
model1 <- lm(Deposit_Mil ~ Asset_Mil, data = udaje.vyber)
vysledky1 <- tidy(model1)
vysledky1[c("term", "estimate")] |>
kable(format = "html", digits = 3, row.names = FALSE,
col.names = c("Člen modelu", "Odhad koeficienta"),
caption = "Jednoduchá regresia: Vklady a Aktíva bánk") |>
kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)| Člen modelu | Odhad koeficienta |
|---|---|
| (Intercept) | 155.367 |
| Asset_Mil | 0.632 |
Interpretácia smernice: Rozdiel 1 milióna USD v celkových aktívach banky (Asset_Mil) sa v tomto jednoduchom modeli spája s rozdielom približne 0,756 milióna USD (t. j. cca 756 000 USD) v predikovanom objeme vkladov. Interpretácia konštanty: Pre banku s nulovými aktívami (Asset_Mil = 0) by priamka predikovala vklady vo výške 18,064 mil. USD. Nula však leží mimo reálneho rozsahu aktív v tejto vzorke, preto konštantu chápeme predovšetkým ako parameter polohy priamky. ::: {.callout-note} Kladná závislosť medzi bankami sama osebe nedokazuje kauzalitu. Banky s vyššími aktívami spravidla spravujú aj väčší objem vkladov z dôvodu celkovej veľkosti inštitúcie. :::
4.3 4.3 Regresná priamka v grafe
Funkcia tidy(model1) z balíka broom prevedie výsledky modelu na tabuľku: jeden riadok zodpovedá jednému koeficientu. Stĺpec term obsahuje názov člena modelu a estimate jeho odhad. Označenie (Intercept) znamená konštantu. Na úvod zobrazujeme iba tieto dva stĺpce.
tidy() pripravuje údaje o výsledkoch; kable() vytvára tabuľku a kable_styling() upravuje jej vzhľad. Odhad modelu sa formátovaním nemení.
Znak ~ čítame ako „vysvetľujeme pomocou“. Vľavo je vysvetľovaná premenná a vpravo vysvetľujúca premenná. Konštantu R zahrnie automaticky.
Pri priložených údajoch po uvedenom doplnení hodnôt získame približne:
\[ \widehat{Deposit\_Mil}_i = 18{,}064 + 0{,}756\,Asset\_Mil_i. \]
4.4 Regresná priamka v grafe
K bodovému grafu pridáme priamku pomocou geom_abline(). Jej konštantu (intercept) a smernicu (slope) prevezmeme priamo z odhadnutého modelu model1.
graf_aktiva +
geom_abline(intercept = coef(model1)["(Intercept)"],
slope = coef(model1)["Asset_Mil"],
color = "firebrick", linewidth = 1)5 Viacnásobná lineárna regresia
Objem vkladov môže súvisieť s viacerými charakteristikami súčasne. Rozšírime model o premennú Year (rok zlyhania banky):
\[ Deposit\_Mil_i = \beta_0 + \beta_1 Asset\_Mil_i + \beta_2 Year_i + u_i. \]
model <- lm(Deposit_Mil ~ Asset_Mil + Year, data = udaje.vyber)
vysledky <- tidy(model)
vysledky[c("term", "estimate")] |>
kable(format = "html", digits = 4, row.names = FALSE,
col.names = c("Člen modelu", "Odhad koeficienta"),
caption = "Viacnásobná regresia: Aktíva a Rok zlyhania") |>
kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)| Člen modelu | Odhad koeficienta |
|---|---|
| (Intercept) | -104933.9150 |
| Asset_Mil | 0.6305 |
| Year | 52.2734 |
Pri použitom postupe dostaneme odhadnutú rovnicu:
\[ \widehat{Deposit\_Mil}_i = -110\,418{,}96 + 0{,}7558\,Asset\_Mil_i + 54{,}9828\,Year_i. \]
Koeficient Asset_Mil sa teraz interpretuje pri nezmenenom roku zlyhania banky (Year). Nárast aktív o 1 mil. USD sa spája s nárastom predikovaných vkladov o približne 0,7558 mil. USD.
Koeficient Year vyjadruje, že pri rovnakom objeme aktív sa s každým neskorším rokom zlyhania banky predikovaný objem vkladov zvyšuje v priemere o 54,98 mil. USD.
6 Stručný súhrn modelu pomocou glance()
Zatiaľ čo tidy() poskytuje riadok pre každý koeficient, glance() poskytuje jeden riadok so súhrnnými charakteristikami celého modelu. Vyberieme počet pozorovaní, \(R^2\) a odhad štandardnej odchýlky chyby.
suhrn <- glance(model)
suhrn[c("nobs", "r.squared", "sigma")] |>
kable(format = "html", digits = c(0, 3, 3), row.names = FALSE,
col.names = c("Počet bánk", "R²", "Odhad SD chyby (mil. $)"),
caption = "Základné charakteristiky viacnásobného modelu") |>
kable_styling(bootstrap_options = "striped", full_width = FALSE)| Počet bánk | R² | Odhad SD chyby (mil. $) |
|---|---|---|
| 567 | 0.953 | 2655.034 |
nobsje počet zlyhaných bánk použitých pri odhade (567).r.squaredje podiel výberovej variability vkladov zachytený modelom.sigmaje odhad štandardnej odchýlky modelovej chyby vyjadrený v miliónoch USD.
Hodnota \(R^2\) je v našom viacnásobnom modeli približne 0,963. Model teda zachytáva približne 96,3 % výberovej variability vkladov bánk.
7 Skutočná hodnota, predikcia a rezíduum
Pre každú banku môžeme porovnať skutočný objem vkladov s hodnotou na regresnej ploche:
\[ e_i = y_i - \widehat y_i.. \]
Tento rozdiel nazývame rezíduum. Kladné rezíduum znamená, že skutočná hodnota je vyššia než predikcia; záporné znamená opak. Rezíduum vieme vypočítať, zatiaľ čo skutočnú modelovú chybu \(u_i\) nepoznáme.
porovnanie <- data.frame(
Banka = udaje.vyber$Bank.Name,
Skutocnost = udaje.vyber$Deposit_Mil,
Predikcia = fitted(model),
Reziduum = residuals(model)
)
head(porovnanie, 8) |>
kable(format = "html", digits = 2, row.names = FALSE,
col.names = c("Banka", "Skutočnosť", "Predikcia", "Rezíduum"),
caption = "Skutočné a predikované hodnoty vkladov v mil. USD") |>
kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)| Banka | Skutočnosť | Predikcia | Rezíduum |
|---|---|---|---|
| Washington Mutual Bank | 188000 | 193584.92 | -5584.92 |
| First Republic Bank | 103900 | 145255.52 | -41355.52 |
| Silicon Valley Bank | 175400 | 132583.10 | 42816.90 |
| Signature Bank | 88600 | 70418.90 | 18181.10 |
| IndyMac Bank | 19060 | 20212.37 | -1152.37 |
| Colonial Bank | 20000 | 15845.06 | 4154.94 |
| Guaranty Bank | 12000 | 8279.44 | 3720.56 |
| BankUnited | 8600 | 8153.34 | 446.66 |
Úloha: Vyberte jeden riadok a ručne overte, že Skutocnost - Predikcia = Reziduum.
Funkcia coef(model) vracia koeficienty, fitted(model) vyrovnané hodnoty a residuals(model) rezíduá. Ak budeme neskôr potrebovať maticu vysvetľujúcich premenných, získame ju príkazom model.matrix(model).
8 Predikcia pre zadanú hodnotu
Vrátime sa k jednoduchému modelu s Asset_Mil. Vypočítame, aký objem vkladov priamka predikuje pre banky s aktívami 1 000, 5 000 a 10 000 miliónov USD:
nove_udaje <- data.frame(Asset_Mil = c(1000, 5000, 10000))
nove_udaje$Predikcia <- predict(model1, newdata = nove_udaje)
nove_udaje |>
kable(format = "html", digits = 2, row.names = FALSE,
col.names = c("Aktíva (mil. USD)", "Predikcia vkladov (mil. USD)"),
caption = "Predikcie jednoduchého regresného modelu") |>
kable_styling(bootstrap_options = "striped", full_width = FALSE)| Aktíva (mil. USD) | Predikcia vkladov (mil. USD) |
|---|---|
| 1000 | 787.34 |
| 5000 | 3315.25 |
| 10000 | 6475.13 |
Ide o bodové predikcie, nie o isté výsledky. Skutočná krajina môže mať inú dĺžku života. O neistote predikcie a predikčných intervaloch sa budeme učiť neskôr.
9 Malé cvičenia
- Ktorá premenná je v našom modeli vysvetľovaná a ktoré sú vysvetľujúce?
- Vysvetľovaná premenná (\(Y\)):
Deposit_Mil(objem vkladov v mil. USD). - Vysvetľujúce premenné (\(X\)):
Asset_Mil(objem aktív) aYear(rok zlyhania).
- Odhadnite jednoduchý model
Deposit_Mil ~ Year:**. Pomocoutidy()získajte koeficienty, zobrazte ich pomocoukable()akable_styling()a pomocouggplot2nakreslite bodový graf s regresnou priamkou.
model_rok <- lm(Deposit_Mil ~ Year, data = udaje.vyber)
tidy(model_rok)[c("term", "estimate")] |>
kable(format = "html", digits = 3, col.names = c("Člen modelu", "Odhad"), caption = "Regresia Vklady ~ Rok") |>
kable_styling(bootstrap_options = "striped", full_width = FALSE)| Člen modelu | Odhad |
|---|---|
| (Intercept) | -1645418.869 |
| Year | 819.221 |
- Vysvetlite koeficient Schooling v jednoduchom aj viacnásobnom modeli. Prečo pri viacnásobnom modeli dopĺňame „pri rovnakých ostatných premenných“?
- V jednoduchom modeli vyjadruje koeficient celkový nekorigovaný vzťah.
- Vo viacnásobnom modeli doplňujeme “pri rovnakých ostatných premenných” preto, lebo izolujeme čistý vplyv aktív pri očistení od časového vývoja (roku zlyhania).
- Pomocou
model1vypočítajte predikciu pre Schooling = 13. Overte ju ručným dosadením do odhadnutej rovnice.
- Ručný výpočet: \(\widehat{\text{Deposit\_Mil}} = 18{,}064 + 0{,}7562 \times 1000 = 774{,}264 \text{ mil. USD}\).
- Nájdite v tabuľke
porovnaniejednu krajinu s kladným a jednu so záporným rezíduom. Vysvetlite význam znamienka.
- Kladné rezíduum: Skutočný objem vkladov banky je vyšší, než predikoval model na základe jej aktív.
- Záporné rezíduum: Skutočný objem vkladov banky je nižší, než predikoval model.
- Vysvetlite, prečo z kladného koeficienta nemožno automaticky vyvodiť kauzálny účinok.
- Kladný koeficient indikuje len štatistickú asociáciu (koreláciu), nie kauzalitu. Vyššie aktíva a vklady sa vyskytujú spoločne z dôvodu celkovej veľkosti banky, nie preto, že by jedno automaticky vytváralo druhé.