% This graph creates Figure 4 for the paper "Only time will tell" using
% LaTeX font
clear all
clc

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                              Example 1
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
rho=2/3;
pL=0.5;
tstar=log(pL)/log(rho);
LupperH=1;
LlowerH=1;
kappaIC1=4;
kappaPC1=10;
Deltar=0.13;

T1=1:1:5;
for t=T1
    pH=1-rho^t;
    if t<=tstar
        LupperH=[LupperH;LupperH(end)*(1-pH)];
        LlowerH=[LlowerH;LlowerH(end)*pH];
    else
        LupperH=[LupperH;LupperH(end)*pH];
        LlowerH=[LlowerH;LlowerH(end)*(1-pH)];
    end
end
LupperH=LupperH(2:end);
LlowerH=LlowerH(2:end);

Lupper1=[];%failures until tstar and then successes
Llower1=[];%successes until tstar and then failures
for t=T1
    Lupper1=[Lupper1;1-(pL^t/LupperH(t))];
    Llower1=[Llower1;1-(pL^t/LlowerH(t))];
end    
Llower1=[0;Llower1];
Lupper1=[0;Lupper1];


L11=1-0.5/rho;
L12=1-0.5^2/(rho*rho^2);
L13=1-0.5^3/(rho*rho^2*rho^3);
L14=1-0.5^4/(rho*rho^2*rho^3*(1-rho^4));
L15=1-0.5^5/(rho*rho^2*rho^3*(1-rho^4)*(1-rho^5));
L1=[0;L11;L12;L13;L14];
%failure, failure, failure, success, success

steps=5000;
t_min=0;
t_max=5;
T1=t_min:(t_max-t_min)/steps:t_max;

LR1=zeros(steps+1,steps+1);

for i = 1:steps+1
    if T1(i)>=1&&T1(i)<2
        LR_min=Llower1(2);
        LR_max=Lupper1(2);
        LR1(:,i)=LR_min:(LR_max-LR_min)/steps:LR_max;
    elseif T1(i)>=2&&T1(i)<3
        LR_min=Llower1(3);
        LR_max=Lupper1(3);
        LR1(:,i)=LR_min:(LR_max-LR_min)/steps:LR_max;
    elseif T1(i)>=3&&T1(i)<4
        LR_min=Llower1(4);
        LR_max=Lupper1(4);
        LR1(:,i)=LR_min:(LR_max-LR_min)/steps:LR_max;
    elseif T1(i)>=4
        LR_min=Llower1(5);
        LR_max=Lupper1(5);
        LR1(:,i)=LR_min:(LR_max-LR_min)/steps:LR_max;
    end
end



T1=kron(ones(steps+1,1),T1);
C1=real(exp(-2.*Deltar.*T1).*(kappaIC1.*LR1+kappaPC1).^2);
    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                              Example 2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Single loan with exponential default time distribution
lambdaH=1;
lambdaL=2;
step=0.01;
T2=0:step:4;
T22=T2;
Lupper2=[];%survival up to t
Llower2=[];
for t=T2
    Lupper2=[Lupper2;1-exp(-(lambdaL-lambdaH)*t)];
    Llower2=[Llower2;1-exp(-(lambdaL-lambdaH)*0)*lambdaL/lambdaH]; %default at t=0
end
% Llower2=1-exp(-(lambdaL-lambdaH)*0)*lambdaL/lambdaH.*ones(length(T2)); %default at t=0

L2=[];
i=0;
for t=T2
    if t<=1.5
        i=i+1;
        L2=[L2;Lupper2(i)];
    else 
        L2=[L2;Lupper2(i)-(1-Lupper2(i))*(lambdaL/lambdaH-1)]; %1-exp(-(lambdaL-lambdaH)*2.5)*lambdaL/lambdaH
    end
end    

lambdaL=2;
lambdaH=1;
Deltal=lambdaL-lambdaH;
Deltar=0.2;
rA=2;
cH=1;
Deltac=1;
v=7;
kappaIC2=Deltac/(2*Deltal^2/(lambdaH*(Deltal+lambdaL))*(1/(Deltar+rA)-1/(Deltal+lambdaL+Deltar+rA)));
kappaPC2=(Deltar+rA)*(v+cH)/2;


LR2=[]; %ones(T2,T2).*(1-lambdaL/lambdaH);
b2=[];
for t=T2
    LL=[];
    bb=[];
    for s=step:step:t
        ll=1-exp(-Deltal*s)*lambdaL/lambdaH;      
        LL=[LL;ll];
        bb=[bb;exp(-2*Deltar*s)*(kappaIC2*ll+kappaPC2)^2];
    end
    for s=t+step:step:T2(end)
        ll=1-exp(-Deltal*t)*lambdaL/lambdaH;      
        LL=[LL;ll];
        bb=[bb;exp(-2*Deltar*t)*(kappaIC2*ll+kappaPC2)^2];
    end
    I=1-exp(-Deltal*t);
    LL=[LL;I];
    LR2=[LR2,LL];
    bb=[bb;exp(-2*Deltar*t)*(kappaIC2*I+kappaPC2)^2];
    b2=[b2,bb];
end

LR_max  =       Lupper2;
LR_min  =       Llower2;
LR2      =      zeros(length(T2),length(T2));

for i = 1:1:length(T22)
LR2(:,i) =       (LR_min(i):(LR_max(i)-LR_min(i))/(length(T2)-1):LR_max(i))';
end



T2       =       kron(ones(length(T2),1),T2);
C2       =       real(exp(-2.*Deltar.*T2).*(kappaIC2.*LR2+kappaPC2).^2);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                              Example 3
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
muH=1.5;
muL=-1.5;
sigma=2;
Lupper3=[1,1];
T3=0:0.01:4;

N=length(T3)-1;
T=T3(end);
t = (0:1:N)'/N;                     % t is the column vector [0 1/N 2/N ... 1]
W = [0; cumsum(randn(N,1))]/sqrt(N); % S is running sum of N(0,1/N) variables
t = t*T;
W = W*sqrt(T);
X = muH*t + sigma*W;
L3=1-exp(-((muH-muL)/sigma^2).*X);

kappaIC3=4.2416;
kappaPC3=0.2651;
aH=1/2;
aL=-aH;
rA=0.6;
rP=0.1;
dr=rA-rP;

steps   =       500;
t_min   =       0;
t_max   =       1.25;

T       =       t_min:(t_max-t_min)/steps:t_max;
LR_max  =       1.*ones(size(T,2),1);
xstar   =       -sigma.^2./(aH-aL).*log(1+kappaPC3./kappaIC3);
LR_min  =       (1-exp(-(aH-aL)./sigma.^2.*xstar)).*ones(size(T,2),1);

T3       =       kron(ones(steps+1,1),T);
LR3      =       zeros(steps+1,steps+1);

for i = 1:steps+1
LR3(:,i) =       LR_min(i):(LR_max(i)-LR_min(i))/steps:LR_max(i);
end

C3       =       real(exp(-2.*dr.*T3).*(kappaIC3.*LR3+kappaPC3).^2);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                              Figure
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
figure;
set(0,'DefaultTextInterpreter', 'latex')

subplot(1,3,1);
hold on

contourf(T1,LR1,C1,20)

T1=0:1:5;
plot([T1(5),T1(5)],[Lupper1(4),1],'color', 'white','Linewidth',2);
plot([T1(1),T1(2)],[Lupper1(1),Lupper1(1)],'color', 'black','Linewidth',2);
plot([T1(2),T1(3)],[Lupper1(2),Lupper1(2)],'color', 'black','Linewidth',2);
plot([T1(3),T1(4)],[Lupper1(3),Lupper1(3)],'color', 'black','Linewidth',2);
plot([T1(4),T1(5)],[Lupper1(4),Lupper1(4)],'color', 'black','Linewidth',2);
% plot([T1(5),T1(6)],[Lupper1(5),Lupper1(5)],'color', 'black','Linewidth',2);

% plot([T1(2),T1(2)],[Lupper1(1),Lupper1(2)],'color', 'black','Linewidth',2,'Linestyle',':');
% plot([T1(3),T1(3)],[Lupper1(2),Lupper1(3)],'color', 'black','Linewidth',2,'Linestyle',':');
% plot([T1(4),T1(4)],[Lupper1(3),Lupper1(4)],'color', 'black','Linewidth',2,'Linestyle',':');
% plot([T1(5),T1(5)],[Lupper1(4),Lupper1(5)],'color', 'black','Linewidth',2,'Linestyle',':');
% plot([T1(6),T1(6)],[Lupper1(5),Lupper1(6)],'color', [.99 .6 0],'Linewidth',2,'Linestyle',':');

plot([T1(1),T1(2)],[Llower1(1),Llower1(1)],'color', 'black','Linewidth',2);
plot([T1(2),T1(3)],[Llower1(2),Llower1(2)],'color', 'black','Linewidth',2);
plot([T1(3),T1(4)],[Llower1(3),Llower1(3)],'color', 'black','Linewidth',2);
plot([T1(4),T1(5)],[Llower1(4),Llower1(4)],'color', 'black','Linewidth',2);




xlabel('\bf{Time $$t$$}','fontsize',16,'color','black');
ylabel('\bf{Likelihood Ratio $$L_t$$}','fontsize',16,'color','black');
ax = gca;
ax.XTick = [0 1 2 3 4];
ax.YTick =  [-2 -1 0 1];
ax.TickLabelInterpreter='latex';
ax.XTickLabel = {'$$0$$','$$1$$','$$2$$','$$3$$','$$4$$'};
ax.YTickLabel = {'$$-2$$','$$-1$$','$$0$$','$$1$$',};
ax.FontSize = 16;
axis([0 4 Llower1(4)-0.1 1.01]);
% legendh = legend('\bf{Support $$\left[\underline{L}_t,\overline{L}_t \right]$$}', '\bf{Sample path $$l_t$$}', '\bf{Cutoff $$\widehat{LR_t}$$}','Location','NW');
% set(legendh, 'Interpreter','latex','FontSize',16);
% legend('boxoff');
title('\bf{Example 1}','fontsize',16,'color','black');

% plot(T1,[Lupper1,Llower1]);
hold off

subplot(1,3,2);
hold on

contourf(T2,LR2,C2,20)


plot(T22,Lupper2,'Linewidth',2,'color','black');
plot(T22,Llower2,'Linewidth',2,'color','black');


xlabel('\bf{Time $$t$$}','fontsize',16,'color','black');
ylabel('\bf{Likelihood Ratio $$L_t$$}','fontsize',16,'color','black');
ax = gca;
ax.XTick = [0 1 2 3 4];
ax.YTick =  [-2 -1 0 1];
ax.TickLabelInterpreter='latex';
ax.XTickLabel = {'$$0$$','$$1$$','$$2$$','$$3$$','$$4$$'};
ax.YTickLabel = {'$$-2$$','$$-1$$','$$0$$','$$1$$',};
ax.FontSize = 16;
axis([0 4 -1.05 1.01]);
title('\bf{Example 2}','fontsize',16,'color','black');

hold off

subplot(1,3,3);
hold on
contourf(T3,LR3,C3,20);

colormap(jet);
colorbar('TickLabelInterpreter','latex');

plot([T3(1),T3(end)],Lupper3,'Linewidth',2,'color','black');

% plot(T3,L3,'Linewidth',2,'color',[.99 .6 0]);
xlabel('\bf{Time $$t$$}','fontsize',16,'color','black');
ylabel('\bf{Likelihood Ratio $$L_t$$}','fontsize',16,'color','black');
ax = gca;
ax.XTick = [0 1 2 3 4];
ax.YTick =  [-2 -1 0 1];
ax.TickLabelInterpreter='latex';
ax.XTickLabel = {'$$0$$','$$1$$','$$2$$','$$3$$','$$4$$'};
ax.YTickLabel = {'$$-2$$','$$-1$$','$$0$$','$$1$$',};
ax.FontSize = 16;
axis([0 T3(end,end) 0 1.05]);
% legendh = legend('\bf{Support $$\left(- \infty,\overline{L}_t \right]$$}', '\bf{Sample path $$l_t$$}','Location','NW');
% set(legendh, 'Interpreter','latex','FontSize',16);
% legend('boxoff');
title('\bf{Example 3}','fontsize',16,'color','black');

hold off
% 
set(gcf,'renderer','painters')

% RESIZE GRAPH FIRST
print(gcf,'-depsc2','MATLAB_graph.eps')
system('epstopdf MATLAB_graph.eps');

