
%% prepare data
r=length(fblks);
mqp=max(qp);
[yyy,~,~,~,~,~,T_est,xxx]=prep_data_TV_VAR(Y_est,mqp);  %y is nxT
x_bul=reshape(xxx(2:end,:),[n,mqp,T_est]);

Y_actual=Y_est(:,mqp+1:length(Y_est));

lagsMA=30;

%%

structdecomp=struct('totresids',[],'Omegasdiag',[],'CHISdiag',[],'comresids',[],...
    'diagcomresidsVar',[],'idioresids',[],'diagidioresidsVar',[],...
    'Y_ComComp' ,[], 'RsqCom'  ,[] ,'Y_Idio'  ,[],'LocalProjComm',[]);

tic
parfor js=1:l_c
    %for js=1:l_c  % if you do not have the Parallel toolbox
    
    structdecomp(js)= ...
        make_resdecomp_MAIARSV(Y_actual,x_bul,b0_tildes(:,js),as(:,js),gammas(:,js),gs(:,js),sigmas(:,:,js),fblks,lagsMA);
end
toc

%%

totresids=cat(3,structdecomp.totresids);
Omegasdiag=cat(3,structdecomp.Omegasdiag);
CHISdiag=cat(3,structdecomp.CHISdiag);

comresids=cat(3,structdecomp.comresids);
diagcomresidsVar=cat(3,structdecomp.diagcomresidsVar);

idioresids=cat(3,structdecomp.idioresids);
diagidioresidsVar=cat(3,structdecomp.diagidioresidsVar);

share_aggOmegaVar=permute(cat(4,diagcomresidsVar,diagidioresidsVar)./repmat(Omegasdiag,[1,1,1,2]),[1,2,4,3]);

Y_ComComp=cat(3, structdecomp.Y_ComComp);
Y_Idio=cat(3, structdecomp.Y_Idio);
RsqCom=cat(2, structdecomp.RsqCom);
LocalProjComm=cat(3,structdecomp.LocalProjComm);

mTresids=median(totresids,3);
mcomresids=median(comresids,3);
mIDIOresids=median(idioresids,3);
mOmegasdiag=median(Omegasdiag,3);
mCHISdiag=median(CHISdiag,3);
mdiagcomresidsVar=median(diagcomresidsVar,3);
mdiagidiovar=median(diagidioresidsVar,3);
mshare_aggOmegaVar=mean(share_aggOmegaVar,4);

mY_ComComp=median(Y_ComComp,3);
mY_Idio=median(Y_Idio,3);
mRsqCom=median(RsqCom,2);

qOmegasdiag=quantile(Omegasdiag,quants,3);
qCHISdiag=quantile(CHISdiag,quants,3);
qLocalProjComm=quantile(LocalProjComm,quants,3);

Test_minusp=length(Y_actual);
dates_Testmp=(dateslab(end-Test_minusp+1:end));


clear structdecomp



