Naloga
Julija Gorenc
2024-02-13
Raziskovalno vprašanje: Kako se telesna masa in dolžina kljuna pingvinov medsebojno povezujejo in kakšen vpliv imata ti dve značilnosti na razlikovanje med posameznimi pingvini?
Izbrala sem si to vprašanje, saj je primerno za MGK analizo, ker želim raziskati, kako se različne merjene značilnosti medsebojno povezujejo. MGK mi bo omogočila zmanjšanje dimenzionalnosti podatkov in identifikacijo glavnih komponent, ki pojasnujejo največji del variabilnosti v podatkih. Potrebno je poudariti, da je za glavne komponente značilno nerastoče zaporedje in da pri tem prva glavna komponenta vsebuje največji del variabilnosti začetnih spremenljivk, vsaka nadaljna glavna komponenta pa vsebuje manj variabilnosti.
Podatki in analiza
Podatke sem pridobila na spletni strani Kaggle, izbrala se si nabor podatkov o pingvinih, saj imam rada živali ter lahko s temi podatki naredim nekaj zanimivega. Iskala sem številske podatke, saj je MGK analiza za njih ustrezna.
Vir podatkov: https://www.kaggle.com/datasets/youssefaboelwafa/clustering-penguins-species/data
Sedaj bom uvozila podatke v r studio za opis le teh. Ločeni so z vejico, imajo glavo, ter decimalno piko. S funkcijo read.table preberemo podatke, s funkcijo head pa prikažemo prvih 6 vrstic:
podatki <- read.table("./penguins.csv", header=TRUE, sep=",", dec=".")
head(podatki)## culmen_length_mm culmen_depth_mm flipper_length_mm body_mass_g sex
## 1 39.1 18.7 181 3750 MALE
## 2 39.5 17.4 186 3800 FEMALE
## 3 40.3 18.0 195 3250 FEMALE
## 4 NA NA NA NA <NA>
## 5 36.7 19.3 193 3450 FEMALE
## 6 39.3 20.6 190 3650 MALEOpis podatkov:
culmen_length_mm: dolžina kljuna (mm)
culmen_depth_mm: globina kljuna (mm)
flipper_length_mm: dolžina plavuti (mm)
body_mass_g: telesna masa (g)
sex: spol pingvina (MALE, FEMALE)
Opazimo lahko, da nekaj podatkov manjka, zato jih je potrebno ustrezno prečistiti.
Z na.omit funkcijo, bomo izbrisali prazne vrednosti iz nabora podatkov:
podatki <- na.omit(podatki)
head(podatki)## culmen_length_mm culmen_depth_mm flipper_length_mm body_mass_g sex
## 1 39.1 18.7 181 3750 MALE
## 2 39.5 17.4 186 3800 FEMALE
## 3 40.3 18.0 195 3250 FEMALE
## 5 36.7 19.3 193 3450 FEMALE
## 6 39.3 20.6 190 3650 MALE
## 7 38.9 17.8 181 3625 FEMALEVidimo lahko, da je npr. ID pingvina 4 bil izbrisan, saj je imel prazne vrednosti.
Spremenimo še spol pingvinov v numerične vrednosti. V končni fazi bi lahko kot smo delali na predavanjih tu uporabili funkcijo factor, vendar sem zasledila tudi za funkcijo recode:
library(dplyr)##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionpodatki$sex <- recode(podatki$sex, "MALE"=1, "FEMALE"=0)## Warning: Unreplaced values treated as NA as `.x` is not compatible.
## Please specify replacements exhaustively or supply `.default`.Napaka pravi, da niso vsi podatki v numerični obliki in jih je pretvorilo kot NA, zato jih odstranimo, kot smo to storili prej:
podatki <- na.omit(podatki)Lahko nadaljujemo z opisno statistiko:
library(psych)
describe(podatki[ , -1])## vars n mean sd median trimmed mad min max
## culmen_depth_mm 1 334 17.17 1.97 17.3 17.20 2.22 13.1 21.5
## flipper_length_mm 2 334 214.35 263.65 197.0 200.39 16.31 -132.0 5000.0
## body_mass_g 3 334 4207.19 804.01 4050.0 4159.79 889.56 2700.0 6300.0
## sex 4 334 0.51 0.50 1.0 0.51 0.00 0.0 1.0
## range skew kurtosis se
## culmen_depth_mm 8.4 -0.15 -0.91 0.11
## flipper_length_mm 5132.0 17.90 322.03 14.43
## body_mass_g 3600.0 0.47 -0.75 43.99
## sex 1.0 -0.02 -2.01 0.03Naredimo LM model. VIF je namenjen številskim spremenljivkam, katere sedaj imamo:
fit <- lm(culmen_depth_mm ~ flipper_length_mm + body_mass_g + culmen_length_mm,
data = podatki)
library(car)## Loading required package: carData##
## Attaching package: 'car'## The following object is masked from 'package:psych':
##
## logit## The following object is masked from 'package:dplyr':
##
## recodevif(fit)## flipper_length_mm body_mass_g culmen_length_mm
## 1.002424 1.534983 1.532086Visoke vrednosti VIF (tipično vrednosti nad 5 ali 10) kažejo na visoko stopnjo multikolinearnosti, kar lahko poslabša zanesljivost statističnih sklepov modela.
V naše primeru kaže VIF na nizko stopnjo multikolinearnosti, kar je idealno stanje, saj zagotavlja, da so koeficienti regresije zanesljivi in natančno odražajo vpliv vsake spremenljivke na odvisno spremenljivko. Želimo si, da je povprečen VIF čim bližje 1.
Preverimo še povprečen VIF:
mean(vif(fit))## [1] 1.356498Tudi tukaj lahko vidimo, da je vrednost VIF-a zadovoljiva. Multikolinearnost vpliva na standardne napake, regresijski koeficienti so napačni. Imamo 2 načina reševanja:
Vržemo ven spremenljivke, ampak s tem izgobimo njene informacije
Združimo v Glavne Komponente.
Uporabimo GK, ki jo delamo kot vmesno metodo:
podatki_MGK <- podatki[, c("flipper_length_mm", "body_mass_g", "culmen_length_mm")]
library(pastecs)##
## Attaching package: 'pastecs'## The following objects are masked from 'package:dplyr':
##
## first, lastround(stat.desc(podatki_MGK, basic=FALSE), 2)## flipper_length_mm body_mass_g culmen_length_mm
## median 197.00 4050.00 44.45
## mean 214.35 4207.19 43.99
## SE.mean 14.43 43.99 0.30
## CI.mean.0.95 28.38 86.54 0.59
## var 69510.50 646430.94 29.83
## std.dev 263.65 804.01 5.46
## coef.var 1.23 0.19 0.12Pogledamo koeficiente in vidimo, da ima Flipper_length_mm največji koeficient variacije. Želimo, da so vse spremenljivke enako zastopane, ko bomo razlagali rezultate.
Naredimo korelacijsko matriko, ki omogoča razumevanje linearne povezanosti med spremenljivkami:
R <- cor(podatki_MGK)
round(R, 3)## flipper_length_mm body_mass_g culmen_length_mm
## flipper_length_mm 1.000 0.048 0.021
## body_mass_g 0.048 1.000 0.589
## culmen_length_mm 0.021 0.589 1.000Visoke vrednosti korelacij med spremenljivkami nakazujejo, da PCA lahko učinkovito zmanjša dimenzionalnost podatkov, saj obstajajo skupni vzorci variabilnosti, ki jih PCA lahko izkoristi.
Tukaj vidimo Pearsonove korelacijske koeficiente (izven diagonale). 0.589 interpretiramo kot: linearna povezava med culmen_length_mm in body_mass_g je pozitivna in srednje močna.
Naredimo Bartletov preizkus sferičnosti, preverjamo ali je populacijska korelacijska matrika enaka enotski matrik:
library(psych)
cortest.bartlett(R, n = nrow(podatki_MGK))## $chisq
## [1] 142.0607
##
## $p.value
## [1] 1.358657e-30
##
## $df
## [1] 3Postavimo H0 in H1. V H0 predpostavljamo, da je populacijska korelacijska matrika P enaka enotski matriki I (H0: P = I). To pomeni, da spremenjlivke nimajo nič skupnega. V H1 predpostavljamo, da populacijska korelacijska matrika ni enaka enotski matriki (H1: P ni enak I). To pa pomeni, da imajo spremenljivke vsaj nekaj skupnega.
Ker je p vrednost izjemno nizka zavrnemo H0 pri p < 0,001 in sprejmemo H1, kar pomeni, da so podatki primerni za PCA analizo.
Preverimo še determinanto korelacijske matrike, ki odraža skupno variabilnost med spremenljivkami:
det(R)## [1] 0.6511792Nizka vrednost determinante (bližje 0) nakazuje, da spremenljivke delijo nekaj skupne variabilnosti in so potencialno povezane, kar je ugodno za PCA.
Nato se lotimo KMO, ki je splošna mera vzorčne ustreznosti in ocenjuje, kako primerni so podatki za MGK:
library(psych)
KMO(R)## Kaiser-Meyer-Olkin factor adequacy
## Call: KMO(r = R)
## Overall MSA = 0.5
## MSA for each item =
## flipper_length_mm body_mass_g culmen_length_mm
## 0.58 0.50 0.50Višje vrednosti (bližje 1) kažejo, da je vzorec primeren. Vrednost nad 0.5 se šteje za sprejemljivo, kar pomeni, da so podatki dovolj dobri za izvedbo PCA. KMO je splošna mera vzorčne ustreznosti, medtem ko MSA ovrednoti posamezne spremenljivke posebej.
MSA statistika ovrednoti vsako posamezno spremenljivko. flipper_length_mm iam največ skupnega z ostalimi.
Naredimo sedaj PCA - Princible component analysis ali MGK:
library(FactoMineR)
mgk <- PCA(podatki_MGK,
scale.unit = TRUE,
graph = FALSE)
library(factoextra)## Loading required package: ggplot2##
## Attaching package: 'ggplot2'## The following objects are masked from 'package:psych':
##
## %+%, alpha## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBaget_eigenvalue(mgk)## eigenvalue variance.percent cumulative.variance.percent
## Dim.1 1.5933549 53.11183 53.11183
## Dim.2 0.9965431 33.21810 86.32993
## Dim.3 0.4101020 13.67007 100.00000S scale.unit = TRUE smo opravili standardizacijo, ki je potrebna.
Glavni cilj MGK je poenostavitev podatkov, torej lahko rečemo, da gre za metodo poenostavitve kompleknosti podatkov. Informacije pri MGK se meri z varianco, torej bolj kot spremenljivka variira, več informacije nosi v sebi. Vsaka glavna komponenta je linearna kombinacija vseh p-tih začetnih spremenljivk. Prva glavna komponenta (GK1) je narejena na način, da ujame čim večji del začetne variabilnosti, potem pa apliciramo npr. še drugo glavno komponento, na katero pa ne prenašamo informacij, ki so že aplicirane na GK1 in iz tega razloga mora biti druga glavna komponenta obvezno PRAVOKOTNA (ortogonalna) na prvo, da lahko poberemo še preostanek informacij. Vprašanje se mogoče pojavlja tudi glede standardizacije. Torej, zakaj spremenljivke standardiziramo? Zato, da zagotovimo, da imajo vse spremenljivke isto količino informacij, da so med seboj primerljive. Se pravi, če začetne spremenljivke ne standardiziramo, bo spremenljivka, ki bo imela največjo varianco imela tudi največ informacij in bo zato najbolj prispevala k ustvarjanju glavne komponente. Pomembno je tudi, da med spremenljivkami obstaja KORELACIJA. Z KMO statistiko preverimo splošno mero vzorčne ustreznosti; z njo dejansko ovrednotimo, ali so podatki kot celota primerni za združevanje v GK. MSA statistika pa služi za ovrednotenje posamezne spremenljivke. Enostavno povedano, več skupne variabilnosti imajo spremenljivke, večji bosta KMO in MSA statistika.
Pomembno je tudi, da vemo, da varianco vsake posamezne komponente merimo z LASTNO VREDNOSTJO. Se pravi, če je lastna vrednost neke glavne komponente enaka 3, to pomeni, da je varianca te glavne komponente enaka 3. Če so podatki primerni, lahko preverimo tudi z Bartlettovim preizkusom sferičnosti. Tukaj dejansko preverjamo, ali je populacijska korelacijska matrika (P) enaka enotski matriki (I). Želimo si zavračati H0, ker potem vemo, da imajo spremenljivke vsaj nekaj skupnega. Torej, da ugotovimo, ali so podatki primerni moramo preveriti korelacije (morajo biti dovolj visoke), preveriti moramo tudi KMO in MSA statistiko in preveriti moramo Bartlettov preizkus sferičnosti.
Če pogledamo lastne vrednosti iz izpisa, vidimo, da gre za padajoče zaporedje. Pojavi se nam maksimalno št. komponent (Dim1, Dim2, Dim3). Kaiserjevo pravilo pravi, da moramo vzeti toliko glavnih komponent kot ima lastna vrednost višjo od 1. Prva glavna komponenta oz. GK1 (Dim.1) pojasnjuje 53.11% celotne variabilnosti, kar pomeni, da velik del informacij v podatkih o pingvinih lahko dobimo iz kombinacije telesne mase pingvinov in njihove dolžine kljuna. Če delimo lastno vrednost Dim1 (1.59, tukaj sem zaokrožila) s številom vseh komponent, ki so prikazane (3), torej 1.59/3 x 100 dobimo zaokroženih 53%, ki ga lahko vidimo pod variance.percent v izpisu. “variance.percent” v % izrazi koliko informacije je padlo na posamezno komponento. Vidimo torej, da je na prvo komponento padlo 53.11% informacij, na drugo komponento je padlo 33.22% informacij in na tretjo, ki nosi najmanj informacij, je padlo 13.67% informacije. Na podlagi Kaiserjevega pravila bi torej vzeli le eno komponento (dim1).
Naredimo diagram lastnih vrednosti, ki se uporabi tudi za določanje števila glavnih komponent:
library(factoextra)
fviz_eig(mgk,
choice = "eigenvalue",
main = "Diagram lastnih vrednosti",
ylab = "Lastna vrednost",
xlab = "Glavna komponenta",
addlabels = TRUE)
Diagram lastnih vrednosti se torej uporablja za določanje števila glavnih komponent. Iščemo prelom oz. lahko si pomagamo z vprašanjem “kje se naklon pomembno spremeni”. Splošno priporočilo pravi, da moramo vzeti eno glavno komponento manj kot tam, kjer je prelom. Problem tukaj lahko nastane zaradi subjektivnosti glede preloma. Glede na sliko, lahko rečemo, da tukaj preloma ni.
Naredimo še paralelno analizo, ki je najbolj zanesljiva in najbolj uporabna:
library(psych)
fa.parallel(podatki_MGK,
sim = FALSE,
fa = "pc")
## Parallel analysis suggests that the number of factors = NA and the number of components = 1Glavno pravilo tukaj je, Lastna vrednostEMP > Lastna vrednostTEOR. Teoretične so prikazane s črtkano rdečo črto, empirične oz. dejanske lastne vrednosti pa so prikazane z modro barvo. Vidimo, da je pri prvi glavni komponenti lastna vrednost (empirična) višja od lastne vrednosti (teoretične), zato vzamemo 1 glavno komponento. Tudi v izpisanem tekstu spodaj je nakazano, naj vzamemo 1 glavno komponento.
Rezultati MGK:
library(FactoMineR)
mgk <- PCA(podatki_MGK,
scale.unit = TRUE,
graph = FALSE,
ncp = 1)
mgk## **Results for the Principal Component Analysis (PCA)**
## The analysis was performed on 334 individuals, described by 3 variables
## *The results are available in the following objects:
##
## name description
## 1 "$eig" "eigenvalues"
## 2 "$var" "results for the variables"
## 3 "$var$coord" "coord. for the variables"
## 4 "$var$cor" "correlations variables - dimensions"
## 5 "$var$cos2" "cos2 for the variables"
## 6 "$var$contrib" "contributions of the variables"
## 7 "$ind" "results for the individuals"
## 8 "$ind$coord" "coord. for the individuals"
## 9 "$ind$cos2" "cos2 for the individuals"
## 10 "$ind$contrib" "contributions of the individuals"
## 11 "$call" "summary statistics"
## 12 "$call$centre" "mean of the variables"
## 13 "$call$ecart.type" "standard error of the variables"
## 14 "$call$row.w" "weights for the individuals"
## 15 "$call$col.w" "weights for the variables"Zanima nas “\(var\)cor”, “\(var\)contrib”, “\(ind\)coord” in “\(ind\)cos2”.
print(mgk$var$cor)## Dim.1
## flipper_length_mm -0.1043677
## body_mass_g -0.8907073
## culmen_length_mm -0.8883146Te uteži predstavljajo Perasonove korelacijske koeficiente med merjenimi spremenljivkami in glavnimi komponentami.
Preverimo kakšten % celotne informacije GK je prispevala posamezna merjena spremenljivka:
print(mgk$var$contrib)## Dim.1
## flipper_length_mm 0.6836283
## body_mass_g 49.7917620
## culmen_length_mm 49.5246097Prikažemo koordinate, v katerih so shranjene pozicije posameznih enot:
head(mgk$ind$coord)## Dim.1
## 1 1.0429495
## 2 0.9458125
## 3 1.3231762
## 5 1.6125802
## 6 1.1022094
## 7 1.1786288Ugotovitve:
Lahko naredimo interpretacijo glavne komponente iz MGK metode o pingvinih, ki temelji na utežeh (koeficientih),ki smo jih dobili, ti pa so pomembni za določanje, kako močno vsaka spremenljivka prispeva k tej komponenti. Uteži, ki smo jih dobili iz izpisa, so bile:
flipper_length_mm: - 0.104
body_mass_g: - 0.891
culmen_length_mm: - 0.888
Na podlagi tega lahko vidimo, da ima najmanj skupnega spremenljivka flipper_length oz. lahko bi rekli, da je najbolj drugačna od ostalih. Pogledamo lahko tudi, kako prispevajo predstavljene posamezne spremenljivke k variabilnosti, ki jo pojasnjuje glavna komponenta (iz izpisa): flipper_length_mm: 0.68%body_mass_g: 49.79%culmen_length_mm: 49.52%
Ključne ugotovitve: Na podlagi analiz, ki sem jih naredila v programu R, ugotavljam in opažam, da sta h glavni komponenti pomembneje prispevali spremenljivki body_mass_g (telesna masa pingvina) in culmen_length_mm (dolžina kljuna). Obe spremenljivki imata pomemben prispevek k variabilnosti, ki jo pojasnjuje glavna komponenta zato lahko rečem, da sta spremenljivki tesno povezani. Manjši prispevek h glavni komponenti doprinaša flipper_length_mm (dolžina plavuti) kar lahko vidimo tudi iz izpisa, zato lahko sklepam, da ta spremenljivka nima tako velika pomena pri pojasnjevanju variabilnosti, če jo primerjamo s telesno maso (body_mass_g) in z dolžino kljuna (culmen_length_mm)
Interpretacija rezultatov:
Kot sem omenila sem na podlagi podatkov ugotovila, da telesna masa pingvinov in dolžina kljuna pingvinov predstavljata dve pomembne oz. lahko rečemo tudi ključni spremenljivki, ki pomembno prispevata h glavni komponenti. Ti dve spremenljivki predstavljata pomembne lastnosti na podlagi katerih se pingvini razlikujejo.
Zaključek
V dani domači nalogi sem se osredotočila na raziskovalno vprašanje “kako se telesna masa in dolžina kljuna pingvinov medsebojno povezujeta in kakšen vpliv imata ti dve značilnosti na razlikovanje med posameznimi pingvini”. S pomočjo tehnične metode MGK sem kompleksnost podatkov poenostavila in ugotovila, katere so ključne spremenljivke, ki pojasnjujejo največji del variabilnosti med pingvini.
Ugotovila sem, da sta telesna masa in dolžina kljuna dve zelo pomembni spremenljivki, ki skupaj lahko pojasneta velik delež variabilnosti v razlikovanju pingvinov, torej gre za neke fizične lastnosti pingvinov. Na podlagi izpisa sem videla, da prva glavna komponenta pojasnjuje več kot 53% celotne variabilnosti, kar pomeni, da obstaja močna povezanost med telesno težo pingvinov in njihovo dolžino kljunov.