1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151
| function results=varin(data)
% data is a column vector of the returns defined as ln(p(t)/p(t-1))*100
% or (p(t)-p(t-1))/p(t-1)*100
% the model is,
% r(t)=m+e(t)
% s(t)^2=v+a*e(t-1)^2+b*s(t-1)^2
% results, is a structure of the constant in mean (m), constant in variance
% (v), arch term (a), garch term (b), the maximum likelihood estimation
% (mle), and a column vector of the conditional standard deviation (sigma).
alpha=[.95;.975;.99;.995;.9975;.999]; % confidence intervals
% matrices for the optimization routine. All the inequalities are
% A*x<=b
Coeff.A=[0 -1 0 0;0 0 -1 0;0 0 0 -1;0 0 1 1;0 0 1 0;0 0 0 1];
Coeff.b=[0;0;0;1;1;1];
% Initialization of the parameter vector
Coeff.c=[0.1;0.1;0.1;0.9];
%mean equation
m=a*(data(t-1)/data(t-1)^2)-d*data(t-1)
% maximum likelihood estimation of Garch(1,1) and Gaussian distribution
function f=garch11nmle(x)
r2=(data-x(1)).^2;
sigma2=zeros(length(data),1);
sigma2(1)=var(data);
for index1=2:length(data)
sigma2(index1)=x(2)+x(3)*r2(index1-1)+x(4)*sigma2(index1-1);
end
f=zeros(length(data),1);
for index2 =1:length(data)
% matlab minimizes the mle so we need an extra minus
f(index2)=-(-1/2*log(2*pi)-1/2.*log(x(2)+x(3)*r2(index2)+x(4)*sigma2(index2))-1/2*r2(index2)/sigma2(index2));
end
f=sum(f);
end
% Optimization
options = optimset('Algorithm','interior-point','Display','off');
[x,fval] = fmincon(@garch11nmle,Coeff.c,Coeff.A,Coeff.b,a,d,[],[],[],[],[],options);
Sigmas=sqrt(sigma2); % standard deviations
% Output structure
results.m=x(1);
results.v=x(2);
results.arch=x(3);
results.garch=x(4);
results.mle=-fval;
results.sigmas=Sigmas;
% Plot of the Value at Risk levels
for index3=1:length(alpha);
varL=x(1)-norminv(alpha(index3)).*Sigmas;
varS=x(1)-norminv(1-alpha(index3)).*Sigmas;
figure
plot(data,'or');
hold on
plot(varL,'-g');
plot(varS,'-b');
hold off
xlabel('Time');
ylabel('VaR and returns');
title(['VaR in sample for alpha=',num2str(alpha(index3),'%6.4f')]);
end
% In sample VaR backtesting
disp([' ','IN SAMPLE VALUE AT RISK BACKTESTING']);
disp('-------------------------------------------------------------------------------------------------------');
disp([' ','SHORT POSITION' ]);
disp('-------------------------------------------------------------------------------------------------------');
disp(['alpha',' ','Success rate',' ','LR value',' ','Kupiec-p',' ','LR Christ',' ','Christ']);
% Short position
for index4=1:length(alpha);
varS=x(1)-norminv(1-alpha(index4)).*Sigmas;
% Kupiec PF proportion of failures test
pval=1-alpha(index4);
Tval=length(data);
xvalS1=(data>varS);
xvalS=sum(xvalS1);
LR_S=-2*(xvalS*log(pval)+(Tval-xvalS)*log(1-pval)-xvalS*log(xvalS/Tval)-(Tval-xvalS)*log(1-xvalS/Tval));
Kupiec_PFS=1-chi2cdf(LR_S,1);
% Christoffersen independence test for conditional coverage
n00=0; n01=0; n10=0; n11=0;
for index5=1:length(xvalS1)-1
if xvalS1(index5)==0 && xvalS1(index5+1)==0
n00=n00+1;
elseif xvalS1(index5)==0 && xvalS1(index5+1)==1
n01=n01+1;
elseif xvalS1(index5)==1 && xvalS1(index5+1)==0
n10=n10+1;
elseif xvalS1(index5)==1 && xvalS1(index5+1)==1;
n11=n11+1;
end
end
p0=n01/(n00+n01);
p1=n11/(n10+n11);
p=(n01+n11)/(n00+n01+n10+n11);
LR_ITS=-2*(n00+n10)*log(1-p)-2*(n01+n11)*log(p)+2*n00*log(1-p0)+2*n01*log(p0)+2*n10*log(1-p1)+2*n11*log(p1);
% LR_ITS=-2*log((1-p)^(n00+n10)*p^(n01+n11)/((1-p0)^n00*p0^n01*(1-p1)^n10*p1^n11));
CIT_S=1-chi2cdf(LR_ITS,1);
fprintf('%6.4f %10.4f %12.4f %14.4e %11.4f %13.4f',alpha(index4),1-xvalS/Tval,LR_S,Kupiec_PFS,LR_ITS,CIT_S);
fprintf('\n');
% disp([num2str(alpha(index4),'%10.4f'),' ',num2str(1-xvalS/Tval,'%10.6f'),' ',num2str(LR_S,'%10.6f'),...
% ' ',num2str(Kupiec_PFS,'%10.6e'),' ',num2str(LR_ITS,'%10.6f'),' ',num2str(CIT_S,'%10.6f')]);
end
disp('-------------------------------------------------------------------------------------------------------');
disp([' ','LONG POSITION' ]);
disp('-------------------------------------------------------------------------------------------------------');
disp(['alpha',' ','Failure rate',' ','LR value',' ','Kupiec-p',' ','LR Christ',' ','Christ']);
% Long position
for index6=1:length(alpha);
varL=x(1)-norminv(alpha(index6)).*Sigmas;
% Kupiec PF proportion of failures test
pval=1-alpha(index6);
Tval=length(data);
xvalL1=(data<varL);
xvalL=sum(xvalL1);
%LR_L=-2*log((pval^xvalL*(1-pval)^(Tval-xvalL))/((xvalL/Tval)^xvalL*(1-xvalL/Tval)^(Tval-xvalL)));
LR_L=-2*(xvalL*log(pval)+(Tval-xvalL)*log(1-pval)-xvalL*log(xvalL/Tval)-(Tval-xvalL)*log(1-xvalL/Tval));
Kupiec_PFL=1-chi2cdf(LR_L,1);
% Christoffersen independence test for conditional coverage
n00=0; n01=0; n10=0; n11=0;
for index7=1:length(xvalL1)-1
if xvalL1(index7)==0 && xvalL1(index7+1)==0
n00=n00+1;
elseif xvalL1(index7)==0 && xvalL1(index7+1)==1
n01=n01+1;
elseif xvalL1(index7)==1 && xvalL1(index7+1)==0
n10=n10+1;
elseif xvalL1(index7)==1 && xvalL1(index7+1)==1;
n11=n11+1;
end
end
p0=n01/(n00+n01);
p1=n11/(n10+n11);
p=(n01+n11)/(n00+n01+n10+n11);
LR_ITL=-2*(n00+n10)*log(1-p)-2*(n01+n11)*log(p)+2*n00*log(1-p0)+2*n01*log(p0)+2*n10*log(1-p1)+2*n11*log(p1);
% LR_ITL=-2*log((1-p)^(n00+n10)*p^(n01+n11)/((1-p0)^n00*p0^n01*(1-p1)^n10*p1^n11));
CIT_L=1-chi2cdf(LR_ITL,1);
fprintf('%6.4f %10.4f %12.4f %14.4e %11.4f %13.4f',alpha(index6),1-xvalL/Tval,LR_L,Kupiec_PFL,LR_ITL,CIT_L);
fprintf('\n');
end
end |
Partager