##########################################################################################
################################## SOFTWARE R - CAST 1 ###################################
##########################################################################################

# software R volne stazitelny z http://www.r-project.org/

# Nastaveni pracovniho adresare: File -> Change Dir...
# nebo pomoci prikazu setwd(), napr.:
setwd("f:/Prace/Vyuka/FSTA_2015")
getwd()

############### Nacteni dat a jejich zobrazeni ###############

readLines("Data_pro_R.txt",3) # nacteni prvnich 3 radku textoveho souboru
# prikaz readLines vhodny na zkontrolovani:
# 1) typu oddelovace sloupcu - zde tabulator -> v read.table bude: sep="\t"
# 2) zda jsou v souboru pouzity desetinne carky ci tecky - zde tecky -> v read.table bude: dec="."
# 3) zda ma soubor hlavicku - zde ma soubor hlavicku -> v read.table bude: header=T

kardio <- read.table(file="Data_pro_R.txt", sep="\t", dec=".", header=T) # nacteni dat

dim(kardio)                       # pocet radku a sloupcu datove matice
nrow(kardio)                      # pocet radku datove matice
ncol(kardio)                      # pocet sloupcu datove matice
names(kardio)                     # nazvy promennych
head(kardio)                      # prvnich 10 radku tabulky
tail(kardio)                      # poslednich 10 radku tabulky
str(kardio)                       # sumarni informace o jednotlivych promennych

# zobrazeni a uprava datove tabulky:
# (pro zobrazeni je nutne zmenit textove promenne z "factor" na "character") 
for (i in 1:ncol(kardio)) {
  if (is.factor(kardio[,i])) {
    kardio[,i]<-as.character(kardio[,i])
  }
}
data.entry(kardio)                # zobrazeni datove tabulky; hodnoty je mozne editovat, ale moc to nedoporucuji
class(kardio)                     # trida objektu kardio je "list"
kardio <- as.data.frame(kardio)   # prikaz data.entry zmeni tridu objektu kardio na "list"; je nutne zmenit tridu zpet na "data.frame",
# aby se s daty dalo dále pracovat



############### Vyhledavani funkci a napoveda k funkcim ###############

# vyhledavani funkci - 2 mozne zpusoby:
help.search("read")
??"read"

# napoveda k funkcim - 2 mozne zpusoby:
help(read.table)
?read.table


## Ukol 1: Pomoci napovedy zjistete, jak je mozne vypsat prvnich 20 radku tabulky.
# Reseni:
?head
head(kardio,n=20L) # tohle je z napovedy 
head(kardio,20) # tohle je jednodussi zpusob



############### Indexovani ###############

# vypsani hodnot veku u vsech pacientu (nasledujici 3 prikazy jsou ekvivalentni):
kardio$vek 
kardio[,"vek"]
kardio[,2]

# vypsani hodnot veku u prvnich 10 pacientu (nasledujici 3 prikazy jsou ekvivalentni):
kardio$vek[1:10]
kardio[1:10,"vek"]
kardio[1:10,2]

# vyber vsech sloupcu krome prvniho:
kardio.pom <- kardio[,-1]
dim(kardio.pom)
names(kardio.pom)

# vyber vsech sloupcu krome druheho a pateho (tzn. krome mesic_narozeni a pohlavi)
kardio.pom <- kardio[,-c(2,5)]
dim(kardio.pom)
names(kardio.pom)

# vytvoreni nove datove tabulky, kde budou pouze zeny:
kardio.zeny <- kardio[kardio$pohlavi_rek=="Z",]
dim(kardio.zeny)



## Ukol 2:
# 1. Vytvorte datovou tabulku, ktera bude obsahovat pouze muze. Kolik ma tato nova tabulka radku (tzn. kolik je muzu v souboru)?
# 2. Vytvorte datovou tabulku, ktera bude obsahovat nekuraky (nekuraci jsou oznaceni hodnotou "2"). Kolik ma tato nova tabulka radku?
# 3*. Vytvorte datovou tabulku, ktera bude obsahovat parametry: id, vek, pohlavi_rek
# 4. Vypiste hodnoty vahy u 50. az 70. pacienta.
# 5*. Vypiste hodnoty pohlavi u 5., 12. a 40. pacienta (pouzijte promennou pohlavi_rek).
# 6**. vypiste hodnoty vysky u kazdeho 10. pacienta:

## Reseni:
# ad 1.
kardio.muzi <- subset(kardio,pohlavi_rek=="M")
dim(kardio.muzi) # 30
# ad 2.
kardio.nekurak <- subset(kardio,koureni==2)
dim(kardio.nekurak) # 63
# ad 3.
names(kardio) # na zkontrolovani, kolikate sloupce jsou id, vek a pohlavi_rek
kardio.pom <- kardio[,c(1,2,5)]
head(kardio.pom)
# ad 4.
kardio$vaha[50:70]
# [1]  68 112  78  60  58  60  84  71  52  89  70  64  45  73  85  85  64  70  86
#[20]  68  99
# ad 5.
kardio$pohlavi_rek[c(5,12,40)]
# [1] M Z Z
# ad 6.
kardio$vyska[seq(10,nrow(kardio),by=10)]
# [1] 164 176 173 165 175 175 179 170 172



############### Pridani a odstraneni promenne ###############

kardio$log.vaha <- log(kardio$vaha) # pridani promenne
dim(kardio)
names(kardio)
kardio$log.vaha <- NULL            # odstraneni promenne 
dim(kardio)
names(kardio)


## Ukol 3: 
# 1. Pridejte do datove tabulky promennou se zlogaritmovanou vyskou.
# 2. Pridejte do datove tabulky promennou s vyskou v metrech.

## Reseni:
# ad 1.
kardio$log.vyska <- log(kardio$vyska) # 1. zpusob
kardio <- transform (kardio, log.vyska=log(vyska)) # 2. zpusob
# ad 2.
kardio$vyska.m <- kardio$vyska/100 # 1. zpusob
kardio <- transform(kardio, vyska.m = vyska/100) # 2. zpusob



############### Sumarizace promennych - prvni pohled na data ###############

# sumarizace kategorialnich promennych:
summary(kardio$pohlavi)           # je patrne, ze u jedne zeny je misto "Z" zadano "F" - vytvorila jsem proto novou promennou "pohlavi_rek", ktera obsahuje jen "Z" a "M" (pozdeji si ukazeme, jak promennou "pohlavi_rek" vytvorit)
summary(kardio$pohlavi_rek)
levels(kardio$pohlavi_rek)        # vypise jen kategorie

# sumarizace spojitych promennych:
summary(kardio$vek)               # je patrne, ze funkce summary pozna typ promenne a podle toho vypise vysledek sumarizace
summary(kardio$vek_chybny)        # co je na promenne vek_chybny spatne? Je tam odlehla hodnota - vsimnete si, ze median zustava stejny, ale prumer se hodne zmeni


## Ukol 4: Sumarizujte vysku.
## Reseni:
summary(kardio$vyska) 


############### Vizualizace promennych pomoci grafu ###############

# vizualizace kategorialnich promennych:
pie(summary(kardio$pohlavi_rek))  # kolacovy graf
pie(summary(kardio$pohlavi))
# zapnuti nahravani historie vykreslenych grafu: kliknout na graficke okno -> History -> Recoding

# vizualizace spojitych promennych:
boxplot(kardio$vaha)              # krabicovy graf
boxplot(kardio$vek_chybny)        # je patrne, ze je tu odlehla hodnota
boxplot(kardio$vaha~kardio$pohlavi_rek) # krabicovy graf vahy podle pohlavi
hist(kardio$vaha)                 # histogram (pocet sloupcu je vypocitan pomoci Sturgesova pravidla)
hist(kardio$vaha,20,col="blue")   # histogram s 20 sloupci a modrou barvou
hist(log(kardio$vaha))            # histogram pro zlogaritmovanou vahu
hist(kardio$vek_chybny)

# vykresleni histogramu s krivkou hustoty normalniho rozdeleni:
x <- kardio$vyska
h<-hist(x, breaks=10, col="yellow", xlab="vyska",
        main="Histogram s krivkou hustoty normalniho rozdeleni")
xfit<-seq(min(x),max(x),length=40)
yfit<-dnorm(xfit,mean=mean(x),sd=sd(x))
yfit <- yfit*diff(h$mids[1:2])*length(x)
lines(xfit, yfit, col="blue", lwd=2) 

#ulozeni grafu
windows()

x <- kardio$vyska
h<-hist(x, breaks=10, col="yellow", xlab="vyska",
        main="Histogram s krivkou hustoty normalniho rozdeleni")
xfit<-seq(min(x),max(x),length=40)
yfit<-dnorm(xfit,mean=mean(x),sd=sd(x))
yfit <- yfit*diff(h$mids[1:2])*length(x)
lines(xfit, yfit, col="blue", lwd=2) 

savePlot("Graf_01", type="jpg", device = dev.cur(),restoreConsole = TRUE)
dev.off()


## Ukol 5.
# 1. Vykreslete histogram vysky se 4 kategoriemi a hnedou barvou.
# 2. Vykreslete krabicovy graf veku podle pohlavi, ktery bude mit cervenou barvu.
# 3*. Vykreslete kolacovy graf pro pohlavi, kde zeny budou mit cervenou barvu a muzi modrou barvu.
# 4*. Vykreslete kolacovy graf pro koureni (tip: je nutne prevest promennou koureni na faktor pomoci funkce as.factor).

## Reseni:
# ad 1.
hist(kardio$vyska,4,col="brown")
# ad 2.
boxplot(kardio$vek~kardio$pohlavi_rek,col="red")
# ad 3.
pie(summary(kardio$pohlavi_rek),col=c("blue","red"))
# ad 4.
pie(summary(as.factor(kardio$koureni)))



############### Rekodovani promennych ###############
# vytvoreni nove promenne "pohlavi_rek2", ktera bude obsahovat pouze hodnoty "F" a "M":
kardio$pohlavi_rek2 <- ifelse(kardio$pohlavi_rek=="Z","F","M") # rekodovani nove promenne
table(kardio$pohlavi_rek2) # overeni, ze promenna pohlavi_rek obsahuje pouze hodnoty "F" a "M"

# vytvoreni kategorii ze spojite promenne:
kardio$vek.2kat <- kardio$vek>=50 # vytvoreni dvou kategorii: 0(FALSE)...<50, 1(TRUE)...>=50
table(kardio$vek.2kat)
kardio$vek.3kat <- 1+(kardio$vek>=30)+(kardio$vek>=50) # vytvoreni trech kategorii: 1...<30, 2...30-50, 3...>=50 
table(kardio$vek.3kat) #overeni
summary(as.factor(kardio$vek.3kat)) #jiny zpusob overeni


# slouceni kategorii u kategorialni promenne:
install.packages("car")
library(car) 

table(kardio$cel_chol_norma,useNA="ifany") # 0...pod normou, 1...v norme, 2...nad normou, <NA>...chybejici hodnota
kardio$cel.chol.mimo.normu <- recode(kardio$cel_chol_norma,'c(0,2)=1;1=0;else=99') # 0...v norme, 1...mimo normu, 99...chybejici hodnoty
table(kardio$cel.chol.mimo.normu)


## Ukol 6. 
# 1. Vykreslete histogram vahy a navrhnete nejlepsi rozdeleni vahy do trech skupin.
# 2. Rekodujte promennou "pohlavi" tak, aby muzi meli hodnotu 1 a zeny hodnotu 0.
# 3. Rekodujte promennou "koureni", ktera obsahuje kategorie 0..kouri, 1..kouril(a) drive, 2..nikdy nekouril(a), tak, aby nova promenna obsahovala dve kategorie: 1..kurak x 0..nekurak.

## Reseni:
# ad 1.
hist(kardio$vaha,20)
kardio$vaha.4kat <- 1+(kardio$vaha>=60)+(kardio$vaha>=80) # vytvoreni ctyr kategorii: 1...<60, 2...60-80, 3...>=80 
summary(as.factor(kardio$vaha.4kat))
# ad 2.
kardio$pohlavi_01 <- recode(kardio$pohlavi,'c("F","Z")=0;"M"=1;else=99')
summary(kardio$pohlavi_01)
# ad 3.
kardio$koureni_01 <- recode(kardio$koureni,'c(0,1)=1;2=0;else=99')
summary(as.factor(kardio$koureni_01))



#  ------------------------------------------------------------------------

# Sumarizace dat ----------------------------------------------------------

#  ------------------------------------------------------------------------


############### Popisna statistika kategorialnich promennych ###############

# frekvencni tabulka pro jednu promennou:
table(kardio$pohlavi_rek,useNA="ifany") # frekvencni tabulka

# kontingencni tabulka pro dve promenne:
kt <- table(kardio$pohlavi_rek,kardio$koureni,useNA="ifany") # pohlavi: M..muz,Z..zena; koureni:0..kouri,1..drive,2..nikdy
kt
addmargins(kt)                    # pridani souctu radku a sloupcu

# kontingencni tabulka procent:
kt.procenta <- prop.table(kt)     # kontingencni tabulka procent (procenta z celkoveho poctu)
addmargins(kt.procenta)
kt.procenta <- prop.table(kt,1)   # kontingencni tabulka procent (radkova procenta)
addmargins(kt.procenta)
kt.procenta <- prop.table(kt,2)   # kontingencni tabulka procent (sloupcova procenta)
addmargins(kt.procenta)

cbind(kt,kt.procenta) # slouceni absolutnich hornot se sloupcovami procenty

write.table(cbind(kt,kt.procenta),file="Tab1.csv",sep=";",row.names=TRUE) # ulozeni tabulky do textoveho souboru
# Pozor! Pri otevreni souboru v Excelu je horni radek o jednu bunku posunuty!

## Ukol 1. Sumarizujte vysku podle pohlavi (vysku vhodne kategorizujte: 1...<170, 2...170-180, 3...>=180). Dale vytvorte kontingencni tabulku 
##         procent, z niz bude patrne, zda se v jednotlivych vyskovych kategoriich zeny a muzi lisi.

## Reseni:
hist(kardio$vyska,20)
kardio$vyska.3kat <- 1+(kardio$vyska>=170)+(kardio$vyska>=180) # vytvoreni trech kategorii: 1...<170, 2...170-180, 3...>=180 
summary(as.factor(kardio$vyska.3kat))
kt <- table(kardio$pohlavi_rek,kardio$vyska.3kat,useNA="ifany")
addmargins(kt)                    # pridani souctu radku a sloupcu
kt.procenta <- prop.table(kt,2)   # kontingencni tabulka procent (sloupcova procenta)
addmargins(kt.procenta)

############### Popisna statistika spojitych promennych ###############

# popisna statistika jedne promenne:
summary(kardio$vek)
summary(kardio$vek_chybny)          # pokud je v souboru odlehla hodnota, prumer je touto odlehlou hodnotou ovlivneny (je to videt ze srovnani sumarizace promennych vek a vek_chybny)

# smerodatna odchylka:
sd(kardio$vek)
sd(kardio$vek_chybny)               # pokud je v souboru odlehla hodnota, i smerodatna odchylka je touto odlehlou hodnotou ovlivnena

# popisna statistika spojite promenne podle kategorialni promenne:
summary(kardio$vek[kardio$pohlavi=="M"])
summary(kardio$vek[kardio$pohlavi=="Z"])

# popisna statistika spojite promenne podle kategorialni promenne v jednom vystupu:
aggregate(kardio$vek ~ kardio$pohlavi_rek, FUN = function(x) c(M=mean(x),Med=median(x), SD=sd(x)))

# vytvareni hezcich tabulek:
install.packages("Epi")
library(Epi)
pp<-stat.table(pohlavi_rek,contents=list(count(),mean(vek),median(vek),min(vek),max(vek)),data=kardio) # N, prumer, median, minimum a maximum veku podle pohlavi
pp

# uprava popisku sloupcu, aby tabulka vypadala lepe:
pp<-stat.table(index=list("Pohlavi"=pohlavi_rek),contents=list(N=count(),"Prumer"=mean(vek),"Median"=median(vek),"Minimum"=min(vek),"Maximum"=max(vek)),data=kardio)
pp

# ulozeni tabulky do promenne, prevod na typ "data frame" a ulozeni do textoveho souboru:
tab1 <- stat.table(index=list("Pohlavi"=pohlavi_rek),contents=list(N=count(),"Prumer"=mean(vek),"Median"=median(vek),"Minimum"=min(vek),"Maximum"=max(vek)),data=kardio) # ulozeni do promenne tab1
tab1 <- as.data.frame(t(tab1))    # prevod na typ "data frame"; je nutne tab1 transponovat pomoci t(tab1)
tab1[,2]<-round(tab1[,2],1)       # zaokrouhleni prumeru na jedno desetinne misto
write.table(tab1,file="Tab2.csv",sep=";",row.names=TRUE) # ulozeni tabulky do textoveho souboru
# Pozor! Pri otevreni souboru v Excelu je horni radek o jednu bunku posunuty!


seznam funkci pro popisnou statistiku:
  x <- kardio$vek
length(x)                         # pocet prvku promenne x
mean(x)                           # prumer
median(x)                         # median
min(x)                            # minimum
max(x)                            # maximum
range(x)                          # rozsah hodnot promenne x
c(min(x),max(x))                  # stejny vysledek jako range(x)
var(x)                            # rozptyl
sd(x)                             # smerodatna odchylka
quantile(x,0.25)                  # 1. kvartil (tzn. 25. percentil)
quantile(x,0.5)                   # 2. kvartil (tzn. 50. percentil) - totez co median(x) 
quantile(x,0.75)                  # 3. kvartil (tzn. 75. percentil)
sum(x)                            # suma hodnot x
prod(x)                           # soucin hodnot x

x <- kardio$hdl_cholesterol
mean(x)                           # prumer se nespocita, protoze promenna "hdl_cholesterol" obsahuje alespon jednu hodnotu NA
x <- kardio$hdl_cholesterol[!is.na(kardio$hdl_cholesterol)] # vyber vsech hodnot promenne "hdl_cholesterol" krome hodnot NA
mean(x)

