%% Old plots from MasterOneNorm
%% Table Rankings
clear bumpRanks sortings NewSortings outDeg sortings outDegOff outDegDiag namestring
k=0;
TITLE = {'1972 - 1983','1984 - 1995', '1996 - 2007'};
for t=ts
    k=k+1;
    weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));
    FEVD = DfinalReg(:,:,t)-diag(diag(DfinalReg(:,:,t)));
    FEV = DfinalRegIRF(:,:,t)-diag(diag(DfinalRegIRF(:,:,t)));
    FEV2 = diag(weightComponent)*DfinalRegFEV(:,:,t)-diag(diag(DfinalRegFEV(:,:,t)))*diag(weightComponent);
    
    sorted = sort(FEVD(:,2),'descend');
    sorted2 = sort(FEV(:,2),'descend');
    sorted3 = sort(FEV2(:,2),'descend');
    
    for i = 1:7
        [row, col] = find(sorted(i)==FEVD);
        FromName = convertStringsToChars(names{col}); FromName = FromName(5:end-16);
        ToName = convertStringsToChars(names{row});ToName = ToName(5:end-16);
        namestring{k,i} = {['From ' FromName ' to ' ToName]};
        topvalues(k,i) = FEVD(row,col);
        [row, col] = find(sorted2(i)==FEV);
        row = row(1); col = col(1);
        FromName = convertStringsToChars(names{col}); FromName = FromName(5:end-16);
        ToName = convertStringsToChars(names{row});ToName = ToName(5:end-16);
        topvalues2(k,i) = FEV(row,col);
        namestring2{k,i} = {['From ' FromName ' to ' ToName]};
        namestring2{k,i} = ToName;
        [row, col] = find(sorted3(i)==FEV2);
        row = row(1); col = col(1);
        FromName = convertStringsToChars(names{col}); FromName = FromName(5:end-16);
        ToName = convertStringsToChars(names{row});ToName = ToName(5:end-16);
        topvalues3(k,i) = FEV2(row,col);
        namestring3{k,i} = {['From ' FromName ' to ' ToName]};
    end
end
TableString= ['\begin{tabular}{l r r}'  ' & & From To &  Value \\' '\cline{2-3} '];
for k=1:3
    for i = 1:3
        if i==1
            TableString = [TableString TITLE{k}];
        end
        TableString= [TableString ' & '  char(namestring{k,i}) '& ' num2str(round(topvalues(k,i),2)) ' \\'];
        if i==3
            TableString = [TableString '[0.4cm]'];
        end
    end
end
TableString = [TableString ' \bottomrule \end{tabular}'];
disp(TableString)
TableString2= ['\begin{tabular}{l r r}'  ' & From ' FromName ' To &  Value \\' '\cline{2-3} '];
for k=1:3
    for i = 1:3
        if i==1
            TableString2 = [TableString2 TITLE{k}];
        end
        TableString2= [TableString2 ' & '  char(namestring2{k,i}) '& ' num2str(round(topvalues2(k,i),2)) ' \\'];
        if i==3 && k<3
            TableString2 = [TableString2 '[0.4cm]'];
        end
    end
end
TableString2 = [TableString2 ' \bottomrule \end{tabular}'];
disp(TableString2)
TableString3= ['\begin{tabular}{l l r r}'  ' & & From To &  Value \\' '\cline{3-4}' ...
    '1972 - 1983 & Top 1 & ' char(namestring3{1,1}) '& ' num2str(topvalues3(1,1)) ' \\ ' ...
    ' & Top 2& ' char(namestring3{1,2})  '& ' num2str(topvalues3(1,2)) ' \\' ...
    ' & Top 3& ' char(namestring3{1,3})   '& ' num2str(topvalues3(1,3)) ' \\[0.4cm]'...
    '1984-1995 & Top 1 & ' char(namestring3{2,1})  '& ' num2str(topvalues3(2,1)) ' \\' ...
    '  & Top 2&  ' char(namestring3{2,2}) '& ' num2str(topvalues3(2,2)) ' \\ '...
    '  & Top 3& ' char(namestring3{2,3})   '& ' num2str(topvalues3(2,3)) ' \\[0.4cm]'...
    '1996-2007 & Top 1 & ' char(namestring3{3,1})  '& ' num2str(topvalues3(3,1)) ' \\'...
    '  & Top 2 & '  char(namestring3{3,2}) '& ' num2str(topvalues3(3,2)) ' \\'...
    ' & Top 3 & ' char(namestring3{3,3})  '& ' num2str(topvalues3(3,3)) ' \\ \bottomrule \end{tabular}'];
disp(TableString3)
%% TableIRV long
clear bumpRanks sortings NewSortings outDeg sortings outDegOff outDegDiag outDegSingleFact
k=0;
for t=ts
    k=k+1;
    weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));
    weightComponent2 = weightComponent; weightComponent2(:) = 1/N;
    weightComponent3 = weightComponent-1/N;
    weightInnovations = diag(FinalSigmaReg(:,:,t)).^(-0.5);ones(size(weightComponent));
%     weightInnovations = weightInnovations./sum(weightInnovations);
    weightsFactors = diag(kron(eye(13),covFactor(:,:,t))).^(-0.5);

    Table = (1/9)*DfinalRegIRF(:,:,t);
    TableFact = (1/9)*DfinalRegIRFfactor(:,:,t);
    TableFactSpillovers = (1/9)*DfinalRegIRFfactorSpillovers(:,:,t);
        
    outDeg(k) = weightComponent'*Table*weightInnovations;
    outDegOff(k)= weightComponent'*(Table-diag(diag(Table)))*weightInnovations;
    outDegOff2(k)=weightComponent2'*(Table-diag(diag(Table)))*weightInnovations;
    outDegOff3(k)=outDegOff(k)-outDegOff2(k);
    

    outDegFact(k) = weightComponent'*TableFact*weightsFactors;
    outDegSingleFact(:,k) = weightComponent'*TableFact*kron(ones(FactNumLags+1,1),diag(weightsFactors(1:3)));
    outDegFactSpillover(k) = weightComponent'*TableFactSpillovers*weightsFactors;
    
    outDegDiag(k)=weightComponent'*(diag(diag(Table)))*weightInnovations;
    outDegDiag2(k)= weightComponent2'*(diag(diag(Table)))*weightInnovations;
    outDegDiag3(k)=outDegDiag(k)-outDegDiag2(k);
    
    SymmetricTable = sign(sign(Table')+sign(Table)).*(min(abs(Table'),abs(Table)))-diag(diag(Table));%DfinalRegFEV(:,:,t)-AsymmetricTable-diag(diag(DfinalRegFEV(:,:,t)));
    AsymmetricTable = Table-SymmetricTable-diag(diag(Table));min(DfinalRegIRF(:,:,t)',DfinalRegIRF(:,:,t));
    AsymetriCon(k)=  weightComponent'*AsymmetricTable*weightInnovations;
    SymetriCon(k)=  weightComponent'*SymmetricTable*weightInnovations;
    AsymetriCon2(k)=weightComponent2'*AsymmetricTable*weightInnovations;
    SymetriCon2(k)=weightComponent2'*SymmetricTable*weightInnovations;
    AsymetriCon3(k)=AsymetriCon(k)-AsymetriCon2(k);
    SymetriCon3(k)=SymetriCon(k)-SymetriCon2(k);
    laggedCon(k) = weightComponent'*(Table - (1/3)*FinalSigmaReg(:,:,t))*weightInnovations;
    contempCon(k) = weightComponent'*((1/3)*FinalSigmaReg(:,:,t))*weightInnovations;

end
TableString= ['\begin{tabular}{l l r r r}'  ' & & {1972-1983} &  {{1984-1995}} & {{1996-2007}} \\' '\cline{3-5}' ' \midrule ' ...
    'Data variance & $v^{IP}=var(w' char(39) 'y_t^q)$ &' num2str(round(sampleVar(1),1)) '& ' num2str(round(sampleVar(2),1))  '& ' num2str(round(sampleVar(3),1)) ' \\ ' ...%  'Innovation Variance & $v^{inn}=w' char(39) 'diag(FEV)w$ & ' num2str(round(InnovDiag(1),1)) '& ' num2str(round(InnovDiag(2),1))  '& ' num2str(round(InnovDiag(3),1)) ' \\[0.4cm] '...   '$\quad$  of which  Contemporaenous  & $v^{Cont}=w' char(39) '3*\Sigma_uw$ & ' num2str(round(contCon(1),1)) '& ' num2str(round(contCon(2),1))  '& ' num2str(round(contCon(3),3)) ' \\' ... '$\quad$  of which  Lagged & $v^{lagged}=v^{agg}-v^{Cont}$ &  ' num2str(round(laggedCon(1),1)) '& ' num2str(round(laggedCon(2),1))  '& ' num2str(round(laggedCon(3),1)) ' \\[0.4cm] '...
    'Aggregate variance & $v^{agg}=w' char(39) 'IRVw$ &' num2str(round(outDeg(1),1)) '& ' num2str(round(outDeg(2),1))  '& ' num2str(round(outDeg(3),1)) ' \\ ' ...%  'Innovation Variance & $v^{inn}=w' char(39) 'diag(FEV)w$ & ' num2str(round(InnovDiag(1),1)) '& ' num2str(round(InnovDiag(2),1))  '& ' num2str(round(InnovDiag(3),1)) ' \\[0.4cm] '...   '$\quad$  of which  Contemporaenous  & $v^{Cont}=w' char(39) '3*\Sigma_uw$ & ' num2str(round(contCon(1),1)) '& ' num2str(round(contCon(2),1))  '& ' num2str(round(contCon(3),3)) ' \\' ... '$\quad$  of which  Lagged & $v^{lagged}=v^{agg}-v^{Cont}$ &  ' num2str(round(laggedCon(1),1)) '& ' num2str(round(laggedCon(2),1))  '& ' num2str(round(laggedCon(3),1)) ' \\[0.4cm] '...
    'Factor Contribution & $v^{factIP}=w' char(39) 'B\Sigma_{fact}B' char(39) 'w$ &' num2str(round(outDegFact(1),1)) '& ' num2str(round(outDegFact(2),1))  '& ' num2str(round(outDegFact(3),1)) ' \\ ' ...%  'Innovation Variance & $v^{inn}=w' char(39) 'diag(FEV)w$ & ' num2str(round(InnovDiag(1),1)) '& ' num2str(round(InnovDiag(2),1))  '& ' num2str(round(InnovDiag(3),1)) ' \\[0.4cm] '...   '$\quad$  of which  Contemporaenous  & $v^{Cont}=w' char(39) '3*\Sigma_uw$ & ' num2str(round(contCon(1),1)) '& ' num2str(round(contCon(2),1))  '& ' num2str(round(contCon(3),3)) ' \\' ... '$\quad$  of which  Lagged & $v^{lagged}=v^{agg}-v^{Cont}$ &  ' num2str(round(laggedCon(1),1)) '& ' num2str(round(laggedCon(2),1))  '& ' num2str(round(laggedCon(3),1)) ' \\[0.4cm] '...
    '$\quad$  of which Commodities  & $v^{Com}=w' char(39) '$  & ' num2str(round(outDegSingleFact(1,1),1)) '& ' num2str(round(outDegSingleFact(1,2),1))  '& ' num2str(round(outDegSingleFact(1,3),1)) ' \\' ...
    '$\quad$  of which Mon. policy & $v^{MP}= w' char(39) ' $ &  ' num2str(round(outDegSingleFact(2,1),1)) '& ' num2str(round(outDegSingleFact(2,2),1))  '& ' num2str(round(outDegSingleFact(2,3),1)) ' \\ '...
    '$\quad$  of which  TFP  & $v^{TFP}=w' char(39) '$ & ' num2str(round(outDegSingleFact(3,1),1)) '& ' num2str(round(outDegSingleFact(3,2),1))  '& ' num2str(round(outDegSingleFact(3,3),1)) ' \\[0.4cm]' ...
    ' Intrinsic component &$v^{IC}=w' char(39) 'diag(IRV)w$ &  ' num2str(round(outDegDiag(1),1)) '& ' num2str(round(outDegDiag(2),1))  '& ' num2str(round(outDegDiag(3),1)) ' \\'...
    '$\quad$  of which  equal weights  & $v^{ICeq}=v^{IC}|_{w^*_i=1/N}$  & ' num2str(round(outDegDiag2(1),1)) '& ' num2str(round(outDegDiag2(2),1))  '& ' num2str(round(outDegDiag2(3),1)) ' \\' ...
    '$\quad$  of which  unequal weights & $v^{ICuw}=v^{IC}|_{w^*_i=w_i-1/N} $ &  ' num2str(round(outDegDiag3(1),1)) '& ' num2str(round(outDegDiag3(2),1))  '& ' num2str(round(outDegDiag3(3),1)) ' \\ '...
    ' Extrinsic component & $v^{EC}=v^{agg}-v^{IC}$ & ' num2str(round(outDegOff(1),1)) '& ' num2str(round(outDegOff(2),1))  '& ' num2str(round(outDegOff(3),1)) ' \\' ...
    '$\quad$  of which  equal weights  & $v^{ECeq}=v^{EC}|_{w^*_i=1/N}$  & ' num2str(round(outDegOff2(1),1)) '& ' num2str(round(outDegOff2(2),1))  '& ' num2str(round(outDegOff2(3),1)) ' \\' ...
    '$\quad$  of which  unequal weights & $v^{ECuw}=v^{EC}|_{w^*_i=w_i-1/N} $ &  ' num2str(round(outDegOff3(1),1)) '& ' num2str(round(outDegOff3(2),1))  '& ' num2str(round(outDegOff3(3),1)) ' \\ '...
    '$\quad$  of which  symmetric  & $v^{Sy}=w' char(39) '(IRV^s-diag(IRV^s))w$ & ' num2str(round(SymetriCon(1),1)) '& ' num2str(round(SymetriCon(2),1))  '& ' num2str(round(SymetriCon(3),1)) ' \\' ...
    '$\qquad$  of which  equal weights  & $v^{Syeq}=v^{Sy}|_{w^*_i=1/N}$  & ' num2str(round(SymetriCon2(1),1)) '& ' num2str(round(SymetriCon2(2),1))  '& ' num2str(round(SymetriCon2(3),1)) ' \\' ...
    '$\qquad$  of which  unequal weights & $v^{Syuw}=v^{Sy}|_{w^*_i=w_i-1/N} $ &  ' num2str(round(SymetriCon3(1),1)) '& ' num2str(round(SymetriCon3(2),1))  '& ' num2str(round(SymetriCon3(3),1)) ' \\ '...
    '$\quad$  of which  asymmetric & $v^{Asy}=v^{EC}-v^{Sy}$ &  ' num2str(round(AsymetriCon(1),1)) '& ' num2str(round(AsymetriCon(2),1))  '& ' num2str(round(AsymetriCon(3),1)) ' \\ '...
    '$\qquad$  of which  equal weights  & $v^{Asyeq}=v^{Asy}|_{w^*_i=1/N}$  & ' num2str(round(AsymetriCon2(1),1)) '& ' num2str(round(AsymetriCon2(2),1))  '& ' num2str(round(AsymetriCon2(3),1)) ' \\' ...
    '$\qquad$  of which  unequal weights & $v^{Asyuw}=v^{Asy}|_{w^*_i=w_i-1/N} $ &  ' num2str(round(AsymetriCon3(1),1)) '& ' num2str(round(AsymetriCon3(2),1))  '& ' num2str(round(AsymetriCon3(3),1)) ' \\[0.4cm] '...
    'Extrinsic share & $ v^{EC}/v^{agg}$& '  num2str(round(outDegOff(1)/outDeg(1),2)) '& ' num2str(round(outDegOff(2)/outDeg(2),2))  '& ' num2str(round(outDegOff(3)/outDeg(3),2)) ' \\'...
    'Asymmetric share & $v^{Asy}/v^{agg}$ & ' num2str(round(AsymetriCon(1)/outDeg(1),2)) '& ' num2str(round(AsymetriCon(2)/outDeg(2),2))  '& ' num2str(round(AsymetriCon(3)/outDeg(3),2)) ' \\'...
    'Asymmetric/extrinsic & $ v^{Asy}/v^{EC}$& '  num2str(round(AsymetriCon(1)/outDegOff(1),2)) '& ' num2str(round(AsymetriCon(2)/outDegOff(2),2))  '& ' num2str(round(AsymetriCon(3)/outDegOff(3),2)) ' \\ \bottomrule \end{tabular}'];
disp(TableString)
%% Bias Variance
close all

fig1=figure;
fig1.PaperPositionMode = 'manual';
fig1.PaperOrientation = 'landscape';
fig1.PaperSize = [ 10  5 ];
fig1.PaperPosition = [ 0 0 1 1];
Tile = tiledlayout(1,2);
tile = tiledlayout(Tile,2,1,'TileSpacing','none');

Titles={'Bias','Frob-Norm','Frobenius Norm to $1_N$','','Variance'};
yTitles={'Average difference to true FEVD (log-scale)','','','','Average Variance of an entry (log-scale)'};
Ylims=[ min(abs([regCurve(:,1);regCurve2(:,1);regCurve3(:,1)])) 1; 0 105; 0 105;0 105; min(abs([regCurve(:,5);regCurve2(:,5);regCurve3(:,5)])) 10^-1 ];
order=[1 1 1 1 2];
for k = [1 5]
    subplot(1,2,order(k))
    semilogy(T,olsCurve(1:end,k),T,abs(regCurve(:,k)),':',T,abs(regCurve2(:,k)),'-.',T,abs(regCurve3(:,k)),'--','LineWidth',1.5)
    if mean(sign(regCurve3(:,k))<0)==1
        hold on
        semilogy(T(floor(linspace(1,length(T),10))),abs(regCurve3(floor(linspace(1,length(T),10)),k)),'xk','LineWidth',2,'MarkerSize',6)
    end
    if mean(sign(regCurve2(:,k))<0)==1
        hold on
        semilogy(T(floor(linspace(1,length(T),10))),abs(regCurve2(floor(linspace(1,length(T),10)),k)),'xk','LineWidth',2,'MarkerSize',6)
    end
%     if mean(regCurve2(:,k))<0.0001*mean(regCurve3(:,k))
%             semilogy(T,olsCurve(1:end,k),T,abs(regCurve(:,k)),':',T,10000*abs(regCurve2(:,k)),'-.',T,abs(regCurve3(:,k)),'--','LineWidth',1.5)
%     end
    ax = gca;
    ax.FontSize = 10;
    xticks([100:100:max(T)]);
%     axis([min(T)-5 max(T)+5 Ylims(k,:)]);
    title(Titles{k})
    ylabel(yTitles{k})
    xlabel('T')
    box off
    if k==1
        legend('OLS + Sample Cov',[Regularization.A ' + ' Regularization.Sigma], [Regularization2.A ' + ' Regularization2.Sigma],...
            [Regularization3.A ' + ' Regularization3.Sigma] , ['indicates that the' newline 'true line is negative' ],'FontSize',8,'location','northeast')
        legend boxoff
%         hold on
%         px=[0 10];
%         py1=[min(olsCurve(:,k))]*ones(1,2);
%         py2=[max(regCurve2(:,k))*10000]*ones(1,2);
%         plot(px,py1,'k','LineWidth',2);hold all;
%         plot(px,py2,'k','LineWidth',2);hold all;
%         px = [0 max(T)+90];
%         fill([px flip(px)],[py1 flip(py2)],'w','EdgeColor','none');
    end
    ax.YGrid='on';
end

print('-dpdf', './Graphs/BiasOLS.pdf','-fillpage')

%% Old graphs for running after the application file
%%
k=0;
for t=ts
    k= k+1;
    Table = diag(diag(FinalSigmaReg(:,:,t)))\FinalSigmaReg(:,:,t)/diag(diag(FinalSigmaReg(:,:,t)));
    Table= Table - diag(diag(Table));
    subplot(1,3,k)
    bins=linspace(-3,10,50)*10^-4;
    histogram(Table(Table~=0),bins)
    title(num2str(sum(sum(Table<0))))
end
%% scatter OutgoingFEV vs OutgoingIRF
scatter(sum(DfinalRegFEV(:,:,ts(1))-diag(diag(DfinalRegFEV(:,:,ts(1)))),2),sum(DfinalRegIRF(:,:,ts(1))-diag(diag(DfinalRegIRF(:,:,ts(1)))),2))
sum(DfinalRegIRF(:,:,ts(1))-diag(diag(DfinalRegIRF(:,:,ts(1)))),1)
fig=figure;
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [ 10  4];
fig.PaperPosition = [ 0 0 1 1];

k=0;
TITLE={'1972 to 1983','1984 to 1995','1996 to 2007'};

for t=ts
    k=k+1;
    weightComponent = meanWeight(:,k)/sum(meanWeight(:,k));
    Table = DfinalReg(:,:,t)-diag(diag(DfinalReg(:,:,t)));
    Table2 = DfinalRegIRF(:,:,t)-diag(diag(DfinalRegIRF(:,:,t)));
    
    deg2 = sum(Table,1)';
    degIRF = sum(Table2,1)';
    
    H1 = subplot(1,3,k);
    sz = 100*meanWeight(:,k);
    H2 = scatter(deg2,degIRF,sz,'filled','MarkerFaceAlpha',0.85,'LineWidth',1.5);
    colormap(jet)
    ls=lsline;
    ls.LineWidth=1.5; ls.Color=0.25*[1 1 1];
    title(TITLE{k})
    hold off
    if k>1
        yticklabels({})
        H1.YLim = limsY;
    else
        H1.YLim(2) = H1.YLim(2)*1.05;
        limsY = H1.YLim;
        ylabel('IRF-Out')
    end
    if k==2
        xlabel('FEV-Out')
    end
    %     H1.XLim(2) = 0.25;
    %         H1.YLim = [0 0.2];
    box on
    grid on
end

print('-dpdf', './Graphs/ScatteroutFEVIRF.pdf','-fillpage')


%% scatter ingoing  vs vola
fig=figure;
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [ 10  4];
fig.PaperPosition = [ 0 0 1 1];

k=0;
TITLE={'1972 to 1983','1984 to 1995','1996 to 2007'};

for t=ts
    k=k+1;
    weightComponent = meanWeight(:,k)/sum(meanWeight(:,k));
    Table=DfinalReg(:,:,t);%-diag(diag(DfinalReg(:,:,t))); %
    deg2 = sum(Table,2)';
    for i=1:N
        CI=concentrationIndices(Table(i,:));
        concentration(i) = CI.HH;
    end
    
    H1 = subplot(1,3,k);
    sz = 100*meanWeight(:,k);
    H2 = scatter(concentration,sz,sz,'filled','MarkerFaceAlpha',0.85,'LineWidth',1.5);
    colormap(jet)
    ls=lsline;
    ls.LineWidth=1.5; ls.Color=0.25*[1 1 1];
    title(TITLE{k})
    hold off
    if k>1
        yticklabels({})
        H1.YLim = limsY;
    else
        H1.YLim(2) = H1.YLim(2)*1.05;
        limsY = H1.YLim;
        ylabel('Vola')
    end
    if k==2
        xlabel('In-Connectedness')
    end
    %     H1.XLim(2) = 0.25;
    %         H1.YLim = [0 0.2];
    box on
    grid on
end

print('-dpdf', './Graphs/ScatterInVola.pdf','-fillpage')

%% Scatter Std of outgoing vs. sector weights
fig=figure;
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [ 10  4];
fig.PaperPosition = [ 0 0 1 1];
k=0;
TITLE={'1972 to 1983','1984 to 1995','1996 to 2007'};

for t=ts
    k=k+1;
    weightComponent = meanWeight(:,k)/sum(meanWeight(:,k));
    Table=DfinalReg(:,:,t)-diag(diag(DfinalReg(:,:,t))); %
    %     Table=round( diag(weightComponent)*(DfinalRegFEV(:,:,t)-diag(diag(DfinalRegFEV(:,:,t))))  +diag(diag(DfinalRegFEV(:,:,t)))*diag(weightComponent) ,3);
    deg2=std(Table,[],1);
    [totest, rankings]=sort(meanWeight(:,k),'descend');
    deg2ranked=deg2(rankings);
    deg2stairs=deg2ranked(1:floor(N/2));
    
    [~, maxI]=max(meanWeight(:,k));
    
    H1=subplot(1,3,k);
    sz = 100*meanWeight(:,k);
    H2 = scatter(deg2,weightComponent,sz,'filled','MarkerFaceAlpha',0.85,'LineWidth',1.5);
    colormap(jet)
    ls=lsline;
    ls.LineWidth=1.5; ls.Color=0.25*[1 1 1];
    title(TITLE{k})
    hold off
    if k>1
        yticklabels({})
        H1.YLim = limsY;
    else
        H1.YLim(2) = H1.YLim(2)*1.05;
        
        limsY = H1.YLim;
        ylabel('sector weight')
    end
    if k==2
        xlabel('In-Connectedness')
    end
    %     H1.XLim = [0 1];
    box on
    grid on
end

print('-dpdf', './Graphs/ScatterInCounting.pdf','-fillpage')

%% Scatter Out off and diag
fig=figure;
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [ 10  4];
fig.PaperPosition = [ 0 0 1 1];
k=0;
TITLE={'1972 to 1983','1984 to 1995','1996 to 2007'};
for t=ts
    k=k+1;
    weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));
    outDeg=weightComponent'*(DfinalRegFEV(:,:,t))*diag(weightComponent);%/(N*averageVAR(t));-diag(diag(DfinalRegFEV(:,:,t)))
    outDegOff=weightComponent'*(DfinalRegFEV(:,:,t)-diag(diag(DfinalRegFEV(:,:,t))))*diag(weightComponent);
    outDegDiag=weightComponent'*(diag(diag(DfinalRegFEV(:,:,t))))*diag(weightComponent);
    H1=subplot(1,3,k);
    sz = 100*meanWeight(:,k);
    c = outDegDiag./outDeg;
    H2 = scatter(outDegDiag,outDegOff,sz,outDeg,'filled','MarkerFaceAlpha',0.85,'LineWidth',1.5);
    if k==1
        
        yLim = H1.YLim ;[0 max(outDegOff)*1.05];
        xLim = [0 max(outDegDiag)*1.05];
    end
    H1.YLim =yLim;
    H1.XLim = xLim;
    H1.XScale ='log';
    H1.YScale ='log';
    
    colormap(0.7*jet)
    grid on
    box on
    ylabel(['Contr. via Spillovers ' num2str(round(sum(outDegOff),1))])
    xlabel(['Contr. via Idiosyncratic ' num2str(round(sum(outDegDiag),1))])
    title(TITLE{k})
end
print('-dpdf', './Graphs/ContributionViaSpillovers.pdf','-fillpage')

%% Scatter Out FEV changes
fig=figure;
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [ 10  4];
fig.PaperPosition = [ 0 0 1 1];
k=1;
weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));

outDegRef=weightComponent'*(DfinalRegFEV(:,:,t(1)))*diag(weightComponent);%/(N*averageVAR(t));-diag(diag(DfinalRegFEV(:,:,t)))

TITLE={'1972 to 1983','1984 to 1995','1996 to 2007'};
for t=ts(2:end)
    k=k+1;
    weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));
    outDeg=weightComponent'*(DfinalRegFEV(:,:,t))*diag(weightComponent);%/(N*averageVAR(t));-diag(diag(DfinalRegFEV(:,:,t)))
    
    H1=subplot(1,2,k-1);
    sz = 100*meanWeight(:,k);
    c = outDeg-outDegRef;
    H2 = scatter(outDeg,outDegRef,sz,c,'filled','MarkerFaceAlpha',0.85,'LineWidth',1.5);
    
    H1.XScale ='log';
    H1.YScale ='log';
    hold on
    plot(linspace(eps,max(outDegRef)*1.05,10),linspace(eps,max(outDegRef)*1.05,10),'k')
    H1.YLim =[min(outDegRef) max(outDegRef)*1.05];
    H1.XLim = [min(outDeg) max(outDeg)*1.05];
    colormap(0.7*jet)
    grid on
    box on
    lsline
    
    ylabel(['Contr. ' TITLE{1}])
    xlabel(['Contr. ' TITLE{k}])
end
print('-dpdf', './Graphs/ContributionScatterChange.pdf','-fillpage')

%% CDF outgoing
fig=figure;
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [ 10  4];
fig.PaperPosition = [ 0 0 1 1];

k=0;
TITLE={'1972 to 1983','1984 to 1995','1996 to 2007'};

for t=ts
    k=k+1;
    weightComponent = N*meanWeight(:,k)/sum(meanWeight(:,k));
    Table=diag(weightComponent)*(DfinalRegFEV(:,:,t)-diag(diag(DfinalRegFEV(:,:,t))))*diag(weightComponent); %
    %     Table=round( diag(weightComponent)*(DfinalRegFEV(:,:,t)-diag(diag(DfinalRegFEV(:,:,t))))  +diag(diag(DfinalRegFEV(:,:,t)))*diag(weightComponent) ,3);
    
    H2 = subplot(1,3,k);
    sectorSelect = 3;
    refSector = Table(sectorSelect,:);
    otherSectors = Table([1:sectorSelect-1 sectorSelect+1:N],:);
    bins=linspace(0,max(max(otherSectors)),N^2);
    cdfplot(refSector);
    hold on
    cdfplot(otherSectors(:));
    H2.XScale ='log';
    title(TITLE{k})
    hold off
    if k>1
        yticklabels({})
    end
end
print('-dpdf', './Graphs/HistOutCoalvsRest.pdf','-fillpage')

%% BumpChart
clear bumpRanks sortings NewSortings sortingsNew
fig=figure;
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [20 9];
fig.PaperPosition = [ 0 0 1 1];

k=0;
for t=ts
    k=k+1;
    weightComponent = N*meanWeight(:,k)/sum(meanWeight(:,k));
    Table=diag(weightComponent)*(DfinalReg(:,:,t)-diag(diag(DfinalReg(:,:,t)))); %
    deg2=sum(Table,1);
    [bumpRanks(k,:),sortings(k,:)]=sort(deg2,'descend');
    sortings2(k,sortings(k,:))=1:length(deg2);
end
for i=1:10
    for j=1:length(ts)
        sortingsNew(j,sortings(j,i))=11-i;
    end
end

k=0;
Xs=1:0.05:length(ts);
for x=Xs
    k=k+1;
    weightings=zeros(length(ts),1);
    weightings(floor(x):ceil(x))=pdf(makedist('Normal','mu',x,'sigma',0.2),floor(x):ceil(x));
    weightings=weightings/sum(weightings);
    NewSortings(k,:)=weightings'*sortingsNew;
end
plot(Xs,NewSortings,'LineWidth',2.5)
for i=1:10
    for j=linspace(1,length(ts),3)
        txt=convertStringsToChars(names{sortings(j,i)});
        txt=txt(5:end-16);
        text(j,11.1-i,txt)
    end
end
xlim([0.7 length(ts)+0.25])
ylim([-0.3 10.3])
xticks(1:3)
xticklabels({'1972 to 1983','1984 to 1995','1996 to 2007'});


%%
figure
startDate = datenum('01-01-1972');
endDate = datenum('01-04-2008');
[~,order] = sort(std(Y(:,1:147),[],2));
sectorPlot = [77 3];order(end-3:end);
plotData =  Y;
timeDate = linspace(startDate,endDate,size(plotData,2));
plot(timeDate,plotData(sectorPlot,:)','LineWidth',2)
title('IP growth by sector')
datetick('x','yyyy','keeplimits')
legend(names(sectorPlot))
figure
k = 0;
for t=ts
    k=k+1;
    subplot(1,3,k,0.03)
    bins = linspace(0,max(max(varBySector.^0.5)),N/2);
    histogram(varBySector(:,t).^0.5,bins)
end
histogram( ((mean(varBySector(:,ts(2:3)).^0.5,2)-varBySector(:,ts(1)).^0.5))./(varBySector(:,ts(1)).^0.5),linspace(-.8,0.8,100))


%% Scatter Out IRV changes
fig = figure('Name','IRV-out changes: scatter');
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [ 10  4];
fig.PaperPosition = [ 0 0 1 1];
k=1;
weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));

outDegRef=weightComponent'*(DfinalRegIRF(:,:,t(1)))*diag(weightComponent);%/(N*averageVAR(t));-diag(diag(DfinalRegFEV(:,:,t)))

TITLE={'1972 to 1983','1984 to 1995','1996 to 2007'};
for t=ts(2:end)
    k=k+1;
    weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));
    outDeg=weightComponent'*(DfinalRegIRF(:,:,t))*diag(weightComponent);%/(N*averageVAR(t));-diag(diag(DfinalRegFEV(:,:,t)))
    
    H1 = subplot(1,2,k-1);
    sz = 100*meanWeight(:,k);
    c = log(outDeg)-log(outDegRef);
    H2 = scatter(outDeg,outDegRef,sz,c,'filled','MarkerFaceAlpha',0.85,'LineWidth',1.5);
    
    H1.XScale ='log';
    H1.YScale ='log';
    hold on
    H1.YLim =[min(outDegRef) max(outDegRef)*15];
    H1.XLim = [min(outDeg) max(outDeg)*15];
    plot(linspace(min(abs(outDegRef)),max(outDegRef)*15,10),linspace(min(abs(outDegRef)),max(outDegRef)*15,10),'k')
    
    colormap(0.7*jet)
    grid on
    box on
    
    ylabel(['Contr. ' TITLE{1}])
    xlabel(['Contr. ' TITLE{k}])
end
print('-dpdf', './Graphs/ContributionIRVScatterChange.pdf','-fillpage')

%% Scatter Weights vs Out
fig=figure('Name','In- vs Outgoing Connectedness');
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [ 10  4];
fig.PaperPosition = [ 0 0 1 1];
k=0;

for t=ts(1:end)
    k=k+1;
    weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));
    outDeg=weightComponent'*(DfinalRegIRF(:,:,t))*diag(weightComponent);
    x = weightComponent;
    y = outDeg/sum(outDeg);
    
    H1 = subplot(1,length(ts),k);
    sz = 100*meanWeight(:,k);
    c = log(x)-log(y);
    c = c + min(c)+eps;
    H2 = scatter(x,y,sz,'filled','MarkerFaceAlpha',0.85,'LineWidth',1.5);
    
%     H1.XScale ='log';
%     H1.YScale ='log';
    hold on
    H1.YLim = [min(min(y),min(x)) max(max(y),max(x))];
    H1.XLim = H1.YLim;
    plot(linspace(min(min(y),min(x)),max(max(y),max(x)),10),linspace(min(min(y),min(x)),max(max(y),max(x)),10),'k')
    hold off
    colormap(0.7*jet)
    grid on
    box on
    
    ylabel(['Out-Con. IRV (normalized) ' TITLE{1}])
    xlabel(['Weights ' TITLE{k}])
end

%% Table MSFE

% calculate FEVDs
clear FEVD DfinalRegIRFfactor
for t = ts
    Factors = FactorsSave{t};
    X = Y(:,t-RW+1:t);
    Z = [ones(1,size(X,2));  Factors; Y(:,t-RW:t-1)];%;  Y(:,t-RW-1:t-2)]; %
    Est2 =  EstSave2(t) ;
    Est = EstSave1(t);
    factors = size(Z,1) - size(X,1);
    % Mean squared forecast error
    D = FEVD(Est2,0,12,'MSFE');
    covFactor = cov(Factors');

    DfinalRegMSFE(:,:,t) = D.Reg(1:N,1:N);
    DfinalRegMSFEfactor(:,:,t) =  Est.Areg(1:end,(2:factors))*covFactor*Est.Areg(1:end,(2:factors))';
    Est3 = Est2; Est3.SigReg = DfinalRegMSFEfactor(:,:,t);   D = FEVD(Est3,0,12,'MSFE');
    DfinalRegMSFEfactorPlusSpill = D.Reg(1:N,1:N);   
    DfinalRegMSFEfactorSpillovers(:,:,t) = DfinalRegMSFEfactorPlusSpill - DfinalRegMSFEfactor(:,:,t);

end
clear bumpRanks sortings NewSortings outDeg sortings outDegOff outDegDiag
k=0;
for t=ts
    k=k+1;
    weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));
    outDeg(k) = weightComponent'*(DfinalRegMSFE(:,:,t))*weightComponent;
    outDegOff(k) = weightComponent'*(DfinalRegMSFE(:,:,t)-diag(diag(DfinalRegMSFE(:,:,t))))*weightComponent;
    outDegDiag(k) = weightComponent'*(diag(diag(DfinalRegMSFE(:,:,t))))*weightComponent;
    
    outDegFactor(k) = weightComponent'*(DfinalRegMSFEfactor(:,:,t))*weightComponent;

    outDegFactorSpill(k) = weightComponent'*(DfinalRegMSFEfactorSpillovers(:,:,t))*weightComponent;

    laggedCon(k) = weightComponent'*(DfinalRegMSFE(:,:,t) - FinalSigmaReg(:,:,t))*weightComponent;
    contempCon(k) = weightComponent'*(FinalSigmaReg(:,:,t))*weightComponent;
end
TableString= ['\begin{tabular}{l l r r r}'  ' & & {1972-1983} &  {{1984-1995}} & {{1996-2007}} \\' '\cline{3-5}' ' \midrule ' ...
    'Data variance & $v^{IP}=var(w' char(39) 'y_t^q)$ &' num2str(round(sampleVar(1),1)) '& ' num2str(round(sampleVar(2),1))  '& ' num2str(round(sampleVar(3),1)) ' \\[0.2cm] ' ...%  'Innovation Variance & $v^{inn}=w' char(39) 'diag(FEV)w$ & ' num2str(round(InnovDiag(1),1)) '& ' num2str(round(InnovDiag(2),1))  '& ' num2str(round(InnovDiag(3),1)) ' \\[0.4cm] '...   '$\quad$  of which  Contemporaenous  & $v^{Cont}=w' char(39) '3*\Sigma_uw$ & ' num2str(round(contCon(1),1)) '& ' num2str(round(contCon(2),1))  '& ' num2str(round(contCon(3),3)) ' \\' ... '$\quad$  of which  Lagged & $v^{lagged}=v^{agg}-v^{Cont}$ &  ' num2str(round(laggedCon(1),1)) '& ' num2str(round(laggedCon(2),1))  '& ' num2str(round(laggedCon(3),1)) ' \\[0.4cm] '...
    'Factor variance & $v^{agg}=w' char(39) 'B\Sigma_{f}B' char(39) 'w$ &' num2str(round(outDegFactor(1),1)) '& ' num2str(round(outDegFactor(2),1))  '& ' num2str(round(outDegFactor(3),1)) ' \\[0.2cm] '...
    'Spillover variance & $v^{agg}=w' char(39) 'MSFEw$ &' num2str(round(outDeg(1),1)) '& ' num2str(round(outDeg(2),1))  '& ' num2str(round(outDeg(3),1)) ' \\ ' ...
    ' $\quad$ Intrinsic Component &$v^{IC}=w' char(39) 'diag(MSFE)w$ &  ' num2str(round(outDegDiag(1),1)) '& ' num2str(round(outDegDiag(2),1))  '& ' num2str(round(outDegDiag(3),1)) ' \\'...
    ' $\quad$ Extrinsic Component & $v^{EC}=v^{agg}-v^{IC}$ & ' num2str(round(outDegOff(1),1)) '& ' num2str(round(outDegOff(2),1))  '& ' num2str(round(outDegOff(3),1)) ' \\[0.2cm]' ...
    ' $\quad$ Contemporaneous Component & $v^{cont.}=w' char(39) '\Sigma w$& ' num2str(round(contempCon(1),1)) '& ' num2str(round(contempCon(2),1))  '& ' num2str(round(contempCon(3),1)) ' \\' ...
    ' $\quad$ Lagged Component &$v^{lag}= v^{agg}-v^{cont.}$ &  ' num2str(round(laggedCon(1),1)) '& ' num2str(round(laggedCon(2),1))  '& ' num2str(round(laggedCon(3),1)) ' \\[0.2cm]'...
    'Extrinsic Share & $ v^{EC}/v^{agg}$& '  num2str(round(outDegOff(1)/outDeg(1),2)) '& ' num2str(round(outDegOff(2)/outDeg(2),2))  '& ' num2str(round(outDegOff(3)/outDeg(3),2)) ' \\ '...
    'Lagged Share & $ v^{lag}/v^{agg}$& '  num2str(round(laggedCon(1)/outDeg(1),2)) '& ' num2str(round(laggedCon(2)/outDeg(2),2))  '& ' num2str(round(laggedCon(3)/outDeg(3),2)) ' \\ \bottomrule \end{tabular}'];
disp(TableString)

%% Joyplot
figure('Name','Joyplot')
JoyData = DfinalOLS(:,:,ts(1));
JoyData = JoyData - diag(diag(JoyData));
xLin = linspace(0,1.5*max(mean(JoyData,2)),size(JoyData,2));
tiledlayout(1,length(ts))
for t=ts
    nexttile
    JoyData = DfinalReg(:,:,t);
    JoyData = JoyData - diag(diag(JoyData));
    
    for j=1:size(JoyData,2)
        pd = fitdist(JoyData(:,j),'kernel','Kernel','normal','Width',0.00075);
        JoyPLotData(:,j) = pdf(pd,xLin);
    end
    if t == ts(1)
        offset = 0.07*max(max(JoyPLotData));
    end
    joyPlot(JoyPLotData,xLin,offset,'FaceColor',mean(JoyData,1),'FaceAlpha',0.6,'LineColor','k','StrokeColor','w')
    set(gcf,'position',[500,100,560,680])
    cbar = colorbar;
    
    if t == ts(1)
        axCur = gca;
        YlimSave = axCur.YLim;
        cLims = cbar.Limits;
    end
        caxis(cLims)
        cbar.Limits = cLims;

    set(gca,'box','off','YLim',YlimSave,'YTick',[],'Color',[0.93 0.93 0.93])
    set(gca,'YGrid','on')
    colormap(jet)
    title(TITLE{k}, 'FontSize',titleFontSize )
    
end

%%
k=0;
for t=ts
    k=k+1;
    weightComponent = N*meanWeight(:,k)/sum(meanWeight(:,k));
    Table=round( diag(weightComponent)*(DfinalReg(:,:,t)-diag(diag(DfinalReg(:,:,t))))*diag(weightComponent) ,3);
    deg2=sum(Table,1);
    [totest, rankings]=sort(deg2,'descend');
    
    for i=1:sum(deg2>1)
        top=names{rankings(i)};
        top=char(string(top));   top=top(5:end-16);
        HighDegree{i,k}=top;
    end
end

%% Average Pairwise Correlation
IndexIP = zeros(size(Y,2),1);
IndexIP(1+DataOffset:end) =  sum(weights(:,2:end).*Y(:,1+DataOffset:end),1);
startIndex = DataOffset;
d = datetime(year_date(startIndex),month_date(startIndex):(size(Y,2)- (startIndex)+ month_date(startIndex)),1);
for t = startIndex:size(Y,2)
    singleWeights = false(size(Y,2),1);
    singleWeights(max(t-RW+1,1):t) = 1;
    
    avPcorr(t) = mean(mean(corr(Y(:,singleWeights)'))-1/N);
    avVar(t) = var(IndexIP(singleWeights));
end
figure('Name','Detect the great Moderation')
tiledlayout(1,2);
nexttile
plot(d,avPcorr(startIndex:size(Y,2)))
title('Average Pairwise Correlation')
nexttile
plot(d,avVar(startIndex:size(Y,2)))
title('Average Variance')

%% BumpChart FEV
clear bumpRanks sortings NewSortings outDeg sortings outDegOff outDegDiag
fig = figure('Name','Bumpchart IRF');
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [12  7];
fig.PaperPosition = [ 0 0 1 1];
dimensionSum = 2;
k=0;
for t=ts
    k=k+1;
    weightComponent = (meanWeight(:,k)/sum(meanWeight(:,k)));
    weightInnovations = ones(size(weightComponent)); diag(FinalSigmaReg(:,:,t)).^(-0.5);

    outDeg(:,k)=sum(diag(weightComponent)'*(DfinalRegFEV(:,:,t))*diag(weightInnovations),dimensionSum);
    [bumpRanks(k,:), sortings(k,:)] = sort(outDeg(:,k),'descend');
end
sum(outDeg,1)
outDegSave= outDeg;
for i=1:length(outDeg)
    for j=1:length(ts)
        if bumpRanks(j,:)>0
            EndPositiv=sum(bumpRanks(j,:)>0);
            sortingsNew(j,sortings(j,i),1) = sum(bumpRanks(j,:)) - sum(bumpRanks(j,i:EndPositiv) )+distWhite(i);

            sortingsNew(j,sortings(j,i),2) = sum(bumpRanks(j,:)) - sum(bumpRanks(j,i+1:EndPositiv) ); %
        else
            sortingsNew(j,sortings(j,i),1) = sum(bumpRanks(j,:)) - sum(bumpRanks(j,i:end) );

            sortingsNew(j,sortings(j,i),2) = sum(bumpRanks(j,:)) - sum(bumpRanks(j,i+1:end)); %
        end
    end
end

k=0;
Xs=0.5:0.025:(length(ts)+0.5);
for x=Xs
    k=k+1;
    weightings = zeros(length(ts),1);
    if x<1
        firstWeight=1;
        lastWeight=ceil(x);
    elseif x>length(ts)
        firstWeight=floor(x);
        lastWeight=length(ts);
    else
        firstWeight=floor(x);
        lastWeight=ceil(x);
    end
    weightings(firstWeight:lastWeight) = pdf(makedist('Normal','mu',x,'sigma',0.2),firstWeight:lastWeight);
    weightings = weightings/sum(weightings);
    NewSortings(k,:,1) = weightings'*sortingsNew(:,:,1);
    NewSortings(k,:,2) = weightings'*sortingsNew(:,:,2);
end
hold on

for i=1:length(outDeg)-0
    fill([Xs, fliplr(Xs)],[NewSortings(:,length(outDeg)+1-i,1)' fliplr( NewSortings(:,length(outDeg)+1-i,2)' ) ],ColorsScale(i,:),'EdgeColor',0.125*[1 1 1])
end

for i=1:5
    for j=linspace(1,length(ts),3)
        txt=convertStringsToChars(names{sortings(j,i)});
        txt=txt(5:end-16);
        if sortings(j,i)==7
            txt='Electric Power Gen.';
        end
        %         txt =[txt ' (' num2str( round( 100*outDegOff(sortings(j,i),j)/outDeg(sortings(j,i),j) ,0) ) '%)'];
        text(j,0.5*sortingsNew(j,sortings(j,i),1)+0.5*sortingsNew(j,sortings(j,i),2),txt,'HorizontalAlignment','center','Color',[ 1 1 1])
    end
end
xlim([0.7 length(ts)+0.25])
ylim([0 max(sum(outDeg))])
box off
% axis off
xticks(1:3)
ylabel('Contribution to IP index variance (FEV)')
xticklabels({'1972 to 1983','1984 to 1995','1996 to 2007','2008 to 2019'});
print('-dpdf', './Graphs/BumpChartIRF.pdf','-fillpage')

%% Scatter In and Out
fig=figure('Name','In- vs Outgoing Connectedness');
fig.PaperPositionMode = 'manual';
fig.PaperOrientation = 'landscape';
fig.PaperSize = [ 10  4];
fig.PaperPosition = [ 0 0 1 1];
k=0;

for t=ts
    k=k+1;
    Table=DfinalReg(:,:,t)-diag(diag(DfinalReg(:,:,t))); %
    
    H1 = subplot(1,length(ts),k);
    deg2 = sum(Table,1);
    deg = sum(Table,2)';
    sz = 100*meanWeight(:,k);
    c = deg2-deg;
    H2 = scatter(deg2,deg,sz,c,'filled','MarkerFaceAlpha',0.85,'LineWidth',1.5);
    colormap(0.7*jet)
    [totest, rankings]=sort(c);
    colormap(0.7*jet)
    Alignment={'left','right'};
    for i=1:length(deg2)
        if i<=4 || i >= length(deg2) -4
             top = names(rankings(i)); top=char(string(top));  
            [~,ranks2]=sort(deg2ref,'descend');
            I = find(rankings(i)==ranks2);
            if i<=4
                %                     line([deg2(rankings(i)) 0.6],[deg2ref(rankings(i)) max(deg2ref)-i*0.075] ,'LineWidth',1,'Color',[0 0 0])
                text(0,max(deg2ref)-i*0.075,top); %+0.05*(-floor(rem(I,2))+1)
            elseif i >= length(deg2) -3
                %                     line([deg2(rankings(i)) 1],[deg2ref(rankings(i)) min(deg2ref)-(i-length(deg2ref))*0.075] ,'LineWidth',1,'Color',[0 0 0])
                text(0.6,min(deg2ref)-(i-length(deg2ref))*0.075 + 0.05,top); %+0.05*(-floor(rem(I,2))+1)
            end
        end
    end
    hold on
    plot(linspace(0,max(deg2ref)*1.05,10),linspace(0,max(deg2ref)*1.05,10),'k')
    H1.YLim = [0 max(deg2ref)*1.05];
    H1.XLim = [0 max(deg2ref)*1.05];
    grid on
    box on
    ylabel(['In-Con. ' TITLE{k}])
    xlabel(['Out-Con. ' TITLE{k}])
end
print('-dpdf', './Graphs/ScatterOutIn.pdf','-fillpage')
