%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Estimate the Hedge Ratio of MRS models for Four pairs spot &
%% futures, with 500-winodw, one-step out-sample forecast, weekly data
%% MRS-OLS, MRS-VAR, MRS-ECM, MRS-SUR models with MLE method
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear
future=xlsread('C:\Research_Weijun_Xu\Research\Hedging_MRS\data\future_week.xls');  % Weekly data sets
spot=xlsread('C:\Research_Weijun_Xu\Research\Hedging_MRS\data\spot_week.xls');
forecast_step=1;
T=size(spot,1);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Let's first orginize the return data sets
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for index=1:4
    for i=1:T-500
        ret=[100*diff(log(spot(:,index))), 100*diff(log(future(:,index)))]; % return of spot and future
        sample_ret{index,i}=ret(i:i+500-1,:);
    end
end

% Let's get the Cointegration term for MRS-ECM model use
for index=1:4
    for i=1:T-500
        log_spot=100*log(spot(i:i+500-1,index));         
        log_future=100*log(future(i:i+500-1,index)); 
        Cointegrated_coeff=[ones(size(log_future)),log_future]\[log_spot];
        Cointegrated_error{index,i}=log_spot-[ones(size(log_future)),log_future]*Cointegrated_coeff;                
    end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Model 4: MRS-ECM
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for index=1:4 % index=1~4,represents: SP500, FTSE100, Nikk225, Hengsheng
    for i=1:ceil((T-500)/forecast_step)
        y=sample_ret{index,(i-1)*forecast_step+1}(2:end,:);
        x=[ones(size(y,1),1), Cointegrated_error{index,(i-1)*forecast_step+1}(1:end-1,:),...
            sample_ret{index,(i-1)*forecast_step+1}(1:end-1,:)];

        % y=[Rs_t, Rf_t], x=[1 Cointegrated_error_t-1, Rs_t-1, Rf_t-1]
        tic
        [parameters, LLF, smooth_t,filter_t, cov_matrix]=MLE_RS_Multi_intercept_VAR_fit(y, x);
        % Since we use two stage to estimate the model, first calculate the
        % cointegrated error term first, then use MLE, thus just treat
        % error correction term as new independent variable, we still can  use 'VAR_fit fucntion'
        hedge0=cov_matrix(1,2,1)/cov_matrix(2,2,1);
        hedge1=cov_matrix(1,2,2)/cov_matrix(2,2,2);

        P=[parameters(end-1), 1-parameters(end-1); 1-parameters(end), parameters(end)]; % Transistion probability Matrix
        state_probability(1,:)=filter_t(end,:)*P;
        if forecast_step<=1
            % when forecast=1

            hedge_MRS_ECM_MLE(i,index)=state_probability(1,:)*[hedge1,hedge0]';
        elseif forecast_step>1
            % when forecast>1, we use state_probability(1,:) as new filter probability
            for j=2:forecast_step
                state_probability(j,:)=state_probability(j-1,:)*P;
            end
            hedge_MRS_ECM_MLE(i,index)=state_probability(end,:)*[hedge1,hedge0]';
        end
        fprintf('\n *********************************************************************** \n');
        save('C:\Research_Weijun_Xu\Research\Hedging_MRS\results\hedge_MRS\hedge_MRS_ECM_MLE', 'hedge_MRS_ECM_MLE');
        fprintf('\n Time cost for the Multivarate MRS-ECM model in sub_sample %3g  is: %3g \n', index, toc)
    end
end

clear parameters  filter_t  state_probability y x hedge1 hedge0 P

forecast_step=4;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Model 4: MRS-ECM
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for index=1:4 % index=1~4,represents: SP500, FTSE100, Nikk225, Hengsheng
    for i=1:ceil((T-500)/forecast_step)
        y=sample_ret{index,(i-1)*forecast_step+1}(2:end,:);
        x=[ones(size(y,1),1), Cointegrated_error{index,(i-1)*forecast_step+1}(1:end-1,:),...
            sample_ret{index,(i-1)*forecast_step+1}(1:end-1,:)];

        % y=[Rs_t, Rf_t], x=[1 Cointegrated_error_t-1, Rs_t-1, Rf_t-1]
        tic
        [parameters, LLF, smooth_t,filter_t, cov_matrix]=MLE_RS_Multi_intercept_VAR_fit(y, x);
        % Since we use two stage to estimate the model, first calculate the
        % cointegrated error term first, then use MLE, thus just treat
        % error correction term as new independent variable, we still can  use 'VAR_fit fucntion'
        hedge0=cov_matrix(1,2,1)/cov_matrix(2,2,1);
        hedge1=cov_matrix(1,2,2)/cov_matrix(2,2,2);

        P=[parameters(end-1), 1-parameters(end-1); 1-parameters(end), parameters(end)]; % Transistion probability Matrix
        state_probability(1,:)=filter_t(end,:)*P;
        if forecast_step<=1
            % when forecast=1

            hedge_MRS_ECM_MLE_4(i,index)=state_probability(1,:)*[hedge1,hedge0]';
        elseif forecast_step>1
            % when forecast>1, we use state_probability(1,:) as new filter probability
            for j=2:forecast_step
                state_probability(j,:)=state_probability(j-1,:)*P;
            end
            hedge_MRS_ECM_MLE_4(i,index)=state_probability(end,:)*[hedge1,hedge0]';
        end
        fprintf('\n *********************************************************************** \n');
        save('C:\Research_Weijun_Xu\Research\Hedging_MRS\results\hedge_MRS\hedge_MRS_ECM_MLE_4', 'hedge_MRS_ECM_MLE_4');
        fprintf('\n Time cost for the Multivarate MRS-ECM model in sub_sample %3g  is: %3g \n', index, toc)
    end
end

clear parameters  filter_t  state_probability y x hedge1 hedge0 P














