%% Multivariate regression between r* and macro-variables

clc
clear
close all

load('USrstar.mat')
load('UKrstar.mat')
load('GMrstar.mat')

load('USdata.mat')
load('UKdata.mat')
load('GMdata.mat')

% usrstar = mean(USrstar,1)';  % For US the sample is from 197110-202204
usrstar = mean(USrstar,1)';  % For US the sample is from 197111-202204
ukrstar = mean(UKrstar,1)';  % For UK the sample is from 198301-202203
gmrstar = mean(GMrstar,1)';  % For GM the sample is from 197210-202203

% n1 = size(usrstar,1)-1;
n1 = size(usrstar,1);
n2 = size(ukrstar,1);
n3 = size(gmrstar,1);

usr = nan(n1/3,1);  % For US the sample is from 1971Q4-2022Q1
ukr = nan(n2/3,1);  % For UK the sample is from 1983Q1-2022Q1
gmr = nan(n3/3,1);  % For GM the sample is from 1972Q4-2022Q1

% for i = 1:n1
%     if mod(i,3)==0
%          % usr(i/3) = mean(usrstar(i-2:i));
%         usr(i/3) = usrstar(i);
%     end
% end

usrstar = [nan;usrstar];  % LW data starts at 197111
n1 = n1+1;
for i = 1:n1
    if mod(i,3)==0
         % usr(i/3) = mean(usrstar(i-2:i));
        usr(i/3) = usrstar(i);
    end
end
n1 = n1-1;

for i = 1:n2
    if mod(i,3)==0
        % ukr(i/3) = mean(ukrstar(i-2:i));
        ukr(i/3) = ukrstar(i);
    end
end

for i = 1:n3
    if mod(i,3)==0
         % gmr(i/3) = mean(gmrstar(i-2:i));
        gmr(i/3) = gmrstar(i);
    end
end


xdate = datetime("1971-10-01") + calmonths(0:3:605)';
rstars = table(xdate);
rstars.US = usr;
rstars.UK = nan(n1/3,1);
rstars.GM = nan(n1/3,1);
rstars.UK(46:end) = ukr;
rstars.GM(5:end) = gmr;

% save("rsatr_Qend","rstars");
% save("rsatr_Qend_LW.mat","rstars");
save("rsatr_Qend_LW_RR2.mat","rstars");

%% Annual mean of r*

% usrstar = usrstar(2:end);
tempus = usrstar(3:end-4);
tempuk = ukrstar(1:end-3);
tempgm = gmrstar(4:end-3);
n1 = length(tempus)/12;
n2 = length(tempuk)/12;
n3 = length(tempgm)/12;
usr = nan(n1,1);  % For US the sample is from 1972-2021
ukr = nan(n2,1);  % For UK the sample is from 1983-2021
gmr = nan(n3,1);  % For GM the sample is from 1973-2021


for i = 1:length(tempus)
    if mod(i,12)==0
        usr(i/12) = mean(tempus(i-11:i));
        % usr(i/12) = tempus(i);
    end
end

for i = 1:length(tempuk)
    if mod(i,12)==0
        ukr(i/12) = mean(tempuk(i-11:i));
        % ukr(i/12) = tempuk(i);
    end
end

for i = 1:length(tempgm)
    if mod(i,12)==0
        gmr(i/12) = mean(tempgm(i-11:i));
        % gmr(i/12) = tempgm(i);
    end
end

xdate = datetime("1972-1-01") + calyears(0:1:49)';
rstars = table(xdate);
rstars.US = usr;
rstars.UK = nan(n1,1);
rstars.GM = nan(n1,1);
rstars.UK(12:end) = ukr;
rstars.GM(2:end) = gmr;
% 
% save("rsatr_Ymean","rstars");
 % save("rsatr_Yend","rstars");
 % save("rsatr_Ymean_LW","rstars");
 save("rsatr_Ymean_LW_rr2","rstars");


%%
% GDP growth  1991Q2-2022Q1
rgdp = [USvars.rgdp(125:end-3) UKvars.rgdp(125:end-3) GMvars.rgdp(125:end-3)];
macro_var = rgdp;
rs = rstars{79:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_rgdp,~,~,CovB_rgdp] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_rgdp(4)/sqrt(CovB_rgdp(4,4)) ;
tboots_rgdp = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_rgdp)  % 0.23

adftest(rgdp(:,1),Model="ARD",Lags=0:4)
adftest(rgdp(:,2),Model="ARD",Lags=0:4)
adftest(rgdp(:,3),Model="ARD",Lags=0:4)

%% consumption growth
consump = [USvars.consump(141:end-3) UKvars.csonsump(141:end-3) GMvars.csonsump(141:end-3)];
macro_var = consump;
rs = rstars{95:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_consump,~,~,CovB_consump] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_consump(4)/sqrt(CovB_consump(4,4)) ;
tboots_consump = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_consump)   % 0.38

%% variance of gdp growth 
var_gdp = [USvars.VarGDP(145:end-3) UKvars.VarGDP(145:end-3) GMvars.VarGDP(145:end-3)];
macro_var = var_gdp;
rs = rstars(99:end,2:4);
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_vargdp,~,~,CovB_vargdp] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_vargdp(4)/sqrt(CovB_vargdp(4,4)) ;
tboots_vargdp = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_vargdp)   % 0.38

adftest(var_gdp(:,1),Model="ARD",Lags=0:4)
adftest(var_gdp(:,2),Model="ARD",Lags=0:4)
adftest(var_gdp(:,3),Model="ARD",Lags=0:4)

%% variance of GDP inflaton 1983Q3 - 2022Q1
var_pi = [USvars.VarPi(92:end-3) UKvars.VarPi(92:end-3) GMvars.VarPi(92:end-3)];
macro_var = var_pi;
rs = rstars{46:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_varpi,~,~,CovB_varpi] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_varpi(4)/sqrt(CovB_varpi(4,4)) ;
tboots_varpi = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_varpi)   % 0.44

%% variance of market return
var_mkt = [USvars.VarMkt(118:end-3) UKvars.VarMkt(118:end-3) GMvars.VarMkt(118:end-3)];
macro_var = var_mkt;
rs = rstars{72:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_varmkt,~,~,CovB_varmkt] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_varmkt(4)/sqrt(CovB_varmkt(4,4)) ;
tboots_varmkt = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_varmkt) % 0.59

%% covariance of pi and GDP
cov_gdp = [USvars.CovGDP(145:end-3) UKvars.CovGDP(145:end-3) GMvars.CovGDP(145:end-3)];
macro_var = cov_gdp*100;
rs = rstars(99:end,2:4);
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_covgdp,~,~,CovB_covgdp] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_covgdp(4)/sqrt(CovB_covgdp(4,4)) ;
tboots_covgdp = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_covgdp) 

%% covariance of pi and market return
cov_mkt = [USvars.CovMkt(118:end-3) UKvars.Covmkt(118:end-3) GMvars.Covmkt(118:end-3)];
macro_var = cov_mkt;
rs = rstars{72:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_covmkt,~,~,CovB_covmkt] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_covmkt(4)/sqrt(CovB_covmkt(4,4)) ;
tboots_covmkt = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_covmkt) % 0.59

%% investment growth   1983Q1-2022Q1
cap_growth = [USvars.Capital(92:end-3) UKvars.Capital(92:end-3) GMvars.Capital(92:end-3)];
macro_var = cap_growth;
rs = rstars{46:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_invest,~,~,CovB_invest] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_invest(4)/sqrt(CovB_invest(4,4)) ;
tboots_invest = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_invest)  % 0.56

%% debt growth 1991Q1-2022Q1
debt_growth = [USvars.Debt(125:end-3) UKvars.Debt(125:end-3) GMvars.Debt(125:end-3)];
macro_var = debt_growth;
rs = rstars{79:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_debt,~,~,CovB_debt] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_debt(4)/sqrt(CovB_debt(4,4)) ;
tboots_debt = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_debt) % 0.67

%% current account   1991Q1-2022Q1
current =  [USvars.Current(125:end-3) UKvars.Current(125:end-3) GMvars.Current(125:end-3)];
macro_var = current;
rs = rstars{79:end,2:4};
macro_var = diff(current);
rs = rstars{80:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_curr,~,~,CovB_curr] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_curr(4)/sqrt(CovB_curr(4,4)) ;
tboots_curr = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_curr)  %0.01

%% realtive price of investment against consumption price   2009Q1-2022Q1
rel_price = [USvars.Price(196:end-3) UKvars.Price(196:end-3) GMvars.Price(196:end-3)];
macro_var = rel_price;
rs = rstars{150:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_price,~,~,CovB_price] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_price(4)/sqrt(CovB_price(4,4)) ;
tboots_price = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_price)  % 0.67

%% yield spread   1986Q3-2022Q1
Spread = [USvars.Spread(104:end-3) UKvars.Spread(104:end-3) GMvars.Spread(104:end-3)];
macro_var = Spread;
rs = rstars{58:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_spread,~,~,CovB_spread] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_spread(4)/sqrt(CovB_spread(4,4)) ;
tboots_spread = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_spread)  % 0.37

%% GDP inflation

%% M1 growth  1987Q2-2022Q1
M1 = [USvars.M1(109:end-3) UKvars.M1(109:end-3) GMvars.M1(109:end-3)];
macro_var = M1;
rs = rstars{63:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_m1,~,~,CovB_m1] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_m1(4)/sqrt(CovB_m1(4,4)) ;
tboots_m1 = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_m1)  % 0.24

%% M3 growth 1987Q2-2022Q1
M3 = [USvars.M3(109:end-3) UKvars.M3(109:end-3) GMvars.M3(109:end-3)];
macro_var = M3;
rs = rstars{63:end,2:4};
n = size(rs,1);
myx = cell(n,1);
for i = 1:n
    myx{i} = [eye(3), macro_var(i,:)'];
end
[b_m3,~,~,CovB_m3] = mvregress(myx,rs, 'algorithm','mvn');
tvalue = b_m3(4)/sqrt(CovB_m3(4,4)) ;
tboots_m3 = panel_bootstrap(rs, macro_var, 5000, 6);
findprctile(tvalue,tboots_m3)



















