GHADAsim1

Download: GHADAsim1.xpl
Code:
setsize(800,400)
proc(MSE,est,step) = simest(n,sig)
i = 1
MSE = matrix(n,1)*(-100)
est = matrix(800,200)
step = matrix(200,6)
while (i<=n)
path = "D:\WIAS talk\XploRe program\simulation\sim1\"
name = "vol_sim"
file = path + name + string("%4.0f",i) + ".txt"
v = read(file)
vola = v[,4]
rr = v[,1]
; first jump
stepindex = 200:230
s1 = rr[stepindex]~(vola[stepindex]-0.05*0.4)
s2 = rr[stepindex]~(vola[stepindex]-0.05*0.5)
s3 = rr[stepindex]~(vola[stepindex]-0.05*0.6)
step[i,1] = paf(s1,s1[,2]>=0)[1,1]-400
step[i,2] = paf(s2,s2[,2]>=0)[1,1]-400
step[i,3] = paf(s3,s3[,2]>=0)[1,1]-400
; second jump
stepindex = 550:700
s1 = rr[stepindex]~(vola[stepindex]-0.05*0.4)
s2 = rr[stepindex]~(vola[stepindex]-0.05*0.5)
s3 = rr[stepindex]~(vola[stepindex]-0.05*0.6)
step[i,4] = paf(s1,s1[,2]<=0)[1,1]-750
step[i,5] = paf(s2,s2[,2]<=0)[1,1]-750
step[i,6] = paf(s3,s3[,2]<=0)[1,1]-750

MSE[i] = mean(sqrt((vola-sig[rr])^2))
est[,i] = vola
sigl = setmask(1:1000~sig,"line","dotted","red")
volal = setmask(rr~vola, "line", "blue")
volaplot = createdisplay(1,1)
show(volaplot,1,1, sigl, volal)
title = "sim"+string("%4.0f",i)
setgopt(volaplot,1,1,"title",title,"border",0)
fig = path ++string("%1.0f",i) + ".ps"
print(volaplot, fig)
i = i+1
endo
MSE = paf(MSE,MSE!=-100)
write(MSE,path+"MSE.txt",path+"MSE.txt")
endp

library("plot")
sig = matrix(400,1)*0.01|matrix(350,1)*0.05|matrix(250,1)*0.01
MSE = simest(200,sig)
step = MSE.step
library("stats")
descriptive(step)
stepbox = step[,1]|step[,2]|step[,3]
theFactor = string("40% rule", 1:200) | string("50% rule", 1:200)| string("60% rule", 1:200)
plotbox(stepbox, theFactor)
setgopt(boxplot,1,1,"title","Boxplots of the steps detecting the first jump","border",0)
/*
stepbox = step[,4]|step[,5]|step[,6]
theFactor = string("40% rule", 1:200) | string("50% rule", 1:200)| string("60% rule", 1:200)
plotbox(stepbox, theFactor)
setgopt(boxplot,1,1,"title","Boxplots of the steps detecting the second jump","border",0)
*/

est = MSE.est'
eststd = sqrt(var(est))'
confip = mean(est)'+2.58*eststd ; 99% confidence
confin = mean(est)'-2.58*eststd
t = 201:1000
estml = setmask(t~mean(est)',"line")
confipl = setmask(t~confip,"line","dashed","red")
confinl = setmask(t~confin,"line","dashed","red")

estplot = createdisplay(1,1)
show(estplot,1,1,estml,confipl,confinl)
setgopt(estplot,1,1,"title","Estimation mean and 99% confidence interval - HYP","border",0)



05.11.2004