paretoDens<-function(a,b,x){

if(x<a){ return(0)}

else{

return( b*(a^b)/(x^(b+1)))
}

}


paretoQuant<-function(a,b,y){

if(y>=1 || y<0){ return(NaN)}

else{
return(a/((1-y)^(1/b)))
}

}

paretoSample<-function(a,b,n){

return( paretoQuant(a,b, runif(n)))

}

ML_est_pareto<-function(sample){
	a<-min(sample)
	n<-length(sample)
	b<--n/(n*log(a)-sum(log(sample)))
	return(c(a,b))
}

empVaR<-c()
ML_VaR<-c()


for(i in 1:10000){
	
	hist(S)
	empVaR<-c(empVaR, quantile(S, 0.995, type=4))

	ab<- ML_est_pareto(S)
	ML_VaR<-c(ML_VaR, paretoQuant(ab[1], ab[2], 0.995))

}

path<-"C:\\Eigene Dateien\\Lehre\\RTII\\"
png(paste(path, "emp_VaR.png"))

mybreaks<-seq(from=0, to=50, by=2)

hist(empVaR, breaks=mybreaks)

dev.off()
png(paste(path, "ML_VaR.png"))

mybreaks<-seq(from=0, to=50, by=2)

hist(ML_VaR, breaks=mybreaks)
dev.off()

png(paste(path, "Diff_emp-ML_VaR.png"))

mybreaks<-seq(from=-20, to=44, by=2)

hist(empVaR-ML_VaR, breaks=mybreaks)
dev.off()