Analýza stravovacích návykov

Author

Marek Holik

1. Úvod

Dataset Estimation of Obesity Levels Based on Eating Habits and Physical Condition pochádza z repozitára UCI Machine Learning Repository. Dataset bol publikovaný v roku 2019 a je určený na analýzu a predikciu úrovne obezity jednotlivcov na základe ich demografických charakteristík, stravovacích návykov a fyzickej aktivity.

Dáta boli získané od respondentov z krajín Mexiko, Peru a Kolumbia prostredníctvom online dotazníka.

Vzhľadom na nevyvážené zastúpenie jednotlivých kategórií obezity boli pôvodné dáta rozšírené pomocou syntetických vzoriek generovaných technikou SMOTE (Synthetic Minority Over-sampling Technique) v prostredí Weka. Teda dáta sú kombináciou reálnych odpovedí respondentov (23 %) a synteticky generovaných dát (77 %).

Dataset obsahuje 2111 pozorovaní a 17 premenných, z ktorých jedna je cieľová (úroveň obezity) a zvyšné predstavujú vysvetľujúce premenné.

Zdroj: archive.ics.uci.edu

2. Význam Jednotlivých Premenných

Demografické a Antropometrické Premenné

Gender : Pohlavie respondenta (Male, Female).

Age : Vek respondenta v rokoch.

Height : Výška respondenta v metroch.

Weight : Hmotnosť respondenta v kilogramoch.

family_history_with_overweight : Informácia o výskyte nadváhy alebo obezity v rodine (yes, no).

Stravovacie Návyky

FAVC : Indikátor častej konzumácie vysokokalorických jedál (yes, no).

FCVC : Frekvencia konzumácie zeleniny vyjadrená na numerickej škále.

NCP : Počet hlavných jedál skonzumovaných počas dňa.

CAEC : Frekvencia konzumácie jedla medzi hlavnými jedlami (no, Sometimes, Frequently, Always).

CH2O : Priemerný denný príjem vody vyjadrený numerickou škálou.

CALC : Frekvencia konzumácie alkoholických nápojov (no, Sometimes, Frequently, Always).

Životný Štýl a Fyzická Aktivita

SMOKE : Informácia o tom, či respondent fajčí (yes, no).

SCC : Indikátor sledovania alebo kontroly kalorického príjmu (yes, no).

FAF : Frekvencia vykonávania fyzickej aktivity vyjadrená numerickou škálou.

TUE : Čas strávený používaním technológií ako televízia, počítač alebo mobilné zariadenia.

MTRANS : Preferovaný spôsob dopravy (Automobile, Bike, Motorbike, Public_Transportation, Walking).

Cieľová Premenná

NObeyesdad : Kategória úrovne obezity respondenta (Insufficient Weight, Normal Weight, Overweight Level I, Overweight Level II, Obesity Type I, Obesity Type II, Obesity Type III).

3. Regresná Úloha

3.1 Číselný a Grafický Súhrn Premenných

V tejto kapitole je uvedený číselný a grafický popis rozdelenia pravdepodobnosti jednotlivých premenných, ako aj ich vzájomných vzťahov.

Načítanie Dát a Príprava Prostredia

Načítanie knižníc a dát 
library(ggplot2)
library(dplyr)
library(tidyr)
library(knitr)
library(gt) 
library(car)
library(MASS)
library(glmnet)

data <- read.csv("ObesityDataSet.csv", sep = ",")
Premenovanie premenných 
data <- data %>%
  rename(
    gender = Gender,
    age_years = Age,
    height_m = Height,
    weight_kg = Weight,
    family_history_overweight = family_history_with_overweight,

    high_calorie_food_frequent = FAVC,
    vegetable_consumption_frequency = FCVC,
    meals_per_day = NCP,
    snacking_between_meals = CAEC,
    daily_water_intake = CH2O,

    smoking = SMOKE,
    calorie_monitoring = SCC,
    physical_activity_frequency = FAF,
    technology_use_time = TUE,
    alcohol_consumption = CALC,
    transportation_mode = MTRANS,

    obesity_level = NObeyesdad
  )
Vypočítanie štatistík číselných premenných 
num_data <- data[, sapply(data, is.numeric)]

num_summary <- num_data %>%
  summarise(
    across(
      everything(),
      list(
        mean = mean,
        sd = sd,
        min = min,
        q1 = ~quantile(.x, 0.25),
        median = median,
        q3 = ~quantile(.x, 0.75),
        max = max
      ),
      .names = "{.col}_{.fn}"
    )
  ) %>%

  pivot_longer(
    everything(),
    names_to = "variable_stat",
    values_to = "value"
  ) %>%
 
  separate(variable_stat, into = c("variable", "stat"), sep = "_(?=[^_]+$)") %>%

  pivot_wider(names_from = stat, values_from = value)


num_summary %>% knitr::kable(caption = "Súhrn numerických premenných", )
Súhrn numerických premenných
variable mean sd min q1 median q3 max
age_years 24.3125999 6.3459683 14.00 19.947192 22.777890 26.000000 61.00
height_m 1.7016774 0.0933048 1.45 1.630000 1.700499 1.768464 1.98
weight_kg 86.5860581 26.1911717 39.00 65.473343 83.000000 107.430682 173.00
vegetable_consumption_frequency 2.4190431 0.5339266 1.00 2.000000 2.385502 3.000000 3.00
meals_per_day 2.6856280 0.7780386 1.00 2.658738 3.000000 3.000000 4.00
daily_water_intake 2.0080114 0.6129535 1.00 1.584812 2.000000 2.477420 3.00
physical_activity_frequency 1.0102977 0.8505924 0.00 0.124505 1.000000 1.666678 3.00
technology_use_time 0.6578659 0.6089273 0.00 0.000000 0.625350 1.000000 2.00

Vizualizácia Štatistík pre Numerické Premenné

Histogramy numerických premenných 
num_data_long <- num_data %>%
  pivot_longer(cols = everything(), names_to = "variable", values_to = "value")

ggplot(num_data_long, aes(x = value)) +
  geom_histogram(aes(y = after_stat(density)), bins = 30, fill = "skyblue", color = "black", alpha = 0.6) +
  geom_density(color = "red", size = 1) +
  facet_wrap(~variable,nrow=2, ncol=4, scales = "free") +
  theme_minimal(base_size = 14) +
     theme(
    legend.position = "none",
    strip.text = element_text(face = "bold", size = 12),
    plot.title = element_text(face = "bold", hjust = 0.5, size = 16),
    axis.text.x = element_blank(),
    axis.text.y = element_text(size = 12),
    axis.title = element_text(size = 14),
    panel.spacing = unit(1, "lines")
  ) +
  labs(title = "Histogramy numerických premenných s hustotou pravdepodobnosti", x = "Hodnota", y = "Hustota")

Boxploty číselných charakteristík 
ggplot(num_data_long, aes(x = "", y = value, fill = variable)) +
  geom_boxplot(alpha = 0.6, outlier.color = "red", outlier.shape = 16) +
  facet_wrap(~variable, nrow = 2, ncol = 4, scales = "free") +
  labs(
    title = "Boxploty vypočítaných číselných premenných",
    x = "",
    y = "Hodnota"
  ) +
  theme_minimal(base_size = 14) +
  theme(
    legend.position = "none",
    strip.text = element_text(face = "bold", size = 12),
    plot.title = element_text(face = "bold", hjust = 0.5, size = 16),
    axis.text.x = element_blank(),
    axis.text.y = element_text(size = 12),
    axis.title = element_text(size = 14),
    panel.spacing = unit(1, "lines")
  )

Korelačná mapa 
corr_matrix <- cor(num_data, use = "pairwise.complete.obs")
corr_long <- as.data.frame(as.table(corr_matrix))
names(corr_long) <- c("Var1", "Var2", "Correlation")

# heatmapa
ggplot(corr_long, aes(x = Var1, y = Var2, fill = Correlation)) +
  geom_tile(color = "white") +
  scale_fill_gradient2(low = "blue", mid = "white", high = "red", midpoint = 0) +
  geom_text(aes(label = round(Correlation, 2)), color = "black", size = 3) +
  theme_minimal(base_size = 14) +
  theme(
    axis.text.x = element_text(angle = 45, hjust = 1, face = "bold"),
    axis.text.y = element_text(face = "bold"),
    plot.title = element_text(face = "bold", hjust = 0.5, size = 16),
    panel.grid = element_blank()
  ) +
  labs(
    title = "Heatmapa korelácií numerických premenných",
    x = "",
    y = ""
  )

Štatistiky pre Kategoriálne Premenné

Barploty pre kategoriálne premenné 
cat_data <- data[, sapply(data, function(x) !is.numeric(x))]
cat_data_long <- cat_data %>%
  pivot_longer(cols = everything(), names_to = "variable", values_to = "value")

ggplot(cat_data_long, aes(x = value, fill = value)) +
  geom_bar() +
  facet_wrap(~variable, scales = "free") +
  theme_minimal(base_size = 14) +
  theme(
    legend.position = "none",
    strip.text = element_text(face = "bold", size = 12),
    plot.title = element_text(face = "bold", hjust = 0.5, size = 16),
    axis.text.x = element_text(angle = 45, hjust = 1),
    axis.text.y = element_text(size = 12),
    axis.title = element_text(size = 14),
    panel.spacing = unit(1, "lines")
  ) +
  labs(
    title = "Barploty kategóriálnych premenných",
    x = "Hodnota",
    y = "Počet"
  )

3.2 Odhady Lineárnych Predpovedných Modelov

One-Hot Encoding - Transformácia Kategoriálnych Premenných 
cat_vars <- names(data)[sapply(data, function(x) is.character(x) | is.factor(x))]
dummy_matrix <- model.matrix(~ . - 1, data = data[, cat_vars])

dummy_df <- as.data.frame(dummy_matrix)
data <- bind_cols(num_data, dummy_df)
Min Max Rescaler 
#min_max <- function(x) {
#  (x - min(x, na.rm = TRUE)) / (max(x, na.rm = TRUE) - min(x, na.rm = TRUE))
#}
#
#num_vars <- names(data)[sapply(data, is.numeric)]
#
#data[num_vars] <- lapply(data[num_vars], min_max)
Výpočet Štatistík 
summ <- data %>%
  summarise(
    across(
      everything(),
      list(
        mean = mean,
        sd = sd,
        min = min,
        q1 = ~quantile(.x, 0.25),
        median = median,
        q3 = ~quantile(.x, 0.75),
        max = max
      ),
      .names = "{.col}_{.fn}"
    )
  ) %>%

  pivot_longer(
    everything(),
    names_to = "variable_stat",
    values_to = "value"
  ) %>%
 
  separate(variable_stat, into = c("variable", "stat"), sep = "_(?=[^_]+$)") %>%

  pivot_wider(names_from = stat, values_from = value)


summ %>% knitr::kable(caption = "Súhrn preškálovaných numerických premenných", )
Súhrn preškálovaných numerických premenných
variable mean sd min q1 median q3 max
age_years 24.3125999 6.3459683 14.00 19.947192 22.777890 26.000000 61.00
height_m 1.7016774 0.0933048 1.45 1.630000 1.700499 1.768464 1.98
weight_kg 86.5860581 26.1911717 39.00 65.473343 83.000000 107.430682 173.00
vegetable_consumption_frequency 2.4190431 0.5339266 1.00 2.000000 2.385502 3.000000 3.00
meals_per_day 2.6856280 0.7780386 1.00 2.658738 3.000000 3.000000 4.00
daily_water_intake 2.0080114 0.6129535 1.00 1.584812 2.000000 2.477420 3.00
physical_activity_frequency 1.0102977 0.8505924 0.00 0.124505 1.000000 1.666678 3.00
technology_use_time 0.6578659 0.6089273 0.00 0.000000 0.625350 1.000000 2.00
genderFemale 0.4940786 0.5000834 0.00 0.000000 0.000000 1.000000 1.00
genderMale 0.5059214 0.5000834 0.00 0.000000 1.000000 1.000000 1.00
family_history_overweightyes 0.8176220 0.3862473 0.00 1.000000 1.000000 1.000000 1.00
high_calorie_food_frequentyes 0.8839413 0.3203712 0.00 1.000000 1.000000 1.000000 1.00
snacking_between_mealsFrequently 0.1146376 0.3186596 0.00 0.000000 0.000000 0.000000 1.00
snacking_between_mealsno 0.0241592 0.1535795 0.00 0.000000 0.000000 0.000000 1.00
snacking_between_mealsSometimes 0.8360966 0.3702756 0.00 1.000000 1.000000 1.000000 1.00
smokingyes 0.0208432 0.1428931 0.00 0.000000 0.000000 0.000000 1.00
calorie_monitoringyes 0.0454761 0.2083952 0.00 0.000000 0.000000 0.000000 1.00
alcohol_consumptionFrequently 0.0331596 0.1790957 0.00 0.000000 0.000000 0.000000 1.00
alcohol_consumptionno 0.3027001 0.4595354 0.00 0.000000 0.000000 1.000000 1.00
alcohol_consumptionSometimes 0.6636665 0.4725665 0.00 0.000000 1.000000 1.000000 1.00
transportation_modeBike 0.0033160 0.0575025 0.00 0.000000 0.000000 0.000000 1.00
transportation_modeMotorbike 0.0052108 0.0720146 0.00 0.000000 0.000000 0.000000 1.00
transportation_modePublic_Transportation 0.7484604 0.4340007 0.00 0.000000 1.000000 1.000000 1.00
transportation_modeWalking 0.0265277 0.1607365 0.00 0.000000 0.000000 0.000000 1.00
obesity_levelNormal_Weight 0.1359545 0.3428215 0.00 0.000000 0.000000 0.000000 1.00
obesity_levelObesity_Type_I 0.1662719 0.3724128 0.00 0.000000 0.000000 0.000000 1.00
obesity_levelObesity_Type_II 0.1406916 0.3477855 0.00 0.000000 0.000000 0.000000 1.00
obesity_levelObesity_Type_III 0.1534818 0.3605367 0.00 0.000000 0.000000 0.000000 1.00
obesity_levelOverweight_Level_I 0.1373757 0.3443251 0.00 0.000000 0.000000 0.000000 1.00
obesity_levelOverweight_Level_II 0.1373757 0.3443251 0.00 0.000000 0.000000 0.000000 1.00

Model 0: Konštantný Model

Konštantný model predikuje jednu hodnotu pre všetky pozorovania – priemer zvolenej premennej.

Naučenie Modelu 
model0 <- lm(weight_kg ~ 1, data = data)
Vizualizácia Modelu 
ggplot(data, aes(x = seq_along(weight_kg), y = weight_kg)) +
  geom_point(alpha = 0.4) +
  geom_hline(
    yintercept = fitted(model0)[1],
    color = "red",
    linewidth = 1
  ) +
  labs(
    x = "Pozorovanie",
    y = "Hmotnosť",
    title = "Model0: konštantný model (iba priemer)"
  ) +
  theme_minimal()

Model 0: Analýza Reziduií

Analýza Reziduií 
par(mfrow = c(1, 3))

plot(
  fitted(model0),
  resid(model0),
  main = "Residuals vs Fitted",
  xlab = "Fitted values",
  ylab = "Residuals"
)
abline(h = 0, col = "red")

hist(
  resid(model0),
  breaks = 20,
  main = "Histogram rezíduí",
  xlab = "Rezíduá"
)

qqnorm(resid(model0))
qqline(resid(model0), col = "red")

Residuals vs Fitted

Graf rezíduí voči odhadnutým hodnotám slúži na kontrolu linearity a homoskedasticity.

V prípade konštantného modelu sú všetky odhadnuté hodnoty totožné, čo vedie k vertikálnemu rozloženiu bodov.

Rezíduá majú veľký rozptyl a nevykazujú žiadnu systematickú štruktúru.

To potvrdzuje, že model nedokáže zachytiť variabilitu závislej premennej a neposkytuje žiadnu vysvetľujúcu informáciu.

Histogram rezíduí

Histogram rezíduí je približne symetrický, avšak s výrazným rozptylom.

To naznačuje, že hoci rozdelenie rezíduí môže byť približne normálne, ich variancia je vysoká.

Tento výsledok je očakávaný, keďže model neobsahuje žiadne vysvetľujúce premenné a všetka variabilita zostáva v rezíduách.

Normal Q–Q Plot

Q–Q graf ukazuje približne lineárny priebeh v strednej časti rozdelenia, avšak s odchýľkami v chvostoch.

To naznačuje mierne porušenie predpokladu normality rezíduií, čo je pri konštantnom modeli bežné a nie je považované za závažný problém.

Odchýlky sú dôsledkom veľkej nevysvetlenej variability v dátach.

Model 1: Jednoduchá Lineárna Regresia s Jednou Premennou

Model 1 je jednoduchý lineárny regresný model, ktorý skúma vzťah medzi telesnou hmotnosťou a denným príjmom vody.

Trénovanie Modelu 
model1 <- lm(weight_kg ~ daily_water_intake, data = data)
Vizualizácia Modelu 
ggplot(data, aes(x = daily_water_intake, y = weight_kg)) +
  geom_point(alpha = 0.4) +
  geom_smooth(method = "lm", formula = y ~ x, color = "red", se = FALSE) +
  labs(
    x = "Denný príjem vody (škálovaný)",
    y = "Hmotnosť (škálovaná)",
    title = "Model1: weight_kg ~ daily_water_intake"
  ) +
  theme_minimal()

V tomto prípade hmotnosť osoby (v kilogramoch) je slabo lineárne závislá od denného množstva prijatej vody.

Model 1: Analýza Reziduií

Analýza Reziduií 
par(mfrow = c(2, 2))
plot(model1)

Residuals vs Fitted

V grafe Residuals vs Fitted sú reziduá rozptýlené okolo nulovej hodnoty bez výrazného systematického trendu, čo naznačuje, že lineárny tvar modelu je v zásade prijateľný.

Q–Q Residuals

Graf Q–Q Residuals ukazuje odchýlky reziduí od normálneho rozdelenia, najmä v chvostoch.

Reziduá v extrémnych kvantiloch sa odchyľujú od referenčnej priamky, čo naznačuje porušenie predpokladu normality.

Toto je však pri jednoduchých behaviorálnych dátach bežné a neznamená nutne neplatnosť modelu, môže však ovplyvniť inferenciu.

Scale–Location

Je viditeľná mierna zmena variability reziduí v závislosti od predikovaných hodnôt.

To naznačuje slabú heteroskedasticitu, hoci nejde o výrazné porušenie predpokladov.

Residuals vs Leverage

Graf Residuals vs Leverage neodhaľuje pozorovania s extrémne vysokou pákou (leverage) ani body presahujúce hranice Cookovej vzdialenosti.

To znamená, že model nie je výrazne ovplyvnený jednotlivými extrémnymi pozorovaniami.

Model 2: Viacnásobná Lineárna Regresia

Model 2 je viacnásobný lineárny regresný model, ktorý rozširuje jednoduchý prístup z Modelu 1 o viacero vysvetľujúcich premenných.

Model predpokladá, že telesná hmotnosť (v kilogramoch) je lineárne ovplyvnená kombináciou:

  • telesnej výšky,
  • veku,
  • frekvencie fyzickej aktivity,
  • počtu jedál za deň,
  • frekvencie konzumácie zeleniny.
Trénovanie Modelu 
model2 <- lm(weight_kg ~ height_m + age_years + physical_activity_frequency + meals_per_day + vegetable_consumption_frequency,
  data = data
)
Vizualizácia Modelu 
ggplot(data, aes(x = fitted(model2), y = weight_kg)) +
  geom_point(alpha = 0.4) +
  geom_abline(intercept = 0 , slope = 1, color = "red", linewidth = 1) +
  labs(
    x = "Predikovaná hmotnosť",
    y = "Skutočná hmotnosť",
    title = "Model2: weight_kg ~ height_m + age_years + physical_activity_frequency + meals_per_day + vegetable_consumption_frequency"
  ) +
  theme_minimal()

Model 2: Kontrola Kolinearity Viacnásobnej Lineárnej Regresie

VIF (Variance Inflation Factor) slúži na diagnostiku multikolinearity vo viacnásobnej lineárnej regresii, teda na zistenie, do akej miery sú vysvetľujúce premenné navzájom lineárne závislé.

VIF pre jednotlivú premennú vyjadruje, koľkokrát je rozptyl (variancia) odhadu jej regresného koeficientu zväčšený v dôsledku korelácie s ostatnými premennými v modeli. Vyššia hodnota VIF znamená menej stabilný odhad koeficientu a horšiu interpretovateľnosť.

Interpretácia hodnôt VIF:

  • VIF ≈ 1 – premenná nie je korelovaná s ostatnými, multikolinearita nie je problém
  • VIF 1–5 – mierna korelácia, zvyčajne akceptovateľná
  • VIF > 5 (niekedy > 10) – silná multikolinearita, potenciálny problém modelu

V kontexte Modelu 2 VIF ukazuje, či sa výška, vek, frekvencia fyzickej aktivity, počet jedál za deň a frekvencia konzumácie zeleniny navzájom významne prekrývajú v informácii, ktorú nesú.

Kolinearita 
vif(model2) %>% knitr::kable(caption="Tabulka VIF koeficientov", )
Tabulka VIF koeficientov
x
height_m 1.153522
age_years 1.023200
physical_activity_frequency 1.123063
meals_per_day 1.071119
vegetable_consumption_frequency 1.005642

Nízke hodnoty VIF znamenajú, že parciálne vplyvy premenných sú spoľahlivo interpretovateľné, zatiaľ čo vysoké hodnoty naznačujú potrebu úpravy modelu (odstránenie premenných, transformácie alebo regularizácia), tieto tu však nie sú.

Model 2: Analýza Reziduií

Analýza Reziduií 
par(mfrow = c(2,2))
plot(model2)

Residuals vs Fitted

Rezíduá sú rozptýlené okolo nuly bez výrazného trendu.

Variancia rezíduí je relatívne konštantná, čo naznačuje, že homoskedasticita je v poriadku.

Model zachytáva viac štruktúry než model1.

Normal Q–Q Plot

Body sa väčšinou nachádzajú na priamke

Menšie odchýlky v chvostoch naznačujú mierne porušenie normality, typické pre reálne dáta.

Scale–Location

Rezíduá majú mierne stúpajúci trend pri vyšších fitted hodnotách.

To naznačuje miernu heteroskedasticitu – rozptyl rezíduí rastie s predikovanou hmotnosťou.

Väčšina bodov je blízko 1 (štandardizované), takže problém nie je extrémny, ale je vhodné ho spomenúť pri interpretácii modelu.

Residuals vs Leverage

Väčšina bodov má nízku hodnotu leverage → žiadne extrémne vplyvné pozorovania.

Niekoľko bodov na pravej strane má vyšší leverage a mierne vyšší Cookov distance (číslované body: 93°, 134°, 233°), ale nie sú natoľko extrémne, aby dominovali modelu.

Model nie je výrazne ovplyvnený jednotlivými pozorovaniami.

Model 3: Model s Interakciami

Trénovanie Modelu 
model3 <- lm(weight_kg ~ height_m * physical_activity_frequency + age_years + meals_per_day, data = data)
Vizualizácia Modelu 
ggplot(data, aes(x = fitted(model3), y = weight_kg)) +
  geom_point(alpha = 0.4) +
  geom_abline(intercept = 0, slope = 1, color = "red", linewidth = 1) +
  labs(
    x = "Predikovaná hmotnosť",
    y = "Skutočná hmotnosť",
    title = "Model3: weight_kg ~ height_m * physical_activity_frequency + age_years + meals_per_day"
  ) +
  theme_minimal()

Model 3: Kontrola Kolinearity Modelu s Interakciami

Kolinearita 
vif(model3) %>% knitr::kable(caption="Tabulka VIF koeficientov", )
Tabulka VIF koeficientov
x
height_m 2.787642
physical_activity_frequency 372.603931
age_years 1.023032
meals_per_day 1.069062
height_m:physical_activity_frequency 388.599200

Na základe výsledkov VIF pre model3 vidíme, že niektoré premenné vykazujú silnú kolinearitu. Vysoký VIF pre interakčný člen je ale bežný, pretože interakcia je priamo odvodená z dvoch pôvodných premenných.

Model 3: Analýza Reziduií

Analýza Reziduií 
par(mfrow = c(2, 2))
plot(model3)

Residuals vs Fitted

Rezíduá sú rozptýlené okolo nuly bez systematického trendu.

Variancia je približne konštantná → homoskedasticita je splnená.

Mierne odchýlky pri extrémnych fitted hodnotách sú zanedbateľné.

Normal Q–Q Plot

Body sledujú približne lineárnu priamku.

Mierne odchýlky v chvostoch, čo je typické pre reálne dáta.

Predpoklad normality je dostatočne splnený.

Scale–Location

Variancia rezíduí je stabilná naprieč fitted hodnotami.

Žiadne významné systematické odchýlky.

Residuals vs Leverage

Nie sú prítomné výrazne vplyvné pozorovania.

Model3 nie je dominovaný jednotlivými bodmi.

Model3, ktorý zahŕňa interakciu medzi výškou a frekvenciou fyzickej aktivity spolu s vekom a počtom jedál denne, poskytuje najlepšiu popisnú schopnosť spomedzi všetkých navrhnutých modelov.

Rezíduá sú rozptýlené okolo nuly, približne normálne rozdelené, variancia je stabilná a nenachádzajú sa výrazne vplyvné pozorovania.

Interakčný člen umožňuje, aby efekt fyzickej aktivity závisel od výšky, čo má biologický a behaviorálny zmysel.

Model3 je preto najvhodnejší pre interpretáciu skutočného regresného vzťahu medzi hmotnosťou a vybranými vysvetľujúcimi premennými.

3.3 Porovnanie Modelov

Porovnanie Modelov z Predpovedného Hľadiska

V tabuľke nižšie sú zhrnuté metriky štyroch lineárnych modelov pre predikciu hmotnosti na základe rôznych premenných.

Porovnanie Modelov a Chýb 
comparison <- data.frame(
  Model = c("model0", "model1", "model2", "model3"),
  R2 = c(
    summary(model0)$r.squared,
    summary(model1)$r.squared,
    summary(model2)$r.squared,
    summary(model3)$r.squared
  ),
  Adj_R2 = c(
    summary(model0)$adj.r.squared,
    summary(model1)$adj.r.squared,
    summary(model2)$adj.r.squared,
    summary(model3)$adj.r.squared
  ),
  AIC = c(
    AIC(model0),
    AIC(model1),
    AIC(model2),
    AIC(model3)
  )
)

comparison %>% knitr::kable(caption = "Porovnanie modelov", )
Porovnanie modelov
Model R2 Adj_R2 AIC
model0 0.0000000 0.0000000 19780.37
model1 0.0402305 0.0397754 19695.69
model2 0.3446892 0.3431326 18898.17
model3 0.2893388 0.2876508 19069.34
  • model0: Nevysvetľuje žiadnu variabilitu hmotnosti ((R^2 = 0)) a má veľmi vysoké AIC (19780.37), čo naznačuje nevhodnosť modelu.
  • model1: Pridanie prediktorov prinieslo len minimálne zlepšenie ((R^2 = 0.0402)), AIC sa mierne znížila (19695.69), ale model zostáva slabý.
  • model2: Podstatne vyššie (R^2 = 0.358) a Adjusted (R^2 = 0.356) znamenajú, že model vysvetľuje približne 36 % variability hmotnosti. Hodnota AIC (18857.83) je stále vysoká, model však lepšie zachytáva vzťahy v dátach.
  • model3: Nižšie (R^2 = 0.289) a Adjusted (R^2 = 0.288) znamenajú mierne slabšie vysvetlenie, avšak AIC (19069.34) je nižšia než u model0 a model1, čo naznačuje lepší kompromis medzi kvalitou prispôsobenia a zložitosťou modelu.

3.4 Odhad Testovacej Chyby

Rozdelenie Dát na Tréningovú a Testovaciu Množinu

Najprv rozdelíme dataset, 70 % tréningovú a 30 % testovacú množinu. Potom vyhodnotíme modely pomocou RMSE na testovacej množine.

Rozdelenie Dát Train/Test 
set.seed(123)

train_idx <- sample(seq_len(nrow(data)), size = 0.7 * nrow(data))
train_data <- data[train_idx, ]
test_data <- data[-train_idx, ]


rmse <- function(actual, predicted) {
  sqrt(mean((actual - predicted)^2))
}

model0_train <- lm(weight_kg ~ 1, data = train_data)
model1_train <- lm(weight_kg ~ daily_water_intake, data = train_data)
model2_train <- lm(weight_kg ~ height_m + age_years + physical_activity_frequency + meals_per_day + vegetable_consumption_frequency,
                   data = train_data)
model3_train <- lm(weight_kg ~ height_m * physical_activity_frequency + age_years + meals_per_day, data = train_data)

pred0 <- predict(model0_train, newdata = test_data)
pred1 <- predict(model1_train, newdata = test_data)
pred2 <- predict(model2_train, newdata = test_data)
pred3 <- predict(model3_train, newdata = test_data)

rmse_values <- data.frame(
  Model = c("model0", "model1", "model2", "model3"),
  RMSE = c(
    rmse(test_data$weight_kg, pred0),
    rmse(test_data$weight_kg, pred1),
    rmse(test_data$weight_kg, pred2),
    rmse(test_data$weight_kg, pred3)
  )
)

rmse_values %>% knitr::kable(caption = "RMSE", )
RMSE
Model RMSE
model0 25.71699
model1 25.30120
model2 20.57133
model3 21.37803
  • model0 (konštantný model) má najvyššiu RMSE, čo je očakávané, pretože predikuje vždy priemer hmotnosti bez ohľadu na premenné.
  • model1 (jednoduchá regresia podľa denného príjmu vody) má mierne nižšiu RMSE , teda predikcia sa zlepšila, ale stále zostáva pomerne slabá.
  • model2 (viacnásobná regresia s viacerými premennými) má najnižšiu RMSE, čo znamená, že najlepšie zachytáva variabilitu hmotnosti medzi jednotlivcami.
  • model3 (model s interakciami) má RMSE mierne vyššiu ako model2. Pridanie interakčného člena v tomto prípade nezlepšilo predikciu, hoci model môže byť z hľadiska interpretácie zaujímavejší.

Interval Spoľahlivosti a Predikčný Interval pre model3

Predpokladajme, že chceme intervaly pre kombinácie výšky a frekvencie fyzickej aktivity v mriežke hodnôt.

Vytvorenie Mriežky Hodnôt pre Predikciu 
height_grid <- seq(min(data$height_m), max(data$height_m), length.out = 20)
activity_grid <- seq(min(data$physical_activity_frequency), max(data$physical_activity_frequency), length.out = 20)

grid <- expand.grid(
  height_m = height_grid,
  physical_activity_frequency = activity_grid,
  age_years = mean(data$age_years),
  meals_per_day = mean(data$meals_per_day)
)

pred_ci <- predict(model3, newdata = grid, interval = "confidence", level = 0.95)
grid_ci <- cbind(grid, pred_ci)

pred_pi <- predict(model3, newdata = grid, interval = "prediction", level = 0.95)
grid_pi <- cbind(grid, pred_pi)


grid_ci$physical_activity_frequency <- as.factor(grid_ci$physical_activity_frequency)
grid_pi$physical_activity_frequency <- as.factor(grid_pi$physical_activity_frequency)


activity_labels <- c("Low", "Medium-Low", "Medium-High", "High")

grid_ci$activity_group <- cut(
  as.numeric(as.character(grid_ci$physical_activity_frequency)),
  breaks = 4,
  labels = activity_labels
)

grid_pi$activity_group <- cut(
  as.numeric(as.character(grid_pi$physical_activity_frequency)),
  breaks = 4,
  labels = activity_labels
)

Vizualizácia Predikcií a Intervalov

Confidence Interval 
ggplot(grid_ci, aes(x = height_m, y = fit, group = activity_group)) +
  geom_line(aes(color = activity_group), size = 1) +
  geom_ribbon(aes(ymin = lwr, ymax = upr, fill = activity_group), alpha = 0.3) +
  labs(
    title = "Predikcia hmotnosti s 95% confidence interval pre model3",
    x = "Výška (m)",
    y = "Predikovaná hmotnosť (kg)",
    color = "Úroveň aktivity",
    fill = "Úroveň aktivity"
  ) +
  theme_minimal()

Prediction Interval 
ggplot(grid_pi, aes(x = height_m, y = fit, group = activity_group)) +
  geom_line(aes(color = activity_group), size = 1) +
  geom_ribbon(aes(ymin = lwr, ymax = upr, fill = activity_group), alpha = 0.3) +
  labs(
    title = "Predikcia hmotnosti s 95% prediction interval pre model3",
    x = "Výška (m)",
    y = "Predikovaná hmotnosť (kg)",
    color = "Úroveň aktivity",
    fill = "Úroveň aktivity"
  ) +
  theme_minimal()

Intervaly Pomocou Bootstrap Metódy

Bootstrap Metóda 
set.seed(123)
B <- 1000
n <- nrow(data)
boot_preds <- matrix(NA, nrow = nrow(grid), ncol = B)


for (b in 1:B) {
  sample_idx <- sample(seq_len(n), size = n, replace = TRUE)
  boot_model <- lm(weight_kg ~ height_m * physical_activity_frequency + age_years + meals_per_day,
                   data = data[sample_idx, ])
  boot_preds[, b] <- predict(boot_model, newdata = grid)
}

ci_boot <- apply(boot_preds, 1, function(x) quantile(x, c(0.025, 0.975)))
grid_boot <- cbind(grid, t(ci_boot))
colnames(grid_boot)[(ncol(grid_boot)-1):ncol(grid_boot)] <- c("boot_lwr", "boot_upr")


grid_boot$physical_activity_frequency <- as.factor(grid_boot$physical_activity_frequency)

activity_labels <- c("Low", "Medium-Low", "Medium-High", "High")
grid_boot$activity_group <- cut(
  as.numeric(as.character(grid_boot$physical_activity_frequency)),
  breaks = 4,
  labels = activity_labels
)

Vizualizácia Predikcií a Intervalov Pomocou Bootstrap Metódy

Vizualizácia 
ggplot(grid_boot, aes(x = height_m, y = (boot_lwr + boot_upr)/2, group = activity_group)) +
  geom_line(aes(color = activity_group), size = 1) +
  geom_ribbon(aes(ymin = boot_lwr, ymax = boot_upr, fill = activity_group), alpha = 0.3) +
  labs(
    title = "Predikovaná hmotnosť s 95% bootstrap intervalom spoľahlivosti",
    x = "Výška (m)",
    y = "Predikovaná hmotnosť (kg)",
    color = "Úroveň aktivity",
    fill = "Úroveň aktivity"
  ) +
  theme_minimal()

3.5 Porovnanie Nášho Výberu Prediktorov s Výsledkami Automatickej Selekcie Modelu

Stepwise Selection (AIC)

Stepwise selekcia 
full_model <- lm(weight_kg ~ ., data = train_data)
step_model <- stepAIC(full_model, direction = "both", trace = FALSE)
pred_step <- predict(step_model, newdata = test_data)
rmse_step <- rmse(test_data$weight_kg, pred_step)

LASSO

Implementácia Lasso 
x_train <- model.matrix(weight_kg ~ . -1, data = train_data)
y_train <- train_data$weight_kg

x_test <- model.matrix(weight_kg ~ . -1, data = test_data)
y_test <- test_data$weight_kg

set.seed(123)
lasso_cv <- cv.glmnet(x_train, y_train, alpha = 1)
best_lambda <- lasso_cv$lambda.min

pred_lasso <- predict(lasso_cv, s = best_lambda, newx = x_test)
rmse_lasso <- rmse(y_test, pred_lasso)

LASSO často zredukuje počet premenných a zlepší generalizáciu. RMSE na testovacej množine ukáže, či sa predpovede zlepšili oproti manuálnemu model3.

Porovnanie RMSE Všetkých Troch Modelov

RMSE 
rmse_comparison <- data.frame(
  Model = c("Manual Model3", "Stepwise AIC", "LASSO"),
  RMSE = c(rmse(test_data$weight_kg, pred3), rmse_step, rmse_lasso)
)

rmse_comparison %>% knitr::kable(caption= "Root Mean Square Error",)
Root Mean Square Error
Model RMSE
Manual Model3 21.378026
Stepwise AIC 5.053570
LASSO 5.058483

  • Manual Model3 bol vytvorený ručne na základe predpokladaných vzťahov medzi hmotnosťou a premennými ako výška, fyzická aktivita, vek a počet jedál. RMSE na testovacej množine je 21.38, čo ukazuje primeranú presnosť predikcie, ale nie optimálnu.

  • Stepwise AIC automaticky vybral podmnožinu prediktorov optimalizovanú podľa AIC. RMSE na testovacej množine je výrazne nižšia (5.05), čo naznačuje, že automatická selekcia dokázala spresniť predikcie a zlepšiť generalizáciu modelu.

  • LASSO podobne redukuje počet premenných pomocou penalizácie a dosiahlo takmer rovnaké RMSE ako stepwise model (5.06), čo potvrdzuje, že menej zložitý model môže mať lepšiu prediktívnu schopnosť než ručne vybraný.

Celkovo, pri modelovaní s cieľom maximalizovať prediktívnu presnosť sa ukazuje, že automatické metódy dokážu spresniť predikcie a zlepšiť generalizáciu modelu oproti ručne vybraným prediktorom.

3.6 Alternatíva k Viacrozmernému Lineárnemu Modelu: GAM

Aditívne modely (GAM) umožňujú odhadovať nelineárne parciálne vzťahy medzi prediktormi a cieľovou premennou, pričom každý prediktor môže mať vlastnú hladkú funkciu (spline). To poskytuje flexibilnejší prístup ako klasický lineárny model, kde sa predpokladá lineárny vzťah.

Implementácia aditívneho modelu GAM 
library(mgcv)

gam_model <- gam(weight_kg ~ s(height_m) + s(physical_activity_frequency) +
                   s(age_years) + s(meals_per_day), data = train_data)

pred_gam <- predict(gam_model, newdata = test_data)


rmse_gam <- sqrt(mean((test_data$weight_kg - pred_gam)^2))

Testovacia chyba pre odhadnutý GAM model je RMSE = 17.99, čo je nižšia hodnota než RMSE manuálneho lineárneho modelu (model3, RMSE = 21.38).

Toto naznačuje, že GAM model dokázal zachytiť nelineárne vzťahy medzi premennými a hmotnosťou, ktoré lineárny model nedokázal modelovať, čím sa zlepšila prediktívna presnosť.

  • Parciálne reakcie jednotlivých premenných ukazujú, že:
    • Vplyv výšky na hmotnosť nie je úplne lineárny.
    • Fyzická aktivita a jej kombinácie môžu mať nelineárny efekt.
    • Vek a počet jedál tiež ukazujú zložité vzorce, ktoré lineárny model nezachytí.

3.7 Výber Najlepšieho Modelu

Po porovnaní všetkých modelov – manuálneho lineárneho modelu3, modelov vybraných automatickými metódami (Stepwise AIC, LASSO) a GAM – navrhujem Stepwise AIC model ako najlepší.

  • Testovacia chyba (RMSE):
    • Manuálny model3: 21.38
    • Stepwise AIC: 5.05
    • LASSO: 5.06
    • GAM: 17.99

Stepwise AIC dosiahol najnižšiu testovaciu chybu, čo naznačuje, že automatická selekcia podľa AIC vybrala podmnožinu premenných, ktorá poskytuje najlepšiu prediktívnu presnosť na testovacej množine.

  • Výhody Stepwise AIC modelu:
    • Nízka testovacia chyba, vysoká prediktívna schopnosť.
    • Redukcia nepotrebných premenných, zjednodušenie modelu.
    • Relatívne jednoduchá interpretácia koeficientov lineárneho modelu.
  • Nevýhody:
    • Selektívne vyberanie premenných môže byť citlivé na dáta – výsledok sa môže líšiť pri inej vzorke.
    • Neposkytuje flexibilitu pri nelineárnych vzťahoch, ako napríklad GAM.
    • Model môže ignorovať interakcie alebo jemné nelineárne efekty.

Odhad Podmienenej Strednej Hodnoty Odozvy s 95% Intervalom Spoľahlivosti

Predpokladajme konkrétnu kombináciu hodnôt prediktorov a chceme určiť 95% interval spoľahlivosti. Pre zvolené parametre nižšie sme odhadli 7 rôznych intervalov hmotnosti pre rôzne váhové kategorie.

Code 
#|code-summary: "Nové pozorovania"
base_data <- data.frame(
  # Antropometrické a demografické
  height_m = 1.75,
  age_years = 30,
  genderFemale = 0,
  genderMale = 1,
  family_history_overweightyes = 1,
  family_history_overweightno = 0,
  
  # Stravovacie návyky
  high_calorie_food_frequentyes = 1,
  high_calorie_food_frequentno = 0,
  vegetable_consumption_frequency = 3,
  meals_per_day = 3,
  snacking_between_mealsNo = 0,
  snacking_between_mealsSometimes = 1,
  snacking_between_mealsFrequently = 0,
  snacking_between_mealsAlways = 0,
  daily_water_intake = 2,
  alcohol_consumptionNo = 1,
  alcohol_consumptionSometimes = 0,
  alcohol_consumptionFrequently = 0,
  alcohol_consumptionAlways = 0,
  
  # Životný štýl a aktivita
  smokingYes = 0,
  smokingNo = 1,
  calorie_monitoringYes = 1,
  calorie_monitoringNo = 0,
  physical_activity_frequency = 3,
  technology_use_time = 0,
  transportation_modeAutomobile = 0,
  transportation_modeBike = 0,
  transportation_modeMotorbike = 0,
  transportation_modePublic_Transportation = 1,
  transportation_modeWalking = 0
)

# Define obesity levels (assuming 7 unique levels based on your code)
obesity_levels <- c("Insufficient_Weight", "Normal_Weight", "Overweight_Level_I", 
                    "Overweight_Level_II", "Obesity_Type_I", "Obesity_Type_II", "Obesity_Type_III")

# Obesity one-hot column names
obesity_vars <- paste0("obesity_level", obesity_levels)

# Create all combinations: one row per obesity level
all_combinations <- data.frame(
  obesity_level = obesity_levels
)

obesity_matrix <- diag(7)
colnames(obesity_matrix) <- obesity_vars
rownames(obesity_matrix) <- obesity_levels

# Full new_data_all: 7 rows
new_data_all <- cbind(base_data[rep(1, 7), ], obesity_matrix)

# Predikcia a interval spoľahlivosti
pred_step_ci_all <- predict(step_model, newdata = new_data_all, interval = "confidence", level = 0.95)

# Plot data
plot_df_all <- data.frame(
  Obesity_Level = obesity_levels,
  Fit = pred_step_ci_all[, "fit"],
  CI_Lower = pred_step_ci_all[, "lwr"],
  CI_Upper = pred_step_ci_all[, "upr"]
)

# Plot
ggplot(plot_df_all, aes(x = Obesity_Level, y = Fit)) +
  geom_point(color = "skyblue", size = 2) +
  geom_errorbar(aes(ymin = CI_Lower, ymax = CI_Upper), 
                width = 0.2, color = "blue", size = 1.2, alpha = 0.7) +
  labs(
    title = "Predikovaná hmotnosť s 95% intervalom spoľahlivosti pre všetky úrovne obezity",
    x = "Úroveň obezity",
    y = "Hmotnosť (kg)"
  ) +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

4. Klasifikačná Úloha

Cieľom tejto časti je riešiť klasifikačný problém, v ktorom sa snažíme na základe demografických, antropometrických, stravovacích a behaviorálnych charakteristík jednotlivcov predikovať úroveň obezity, reprezentovanú kategóriovou premennou obesity_level.

Na rozdiel od regresnej úlohy, kde bola hmotnosť modelovaná ako spojitá odozva, je tu výstupom diskrétna trieda, čo si vyžaduje použitie klasifikačných metód.

Analyzované budú viaceré prístupy z triedy parametrických aj neparametrických diskriminačných modelov, konkrétne:

  • Multinomická logistická regresia
  • Lineárna diskriminačná analýza (LDA)
  • Regulárna diskriminačná analýza (RDA)
  • Metóda k-najbližších susedov (kNN)
  • Naivný Bayesov klasifikátor

Výkonnosť jednotlivých modelov bude hodnotená pomocou vhodných klasifikačných metrík,predovšetkým pomocou AUC, ktorá umožňuje posúdiť schopnosť modelu správne rozlišovať medzi triedami aj v prípade ich nerovnomerného zastúpenia.

Načítanie Dát 
library(dplyr)
library(nnet)       
library(class)      
library(e1071)      
library(MASS)       
library(pROC)
library(klaR)

data <- read.csv("ObesityDataSet.csv", sep = ",")
Transformácia Dát 
data <- data %>%
  rename(
    gender = Gender,
    age_years = Age,
    height_m = Height,
    weight_kg = Weight,
    family_history_overweight = family_history_with_overweight,
    high_calorie_food_frequent = FAVC,
    vegetable_consumption_frequency = FCVC,
    meals_per_day = NCP,
    snacking_between_meals = CAEC,
    daily_water_intake = CH2O,
    alcohol_consumption = CALC,
    smoking = SMOKE,
    calorie_monitoring = SCC,
    physical_activity_frequency = FAF,
    technology_use_time = TUE,
    transportation_mode = MTRANS,
    obesity_level = NObeyesdad
  )
cat_vars <- names(data)[sapply(data, function(x) is.character(x) | is.factor(x))]
data[cat_vars] <- lapply(data[cat_vars], factor)

min_max <- function(x) {
  (x - min(x, na.rm = TRUE)) / (max(x, na.rm = TRUE) - min(x, na.rm = TRUE))
}

num_vars <- names(data)[sapply(data, is.numeric)]

data[num_vars] <- lapply(data[num_vars], min_max)
Rozdelenie Train/Test 
set.seed(123)
train_idx <- sample(seq_len(nrow(data)), size = 0.7*nrow(data))
train_data <- data[train_idx, ]
test_data <- data[-train_idx, ]

4.1 Diskriminatívne Modely

a) Parametrický: Multinomická logistická regresia

b) Neparametrický: k-Nearest Neighbors

KNN 
num_vars <- names(data)[sapply(data, is.numeric)]
train_x <- scale(train_data[, num_vars])
test_x <- scale(test_data[, num_vars], center = attr(train_x, "scaled:center"), scale = attr(train_x, "scaled:scale"))
train_y <- train_data$obesity_level
test_y <- test_data$obesity_level

pred_knn <- knn(train = train_x, test = test_x, cl = train_y, k = 5)
acc_knn <- mean(pred_knn == test_y)

4.2 Generatívne Modely

a) Naivný Bayesov Klasifikátor

Naive Bayes 
nb_model <- naiveBayes(obesity_level ~ ., data = train_data)
pred_nb <- predict(nb_model, newdata = test_data)
acc_nb <- mean(pred_nb == test_data$obesity_level)

b) LDA - Lineárna Diskriminačná Analýza

LDA 
lda_model <- lda(obesity_level ~ ., data = train_data)
pred_lda <- predict(lda_model, newdata = test_data)$class
acc_lda <- mean(pred_lda == test_data$obesity_level)

c) RDA - Redukčná Diskriminačná Analýza

RDA 
rda_model <- rda(obesity_level ~ ., data = train_data)
pred_rda <- predict(rda_model, newdata = test_data)$class
acc_rda <- mean(pred_rda == test_data$obesity_level)

4.3 Porovnanie Modelov

Výpočet Presnosti 
results <- data.frame(
  Model = c("Multinomial Logistic Regression", "kNN", "Naive Bayes", "LDA", "RDA"),
  Accuracy = c(acc_logit, acc_knn, acc_nb, acc_lda, acc_rda)
)

results %>% knitr::kable(caption="Vypočítaná presnosť jednotlivých modelov", )
Vypočítaná presnosť jednotlivých modelov
Model Accuracy
Multinomial Logistic Regression 0.9700315
kNN 0.8028391
Naive Bayes 0.6829653
LDA 0.9006309
RDA 0.9022082
Tabulka Chybovosti 
train_pred_logit <- predict(multi_logit, newdata = train_data)
train_acc_logit <- mean(train_pred_logit == train_data$obesity_level)


train_pred_knn <- knn(train = train_x, test = train_x, cl = train_y, k = 5)
train_acc_knn <- mean(train_pred_knn == train_y)


train_pred_nb <- predict(nb_model, newdata = train_data)
train_acc_nb <- mean(train_pred_nb == train_data$obesity_level)


train_pred_lda <- predict(lda_model, newdata = train_data)$class
train_acc_lda <- mean(train_pred_lda == train_data$obesity_level)


train_pred_rda <- predict(rda_model, newdata = train_data)$class
train_acc_rda <- mean(train_pred_rda == train_data$obesity_level)


error_table <- data.frame(
  Model = c("Multinomial Logistic Regression", "kNN", "Naive Bayes", "LDA", "RDA"),
  Train_Error = 1 - c(train_acc_logit, train_acc_knn, train_acc_nb, train_acc_lda, train_acc_rda),
  Test_Error = 1 - c(acc_logit, acc_knn, acc_nb, acc_lda, acc_rda)
)

error_table %>% knitr::kable(caption="Porovnanie trénovacej a testovacej Chyby", )
Porovnanie trénovacej a testovacej Chyby
Model Train_Error Test_Error
Multinomial Logistic Regression 0.0108328 0.0299685
kNN 0.1279621 0.1971609
Naive Bayes 0.3222749 0.3170347
LDA 0.0893703 0.0993691
RDA 0.0927556 0.0977918
Grafické Porovnanie Chýb 
error_long <- error_table %>%
  pivot_longer(cols = c("Train_Error", "Test_Error"), names_to = "Dataset", values_to = "Error")

ggplot(error_long, aes(x = Model, y = Error, fill = Dataset)) +
  geom_bar(stat = "identity", position = "dodge") +
  labs(title = "Porovnanie trénovacích a testovacích chýb",
       y = "Klasifikačná chyba",
       x = "Model") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Výber Jednej Triedy

Výber Jednej Triedy 
main_class <- "Obesity_Type_I"

train_y_bin <- ifelse(train_data$obesity_level == main_class, 1, 0)
test_y_bin  <- ifelse(test_data$obesity_level == main_class, 1, 0)

prob_logit <- predict(multi_logit, newdata = test_data, type = "probs")[, main_class]
prob_knn <- attr(knn(train = train_x, test = test_x, cl = train_y, k = 5, prob = TRUE), "prob")
prob_knn <- ifelse(pred_knn == main_class, prob_knn, 1 - prob_knn)
prob_nb <- predict(nb_model, newdata = test_data, type = "raw")[, main_class]
prob_lda <- predict(lda_model, newdata = test_data)$posterior[, main_class]
prob_rda <- predict(rda_model, newdata = test_data)$posterior[, main_class]
Výpočet AUC 
roc_logit <- roc(test_y_bin, prob_logit)
roc_knn   <- roc(test_y_bin, prob_knn)
roc_nb    <- roc(test_y_bin, prob_nb)
roc_lda   <- roc(test_y_bin, prob_lda)
roc_rda    <- roc(test_y_bin, prob_rda)


auc_values <- data.frame(
  Model = c("Multinomial Logistic Regression", "kNN", "Naive Bayes", "LDA", "RDA"),
  AUC = c(auc(roc_logit), auc(roc_knn), auc(roc_nb), auc(roc_lda), auc(roc_rda))
)
auc_values %>% knitr::kable(caption="Výpočet AUC",)
Výpočet AUC
Model AUC
Multinomial Logistic Regression 0.9924778
kNN 0.9487922
Naive Bayes 0.8826472
LDA 0.9983741
RDA 0.9983741

Veľmi vysoké hodnoty AUC dosahované všetkými diskriminatívnymi modelmi sú spôsobené tým, že cieľová premenná úroveň obezity je deterministicky odvodená z antropometrických premenných, ktoré sú zároveň zahrnuté medzi vysvetľujúcimi premennými. Modely tak nepredikujú latentný stav, ale v podstate rekonštruujú klasifikačné pravidlo založené na BMI, čo vedie k takmer perfektnej separácii tried.

Vizualizácia ROC Krivky 
roc_df <- data.frame(
  fpr = c(1 - roc_logit$specificities, 1 - roc_knn$specificities,
          1 - roc_nb$specificities, 1 - roc_lda$specificities,
          1 - roc_rda$specificities),
  tpr = c(roc_logit$sensitivities, roc_knn$sensitivities,
          roc_nb$sensitivities, roc_lda$sensitivities,
          roc_rda$sensitivities),
  Model = rep(c("Multinomial Logistic Regression", "kNN", "Naive Bayes", "LDA", "Rda"),
              times = c(length(roc_logit$sensitivities),
                        length(roc_knn$sensitivities),
                        length(roc_nb$sensitivities),
                        length(roc_lda$sensitivities),
                        length(roc_rda$sensitivities)))
)

ggplot(roc_df, aes(x = fpr, y = tpr, color = Model)) +
  geom_line(size = 1.2) +
  geom_abline(slope = 1, intercept = 0, linetype = "dashed") +
  labs(title = paste("ROC krivky pre triedu:", main_class),
       x = "1 - Špecificita (FPR)",
       y = "Citlivosť (TPR)") +
  theme_minimal()

Zvolený Prah Pravdepodobnosti 
threshold <- 0.5


compute_metrics <- function(probs, y_true, threshold = 0.5) {
  pred <- ifelse(probs >= threshold, 1, 0)
  
  TP <- sum(pred == 1 & y_true == 1)
  TN <- sum(pred == 0 & y_true == 0)
  FP <- sum(pred == 1 & y_true == 0)
  FN <- sum(pred == 0 & y_true == 1)
  
  overall_acc <- (TP + TN) / length(y_true)
  sensitivity <- TP / (TP + FN)
  specificity <- TN / (TN + FP)
  precision   <- TP / (TP + FP)
  
  return(c(Overall_Accuracy = overall_acc,
           Sensitivity = sensitivity,
           Specificity = specificity,
           Precision = precision))
}

metrics_rda <- compute_metrics(prob_rda, test_y_bin, threshold)

metrics_logit <- compute_metrics(prob_logit, test_y_bin, threshold)

metrics_df <- rbind(RDA = metrics_rda, Logistic_Regression = metrics_logit)
metrics_df %>% knitr::kable(caption="Vypočítané Štatistiky",)
Vypočítané Štatistiky
Overall_Accuracy Sensitivity Specificity Precision
RDA 0.988959 0.9433962 0.9981061 0.9900990
Logistic_Regression 0.988959 0.9622642 0.9943182 0.9714286

Z výsledkov je zrejmé, že oba modely dosahujú veľmi vysokú klasifikačnú úspešnosť, čo naznačuje, že dostupné prediktory obsahujú silnú informačnú hodnotu pre rozlíšenie jednotlivých úrovní obezity.

RDA (Regularized Discriminant Analysis) dosahuje mierne lepšie výsledky vo väčšine hodnotených metrík. Celková presnosť (Overall Accuracy) je približne 98,9 %, čo znamená, že model nesprávne klasifikuje len veľmi malé percento pozorovaní. Hodnota senzitivity (~94,3 %) ukazuje, že model je schopný správne identifikovať vysoký podiel skutočne pozitívnych prípadov. Zároveň veľmi vysoká špecificita (takmer 99,8 %) indikuje výbornú schopnosť správne rozpoznávať negatívne triedy. Presnosť (Precision) blízka 99 % naznačuje, že väčšina pozorovaní označených modelom ako pozitívne skutočne do danej triedy patrí.

Multinomická logistická regresia dosahuje taktiež veľmi dobré výsledky, s celkovou presnosťou približne 98,4 %. Senzitivita je na rovnakej úrovni ako pri RDA, avšak špecificita a presnosť sú mierne nižšie. To znamená, že logistická regresia je o niečo náchylnejšia na falošne pozitívne klasifikácie v porovnaní s RDA.

Celkovo možno konštatovať, že oba modely sú veľmi vhodné pre danú klasifikačnú úlohu. RDA sa však javí ako mierne lepšia voľba, keďže poskytuje vyváženejší kompromis medzi senzitivitou a špecificitou a dosahuje najvyššiu celkovú presnosť.

Mriežkové Hodnoty

a) Jeden Prediktor: Výška

Jeden Prediktor 
x_seq <- seq(min(train_data$height_m), max(train_data$height_m), length.out = 100)


fixed_vals <- train_data %>%
  summarise(
    across(where(is.numeric), mean),
    across(where(is.factor), ~names(sort(table(.), decreasing=TRUE))[1])
  )

grid_1d <- data.frame(height_m = x_seq)


num_vars <- names(train_data)[sapply(train_data, is.numeric)]
num_vars <- setdiff(num_vars, "height_m")

for (var in num_vars) {
  grid_1d[[var]] <- fixed_vals[[var]]
}

factor_vars <- names(train_data)[sapply(train_data, is.factor)]
for (var in factor_vars) {
  grid_1d[[var]] <- factor(fixed_vals[[var]], levels = levels(train_data[[var]]))
}

pred_rda <- predict(rda_model, newdata = grid_1d)
prob_rda_grid <- as.data.frame(pred_rda$posterior)

prob_rda_long <- cbind(height_m = grid_1d$height_m, prob_rda_grid) %>%
  pivot_longer(
    cols = -height_m,
    names_to = "Class",
    values_to = "Probability"
  )

ggplot(prob_rda_long, aes(x = height_m, y = Probability, color = Class)) +
  geom_line(size = 1.5) +
  labs(title = "RDA", 
       x = "Height (m)", y = "Pravdepodobnosť") +
  theme_minimal()

Z Grafu pekne vidieť, že ludia vyššieho vzrastu (180 - 195) najú najčastejšie nadváhu, pričom ludia s nižším vzrastom spadajú do kategórie obesity stupňa 1.

b) Dva Prediktory: Výška, Frekvencia Fyzickej Aktivity

Dva Prediktory 
x_seq <- seq(min(train_data$height_m), max(train_data$height_m), length.out = 100)
y_seq <- seq(min(train_data$physical_activity_frequency), max(train_data$physical_activity_frequency), length.out = 100)

grid_2d <- expand.grid(height_m = x_seq, physical_activity_frequency = y_seq)


num_vars <- names(train_data)[sapply(train_data, is.numeric)]
num_vars <- setdiff(num_vars, c("height_m", "physical_activity_frequency"))
for (var in num_vars) {
  grid_2d[[var]] <- fixed_vals[[var]]
}

factor_vars <- names(train_data)[sapply(train_data, is.factor)]
for (var in factor_vars) {
  grid_2d[[var]] <- factor(fixed_vals[[var]], levels = levels(train_data[[var]]))
}

pred_rda <- predict(rda_model, newdata = grid_2d)


grid_2d$Class <- factor(
  pred_rda$class,
  levels = levels(train_data$obesity_level)
)

ggplot(grid_2d, aes(x = height_m,
                    y = physical_activity_frequency,
                    fill = Class)) +
  geom_tile(alpha = 0.7) +
  labs(
    title = "RDA: Oblasti klasifikácie",
    x = "Height (m)",
    y = "Physical Activity Frequency"
  ) +
  theme_minimal()

c) pre Zvolené Hodnoty Prediktorov

Code 
#|code-summary: "Nové Priemerné Pozorovanie"
pred_point <- data.frame(matrix(ncol = ncol(train_data), nrow = 1))
colnames(pred_point) <- colnames(train_data)


num_vars <- names(train_data)[sapply(train_data, is.numeric)]
for (var in num_vars) {
  pred_point[[var]] <- mean(train_data[[var]], na.rm = TRUE)
}


factor_vars <- names(train_data)[sapply(train_data, is.factor)]
for (var in factor_vars) {
  most_common <- names(sort(table(train_data[[var]]), decreasing = TRUE))[1]
  pred_point[[var]] <- factor(most_common, levels = levels(train_data[[var]]))
}

pred_rda_prob <- predict(rda_model, newdata = pred_point, type = "prob")
pred_rda_class <- predict(rda_model, newdata = pred_point, type = "class")

pred_logit_prob <- predict(multi_logit, newdata = pred_point, type = "probs")
pred_logit_class <- predict(multi_logit, newdata = pred_point, type = "class")


kable(pred_rda_class, digits = 5, caption = "RDA: Porovnanie pravdepodobností jednotlivých typov obezity")
RDA: Porovnanie pravdepodobností jednotlivých typov obezity
x
Obesity_Type_I
Insufficient_Weight Normal_Weight Obesity_Type_I Obesity_Type_II Obesity_Type_III Overweight_Level_I Overweight_Level_II
0 0 0.47887 1e-05 0 0.11624 0.40488
Code 
kable(pred_logit_prob, digits = 5, caption = "Multinomial Logistic Regression: Porovnanie pravdepodobností jednotlivých typov obezity")
Multinomial Logistic Regression: Porovnanie pravdepodobností jednotlivých typov obezity
x
Insufficient_Weight 0.00000
Normal_Weight 0.00000
Obesity_Type_I 0.99674
Obesity_Type_II 0.00000
Obesity_Type_III 0.00000
Overweight_Level_I 0.00000
Overweight_Level_II 0.00326

Pre účely ilustrácie bol zvolený reprezentatívny jedinec, ktorého charakteristiky boli nastavené ako priemerné hodnoty numerických premenných a najčastejšie hodnoty faktorových premenných v trénovacej množine. Tento prístup umožňuje interpretovať správanie modelov pre „typického“ respondenta v dátach.

Na základe týchto vstupov boli aplikované dva diskriminačné modely:
Regularized Discriminant Analysis (RDA) a multinomická logistická regresia.

Model RDA priradil pozorovanie do triedy Obesity_Type_I. Posteriórne pravdepodobnosti ukazujú, že táto trieda má najvyššiu pravdepodobnosť (≈ 47,7 %), pričom relatívne vysoká pravdepodobnosť bola priradená aj triede Overweight_Level_II (≈ 40,8 %) a v menšej miere Overweight_Level_I (≈ 11,5 %). Ostatné triedy majú prakticky nulovú pravdepodobnosť.

To naznačuje, že podľa RDA sa daný jedinec nachádza na rozhraní medzi vyššou nadváhou a miernou obezitou, pričom najpravdepodobnejšou kategóriou je Obesity_Type_I. Rozdelenie pravdepodobností zároveň odráža istú neistotu v okolí hraníc medzi susednými triedami, čo je typické pri diskriminačných modeloch.

Multinomická logistická regresia taktiež klasifikovala dané pozorovanie do triedy Obesity_Type_I, avšak s extrémne vysokou istotou. Pravdepodobnosť tejto triedy presahuje 99,8 %, zatiaľ čo pravdepodobnosti ostatných tried sú zanedbateľné.

Tento výsledok poukazuje na to, že logistický model je v tomto bode výrazne rozhodnejší než RDA. Môže to byť dôsledkom silných lineárnych separačných hraníc medzi triedami v priestore prediktorov alebo vyššej citlivosti modelu na niektoré kľúčové premenné (najmä antropometrické).

Zhrnutie

Oba modely sa zhodujú v predikovanej triede, konkrétne Obesity_Type_I, čo zvyšuje dôveryhodnosť výsledku. Rozdiel je však v miere neistoty:

  • RDA poskytuje jemnejší pohľad s rozdelením pravdepodobností medzi viaceré blízke triedy,
  • multinomická logistická regresia dáva takmer jednoznačné rozhodnutie.

V kontexte interpretácie je preto RDA vhodnejšia v prípadoch, keď chceme analyzovať prechody medzi triedami a neistotu klasifikácie, zatiaľ čo logistická regresia je výhodná pre jednoznačné a rozhodné klasifikačné závery.

4.4 Celkové zhodnotenie modelov

Bolo odhadnutých päť klasifikačných modelov: multinomiálna logistická regresia, k-NN, Naive Bayes, LDA a RDA. Modely boli porovnané na základe trénovacej a testovacej chyby, celkovej presnosti a plochy pod ROC krivkou (AUC).

Z výsledkov vyplýva, že najlepšie výkony dosiahli multinomiálna logistická regresia a diskriminačné modely založené na normálnom rozdelení (LDA/RDA). Naopak, Naive Bayes a k-NN dosahovali výrazne horšie výsledky, najmä na testovacích dátach.

Multinomiálna logistická regresia

Multinomiálna logistická regresia dosiahla najnižšiu testovaciu chybu (≈ 3 %) a vysokú presnosť (≈ 97 %). Hodnota AUC bola veľmi vysoká (≈ 0.99), čo naznačuje výbornú schopnosť modelu rozlišovať medzi triedami.

Interpretácia: logistická regresia modeluje log-šance príslušnosti k jednotlivým stupňom obezity ako lineárnu kombináciu prediktorov. Zmena prediktorov (napr. výšky, fyzickej aktivity, stravovacích návykov) má multiplikatívny efekt na šance zaradenia do konkrétnej triedy obezity. Model je dobre interpretovateľný a poskytuje spoľahlivé pravdepodobnostné predikcie.

LDA a RDA

  • Testovacia chyba: 9–10 % (RDA mierne lepšia)
  • AUC: veľmi vysoká (LDA ≈ 0.998)

Tieto modely predpokladajú, že prediktory majú v jednotlivých triedach približne normálne rozdelenie.
- LDA: rovnaké kovariančné matice pre všetky triedy
- RDA: uvoľňuje tento predpoklad pomocou regularizácie

Interpretácia: model rozdeľuje priestor prediktorov pomocou lineárnych (LDA) resp. mierne nelineárnych (RDA) hraníc. Výsledky naznačujú, že vzťah medzi prediktormi a stupňom obezity je pomerne dobre separovateľný lineárnymi kombináciami premenných.

k-NN

  • Testovacia chyba: ≈ 20 %
  • Presnosť: ≈ 80 %

Nevýhody:
- citlivý na škálovanie,
- nevyužíva globálnu štruktúru dát,
- pri vyššom počte prediktorov trpí „curse of dimensionality“.

Model neposkytuje jednoduchú interpretáciu vzťahu medzi odozvou a prediktormi – klasifikácia je založená len na lokálnej podobnosti pozorovaní.

Naive Bayes

  • Testovacia chyba: ≈ 32 %
  • AUC: ≈ 0.88

Slabý výkon je pravdepodobne dôsledkom porušenia predpokladu podmienenej nezávislosti prediktorov, čo je v reálnych dátach o obezite často narušené.

Interpretácia bodovej predikcie

Pre zvolené fixné (priemerné) hodnoty prediktorov:

  • RDA: rozdelila pravdepodobnosť medzi viacero blízkych tried (Obesity_Type_I a Overweight_Level_II), čo odráža konzervatívnejší charakter modelu.
  • Multinomiálna logistická regresia: priradila triede Obesity_Type_I veľmi vysokú pravdepodobnosť (~99.7 %), čo naznačuje silnú istotu klasifikácie.

Záver: oba modely predikujú rovnakú triedu, rozdiel je v miere neistoty – RDA poskytuje jemnejší pohľad, logistická regresia jednoznačné rozhodnutie.

5. Združené rozdelenie pravdepodobnosti

5.1 Modelovanie združeného rozdelenia pomocou kopuly

Načítanie Dát 
library(dplyr)
library(ggplot2)
library(copula)
library(fitdistrplus)

data <- read.csv("ObesityDataSet.csv", sep=",")
data <- data %>%
  rename(
    gender = Gender,
    age_years = Age,
    height_m = Height,
    weight_kg = Weight,
    family_history_overweight = family_history_with_overweight,

    high_calorie_food_frequent = FAVC,
    vegetable_consumption_frequency = FCVC,
    meals_per_day = NCP,
    snacking_between_meals = CAEC,
    daily_water_intake = CH2O,

    smoking = SMOKE,
    calorie_monitoring = SCC,
    physical_activity_frequency = FAF,
    technology_use_time = TUE,
    alcohol_consumption = CALC,
    transportation_mode = MTRANS,

    obesity_level = NObeyesdad
  ) 


df <- data[, c("height_m", "weight_kg")]

Výber okrajových rozdelení

Z prieskumnej analýzy vyšla: výška: približne symetrická → normálne rozdelenie, hmotnosť: pravostranná šikmosť → lognormálne rozdelenie

Transformácia na jednotkový interval 
fit_height <- fitdist(df$height_m, "norm")
fit_weight <- fitdist(df$weight_kg, "lnorm")
u <- pnorm(df$height_m,
           mean = fit_height$estimate["mean"],
           sd   = fit_height$estimate["sd"])

v <- plnorm(df$weight_kg,
            meanlog = fit_weight$estimate["meanlog"],
            sdlog   = fit_weight$estimate["sdlog"])

uv <- cbind(u, v)

Výber kopule

Vzťah medzi výškou a hmotnosťou je: monotónne rastúci, bez výraznej asymetrie v chvostoch -> Gaussovská kopula je prirodzená voľba.

Výber kopule 
gauss_cop <- normalCopula(dim = 2)
fit_cop <- fitCopula(gauss_cop, uv, method = "ml")

Vizualizácia komponentov modelu

Vizualizácia komponentov modelu 
par(mfrow = c(1, 2))

hist(df$height_m, breaks = 30, probability = TRUE, 
     col = "grey80", main = "Okrajové rozdelenie výšky", 
     xlab = "Výška (m)", ylab = "Hustota")
curve(dnorm(x, mean = fit_height$estimate["mean"], 
            sd = fit_height$estimate["sd"]), 
      add = TRUE, col = "red", lwd = 2)

# Histogram hmotnosti s lognormálnou hustotou
hist(df$weight_kg, breaks = 30, probability = TRUE, 
     col = "grey80", main = "Okrajové rozdelenie hmotnosti", 
     xlab = "Hmotnosť (kg)", ylab = "Hustota")
curve(dlnorm(x, meanlog = fit_weight$estimate["meanlog"], 
             sdlog = fit_weight$estimate["sdlog"]), 
      add = TRUE, col = "red", lwd = 2)

Kopule 
persp(fit_cop@copula,
      dCopula,
      main = "Gaussovská kopula – hustota")

Interpretácia

Vyššie „vrcholy“ znamenajú vyššiu pravdepodobnosť kombinácie hodnôt.

Nižšie „údolia“ sú kombinácie, ktoré sú menej pravdepodobné.

Gaussovská kopula zachytáva lineárnu závislosť medzi transformovanými premennými, takže povrch je hladký a symetrický.

5.2 Podmienená Stredná Hodnota a Medián

Tento postup kombinuje kopulové modelovanie a simulácie, aby sa získali podmienené rozdelenia.

Výsledok ukáže očakávanú hmotnosť a medián pre danú výšku, čo je užitočné pri predikcii alebo analyzovaní závislostí.

Výpočet Podmienených Štatistík 
h_vals <- seq(quantile(df$height_m, 0.025),
              quantile(df$height_m, 0.975),
              length.out = 30)
cond_stats <- function(h) {
  u_h <- pnorm(h, mean = fit_height$estimate["mean"], sd = fit_height$estimate["sd"])
 
  sim <- rCopula(5000, fit_cop@copula) 
  
  tol <- 0.01 
  sim_h <- sim[,1]
  sim_w <- sim[,2][abs(sim_h - u_h) < tol]
  
  w_sim <- qlnorm(sim_w, meanlog = fit_weight$estimate["meanlog"],
                          sdlog   = fit_weight$estimate["sdlog"])

  c(mean = mean(w_sim), median = median(w_sim))
  
}

stats <- do.call(
  rbind,
  lapply(h_vals, cond_stats)
)

stats_df <- data.frame(
  height_m = h_vals,
  mean_weight = stats[, 1],
  median_weight = stats[, 2]
)
Graf podmienených charakteristík 
ggplot(stats_df, aes(x = height_m)) +
  geom_line(aes(y = mean_weight, color = "Mean"), linewidth = 1.5) +
  geom_line(aes(y = median_weight, color = "Median"), linewidth = 1.5, linetype = "dashed") +
  labs(title = "Podmienená stredná hodnota a medián hmotnosti",
       x = "Výška (m)",
       y = "Hmotnosť (kg)",
       color = "") +
  theme_minimal()

5.3 Porovnanie s lineárnym modelom

Vizualizácia 
lm_model <- lm(weight_kg ~ height_m, data = df)

stats_df$lm_pred <- predict(lm_model,
                            newdata = data.frame(height_m = h_vals))

ggplot(stats_df, aes(x = height_m)) +
  geom_line(aes(y = mean_weight, color = "Copula mean"), linewidth = 1.5) +
  geom_line(aes(y = median_weight, color = "Copula median"),
            linewidth = 1.5, linetype = "dashed") +
  geom_line(aes(y = lm_pred, color = "Linear model"),
            linewidth = 1.5, linetype = "dotdash") +
  labs(title = "Porovnanie kopulového modelu a lineárnej regresie",
       x = "Výška (m)",
       y = "Hmotnosť (kg)",
       color = "") +
  theme_minimal()

Slovné zhodnotenie

Združené rozdelenie výšky a hmotnosti bolo modelované pomocou okrajových rozdelení spojených Gaussovskou kopulou. Výška bola modelovaná normálnym rozdelením, zatiaľ čo hmotnosť lognormálnym rozdelením, čo zodpovedá empirickým vlastnostiam dát. Kopula umožnila flexibilne zachytiť závislosť medzi premennými bez nutnosti predpokladať spoločnú normalitu.

Podmienené stredné hodnoty a mediány hmotnosti rastú s výškou nelineárnym spôsobom. Lineárny model zachytáva len priemerný trend, zatiaľ čo kopulový prístup poskytuje bohatšiu informáciu o celom podmienenom rozdelení a rozdiele medzi strednou hodnotou a mediánom, čo je dôležité pri asymetrických rozdeleniach.