# Uma Análise dos Modelos de Mensuração de Exposição ao Risco 
# Financeiro no Mercado de Manganês
# Orientando: Lucas Marin Avelleda
# Orientador:Prof.HedibertFreitas Lopes, PhD
#
# Top 10 Manganese Ore Importing Countries in 2017
# China 54%
# India 10%
# Japan 5%
# South Korea 5%
# Norway 4%
# Malaysia 4%
# Russia 3%
# Ukraine 3%
# Spain 2%
# France 1%
# Rest of the world 9%
#
# China is the world’s largest manganese ore importer besides accounting 
# for 9 of the world’s production It absorbs material from the main 
# producers South Africa, Australia, Gabon, Ghana and Brazil.  Other 
# significant manganese ore importers are India, Japan and South Korea
# The most significant manganese ore producers are settled in South Africa 
# (the country produced more than 1/3 of the world’s manganese in 2017), 
# Australia, Gabon, China and Brazil
library("stochvol")

data = read.csv("https://hedibert.org/wp-content/uploads/2024/06/mnprice.txt")

names(data)

attach(data)

n    = nrow(data)

date = data[c(1,seq(50,250,by=50),n),1]

y    = Price[2:n]/Price[1:(n-1)]-1

par(mfrow=c(1,1))
plot(Price,axes=FALSE,xlab="Week",ylab="USD/DMTU (Dry metric tonne unit)",type="l")
axis(2);box()
axis(1,at=c(1,seq(50,250,by=50),n),lab=date)

hist(y,breaks=seq(min(y),max(y),length=30),prob=TRUE,main="",
     xlab="Return")

plot(y,xlab="Week",ylab="Return",type="l",axes=FALSE)
axis(2);box()
axis(1,at=c(1,seq(50,250,by=50),n),lab=date)


fit.sv   = svsample(y,draws=5000,burnin=1000,thin=1)
fit.svt  = svtsample(y,draws=5000,burnin=1000,thin=1)
fit.svl  = svlsample(y,draws=5000,burnin=1000,thin=1)
fit.svtl = svtlsample(y,draws=5000,burnin=1000,thin=1)

q.sv   = t(apply(exp(fit.sv$latent[[1]]/2),2,quantile,c(0.1,0.5,0.9)))
q.svt  = t(apply(exp(fit.svt$latent[[1]]/2),2,quantile,c(0.1,0.5,0.9)))
q.svl  = t(apply(exp(fit.svl$latent[[1]]/2),2,quantile,c(0.1,0.5,0.9)))
q.svtl = t(apply(exp(fit.svtl$latent[[1]]/2),2,quantile,c(0.1,0.5,0.9)))


par(mfrow=c(2,3))

plot(density(fit.sv$para[[1]][,1]),xlim=c(-10,0),xlab=expression(mu),ylim=c(0,1.3),main="",lwd=2)
lines(density(fit.svt$para[[1]][,1]),col=2,lwd=2)
lines(density(fit.svl$para[[1]][,1]),col=3,lwd=2)
lines(density(fit.svtl$para[[1]][,1]),col=4,lwd=2)
legend("topright",legend=c("SV","SV-t","SV-l","SV-lt"),col=1:4,lwd=2,lty=1,bty="n")


plot(density(fit.sv$para[[1]][,2]),xlim=c(0.3,1),xlab=expression(phi),ylim=c(0,10),main="",lwd=2)
lines(density(fit.svt$para[[1]][,2]),col=2,lwd=2)
lines(density(fit.svl$para[[1]][,2]),col=3,lwd=2)
lines(density(fit.svtl$para[[1]][,2]),col=4,lwd=2)

plot(density(fit.sv$para[[1]][,3]),xlim=c(0,2.5),xlab=expression(sigma),
     ylim=c(0,2.5),main="",lwd=2)
lines(density(fit.svt$para[[1]][,3]),col=2,lwd=2)
lines(density(fit.svl$para[[1]][,3]),col=3,lwd=2)
lines(density(fit.svtl$para[[1]][,3]),col=4,lwd=2)

plot(density(fit.svtl$para[[1]][,4]),col=3,xlab=expression(nu),main="",lwd=2,
     xlim=c(2,5),ylim=c(0,8))
lines(density(fit.svt$para[[1]][,4]),col=4,lwd=2)
legend("topright",legend=c("SV-t","SV-lt"),col=3:4,lwd=2,lty=1,bty="n")

plot(density(fit.svl$para[[1]][,5]),col=3,xlab=expression(rho),xlim=c(-0.5,0.5),
     ylim=c(0,5),main="",lwd=2)
lines(density(fit.svtl$para[[1]][,5]),col=4,lwd=2)
legend("topright",legend=c("SV-l","SV-lt"),col=3:4,lwd=2,lty=1,bty="n")


par(mfrow=c(2,2))
plot(q.sv[,2],xlab="Week",main="SV",lwd=2,type="l",ylab="Standard deviation",
ylim=c(0,max(q.sv)))
lines(q.sv[,1],col=2)
lines(q.sv[,3],col=2)


plot(q.svt[,2],xlab="Week",main="SV-t",lwd=2,type="l",ylab="Standard deviation",
ylim=c(0,max(q.svt)))
lines(q.svt[,1],col=2)
lines(q.svt[,3],col=2)

plot(q.svl[,2],xlab="Week",main="SV-l",lwd=2,type="l",ylab="Standard deviation",
ylim=c(0,max(q.svl)))
lines(q.svl[,1],col=2)
lines(q.svl[,3],col=2)

plot(q.svtl[,2],xlab="Week",main="SV-tl",lwd=2,type="l",ylab="Standard deviation",
ylim=c(0,max(q.svtl)))
lines(q.svtl[,1],col=2)
lines(q.svtl[,3],col=2)


par(mfrow=c(1,1))
plot(q.sv[,2],xlab="Week",main="",lwd=2,type="l",ylab="Standard deviation",ylim=c(-0.05,0.32))
lines(q.svt[,2],col=2,lwd=2)
lines(q.svl[,2],col=3,lwd=2)
lines(q.svtl[,2],col=4,lwd=2)
lines(0.2*y,type="h",lwd=2)
legend("topleft",legend=c("SV","SV-t","SV-l","SV-lt"),col=1:4,lwd=2,lty=1,bty="n")