library(actuar)

#1.)
Y=c(99,42,29,28,17,9,3)
c=c(0,7500,17500,32500,67500,125000,300000,Inf)
g=grouped.data(c,Y)

#a)
plot(ogive(g),main='Empiricka distribucni funkce')

#b)

hist(g)    #nefunguje pro nekonecne hodnoty

#upravime nase puvodni data

c2=c(0,7500,17500,32500,67500,125000,300000,1000000)
g2=grouped.data(c2,Y)
hist(g2)

#c)
#logaritmicka verohodnostni funkce pro nase data
f=function(x){
f=0
for (i in 1:7){
f=f+Y[i]*log(pexp(c[i+1],1/x)-pexp(c[i],1/x))
}
return(f)
}

opt=optimize(f,c(0,100000),maximum=T)$max  # hledani maxima v intervalu (0,1 000 000)
opt  # MLE parametru lambda, odhad prumerne vyse plneni jednoho klienta

#d)
n=sum(Y)
#Pearsonova chi^2 statistika pro nase data
chi=function(x){
p=rep(0,7) #teoreticke cetnosti
for (i in 1:7){ 
p[i]=pexp(c[i+1],1/x)-pexp(c[i],1/x)}

chi=sum((Y-n*p)^2/(n*p))
return(chi) 
}

opt2=optimize(chi,c(0,100000))$min    # hledani minima v intervalu (0,1 000 000)
opt2   # odhad parametru lambda metodou minimalniho chi kvadrat



#e)
par(mfrow=c(1,2))
plot(ogive(c,Y))
curve(pexp(x,1/opt),col=2,add=TRUE)
curve(pexp(x,1/opt2),col=3,add=TRUE)

hist(g2)
curve(dexp(x,1/opt),col=2,add=TRUE)
curve(dexp(x,1/opt2),col=3,add=TRUE)

#detailni pohled na nase data
c3=c(0,7500,17500,32500,67500,125000,300000,400000)
g3=grouped.data(c3,Y)
hist(g3)
curve(dexp(x,1/opt),col=2,add=TRUE)
curve(dexp(x,1/opt2),col=3,add=TRUE)


#2.)
#a)
mean(g2)
opt
opt2

quantile(g,0.5)
opt*log(2)
opt2*log(2)


#b)
1-ogive(g)(100000)
1-pexp(100000,1/opt)
1-pexp(100000,1/opt2)

#c)
quantile(g,0.9)
qexp(0.9,1/opt)
qexp(0.9,1/opt2)