Metoda glavnih komponent
podatki <- read.table("/cloud/project/Poglavje 1/Zgled/Olimpijada.csv", header=TRUE, sep=";", dec=",") #Uvozimo podatke
str(podatki) #Prikažemo strukturo spremenljivk## 'data.frame': 200 obs. of 8 variables:
## $ ID : int 1 2 3 4 5 6 7 8 9 10 ...
## $ Ekonomija: num 92.9 103.1 88 92 83.9 ...
## $ Čust : int 96 106 92 102 95 111 94 113 103 83 ...
## $ IQ : int 102 112 91 126 120 84 91 134 102 112 ...
## $ Verb : int 392 593 364 680 616 522 450 740 558 519 ...
## $ Matem : int 487 538 458 550 505 467 457 569 515 467 ...
## $ Starost : int 17 18 18 18 17 17 19 17 17 18 ...
## $ Spol : int 1 1 1 1 1 1 1 1 1 1 ...Opis spremenljivk:
- Ekonomija: Dosežen rezultat na olimpijadi iz ekonomije.
- Čust: Kvocient čustvene inteligence.
- IQ: Inteligenčni kvocient.
- Verb: Verbalne spretnosti.
- Matem: Matematične spretnosti.
- Starost: starost v dopolnjenih letih.
- Spol: 0: Ž, 1: M
library(psych) #Aktiviramo knjižnico
describe(podatki) #Prikažemo opisno statistiko## vars n mean sd median trimmed mad min max
## ID 1 200 100.50 57.88 100.5 100.50 74.13 1.0 200
## Ekonomija 2 200 79.60 17.42 79.5 79.37 16.61 29.5 119
## Čust 3 200 107.52 12.18 108.0 107.74 11.12 75.0 139
## IQ 4 200 100.02 14.56 100.0 99.91 14.08 57.0 134
## Verb 5 200 533.88 105.42 525.0 531.36 109.71 314.0 801
## Matem 6 200 493.22 42.99 487.0 490.62 43.74 364.0 601
## Starost 7 200 17.34 0.65 17.0 17.32 0.00 16.0 19
## Spol 8 200 0.50 0.50 0.0 0.49 0.00 0.0 1
## range skew kurtosis se
## ID 199.0 0.00 -1.22 4.09
## Ekonomija 89.5 0.03 -0.06 1.23
## Čust 64.0 -0.13 -0.22 0.86
## IQ 77.0 0.01 -0.22 1.03
## Verb 487.0 0.22 -0.51 7.45
## Matem 237.0 0.40 -0.25 3.04
## Starost 3.0 0.41 0.14 0.05
## Spol 1.0 0.02 -2.01 0.04LM Model
podatki$SpolF <- factor(podatki$Spol, #Kategorialno spremenljivko spremenimo v faktor
levels = c(0, 1),
labels = c("Ž", "M"))
#Ocenimo linearno regresijsko funkcijo po metodi najmanjših kvadratov
fit <- lm(Ekonomija ~ Čust + IQ + Verb + Matem + Starost + SpolF,
data = podatki)
vif(fit) #Preverimo prisotnost previsoke multikolinearnosti## Čust IQ Verb Matem Starost SpolF
## 2.437805 3.292648 3.323291 2.052659 1.087658 1.444070mean(vif(fit))## [1] 2.273022MGK
podatki_MGK <- podatki[, c("Čust", "IQ", "Verb", "Matem")] #Ustvarimo novo tabelo s podatki, v katero vključimo samo izbrane številske spremenljivke
library(pastecs) #Aktiviramo knjižnico
round(stat.desc(podatki_MGK, basic=FALSE), 2) #Prikažemo opisno statistiko## Čust IQ Verb Matem
## median 108.00 100.00 525.00 487.00
## mean 107.52 100.02 533.88 493.22
## SE.mean 0.86 1.03 7.45 3.04
## CI.mean.0.95 1.70 2.03 14.70 5.99
## var 148.37 211.97 11113.27 1847.97
## std.dev 12.18 14.56 105.42 42.99
## coef.var 0.11 0.15 0.20 0.09R <- cor(podatki_MGK) #Ocena korelacijske matrike
round(R, 3) #Prikaz korelacijske matrike## Čust IQ Verb Matem
## Čust 1.000 0.039 0.550 0.047
## IQ 0.039 1.000 0.639 0.699
## Verb 0.550 0.639 1.000 0.488
## Matem 0.047 0.699 0.488 1.000library(psych) #Aktiviramo knjižnico
cortest.bartlett(R, n=200) #Izvedba Bartlettovega testa sferičnosti## $chisq
## [1] 360.9283
##
## $p.value
## [1] 6.949886e-75
##
## $df
## [1] 6det(R) #Ocena determinante korelacijske matrike## [1] 0.1598252library(psych) #Aktiviramo knjižnico
KMO(R) #KMO statistika## Kaiser-Meyer-Olkin factor adequacy
## Call: KMO(r = R)
## Overall MSA = 0.54
## MSA for each item =
## Čust IQ Verb Matem
## 0.33 0.54 0.53 0.73library(ppcor) #Aktiviramo knjižnico
parc_kor <- pcor(podatki_MGK) #Ocena parcialnih korelacij
round(parc_kor$estimate, 3) #Prikažemo rezultate, zaokrožene na 3 decimalna mesta## Čust IQ Verb Matem
## Čust 1.000 -0.401 0.684 -0.032
## IQ -0.401 1.000 0.594 0.515
## Verb 0.684 0.594 1.000 0.076
## Matem -0.032 0.515 0.076 1.000library(FactoMineR) #Aktiviramo knjižnico
mgk <- PCA(podatki_MGK, #Izvedemo MGK
scale.unit = TRUE, #Standardizacija
graph = FALSE) #Ne prikažemo grafikona lastnih vrednosti
library(factoextra) #Aktiviramo knjižnico
get_eigenvalue(mgk) #Prikažemo lastne vrednosti## eigenvalue variance.percent cumulative.variance.percent
## Dim.1 2.3227044 58.067611 58.06761
## Dim.2 1.1551963 28.879908 86.94752
## Dim.3 0.3536850 8.842125 95.78964
## Dim.4 0.1684142 4.210356 100.00000library(factoextra) #Aktiviramo knjižnico
fviz_eig(mgk, #Narišemo diagram lastnih vrednosti
choice = "eigenvalue", #Na y-osi prikažemo lastne vrednosti
main = "Diagram lastnih vrednosti", #Poimenovanje grafikona
ylab = "Lastna vrednost", #Poimenovanje y-osi
xlab = "Glavna komponenta", #Poimenovanje x-osi
addlabels = TRUE) #Nad stolpci pripišemo lastne vrednosti
library(psych) #Aktiviramo knjižnico
fa.parallel(podatki_MGK, #Izvedemo paralelno analizo
sim = FALSE, #Izključimo simulacijo
fa = "pc") #Izvajamo MGK (pc = principal component)
## Parallel analysis suggests that the number of factors = NA and the number of components = 2library(FactoMineR) #Aktiviramo knjižnico
mgk <- PCA(podatki_MGK, #Izvedemo metodo glavnih komponent
scale.unit = TRUE, #Standardizacija
graph = FALSE, #Izključimo grafičen prikaz, saj ga izvedemo kasneje
ncp=2) #Določimo število obdržanih glavnih komponent
mgk## **Results for the Principal Component Analysis (PCA)**
## The analysis was performed on 200 individuals, described by 4 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"print(mgk$var$cor) #Izpišemo uteži glavnih komponent## Dim.1 Dim.2
## Čust 0.4202514 0.8724243
## IQ 0.8584660 -0.3634767
## Verb 0.8824508 0.3067535
## Matem 0.7939836 -0.4097062print(mgk$var$contrib) #Struktura prispevka merjenih spremenljivk k varianci glavnih komponent## Dim.1 Dim.2
## Čust 7.603691 65.886989
## IQ 31.728696 11.436613
## Verb 33.526411 8.145604
## Matem 27.141203 14.530794library(factoextra) #Aktiviramo knjižnico
fviz_pca_var(mgk, #Grafičen prikaz uteži glavnih komponent
repel = TRUE) #Lažje berljive oznake
library(factoextra) #Aktiviramo knjižnico
fviz_pca_biplot(mgk) #Grafičen prikaz položaja posamezne enote glede na vrednosti prve in druge GK
podatki_MGK_std <- scale(podatki_MGK) #Standardizacija tabele s podatki
podatki_MGK_std[103 , ] #Izpis 103. vrstice tabele s podatki## Čust IQ Verb Matem
## -1.2737270 2.1968702 0.2666964 1.8791274head(mgk$ind$coord) #Prikaz prvih šestih vrstic objekta## Dim.1 Dim.2
## 1 -1.0411731 -1.1452885
## 2 1.3000218 -0.6179274
## 3 -2.0651056 -0.9745598
## 4 2.3770863 -1.0817015
## 5 1.0863805 -1.1832933
## 6 -0.9260634 0.8065728mgk$ind$coord[103, ] #Prikaz 103. vrstice objekta## Dim.1 Dim.2
## 2.024694 -2.423092podatki$GK1 <- mgk$ind$coord[ , 1] #Shranimo vrednosti prve glavne komponente v tabelo s podatki
podatki$GK2 <- mgk$ind$coord[ , 2] #Shranimo vrednosti druge glavne komponente v tabelo s podatki
head(podatki, 3) #Prikaz prvih šestih vrstic tabele s podatki## ID Ekonomija Čust IQ Verb Matem Starost Spol SpolF GK1
## 1 1 92.9 96 102 392 487 17 1 M -1.041173
## 2 2 103.1 106 112 593 538 18 1 M 1.300022
## 3 3 88.0 92 91 364 458 18 1 M -2.065106
## GK2
## 1 -1.1452885
## 2 -0.6179274
## 3 -0.9745598cor(x = podatki$GK1, y = podatki$GK2) #Ocena Pearsonovega korelacijskega koeficienta## [1] 6.066742e-16fit <- lm(Ekonomija ~ GK1 + GK2 + Starost + Spol, #OLS regresija
data = podatki)
vif(fit) #VIF statistika## GK1 GK2 Starost Spol
## 1.066876 1.197743 1.073155 1.282265mean(vif(fit))## [1] 1.15501podatki$StdOstanki <- rstandard(fit) #Izračun stand. ostankov
podatki$StdOcenVred <- scale(fit$fitted.values) #Izračun stand. ocenjenih vrednosti
library(ggplot2) #Aktiviramo knjižnico
ggplot(podatki, aes(y=StdOstanki, x=StdOcenVred)) + #Določimo osnovne parametre
theme_linedraw() + #Določimo temo grafikona
geom_point() + #Določimo način prikaza podatkov
ylab("Standardizirani ostanki") + #Poimenujemo y-os
xlab("Standardizirane ocenjene vrednost") #Poimenujemo x-os
hist(podatki$StdOstanki,
xlab = "Standardizirani ostanki",
ylab = "Frekvenca",
main = "Histogram standardiziranih ostankov")
head(podatki[order(podatki$StdOstanki), ], 3) #Prikaz prvih 6 vrstic tabele s podatki, urejene po standardiziranih ostankih## ID Ekonomija Čust IQ Verb Matem Starost Spol SpolF GK1
## 87 87 62.6 85 114 503 473 18 1 M -0.3842766
## 194 194 29.5 126 107 753 446 17 0 Ž 1.3232520
## 83 83 60.6 93 108 518 448 17 1 M -0.6565740
## GK2 StdOstanki StdOcenVred
## 87 -1.7338595 -3.116408 1.5224539
## 194 2.0867234 -2.617711 -1.5128379
## 83 -0.7967451 -2.287479 0.6183985podatki$CooksD <- cooks.distance(fit) #Izračun Cookovih razdalj
head(podatki[order(-podatki$CooksD), c("ID", "CooksD")], ) #Prikaz prvih 6 vrstic tabele s podatki, urejene po Cookovih razdaljah## ID CooksD
## 87 87 0.04338290
## 194 194 0.03908689
## 75 75 0.02845463
## 165 165 0.02732435
## 163 163 0.02660211
## 109 109 0.02601403podatki <- podatki[c(-87, -194), ] #Odstranitev dveh enotfit <- lm(Ekonomija ~ GK1 + GK2 + Starost + SpolF, #Ocenimo OLS regresijo
data = podatki)
summary(fit) #Prikaz rezultatov##
## Call:
## lm(formula = Ekonomija ~ GK1 + GK2 + Starost + SpolF, data = podatki)
##
## Residuals:
## Min 1Q Median 3Q Max
## -27.1973 -8.8603 0.4681 7.9101 24.2314
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 47.8862 22.2126 2.156 0.03234 *
## GK1 3.3282 0.5482 6.071 6.6e-09 ***
## GK2 -9.8301 0.8322 -11.812 < 2e-16 ***
## Starost 1.6918 1.2884 1.313 0.19071
## SpolFM 5.5028 1.8290 3.009 0.00297 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 11.42 on 193 degrees of freedom
## Multiple R-squared: 0.5629, Adjusted R-squared: 0.5539
## F-statistic: 62.15 on 4 and 193 DF, p-value: < 2.2e-16