clc;
clear;

%Load data:
load('.\problem_var_and_cvar_regression_data.mat')

count = 1;
iargstruc_arr(count) = matrix_pack('matrix_s',H{1},[],c{1});

%Define problem statement for CVaR regression:
problem_statement_cvar_2 = sprintf('%s\n',...
'minimize',...
'  cvar2_err(0.75,matrix_s)',...
'Value:',...
'  L(matrix_s)',...
' ');

%Optimize problem using mpsg_solver function:
[solution_str_2, outargstruc_arr_2] = mpsg_solver(problem_statement_cvar_2, iargstruc_arr);

tail_los_full_cvar = outargstruc_arr_2(1).values(:,2);
intercept_cvar = outargstruc_arr_2(2).values(1);

tail_los_cvar = tail_los_full_cvar+intercept_cvar;

cvar_0_75 = functionvalue('cvar_risk_g', 0.75, tail_los_cvar, [], [], 1);

rep_text = 'CVaR regression with parameter alpha = 0.75';
rep_text = sprintf('%s \n intercept = %f, CVaR = %f',rep_text,intercept_cvar,cvar_0_75);


var_075 = functionvalue('var_risk_g', 0.75, tail_los_cvar, [], [], 1);
loss_null_var = tail_los_cvar - var_075;
loss_positive = loss_null_var(loss_null_var>=0);

%Plot histogram of tail of loss of the CVaR regression
h2 = figure;
histogram(loss_positive,30,'FaceColor',[0.5,0.5,0.5])

%Estimate parameters of GPD
params_cvar = tsallis_harmonik_params_loss(loss_positive);

%Report for parameters estimation of GPD for tail of loss of the CVaR
%regression
rep_text = sprintf('%s\n Mu = %f\n',rep_text,params_cvar(1).mua);
rep_text = sprintf('%s Kappa MLE = %f\n',rep_text,params_cvar(1).kappa);
rep_text = sprintf('%s Kappa Harmonic = %f',rep_text,params_cvar(2).kappa);

disp('    ')
disp(rep_text)