function[quantile, dq] = CAViaROptimisation(y, THETA, C)
% ****************************************************************************************************************************************
% * *
% * Codes for the paper "Measuring Comovements by Regression Quantiles" by Lorenzo Cappiello, *
% * Bruno Gerard and Simone Manganelli *
% * *
% * By SIMONE MANGANELLI, European Central Bank, Frankfurt am Main. *
% * Created on 26 May, 2004. Last modified 13 April 2005. *
% * *
% ****************************************************************************************************************************************
%
%
%************ INPUTS *********************
%
% y = (T,1) time series of returns
% THETA = confidence level of the Value at Risk
% C = Crisis dummy
%
%************ OUTPUT *********************
%
% quantile = (T,1)-vector of estimated quantiles
% dq = vector of quantile derivatives
%
%*****************************************************************************************
% *****************************************************************************************
% Set parameters for optimisation.
% *****************************************************************************************
T = length(y);
%REP = 5; % Number of times the optimization algorithm is repeated.
%nInitialVectors = [1000, 3]; % Number of initial vector fed in the uniform random number generator for AS model.
%nInitialVectors = [1, 5]; % Number of initial vector fed in the uniform random number generator for AS model.
nInitialCond = 5;% Select the number of initial conditions for the optimisation.
MaxFunEvals = 100;% Parameters for the optimisation algorithm. Increase them in case the algorithm does not converge.
MaxIter = 100;
options = optimset('LargeScale','off','HessUpdate','dfp','LineSearchType','quadcubic','MaxFunEvals', MaxFunEvals,...
'display','off','MaxIter', MaxIter,'TolFun', 1e-8,'TolX', 1e-8);
warningoff
% Compute the empirical THETA-quantile for y (the in-sample vector of observations).
%empiricalQuantile = ysort(round(300*THETA));
empiricalQuantile = prctile(y(1:300), THETA*100);
eps = 1e-10;
%
%
%**************************** Optimization Routine ******************************************
initialTargetVectors = [unifrnd(-.1, .1, nInitialCond, 3), unifrnd(.7, .99, nInitialCond, 1), unifrnd(-.1, .1, nInitialCond, 1)];
b=[];
fori = 1:nInitialCond
b0 = initialTargetVectors(i,:)'; b0 = b0; b1 = b0+3*eps;
RQ1 = 1e10;
whilenorm(b1-b0)>eps
b0 = b1;%RQ0 = RQ1;
b1 = fminsearch('CAViaR', b0, [], y, C, THETA, 1);
end
q = CAViaR(b1, y, C, THETA, 0);%q = CAViaR(y, b0, C, .05);
RQ1 = abs(THETA-(y
b = [b;RQ1/100, b1'];
end
b
[aa,a] = min(b(:,1));
b1 = b(a,2:end)';
%************************** Compute variables that enter output *****************************
% Compute VaR and Hit for the estimated parameters of RQ.
quantile = CAViaR(b1, y, C, THETA, 0);
%return
% Compute gradient
dq1 = zeros(T,1); dq2 = dq1; dq3 = dq1; dq4 = dq1; dq5 = dq1; q = dq1;
fori = 3:T
q(i) = b1(1) + b1(2)*C(i) + b1(3)*y(i-1) + b1(4)*q(i-1) - b1(3)*b1(4)*y(i-2) + b1(5)*abs(y(i-1));
dq1(i,1) = 1 + b1(4) * dq1(i-1,1);
dq2(i,1) = C(i) + b1(4) * dq2(i-1,1);
dq3(i,1) = y(i-1) + b1(4) * dq3(i-1,1) - b1(4)*y(i-2);
dq4(i,1) = q(i-1) + b1(4)*dq4(i-1,1) - b1(3)*y(i-2);
dq5(i,1) = abs(y(i-1)) + b1(4) * dq5(i-1,1);
end
dq = [dq1, dq2, dq3, dq4, dq5];