Práca s regresiou

V tejto práci sa budeme zaoberať prvými krokmi v odhade regresnej funkcie. Budeme využívať databázu Life Expectancy (WHO) Fixed dataset, ktorá obsahuje ukazovatele na úrovni krajín, ako je priemerná dĺžka života, HDP, školské vzdelávanie a výdavky na zdravotníctvo.

Využívanie niektorých knižníc

rm(list=ls())
library(lmtest)   #  podpora regresie
library(outliers) # analyza odlahlych hodnot (outliers)
library(gptstudio)
library(kableExtra)
library(knitr)
library(dplyr)
library(broom)
library(corrplot)

Príprava údajov - import z csv súboru

Súbor Life_Expectancy_Data obsahuje databázu determinantov očakávanej dĺžky života. Import údajov urobíme nasledovne

# import the dataset and create a data.frame udaje
#
udaje_svet <- read.csv("udaje/Life-Expectancy-Data-Updated.csv",header=TRUE,sep=",",dec=".",check.names = TRUE)
head(udaje_svet)

#

Databáza obsahuje údaje o 2938 pozorovaniach a 22 premenných. V tejto práci sa budeme zaoberať len časťou z nich, konkrétne tými, ktoré súvisia s dĺžkou dožitia. Na začiatku si vyberieme krajinu, ktorej zdravotný stav chceme analyzovať. V tomto prípade ide o Portugalsko (Slovensko v databáze nemá zastúpenie):

# z databázy udaje_svet si vyberieme len tie pozorovania, ktoré sa týkaj Portugalska 
udaje <- subset(udaje_svet, Country == "Portugal")

Tabuľka uvedená nižšie nám poskuytuje základné popisné štatistiky vybraných kvantitatívnych premenných.

# niektoré štatistiky a ich prehľad v tabuľke KableExtra
library(kableExtra)
udaje %>%
  select(Adult_mortality,Alcohol_consumption, Hepatitis_B, Measles, BMI,Polio,Diphtheria,GDP_per_capita,Population_mln,Life_expectancy,Schooling) %>%
  summary() %>%
  kable() %>%
  kable_styling(bootstrap_options = "striped", full_width = FALSE)

	Adult_mortality	Alcohol_consumption	Hepatitis_B	Measles	BMI	Polio	Diphtheria	GDP_per_capita	Population_mln	Life_expectancy	Schooling
	Min. : 73.28	Min. : 9.52	Min. :58.00	Min. :48.00	Min. :25.50	Min. :93.00	Min. :93.00	Min. :18585	Min. :10.29	Min. :76.30	Min. :6.800
	1st Qu.: 83.07	1st Qu.:10.73	1st Qu.:94.00	1st Qu.:77.25	1st Qu.:25.70	1st Qu.:96.00	1st Qu.:96.75	1st Qu.:18830	1st Qu.:10.41	1st Qu.:77.58	1st Qu.:7.350
	Median : 90.46	Median :11.52	Median :97.00	Median :95.00	Median :25.75	Median :96.50	Median :97.00	Median :19167	Median :10.49	Median :78.45	Median :7.750
	Mean : 91.73	Mean :11.41	Mean :91.56	Mean :85.00	Mean :25.73	Mean :96.50	Mean :97.00	Mean :19199	Mean :10.47	Mean :78.74	Mean :7.831
	3rd Qu.:100.28	3rd Qu.:12.23	3rd Qu.:97.25	3rd Qu.:95.25	3rd Qu.:25.80	3rd Qu.:97.25	3rd Qu.:98.00	3rd Qu.:19399	3rd Qu.:10.54	3rd Qu.:80.42	3rd Qu.:8.350
	Max. :109.63	Max. :12.77	Max. :98.00	Max. :96.00	Max. :25.80	Max. :98.00	Max. :99.00	Max. :19986	Max. :10.57	Max. :81.10	Max. :9.100

Vyššie uvedená tabuľka nám poskytuje prehľad o základných štatistických charakteristikách vybraných premenných, ako sú priemerné hodnoty, rozptyl, minimum a maximum. Tieto informácie nám pomáhajú lepšie pochopiť rozdelenie a rozsah hodnôt v našich dátach. Na druhej strane je zaujímavá aj informácia o vzájomných vzťahoch medzi týmito premennými, čo môžeme merať pomocou korelačnej matice.

# grafický prehľad o korelačných vzťahoch vyjadruje nasledovný obrázok

cor_matrix <- cor(udaje %>% select(Adult_mortality,Alcohol_consumption
,Hepatitis_B, Measles, BMI,Polio,Diphtheria,GDP_per_capita,Population_mln,Life_expectancy,Schooling), use = "complete.obs")
corrplot(cor_matrix, method = "color", type = "upper", tl.col =
"black", tl.srt = 45, title = "Korelačná matica vybraných premenných", mar = c(0,0,1,0))

Uvedený graf nám poskytuje vizuálny prehľad o korelačných vzťahoch medzi vybranými premennými. Farby a intenzita farieb nám umožňujú rýchlo identifikovať silné pozitívne alebo negatívne korelácie. Upozorňujeme, že korelácia neznamená kauzalitu.

V nasledovnom je uvedený graf vývoja očakávanej dĺžky dožitia v Portugalsku v rokoch 2000-2015. Vidíme, že očakávaná dĺžka života sa zvyšovala, čo je pozitívny trend.

#  graf vývoja očakávanej dĺžky dožitia v Portugalsku v rokoch 2000-2015.
library(ggplot2)
ggplot(udaje, aes(x = Year, y = Life_expectancy)) +
  geom_line() +
  geom_point() +
  labs(title = "Vývoj očakávanej dĺžky dožitia v Portugalsku (2000-2015)",
       x = "Rok",
       y = "Očakávaná dĺžka dožitia") +
  theme_minimal()

Na začiatku sa pokúsme o vyrovnanie priebehu tejto premennej v čase pomocou lineárnej regresie, kde nezávislou premennou bude rok a závislou premennou bude očakávaná dĺžka dožitia. Odhadneme koeficienty tejto regresie a posúdime kvalitu vyrovnania pomocou ukazovateľov, ako je R-squared a p-value.

# vyrovnanie priebehu očakávanej dĺžky dožitia v čase
model <- lm(Life_expectancy ~ Year, data = udaje)

library(broom)
library(knitr)
library(kableExtra)

# koeficienty regresie
tidy(model) %>%
  kable(digits = 3, caption = "Odhadnuté koeficienty regresie") %>%
  kable_styling(bootstrap_options = "striped", full_width = FALSE)

Odhadnuté koeficienty regresie

term	estimate	std.error	statistic	p.value
(Intercept)	-575.170	35.209	-16.336	0
Year	0.326	0.018	18.572	0

# kvalita vyrovnania
glance(model) %>%
  kable(digits = 3, caption = "Ukazovatele kvality vyrovnania") %>%
  kable_styling(bootstrap_options = "striped", full_width = FALSE)

Ukazovatele kvality vyrovnania

r.squared	adj.r.squared	sigma	statistic	p.value	df	logLik	AIC	BIC	deviance	df.residual	nobs
0.961	0.958	0.323	344.936	0	1	-3.573	13.145	15.463	1.464	14	16

NA

Výsledky regresie nám ukazujú, že koeficient pre rok je pozitívny a štatisticky významný, čo naznačuje, že očakávaná dĺžka dožitia v Portugalsku sa zvyšovala v priebehu rokov 2000-2015. Jej priemerný ročný nárast dosahoval až štyri roky. Hodnota R-squared hodnota nám hovorí, že model vysvetľuje 96,1 % variability modelu. Podľa hodnoty p-value ( ktorý nám vyšiel 0) môžeme povedať, že model ako celok je štatisticky významný.

# teraz vyššie uvedený obrázok doplníme o regresnú priamku
ggplot(udaje, aes(x = Year, y = Life_expectancy)) +
  geom_line() +
  geom_point() +
  geom_smooth(method = "lm", se = FALSE, color = "red") +
  labs(title = "Vývoj očakávanej dĺžky dožitia v Portugalsku (2000-2015) s regresnou priamkou",
       x = "Rok",
       y = "Očakávaná dĺžka dožitia") +
  theme_minimal()

Vyvoj spotreby alkoholu

Alkohol je predpokladaný negatívny faktor, ktorý vplýva na celkové zdravie a tým aj na dĺžku života.Zvyšuje krvný tlak, má negatívny vplyv na pružnosť ciev.P- value nám vyšla 0, ciže nižšia ako 0,1 a hodnota r-squared vyšla zhruba 0,88 - to znamemná,že v 88% prípadoch mala naozaj konzumácia aloholu veľký vplyv na dĺžku dožitia.

model_alkohol <-lm(Alcohol_consumption ~ Year, data = udaje)
tidy(model_alkohol) %>%
  kable(digits = 3, caption = "Odhadnuté koeficienty regresie pre spotrebu alkoholu") %>%
  kable_styling(bootstrap_options = "striped", full_width = FALSE)

Odhadnuté koeficienty regresie pre spotrebu alkoholu

term	estimate	std.error	statistic	p.value
(Intercept)	416.781	38.734	10.760	0
Year	-0.202	0.019	-10.465	0

# and now the model quality statistics - just R sqared
r2 <- summary(model_alkohol)$r.squared
adj_r2 <- summary(model_alkohol)$adj.r.squared
r2

[1] 0.8866618

adj_r2

[1] 0.8785662

Obrázok o skutočných vyrovnaných hodnotách spotreby alkoholu

fitted_vals <- fitted(model_alkohol)
# vykreslenie skutočných a vyrovnaných hodnôt - Alcohol_consumption a fitted_values
ggplot(udaje, aes(x = Year)) +
  geom_line(aes(y = Alcohol_consumption), color = "blue", size = 1) +
  geom_line(aes(y = fitted_vals), color = "red", size = 1, linetype = "dashed") +
  labs(title = "Skutočné vs. Vyrovnané hodnoty spotreby alkoholu v Rakúsku (2000-2015)",
       x = "Rok",
       y = "Spotreba alkoholu") +
  theme_minimal() +
  scale_y_continuous(limits = c(0, max(udaje$Alcohol_consumption, fitted_vals) * 1.1)) +
  theme(legend.position = "none")