setwd(' ')

xx<-read.table("kostky.csv",header=TRUE,sep=";")
x<-xx$sestky # data ulozime do vektoru x

#1a) prozkoumame data
x<-xx$sestky # data ulozime do vektoru x
table(x)
barplot(table(x),xlab = 'pocet sestek', ylab = 'absolutni cetnost',col=3,density=50)
plot(table(x),xlab = 'pocet sestek', ylab = 'absolutni cetnost')

#1bc)
# zkusime binomicky model; odhadneme jeho parametr p odpovidajici pravdepodobnosti padnuti sestky
n<-length(x)
p<-sum(x)/(10*n) # pocet padnutych sestech/pocet vsech hodu
p # odhadnuta pravdepodobnost padnuti sestky
1/6 # porovnejme s teoretickou hodnotou

#1d)
# muzeme porovnat relativni cetnosti a odhadnute pravdepodobnosti binomickym modelem
prop.table(table(x))
round(dbinom(0:10,10,p),4)

# toto neni prilis prehledne, proto se spise porovnavaji ocekavane o pozorovane cetnosti
table(x) # pozorovane cetnosti
round(n*dbinom(0:10,10,p)) # ocekavane cetnosti

# ulozme tyto hodnoty do vektoru a vykresleme prislusne grafy
poz_cetnosti<-as.numeric(table(x))
oce_cetnosti<-round(n*dbinom(0:10,10,p))

plot(0:10, oce_cetnosti, type = 'h', lty = 2, col = 'black',xlab = 'Pocet sestek', ylab = 'Cetnosti')
points(0:10, oce_cetnosti, pch = 21, col = 'black', bg = 'black')
lines(0:4, poz_cetnosti, type = 'h', col = 'red')
points(0:4, poz_cetnosti, pch = 21, col = 'darkred', bg = 'red')

legend('topright', pch = c(21, 21), col = c('darkred', 'black'), pt.bg = c('red', 'black'),
       legend = c('pozorovane', 'ocekavane'), bty = 'n')

#1e)
# pro popis modelu vykreslime jeste graf odhadnute pravdepodobnostni funkce
plot(0:10,dbinom(0:10,10,p),type='h',xlab = 'Pocet sestek', ylab = 'Pravdepodobnostni funkce',main='Binomicke rozdeleni')
points(0:10,dbinom(0:10,10,p), col = 'red', pch = 19)

# a jeste graf odhadnute distribucni funkce
plot(0:10, pbinom(0:10,10,p), type = 's', main = 'Binomicke rozdeleni',
     xlab = "Pocet sestek", ylab = 'Distribucni funkce', col = 'blue')

# komu se nelibi, muzeme jej malinko vylepsit
plot(0:10, pbinom(0:10,10,p), type = 'n', main = 'Binomicke rozdeleni',
     xlab = "Pocet sestek", ylab = 'Distribucni funkce', col = 'blue')
segments (-1:10, c(0,pbinom(0:10,10,p)) , 0:11, c(0,pbinom(0:10,10,p)))
points (0:10, c(0,pbinom(0:9,10,p)), pch = 21, col=1)
points (0:10, pbinom(0:10,10,p), pch = 19, col=2)

#1f)


########################################################################
########################################################################
#2

x<-xx$neuspechy # data ulozime do vektoru x
table(x)
barplot(table(x),xlab = 'pocet neuspechu', ylab = 'absolutni cetnost',col=3,density=50)
plot(table(x),xlab = 'pocet neuspechu', ylab = 'absolutni cetnost')

#2bc)
# zkusime geometricky model; odhadneme jeho parametr p odpovidajici pravdepodobnosti padnuti sestky
n<-length(x)
p<-1/mean(x+1) # prevracena hodnota prumerneho poctu vsech hodu
p # odhadnuta pravdepodobnost padnuti sestky
1/6 # porovnejme s teoretickou hodnotou

#2d)
# muzeme porovnat relativni cetnosti a odhadnute pravdepodobnosti geometrickym modelem
prop.table(table(x))
round(dgeom(0:15,p),4)

# toto neni prilis prehledne, proto se spise porovnavaji ocekavane o pozorovane cetnosti
table(x) # pozorovane cetnosti
round(n*dgeom(0:15,p)) # ocekavane cetnosti

# ulozme tyto hodnoty do vektoru a vykresleme prislusne grafy
poz_cetnosti<-as.numeric(table(x))
oce_cetnosti<-round(n*dgeom(0:15,p))

plot(0:15, oce_cetnosti, type = 'h', lty = 2, col = 'black',xlab = 'Pocet neuspechu', ylab = 'Cetnosti')
points(0:15, oce_cetnosti, pch = 21, col = 'black', bg = 'black')
lines(c(0,1,2,4,5,6,7,8,9,10,11,14), poz_cetnosti, type = 'h', col = 'red')
points(c(0,1,2,4,5,6,7,8,9,10,11,14), poz_cetnosti, pch = 21, col = 'darkred', bg = 'red')

legend('topright', pch = c(21, 21), col = c('darkred', 'black'), pt.bg = c('red', 'black'),
       legend = c('pozorovane', 'ocekavane'), bty = 'n')

#2e)
# pro popis modelu vykreslime jeste graf odhadnute pravdepodobnostni funkce
plot(0:15,dgeom(0:15,p),type='h',xlab = 'Pocet neuspechu', ylab = 'Pravdepodobnostni funkce',main='
Geometricke rozdeleni',ylim=c(0,0.2))
points(0:15,dgeom(0:15,p), col = 'red', pch = 19)

# a jeste graf odhadnute distribucni funkce
plot(0:15, pgeom(0:15,p), type = 's', main = 'Geometricke rozdeleni',
xlab = "Pocet neuspechu", ylab = 'Distribucni funkce',col = 'blue', ylim=c(0,1))

# komu se nelibi, muzeme jej malinko vylepsit
plot(0:15, pgeom(0:15,p), type = 'n', main = 'Geometricke rozdeleni',
     xlab = "Pocet neuspechu", ylab = 'Distribucni funkce', col = 'blue',ylim=c(0,1))
segments (-1:15, c(0,pgeom(0:15,p)) , 0:16, c(0,pgeom(0:15,p)))
points (0:15, c(0,pgeom(0:14,p)), pch = 21, col=1)
points (0:15, pgeom(0:15,p), pch = 19, col=2)

#2f)


########################################################################
########################################################################
#3

x<-xx$sestky20s # data ulozime do vektoru x
table(x)
barplot(table(x),xlab = 'pocet neuspechu', ylab = 'absolutni cetnost',col=3,density=50)
plot(table(x),xlab = 'pocet neuspechu', ylab = 'absolutni cetnost')

#3bc)
# zkusime Poissonuv model; odhadneme jeho parametr lambda odpovidajici ocekavanemu poctu sestek
n<-length(x)
lambda<-mean(x)
lambda # odhadnuty parametr lambda

#3d)
# muzeme porovnat relativni cetnosti a odhadnute pravdepodobnosti Poissonovym modelem
prop.table(table(x))
round(dpois(0:10,lambda),4)

# toto neni prilis prehledne, proto se spise porovnavaji ocekavane o pozorovane cetnosti
table(x) # pozorovane cetnosti
round(n*dpois(0:10,lambda)) # ocekavane cetnosti

# ulozme tyto hodnoty do vektoru a vykresleme prislusne grafy
poz_cetnosti<-as.numeric(table(x))
oce_cetnosti<-round(n*dpois(0:6,lambda))

plot(0:6, oce_cetnosti, type = 'h', lty = 2, col = 'black',xlab = 'Pocet sestek', ylab = 'Cetnosti')
points(0:6, oce_cetnosti, pch = 21, col = 'black', bg = 'black')
lines(0:5, poz_cetnosti, type = 'h', col = 'red')
points(0:5, poz_cetnosti, pch = 21, col = 'darkred', bg = 'red')

legend('topright', pch = c(21, 21), col = c('darkred', 'black'), pt.bg = c('red', 'black'),
       legend = c('pozorovane', 'ocekavane'), bty = 'n')

#3e)
# pro popis modelu vykreslime jeste graf odhadnute pravdepodobnostni funkce
plot(0:10,dpois(0:10,lambda),type='h',xlab = 'Pocet sestek', ylab = 'Pravdepodobnostni funkce',main='
Poissonovo rozdeleni')
points(0:10,dpois(0:10,lambda), col = 'red', pch = 19)

# a jeste graf odhadnute distribucni funkce
plot(0:10, ppois(0:10,lambda), type = 's', main = 'Poissonovo rozdeleni',
     xlab = "Pocet sestek", ylab = 'Distribucni funkce', col = 'blue',ylim=c(0,1))

# komu se nelibi, muzeme jej malinko vylepsit
plot(0:10, ppois(0:10,lambda), type = 'n', main = 'Poissonovo rozdeleni',
     xlab = "Pocet sestek", ylab = 'Distribucni funkce', col = 'blue',ylim=c(0,1))
segments (-1:10, c(0,ppois(0:10,lambda)) , 0:11, c(0,ppois(0:10,lambda)))
points (0:10, c(0,ppois(0:9,lambda)), pch = 21, col=1)
points (0:10, ppois(0:10,lambda), pch = 19, col=2)

#3f)



########################################################################
########################################################################
#4

data<-read.table("newborns.txt",header=TRUE)
x<-data$prch.N
x<-na.omit(x)