fat <- read.table('DATA/mlrm-fat.txt', header=T)
plot(fat)

boxplot(fat)

R <- cor(fat)
R

n <- nrow(fat)
k <- 7 #pocet promennych
( test.stat <- -n*log(det(R))*(1- (2*k+11)/(6*n)) )
## [1] 407.9905
( kvantil <- qchisq(0.95, df=k*(k-1)/2) )
## [1] 32.67057

fat.PCA <- prcomp(fat, center=T, scale.=T)
fat.PCA

(vl.cisla <- eigen(R)$values)

(vl.cisla.pca <- fat.PCA$sdev^2)

pc <- fat.PCA$rotation

summary(fat.PCA)

plot(fat.PCA, type='l')
abline(h=1, lty=2)

fat.in.pc <- fat.PCA$x #pozorovani v souradnicich hlavnich komponentach
cor(fat,fat.in.pc[,1:2])

biplot(fat.PCA, scale=0)

biplot(fat.PCA, scale=0, ylabs=c('vaha', 'vyska', 'BMI', 'zebro', 'bricho', 'bok', 'stehno'),
       xlab='1. hlavni komponenta', ylab='2. hlavni komponenta' )

(R.reproduced <- pc[,1:2] %*% diag(vl.cisla.pca[1:2]) %*% t(pc[,1:2]))

(R.residual <- R - R.reproduced)