##########################################################################################################
###########################################       CA       ###############################################
##########################################################################################################

# ordinacni metoda pro vyhodnoceni vztahu radku a sloupcu v tabulce
# cil: Redukovat vicerozmerny prostor do maleho poctu novych os, ktere popisuji puvodni vztahy v datech.

# nacteni balicku

install.packages("FactoMineR")
library(FactoMineR)

####################################### Priklad 1

# 1 724 zen vyplnilo dotaznik:

# Jak vypada podle Vas idealni rodina?
# 	pracuji oba - muz i zena ("oba.pracuji")
# 	muz pracuje vice nez zena ("muz.pracuje.vice")
# 	pouze muze pracuje ("pouze.muz.pracuje")
# Co je nejlepsi pro matku, jejiz deti zacaly chodit do skoly?
# 	zustat doma ("zena.v.domacnosti")
# 	prace na castecny uvazek ("zena.castecny.uvazek")
# 	prace na cely uvazek ("zena.cely.uvazek")
# Jaky nazor mate na nasledujici tvrzeni: Nezamestnane zeny se citi "vyrazeny" ze spolecnosti?
# 	naprosto souhlasim ("zena.vyrazena.ze.spolecnosti.souhlasim")
# 	spise souhlasim ("zena.vyrazena.ze.spolecnosti.spise.souhlasim")
# 	spise nesouhlasim ("zena.vyrazena.ze.spolecnosti.spise.nesouhlasim")
# 	zcela nesouhlasim ("zena.vyrazena.ze.spolecnosti.nesouhlasim")
# Data byla sbirana v podobe kontingencni tabulky, kde prvni parametr definuje radky a dalsi dva parametry sloupce tabulky.
# Kazda hodnota v tabulce odkazuje na pocet zen, ktere zapadaji do uvedenych kategorii.

setwd("c:/Users/brozoval/Desktop/5. cvičení - CA, NMDS/")
women_work<-read.table("dotaznik_zeny_francie1974.csv", header=TRUE, row.names=1, sep=";",dec=".")
women_work

sum(women_work[,c(1:3)]) # 1 724 zen
sum(women_work[,c(4:7)]) # 1 724 zen
# prvni 3 sloupce vyuzijeme pro vypocet CA, dalsi 4 sloupce budou v pozici vysvetlujicich promennych ("supplementary variables").

# analyza kontingencni tabulky pred CA
kt <- as.matrix(women_work[,1:3])
addmargins(prop.table(kt)) 
# ze vsech dotazanych by nejvice zen nejradeji pracovalo na castecny uvazek po tom, co deti zacnou chodit do skoly a zastavaji nazor, ze idealni je takova rodina,
# kde pracuje pouze muz nebo muz pracuje vice nez zena

?CA
CAwomen<-CA(women_work,col.sup=4:7)
summary(CAwomen)    # pocet os = min(nrow,ncol)-1

CAwomen$eig			    # inercie, % inercie popsane danou osou (prvni osa popisuje 86% celkove inercie)

CAwomen$row			    # vysledky CA, defaultne pro prvnich 10 radku z puvodni matice
CAwomen$row$coord		# souradnice radkovych kategorii na novych osach 
CAwomen$row$contrib	# prispevek radku k definici nove osy (v sloupecku soucet dava 100%)
CAwomen$row$cos2		# mezi 0-1 (v radku soucet dava 1), blizke jedne pokud je dobre popsano na dane ose 
CAwomen$row$inertia # inercie radkovych kategorii

#analogicke vystupy pro sloupcove kategorie

CAwomen$col			   
CAwomen$col$coord		# prvni osa nejlepe rozlisi kategorie:zena.v.domacnosti a zena.cely.uvazek (tyto budou mit i nejvyssi contrib)
CAwomen$col$contrib
CAwomen$col$cos2	  # zena.castecny.uvazek definuje smer 2. osy
CAwomen$col$inertia # u zena.v.domacnosti je pozorovana nejvyssi inercie (nejvetsi odchylka ocekavanych a pozorovanych cetnosti)
prop.table(kt,2)    # overeni podle kontingencni tabulky - v prvnim sloupci mame nejvetsi rozdily v % zastoupeni radkovych kategorii

# stejne vystupy pro dodatecne promenne (bez $contrib, protoze tyto neprispivaly k definici novych os)
CAwomen$col.sup     # kategorie jsou popsany na 1. ose CA (2. osa popisuje minimum inercie)
CAwomen$row.sup

summary(CAwomen,nbelements=Inf,file="CA results.csv") # vypise vysledky pro vsechny radky/sloupce vstupujici do analyzy,
# ulozi soubor "CA results.csv" do slozky, kam je pomoci setwd nastavena cesta

# vykresleni vysledku CA
plot(CAwomen,invisible=c("col.sup","row.sup"),cex=0.8) # pomoci select.row/select.col muzeme vykreslit jen vybrane klategorie
?plot.CA

# Interpretace
# Zeny, ktere uvedly, ze idealni aktivita po tom, co deti zacnou chodit do skoly, je prace na plny uvazek, uvedly take, ze
# v rodine by meli idealne oba rodice pracovat. Zeny, ktere uvedly, ze idealni aktivita po tom, co deti zacnou chodit do skoly, 
# je zustat doma, nejcasteji take uvedli, ze v rodine by mel pracovat pouze muz (s nizkou frekvenci uvedli, ze maji pracovat oba rodice).

plot(CAwomen,axes=c(1,2),cex=0.8) # vykresleni grafu s doplnkovymi promennymi

# Interpretace
# Zeny, ktere uvedly, ze v rodine by mel pracovat pouze muz, si nemysli, ze by se citily vyrazene ze spolecnosti.
# Zeny, ktere uvedly, ze v rodine by mel pracovat muz vice, castecne nebo zcela souhlasi s tim, ze se zena citi 
# vyrazena ze spolecnosti, pokud bude nezamestnana. 

####################################### Samostatny ukol 1

# nactete soubor "domaci_povinnosti.csv". Radky odkazuji na domaci cinnosti, sloupecky popisuji, kdo tyto cinnosti nejcasteji vykonava, jednotlive bunky
# odkazuji na pocet respondentu, kteri spadaji do uvedenych kategorii

domaci_povinnosti<-read.table("domaci_povinnosti.csv", header=TRUE, row.names=1, sep=";",dec=".")
domaci_povinnosti

# vykresleni - "balonkovy" graf
install.packages("gplots")
library("gplots")
dt <- as.table(as.matrix(domaci_povinnosti)) # pro vykresleni je nutne matici ulozit v podobe tabulky
balloonplot(t(dt), main ="domaci povinnosti", xlab ="", ylab="", label = FALSE, show.margins = FALSE)

# spocitejte CA
CA_povinnosti<-CA(domaci_povinnosti, ncp = 5, graph = TRUE) # ncp = pocet novych dimenzi uchovanych ve vysledku
summary(CA_povinnosti)
?CA

# projdete vypocet inercie v souboru "domaci_povinnosti_inercie.xlsx" a srovnejte s vysledky z R
CA_povinnosti$eig
CA_povinnosti$row$inertia
CA_povinnosti$col$inertia

# Vykreslete a interpretujte biplot:
#		a) ktere cinnosti casteji provadi pouze zeny?
#		b) ktere cinnosti casteji provadi pouze muz?
#		c) ktere cinnosti casto provadi muz a zena spolecne?

plot(CA_povinnosti)

####################################### Samostatny ukol 2

# Nactete soubor "Smoking.csv". Radky odkazuji na pozici zamestnance ve firme, sloupecky popisuji stupen kuractvi 
# a zda zamestnanci konzumuji mimo pracovni dobu alkohol.
# Jednotlive bunky tabulky odkazuji na pocet respondentu, kteri spadaji do uvedenych kategorii.

smoking<-read.table("Smoking.csv", header=TRUE, row.names=1, sep=";",dec=".")
smoking

# Spocitejte CA (konzumaci alkoholu pouzijte jako doplnkovou promennou, ktera nedefinuje smer novych os) a interpretujte:
# - Jaka je hodnota celkove inercie a jak je tato rozlozena do novych os?
# - Jake procento celkove inercie popisuji dve nove osy?
# - Ktere radkove kategorie definuji smer 1. osy?
# - Jaky stupen kuractvi je nejvice vazan s typem zamestnance?
# - Jake pracovni pozice jsou nejcasteji obsazeny silnymi kuraky?
# - Jsou kuraci spise mladsi nebo spise starsi zamestnanci?


##########################################################################################################
###########################################       NMDS     ###############################################
##########################################################################################################

# NMDS - Non-Metric Dimensional Scaling
# aplikuje iterativni proces, pri kterem hleda takove osy, ktere maximalizuji korelaci mezi poradim vzdalenosti v puvodnim
# prostoru a poradim vzdalenosti na predem urcenem poctu novych os.
# Pomoci stresu (odchylek vzdalenosti v puvodnim a novem prostoru) hodnotime, zda nove osy dobre popisuji puvodni vicerozmerny prostor.
# metaMDS provede PCA na novych osach z NMDS - ve vysledku muzeme rici, ze nove osy popisuji max rozptylu v novem prostoru

# nacteni balicku

install.packages("vegan")
library(vegan)

####################################### Priklad 2

# Provedte NMDS na tabulce, ktera popisuje nazory ruzne vzdelanych lidi na to, 
# proc by zena/par mel vahat, zda budou mit deti (budou zakladat rodinu).

data(children)
head(children)
children<-children[-c(15:18),-c(6:8)]
?children

?metaMDS
NMDSchildren<-metaMDS(children, distance="bray",k=2) 
stressplot(NMDSchildren, main="Shepard plot") # vysoka korelace -> vzdalenosti dobre koresponduji puvodnim vzdalenostem
NMDSchildren$stress # hodnoty stresu: idealne pod 0.1 (do 0.2 lze interpretovat, nad 0.2 nelze interpretovat)

plot(NMDSchildren, type="t", main=paste("NMDS stress =", round(NMDSchildren$stress,3))) 

# Interpretujte:
# - Proc se boji nekvalifikovani lide zakladat rodinu?
# - Lide s jakym vzdelanim se boji zakladat rodinu z duvodu, ze prijde valka? 

# domaci ukol: podivejte se na video: https://www.youtube.com/watch?v=Kl49qI3XJKY.