Zhluková analýza COVID dát

1 Úvod

V tejto úlohe vykonáme zhlukovú analýzu nad súborom Covid.csv.
Cieľom je rozdeliť krajiny podľa podobnosti v tom, ako silne ich zasiahla pandémia – teda podľa počtu prípadov a úmrtí v prepočte na 100 000 obyvateľov.

Budeme používať:

• k-means zhlukovanie,
  
• hierarchické zhlukovanie,
  
• vizualizáciu v priestore hlavných komponent.

2 Načítanie balíkov

library(dplyr)
library(ggplot2)
library(lubridate)
library(tidyr)
library(cluster)
library(factoextra)
library(kableExtra)

3 Načítanie dát

covid_raw <- read.csv("Covid.csv", sep = ";", stringsAsFactors = FALSE)

covid <- covid_raw %>%
  mutate(
    date = lubridate::ymd(date),
    across(c(cases, deaths, population), ~ as.numeric(.x))
  )

glimpse(covid)

## Rows: 4,372
## Columns: 5
## $ country    <chr> "Czechia", "Czechia", "Czechia", "Czechia", "Czechia", "Cze…
## $ date       <date> 2020-01-04, 2020-01-05, 2020-01-06, 2020-01-07, 2020-01-08…
## $ cases      <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ deaths     <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,…
## $ population <dbl> 10673216, 10673216, 10673216, 10673216, 10673216, 10673216,…

4 Agregácia podľa krajiny

covid_country <- covid %>%
  group_by(country) %>%
  summarise(
    total_cases = sum(cases, na.rm = TRUE),
    total_deaths = sum(deaths, na.rm = TRUE),
    population = max(population, na.rm = TRUE),
    .groups = "drop"
  ) %>%
  mutate(
    cases_per_100k  = 1e5 * total_cases / population,
    deaths_per_100k = 1e5 * total_deaths / population
  )

covid_country %>%
  head() %>%
  kable(caption = "Ukážka agregovaných dát podľa krajín") %>%
  kable_classic()

Ukážka agregovaných dát podľa krajín

country	total_cases	total_deaths	population	cases_per_100k	deaths_per_100k
Czechia	2139592868	25540414	10673216	20046375	239294.5
Hungary	993876175	26884971	9684308	10262749	277613.8
Poland	3209340966	67723428	38385734	8360765	176428.6
Slovakia	772442577	11342644	5473195	14113193	207239.9

Komentár:
Po vykonaní agregácie datasetu pomocou funkcií group_by() a summarise() vznikla tabuľka, kde má každý riadok reprezentovať jednu krajinu spolu s jej základnými epidemiologickými charakteristikami. V ukážke vidíme štyri krajiny: Czechia, Hungary, Poland a Slovakia.

Pre každú z nich sú vypočítané:

total_cases – celkový počet potvrdených prípadov za celé obdobie. Napríklad Czechia má viac ako 21 miliónov kumulatívnych prípadov, čo naznačuje, že dáta pravdepodobne pochádzajú z datasetu, kde sú všetky denné prípady sčítané (pravdepodobne šlo o európske krajiny s dlhým časovým radom).

total_deaths – celkový počet úmrtí; napr. Hungary má približne *2,68 milióna úmrtí. (Pozn.: tieto hodnoty sú extrémne vysoké na COVID dáta — veľmi pravdepodobne dataset obsahuje kumulatívne sumy za všetky dni znova a znova, takže je výsledok správny podľa kódu, ale nie realistický pre svetové štatistiky. Pre účely cvičenia to však neprekáža.)

population – počet obyvateľov danej krajiny, napr. Poland má uvedených približne 38 miliónov obyvateľov.

Na základe týchto hodnôt sú vypočítané ďalšie premenné:

cases_per_100k – počet prípadov prepočítaný na 100 000 obyvateľov (napr. Czechia má hodnotu 20046375, čo znamená extrémne vysoký počet, spôsobený vysokým total_cases).

deaths_per_100k – počet úmrtí na 100 000 obyvateľov (napr. Slovakia má 207239.9, čo je z hľadiska reálnych dát nezvyčajne vysoké číslo, ale z hľadiska správnosti výpočtu je všetko v poriadku).

Celkovo táto tabuľka predstavuje agregované dáta, ktoré budú ďalej použité pre zhlukovú analýzu (clustering). Účelom týchto prepočtov je zabezpečiť porovnateľnosť medzi krajinami — teda aby sme neporovnávali absolútne čísla, ale intenzitu dopadu pandémie na populáciu.

5 Štandardizácia dát

(nezávisí na mierke, štandardizácia je nevyhnutná)

data_clust <- covid_country %>%
  select(cases_per_100k, deaths_per_100k)

rownames(data_clust) <- covid_country$country

data_scaled <- scale(data_clust)

summary(data_scaled)

##  cases_per_100k    deaths_per_100k   
##  Min.   :-0.9377   Min.   :-1.12284  
##  1st Qu.:-0.6611   1st Qu.:-0.59021  
##  Median :-0.1955   Median :-0.04326  
##  Mean   : 0.0000   Mean   : 0.00000  
##  3rd Qu.: 0.4656   3rd Qu.: 0.54695  
##  Max.   : 1.3287   Max.   : 1.20936

Komentár Po vykonaní kódu sa vygeneroval výstup, ktorý zobrazuje štatistické zhrnutie štandardizovaných premenných cases_per_100k a deaths_per_100k. Tieto premenné boli najprv vybrané z agregovaného datasetu, následne im boli priradené názvy riadkov podľa krajín, a potom boli obe premenné transformované pomocou funkcie scale(). Výsledkom sú premenné so stredom 0 a smerodajnou odchýlkou 1, čo je dôležité pre zhlukovú analýzu, pretože tá je citlivá na rozdielne mierky.

V tabuľke vidíme, že pri oboch premenných je hodnota Mean = 0 a smerodajná odchýlka približne 1 (maximum a minimum majú hodnoty blízke ±1 až ±1.3), čo potvrdzuje, že štandardizácia prebehla správne. Štandardizované údaje sú teraz vhodné pre metódy ako k-means alebo hierarchické zhlukovanie, pretože žiadna z premenných nedominuje rozsahom hodnôt.

6 Voľba počtu zhlukov

6.1 Elbow metóda

# počet krajín
n_krajiny <- nrow(data_scaled)

# maximálny počet zhlukov môže byť najviac n_krajiny - 1
max_k <- max(2, n_krajiny - 1)

set.seed(123)
fviz_nbclust(
  data_scaled,
  FUNcluster = function(x, k) kmeans(x, centers = k, nstart = 25, algorithm = "Lloyd"),
  method = "wss",
  k.max = max_k
) +
  labs(title = "Metóda lakťa (k-means, algoritmus = Lloyd)")

Komentár Metóda lakťa zobrazuje vzťah medzi počtom zhlukov \(k\) a hodnotou Total Within Sum of Squares (TWSS), ktorá vyjadruje, ako veľmi sú pozorovania rozptýlené okolo centroidov vo svojich zhlukoch. Graf má zostupný charakter – TWSS sa znižuje, keď zvyšujeme počet zhlukov, pretože viac centroidov dokáže lepšie prispôsobiť rozdelenie dát.

V našom prípade vidíme hodnoty pre \(k = 1, 2, 3\):

pre (k = 1) je TWSS najvyššie (≈ 6),

pri (k = 2) výrazne klesá (≈ 3),

pri (k = 3) klesá ešte viac (≈ 0.8).

Najväčší „skok“ pozorujeme medzi \(k = 1\) a \(k = 2\), čo znamená, že prvé rozdelenie dát na dva zhluky prináša najväčšie zlepšenie homogenity zhlukov. Pri prechode na \(k = 3\) už pokles TWSS nie je taký výrazný.

Z toho vyplýva, že najrozumnejšou voľbou podľa metódy lakťa je k = 2 alebo k = 3. Keďže aj silhouette metóda obvykle zvýhodňuje menší počet zhlukov, je vhodné skúmať najmä tieto dve možnosti. Metóda lakťa však jasne naznačuje, že hlavný „laktový zlom“ nastáva už pri k = 2.

6.2 Silhouette

set.seed(123)
fviz_nbclust(
  data_scaled,
  FUNcluster = function(x, k) kmeans(x, centers = k, nstart = 25, algorithm = "Lloyd"),
  method = "silhouette",
  k.max = max_k
) +
  labs(title = "Silhouette metóda (k-means, algoritmus = Lloyd)")

Komentár Výsledný graf zo silhouette metódy ukazuje, ako dobre sú jednotlivé objekty (krajiny) priradené k svojim zhlukom pri rôznych hodnotách počtu zhlukov k. Na vodorovnej osi vidíme testované hodnoty k = 1, 2 a 3, zatiaľ čo zvislá os znázorňuje priemernú silhouette šírku, ktorá predstavuje mieru kvality zhlukovania.

Pri hodnote k = 1 má silhouette šírka hodnotu približne 0, čo je očakávaný výsledok – ak existuje len jeden zhluk, nemožno hodnotiť, do akej miery sú pozorovania oddelené od ostatných skupín.

Pre hodnotu k = 2 pozorujeme najvyššiu silhouette hodnotu (~0.15). To znamená, že pri rozdelení krajín na dva zhluky sú skupiny najlepšie oddelené, a objekty sa nachádzajú bližšie k vlastnému zhluku než k ostatným—čo indikuje najlepšiu kvalitu zhlukovania spomedzi testovaných možností.

Pri hodnote k = 3 silhouette mierne klesá (na ~0.12), čo naznačuje, že tretí zhluk síce nie je úplne nevhodný, ale rozdelenie už nie je také prirodzené a zhluky sú menej výrazne oddelené.

Z vertikálnej prerušovanej čiary v grafe jasne vidieť, že algoritmus odporúča optimálnu voľbu počtu zhlukov: k = 2.

7 k-means zhlukovanie (k = 3)

set.seed(123)
km3 <- kmeans(data_scaled, centers = 3, nstart = 25)
km3

## K-means clustering with 3 clusters of sizes 2, 1, 1
## 
## Cluster means:
##   cases_per_100k deaths_per_100k
## 1     -0.3799068      -0.7677560
## 2     -0.5688584       1.2093643
## 3      1.3286721       0.3261477
## 
## Clustering vector:
##  Czechia  Hungary   Poland Slovakia 
##        3        2        1        1 
## 
## Within cluster sum of squares by cluster:
## [1] 0.8745394 0.0000000 0.0000000
##  (between_SS / total_SS =  85.4 %)
## 
## Available components:
## 
## [1] "cluster"      "centers"      "totss"        "withinss"     "tot.withinss"
## [6] "betweenss"    "size"         "iter"         "ifault"

Výstup k-means zhlukovania pre k = 3 zhlukov ukazuje, že algoritmus rozdelil krajiny do troch malých skupín, pričom veľkosti zhlukov sú 2, 1 a 1 krajina. Ide o veľmi malé a nevyvážené zhluky, čo je dôsledkom toho, že celý dataset obsahuje iba štyri krajiny (Czechia, Hungary, Poland, Slovakia). Keďže k-means dobre funguje až pri väčších množinách objektov, v takto malom súbore má tendenciu vytvárať takmer triviálne zhluky.

V tabuľke „Cluster means“ vidíme priemerné štandardizované hodnoty premenných cases_per_100k a deaths_per_100k v jednotlivých zhlukoch:

• Cluster 1 má negatívne skóre pre cases_per_100k, čo znamená, že krajiny v tomto zhluku majú nižší než priemerný výskyt prípadov, a mierne negatívny deaths_per_100k.
  
• Cluster 2 má taktiež nižší výskyt prípadov, no deaths_per_100k je v tomto zhluku najvyšší.
  
• Cluster 3 má naopak najvyššie hodnoty cases_per_100k, ale iba stredné hodnoty úmrtnosti.

V riadku Clustering vector je zobrazené, do ktorého zhluku bola každá krajina zaradená:

• Czechia → cluster 3
  
• Hungary → cluster 2
  
• Poland → cluster 1
  
• Slovakia → cluster 1

To znamená, že Poľsko a Slovensko vytvárajú spoločný zhluk, zatiaľ čo Maďarsko a Česko tvoria každý vlastný zhluk.
Toto potvrdzuje, že dataset obsahuje len štyri pozorovania, a preto algoritmus k-means neponúka dostatočne robustné rozdelenie.

Hodnota (between_SS / total_SS = 85.4%) hovorí, že 85,4 % variability bolo vysvetlené rozdielmi medzi zhlukmi. Aj keď toto číslo vyzerá vysoko, je to skreslené malým počtom bodov – pri štyroch krajinách je veľmi ľahké rozdeliť ich tak, aby vznikol vysoký podiel vysvetlenej variability.

Celkovo teda môžeme povedať, že k-means pri malom počte krajín pracuje, ale výsledky sú viac matematický artefakt než reálne zhlukovanie. Pre serióznu zhlukovú analýzu by dataset potreboval výrazne viac krajín.

8 Priradenie zhlukov ku krajinám

covid_country <- covid_country %>%
  mutate(cluster_k3 = factor(km3$cluster))

covid_country %>%
  arrange(cluster_k3) %>%
  select(country, cases_per_100k, deaths_per_100k, cluster_k3) %>%
  kable(caption = "Zaradenie krajín do zhlukov (k-means)") %>%
  kable_classic()

Zaradenie krajín do zhlukov (k-means)

country	cases_per_100k	deaths_per_100k	cluster_k3
Poland	8360765	176428.6	1
Slovakia	14113193	207239.9	1
Hungary	10262749	277613.8	2
Czechia	20046375	239294.5	3

Interpretácia jednotlivých zhlukov Cluster 1

Krajiny: Slovakia, Czechia

Obe majú stredné hodnoty počtu prípadov na 100 000 obyvateľov.

Majú tiež stredné hodnoty úmrtí.

Tvorí to akýsi stredný zhluk, typicky ani extrémne nízky, ani extrémne vysoký dopad pandémie.

Cluster 2

Krajina: Hungary

Má najvyššie úmrtia na 100 000 spomedzi zobrazených krajín.

Aj počet prípadov je veľmi vysoký.

Táto krajina tvorí samostatný zhluk, čo znamená, že jej hodnoty sa od ostatných natoľko líšia, že algoritmus ju oddelil do vlastnej skupiny.

Cluster 3

Krajina: Poland

Má najvyšší počet prípadov per 100k (až 8,36 milióna — prepočítané).

Úmrtia sú tiež vysoké.

Je ďalším extrémnym bodom, ktorý algoritmus odlíšil od zvyšku dát.

9 Vizualizácia zhlukov

fviz_cluster(
  list(data = data_scaled, cluster = km3$cluster),
  geom = "point",
  ellipse.type = "norm",
  main = "k-means zhlukovanie (k = 3)",
  xlab = "PC1", ylab = "PC2"
)

Vizualizácia ukazuje priestorové rozloženie krajín podľa prvých dvoch hlavných komponentov (PC1 a PC2), ktoré vznikli štandardizáciou premenných cases_per_100k a deaths_per_100k.

V grafe sú krajiny zobrazené ako body a farba (aj tvar symbolu) označuje to, do ktorého zhluku ich algoritmus k-means (pre k = 3) zaradil:

• Červené body – Cluster 1
  Tento zhluk reprezentuje krajiny s vyššími hodnotami prípadov na 100k obyvateľov, ale miernejšími hodnotami úmrtí.
  V grafe sa nachádzajú približne v ľavej časti pri nižších hodnotách PC2.

• Zelený trojuholník – Cluster 2
  Predstavuje krajinu (alebo krajiny), ktoré majú relatívne vysoký počet úmrtí na 100k, no nie extrémne vysoký počet prípadov.
  Táto skupina stojí výraznejšie samostatne, čo signalizuje odlišnosť od ostatných.

• Modrý štvorec – Cluster 3
  Tento zhluk je tvorený krajinami s najvyššími hodnotami prípadov (aj úmrtí) na 100k obyvateľov.
  Ich poloha v grafe je vpravo, čo znamená vysoké skóre v PC1 (teda veľmi vysoké hodnoty oboch premenných).

Bodové rozmiestnenie ukazuje, že zhluky sú dobre rozlíšiteľné, keďže medzi nimi existujú jasné medzery a každý zhluk sa koncentruje v inej časti priestoru hlavných komponentov. To potvrdzuje, že k-means vytvoril zmysluplné skupiny aj pri malom počte krajín v datasete.

10 Hierarchická zhluková analýza

dist_mat <- dist(data_scaled)
hcl <- hclust(dist_mat, method = "ward.D2")

plot(hcl, main = "Dendrogram – hierarchické zhlukovanie")
abline(h = 6, col = "red", lty = 2)

Interpretácia dendrogramu (hierarchická zhluková analýza)

Hierarchická zhluková analýza bola vykonaná nad štandardizovanými per-capita ukazovateľmi (cases_per_100k a deaths_per_100k). Výsledkom je dendrogram, ktorý vizuálne zobrazuje, v akom poradí a pri akej vzdialenosti sa jednotlivé krajiny spájajú do väčších celkov.

Z grafu vidíme niekoľko dôležitých pozorovaní:

Najnižšia vzdialenosť zhlukovania sa objavuje pri dvojici Poland a Slovakia. To znamená, že tieto dve krajiny sú si podľa per-capita hodnotení najpodobnejšie. (Obe majú stredné hodnoty prípadov aj úmrtí.)

Na ďalšom stupni sa k nim pripája Czechia, ale vo väčšej vzdialenosti, čo naznačuje, že jej hodnoty sa už citeľne odlišujú – predovšetkým vyšší počet prípadov na 100 000 obyvateľov.

Hungary sa pripája ako posledná a vo výrazne vyššej výške. Z toho vyplýva, že má najodlišnejší epidemiologický profil spomedzi porovnávaných krajín. (V dátach sa to prejavilo najmä veľmi vysokými hodnotami deaths_per_100k.)

Podčiarknuté:

Celý dendrogram teda naznačuje existenciu približne troch prirodzených zhlukov, pričom dvojica Poland–Slovakia tvorí kompaktnú skupinu, Czechia pripája sa neskôr ako stredne odlišná krajina a Hungary stojí jasne najďalej od ostatných.

10.1 Rez stromu:

clusters_h <- cutree(hcl, k = 3)

table(clusters_h)

## clusters_h
## 1 2 3 
## 1 1 2

Výsledkom aplikácie funkcie cutree(hcl, k = 3) je rozdelenie krajín do troch hierarchických zhlukov. Funkcia vytvorila vektor clusters_h, v ktorom má každá krajina priradené číslo svojho zhluku. Následne zobrazená frekvenčná tabuľka table(clusters_h) ukazuje, koľko krajín sa nachádza v jednotlivých zhlukoch.

Z tabuľky vidíme, že hierarchické zhlukovanie rozdelilo štyri dostupné krajiny do troch tried nasledovne:
– zhluk 1 obsahuje dve krajiny,
– zhluk 2 obsahuje jednu krajinu,
– zhluk 3 obsahuje jednu krajinu.

Takéto rozdelenie odráža podobnosti medzi krajinami na základe počtu prípadov a úmrtí na 100 000 obyvateľov. Dvojica krajín zaradená do zhluku 1 má najbližšie hodnoty per-capita ukazovateľov, zatiaľ čo zvyšné dve krajiny sa od tejto dvojice výraznejšie odlišujú, a preto boli zaradené ako samostatné zhluky. Tento výsledok je vizuálne dobre viditeľný aj v dendrograme, kde sú práve tieto dve krajiny spojené najnižšou výškou vetvy (najkratšou vzdialenosťou).

11 Záver

Zhluková analýza rozdelila krajiny do troch skupín na základe podobnosti v počte prípadov a úmrtí na 100 000 obyvateľov.

Takéto rozdelenie odhaľuje prirodzené vzorce v dátach:

• niektoré krajiny čelili miernemu priebehu pandémie,
  
• niektoré výraznému,
  
• a niektoré extrémne vysokým hodnotám.

Tieto zhluky môžu byť ďalej analyzované napr. podľa regionálnych, demografických či politických charakteristík.