#==========================================================
#                 P R I K L A D    4.12 (4.11)
#==========================================================

1- pnorm(80, 72, 9) # pst, ze jeden student bude mit z testu alespon 80 bodu
1- pnorm(80, 72, sqrt(9^2 / 10)) # pst, ze prumerny pocet bodu u 10 studenti bude alespon 80 bodu

#==========================================================
#                 P R I K L A D    4.13 (4.12)
#==========================================================

# X ~ N(3078.94, 697^2)
xfit <- seq(950, 5300, length = 1000) # 1000 bodu souradnice osy x
fx <- dnorm(xfit, 3078.94, 697) # hustota N(3078.94, 697^2) vypocitana nad body posloupnosti xfit

par(mfrow = c(1, 2)) # rozdeleni grafu na dve okna
plot(xfit, fx, type = 'l', col = 'red', lwd = 2) #graf hustoty

Fx <- pnorm(xfit, 3078.94, 697) # disrt. fce N(3078.94, 697^2) vypocitana nad body posloupnosti xfit
plot(xfit, Fx, type = 'l', col = 'red', lwd = 2) # graf distr. fce
7

#==========================================================
#                 P R I K L A D    8.1
#==========================================================
data <- read.delim('01-one-sample-mean-skull-mf.txt') # nacteni dat
data <- na.omit(data) # odstraneni NA hodnot
head(data)

skull.BM <- data[data$sex == 'm', 'skull.B'] # sirka mozkovny pro muze
skull.BF <- data[data$sex == 'f', 'skull.B'] # sirka mozkovny pro zeny

n1 <- length(skull.BM) # rozsah hodnot pro muze
n2 <- length(skull.BF) # rozsah hodnot pro zeny

# H0: mu1 = mu2 -> H0: mu1 - mu2 = 0
# H1: mu1 != mu2 -> H0: mu1 - mu2 != 0
# alpha <- 0.05
# Test o mu1 - mu2 kdyz sigma1^2 a sigma2^2 nezname.

#------------------------------------------------------
# P R E D P O K L A D Y     T E S T U
#------------------------------------------------------

# Predpoklady: normalita skull.BM
#              normlaita skull.BF
# shoda rozptylu sigma_1^2 = sigma_2^2

# Normalita pro muze (n = 216 > 30)
nortest::lillie.test(skull.BM)$p.val
# p = 0.0766 > alpha -> H0 nezamitame na alpha = 0.05.

# Normalita pro zeny (n = 109 > 30)
nortest::lillie.test(skull.BF)$p.val
# p = 0.06380 > alpha -> H0 nezamitame na alpha = 0.05.

# Shoda rozptylu:
# H0: sigma1^2 = sigma2^2 -> H0: sigma1^2/sigma2^2 = 1
# H1: sigma1^2 != sigma2^2 -> H1: sigma1^2/sigma2^2 != 1
# alpha = 0.05

# Test o podilu rozptylu
# Testovani kritickym oborem:
s1 <- sd(skull.BM)
s2 <- sd(skull.BF)
(t0 <-  s1^2 / s2^2) # testovaci statistika
alpha <- 0.05
qf(alpha / 2, n1 - 1, n2 - 1)
qf(1 - alpha / 2, n1 - 1, n2 - 1)
# W = (-Inf ; 0.7267> u <1.4012 ; Inf)
# t0 = 1.055543 nenalezi W -> H0 nezamitame.

# Testovani IS:
# oboustranna alt -> oboustranny IS
c <- 1
(dh <- (s1^2 / s2^2) / qf(1 - alpha /2, n1 - 1, n2 - 1)) # dolni hranice IS
(fh <- (s1^2 / s2^2) / qf(alpha /2, n1 - 1, n2 - 1)) # horni hranice IS
# IS = (0.7533 ; 1.4526); c = 1 nalezi IS, H0 nezmaitame.

# Testovani p-hodnotou
(p <- 2 * min (1 - pf(t0, n1 - 1, n2 - 1), pf(t0, n1 - 1, n2 - 1)))
# p = 0.7610; p > alpha -> H0 nezamitame na hladine vyznamnosti alpha = 0.05

# Test o rozptylu pomoci funkce var.test()
# chybi akorat hranice kritickeho oboru (nutne dopocitat)
var.test(skull.BM, skull.BF, 
         alt = 'two.sided', # oboustranna alternativa
         conf.level = 0.95) # alpha = 0.05 -> spolehlivost (conf.level) = 0.95

# Zaver: H0 -> nezamitame -> mezi rozptyly sigma1^2 a sigma2^2 neexistuje statisticky vyznamny rozdil
# -> rozptyly sigma1^2 a sigma2^2 jsou shodne

# -> Normalita pro muze OK, normalita pro zeny OK, shoda rozptylu OK
# -> H0 ze zadani testujeme pomoci testu o rozdilu mu1  - mu2, kdyz rozptyly nezname, ale vime, ze jsou shodne

#------------------------------------------------------
# T E S T    H Y P O T E Z Y    Z E    Z A D A N I  (konecne!)
#------------------------------------------------------
# H0: mu1 = mu2  -> H0: mu1 - mu2 = 0
# H1: mu1 != mu2 -> H0: mu1 - mu2 != 0
# alpha <- 0.05
# Test o mu1 - mu2 kdyz sigma1^2 a sigma2^2 nezname (ale vime, ze jsou shodne).

# Test pomoci kritickeho oboru
m1 <- mean(skull.BM)
m2 <- mean(skull.BF)
c <- 0
(sh2 <- ((n1 - 1) * s1^2 + (n2 - 1) * s2^2) / (n1 + n2 - 2)) # vazeny prumer vyberovych rozptylu (s_hvezdicka^2)
(t0 <- ((m1 - m2) - c) / (sqrt(sh2) * sqrt(1/n1 + 1 / n2))) # testovaci statistika
qt(alpha/2, n1 + n2 - 2)
qt(1 - alpha / 2, n1 + n2 - 2)
# t0 = 5.4079; W = (-Inf ; -1.9673> u <1.9673 ; Inf)
# t0 nalezi W -> H0 zamitame

# Testovani IS:
(dh <- m1 - m2 - sqrt(sh2) * sqrt(1 / n1 + 1 / n2) * qt(1 - alpha / 2, n1 + n2 - 1)) # dolni hranice
(hh <- m1 - m2 - sqrt(sh2) * sqrt(1 / n1 + 1 / n2) * qt(alpha / 2, n1 + n2 - 2)) # horni hranice
# IS = (1.9331 ; 4.1437); c = 0 nenalezi IS -> H0 zamitame na hladine vyznamnosti alpha = 0.05

# Testovani p-hodnotou
2 * min (pt(t0, n1 + n2 - 2), 1 - pt(t0, n1 + n2 - 2))
# p = 1.2425 * 10^(-7); p < alpha -> H0 zamitame

# ZAVER: H0 zamitame -> mu1 != mu2 ->
# Mezi nejvetsi sirkou mozkovny muzu a zen staroveke egyptske populace existuje statisticky vyznamny rozdil.

# Test pomoci funkce t.test()
t.test(skull.BM, skull.BF, 
       alt = 'two.sided', # oboustranna alternativa
       conf.level = 0.95, # hl. vyznamnosti alpha = 0.05 
       var.equal = T) # orozptyly jsou shodne


#==========================================================
#                 P R I K L A D    8.2
#==========================================================

data <- read.delim('19-more-samples-correlations-skull.txt')
data <- na.omit(data)
cin <- data[data$pop == 'cin', 'nose.B']
ban <- data[data$pop == 'ban', 'nose.B']

shapiro.test(cin)$p.val
shapiro.test(ban)$p.val

var.test(cin, ban, alt = 'two.sided',
         conf.level = 0.95)

t.test(cin, ban, alt = 'greater',
       conf.level = 0.95, var.equal = F)

# Tip na procviceni: Dopocitat si 8.2 pomoci vzorecku.