%% generate the impluse response to a positive structural monetary policy shock under different capital adequacy ratios


clear;

%% set the type of shocks
sig_tau=1;  % =1 or -1 implies that the economy faces a 1% positive or negative structural monetary policy shock, =0 implies that the economy doesn't experience a structural monetary policy shock
sig_a=0;    % =1 or -1 implies that the economy faces a 1% positive or negative productivity shock, =0 implies that the economy doesn't experience this shock
sig_m=0;    % =1 or -1 implies that the economy faces a 1% positive or negative traditional monetary policy shock, =0 implies that the economy doesn't experience this shock
sig_g=0;    % =1 or -1 implies that the economy faces a 1% positive or negative capital regulatory policy shock, =0 implies that the economy doesn't experience this shock
save shock_type sig_tau sig_a sig_m sig_g

%% set the level of capital adequacy ratio

states_1=0.4; states_2=0.1;         % benchmark case where the strength of structural monetary policy equals 0.4 and the capital adequacy ratio equals 0.1
save regime_type states_1 states_2
dynare   bench_joint.mod 
oo_bench=oo_;                       % define oo_bench: oo_ from benchmark case
save SM_BENCH oo_bench

states_1=0.4; states_2=0.15;       % alternative case where the strength of structural monetary policy equals 0.4 and the capital adequacy ratio equals 0.15
save regime_type states_1 states_2
dynare   bench_joint.mod 
oo_zeta=oo_;                       % define oo_zeta: oo_ from alternative case where tau_p=0.4, GAMMA=0.15;
save SM_GAMMA2 oo_zeta

states_1=0.4; states_2=0.05;       % alternative case where the strength of structural monetary policy equals 0.4 and the capital adequacy ratio equals 0.05
save regime_type states_1 states_2
dynare   bench_joint.mod 
oo_gamma1=oo_;                     % define oo_gamma1: oo_ from alternative case where tau_p=0.4, GAMMA=0.05;
save SM_GAMMA1 oo_gamma1 


%% generate the impluse response to a positive capital regulatory policy shock under different strengths of structural monetary policy
clear;

sig_tau=0;sig_a=0;sig_m=0; sig_g=1;  % set the type of shocks
save shock_type sig_tau sig_a sig_m sig_g

%% set the strength of structural monetary policy

states_1=0.4; states_2=0.1;       % benchmark case where the strength of structural monetary policy equals 0.4 and the capital adequacy ratio equals 0.1
save regime_type states_1 states_2
dynare   bench_joint.mod 
oo_gamma=oo_;                     % define oo_gamma: oo_ from benchmark case where tau_p=0.4, GAMMA=0.1;
save GAMMA_BENCH oo_gamma 


states_1=0.1; states_2=0.1;       % alternative case where the strength of structural monetary policy equals 0.1 and the capital adequacy ratio equals 0.1
save regime_type states_1 states_2
dynare   bench_joint.mod 
oo_tau=oo_;                       % define oo_tau: oo_ from alternative case where tau_p=0.1, GAMMA=0.1;
save GAMMA_TAU1 oo_tau

states_1=0.7; states_2=0.1;       % alternative case where the strength of structural monetary policy equals 0.7 and the capital adequacy ratio equals 0.1
save regime_type states_1 states_2
dynare   bench_joint.mod 
oo_tau1=oo_;                      % define oo_tau1: oo_ from alternative case where tau_p=0.7, GAMMA=0.1;
save GAMMA_TAU2 oo_tau1

%% generate figure 4
clear
load SM_BENCH 
load SM_GAMMA1 
load SM_GAMMA2
load GAMMA_BENCH
load GAMMA_TAU2
load GAMMA_TAU1
figure;
fig_width = 2 * 400;  
fig_height = 400;
set(gcf, 'Position', [100 100 fig_width fig_height]);
subplot(1,2,2)
plot(100*oo_zeta.irfs.pr_exp_tau./oo_zeta.steady_state(38),'k:','LineWidth',1.5);
hold all
plot(100*oo_bench.irfs.pr_exp_tau./oo_bench.steady_state(38),'k-','LineWidth',1.5);
hold all
plot(100*oo_gamma1.irfs.pr_exp_tau./oo_gamma1.steady_state(38),'k--','LineWidth',1.5);
hold all 
axis square; 
title('影子银行违约率','FontSize',18,'interpreter','latex')
legend('$\gamma$=0.15','$\gamma$=0.1','$\gamma$=0.05','FontSize',18,'interpreter','latex')
box off
legend('boxoff')
subplot(1,2,1)
plot(100*oo_tau.irfs.lls_exp_g./oo_tau.steady_state(44),'k:','LineWidth',1.5);
hold all
plot(100*oo_gamma.irfs.lls_exp_g./oo_gamma.steady_state(44),'k-','LineWidth',1.5);
hold all
plot(100*oo_tau1.irfs.lls_exp_g./oo_tau1.steady_state(44),'k--','LineWidth',1.5);
hold all
title('商业银行信贷结构','FontSize',18,'interpreter','latex')
box off
 axis square; 
legend('$\tau_p$=0.1','$\tau_p$=0.4','$\tau_p$=0.7','FontSize',18,'interpreter','latex')
legend('boxoff')

%% Compute steady state values under benchmark economy
clear all;

parameter=calibration;
BETA=parameter(1); DELTA=parameter(2); THETA_S=parameter(3); THETA_P=parameter(4);
ALPHA=parameter(5); As=parameter(6); Ap=parameter(7); ETA=parameter(8); GAMMA=parameter(9);
TAU_S=parameter(10); TAU_P=parameter(11); SIGMA=parameter(12); OMEGA_I=parameter(13); MU=parameter(14);
ZETA =parameter(15);EPSILON=parameter(16);PI=parameter(17); RL=parameter(18); 
KM=parameter(19); BETA_B=parameter(20); BETA_E=parameter(21);PHI=parameter(22); OMEGAM=parameter(32);
ALPHA_S=parameter(23);ALPHA_P =parameter(24);OMEGA_P=parameter(25); ALPHA_SB=parameter(30);
KAPPA=parameter(26); XI=parameter(27);KAPPA_M=parameter(28);BETA_SB=parameter(29); ZETAL=parameter(31);

[lcr,fval,cflagr]=fsolve(@(x) solve_steady(x(1),x(2),x(3),x(4),x(5),PI,DELTA,THETA_S,THETA_P,ALPHA,As,Ap,ETA,GAMMA,TAU_S,TAU_P,KAPPA,...
    SIGMA,MU,ZETA,RL,KM,OMEGAM,BETA_B,PHI,ALPHA_S,ALPHA_P,XI,BETA_E,BETA,EPSILON,BETA_SB,ALPHA_SB,ZETAL),[7.2,6.12,0.924,0.22,0.39],optimset('TolFun',1e-10));
lp=lcr(1);ls=lcr(2);omegas=lcr(3);hp=lcr(4);pp=lcr(5);

[y,yp,ys,ch,h,d,m,ks,wage,DIV_S,hs,kp,lp_bar,DIV_P,ps,Rs,Rp,Rp_bar,rb,Xb,F_xi,PI_CB,PI_SB,b,k,i,DWL,totutility,bg]=solve_star(lp,ls,omegas,hp,pp,PI,DELTA,THETA_S,THETA_P,ALPHA,Ap,ETA,GAMMA,TAU_S,TAU_P,KAPPA,...
    SIGMA,MU,ZETA,RL,KM,OMEGAM,BETA_B,PHI,ALPHA_S,ALPHA_P,XI,BETA,EPSILON,KAPPA_M,BETA_SB,ALPHA_SB,ZETAL);
Y=y;L=lp+ls+lp_bar;LLS=lp/ls;YYS=yp/ys;
save ss_info Y L LLS YYS

%% Compute volatilities of total output, inflation and default rate, as well as welfare gains following a positive specific productivity shocks of POEs
clear all

% set the type of shocks
sig_ap=1;   % =1 or -1 implies that the economy faces a 1% positive or negative specific productivity shocks of POEs, =0 implies that the economy doesn't experience such a shock
sig_a=0;    % =1 or -1 implies that the economy faces a 1% positive or negative productivity shock, =0 implies that the economy doesn't experience this shock
sig_m=0;    % =1 or -1 implies that the economy faces a 1% positive or negative traditional monetary policy shock, =0 implies that the economy doesn't experience this shock
sig_xi=0;   % =1 or -1 implies that the economy faces a 1% positive or negative shadow bank financing efficiency shock, =0 implies that the economy doesn't experience this shock
sig_sigma=0;   % =1 or -1 implies that the economy faces a 1% positive or negative uncertainty shock, =0 implies that the economy doesn't experience such a shock
save shock_types sig_ap sig_a sig_m sig_xi sig_sigma

% Compute volatilities and welfare gains under different response coefficients

vNU1=[-2:0.5:2];
vNU2=[-2:0.5:2];
for mi=1:1:length(vNU1)
    data1=vNU1(mi);
    for mm=1:1:length(vNU2)
        data2=vNU2(mm);
save info_multi data1 data2
dynare bench_mul.mod
totutility_mean=oo_.mean(41);
totutility_ss=oo_.steady_state(41);
var_ch=oo_.var(1,1);var_h=oo_.var(2,2);
ch_ss=oo_.steady_state(1);
vol_y_1(mi,mm)=sqrt(oo_.var(26,26))/oo_.steady_state(26);var_y(mi,mm)=oo_.var(26,26);var_pi(mi,mm)=oo_.var(8,8);
vol_pi_1(mi,mm)=sqrt(oo_.var(8,8))/oo_.steady_state(8);var_gdp(mi,mm)=oo_.var(48,48); var_loan(mi,mm)=oo_.var(57,57);
vol_pr_1(mi,mm)=sqrt(oo_.var(38,38))/oo_.steady_state(38);var_pr(mi,mm)=oo_.var(38,38);
vol_lls(mi,mm)=sqrt(oo_.var(44,44))/oo_.steady_state(44);var_lls(mi,mm)=oo_.var(44,44);
vol_l(mi,mm)=sqrt(oo_.var(57,57))/oo_.steady_state(57);
vol_yys(mi,mm)=sqrt(oo_.var(45,45))/oo_.steady_state(45);var_yys(mi,mm)=oo_.var(45,45);
W_gain_1(mi,mm)=(exp((totutility_mean-totutility_ss)*(1-BETA))-1)*10000;
end
end
save Volatility_TFP 

%% generate figure 5
clear all
load Volatility_TFP 
figure;
subplot(2,2,1)
x=vNU1;y=vNU2;
[xx,yy]=meshgrid(x,y);
surf(xx,yy,vol_y_1);
zlabel('总产出波动率','FontSize',18,'interpreter','latex')
xlabel('$\psi_\gamma$','FontSize',18,'interpreter','latex')
ylabel('$\psi_\tau$','FontSize',18,'interpreter','latex')

subplot(2,2,2)
surf(xx,yy,vol_pi_1);
zlabel('通胀波动率','FontSize',18,'interpreter','latex')
xlabel('$\psi_\gamma$','FontSize',18,'interpreter','latex')
ylabel('$\psi_\tau$','FontSize',18,'interpreter','latex')


subplot(2,2,3)
surf(xx,yy,vol_pr_1);
zlabel('违约率波动率','FontSize',18,'interpreter','latex')
xlabel('$\psi_\gamma$','FontSize',18,'interpreter','latex')
ylabel('$\psi_\tau$','FontSize',18,'interpreter','latex')


subplot(2,2,4)
surf(xx,yy,W_gain_1);
zlabel('福利收益（%%）','FontSize',18,'interpreter','latex')
xlabel('$\psi_\gamma$','FontSize',18,'interpreter','latex')
ylabel('$\psi_\tau$','FontSize',18,'interpreter','latex')

%% Compute volatilities of total output, inflation and default rate, as well as welfare gains following a negative traditional monetary policy shock
clear all

% set the type of shocks
sig_ap=0;   % =1 or -1 implies that the economy faces a 1% positive or negative specific productivity shocks of POEs, =0 implies that the economy doesn't experience such a shock
sig_a=0;    % =1 or -1 implies that the economy faces a 1% positive or negative productivity shock, =0 implies that the economy doesn't experience this shock
sig_m=-1;    % =1 or -1 implies that the economy faces a 1% positive or negative traditional monetary policy shock, =0 implies that the economy doesn't experience this shock
sig_xi=0;   % =1 or -1 implies that the economy faces a 1% positive or negative shadow bank financing efficiency shock, =0 implies that the economy doesn't experience this shock
sig_sigma=0;   % =1 or -1 implies that the economy faces a 1% positive or negative uncertainty shock, =0 implies that the economy doesn't experience such a shock
save shock_types sig_ap sig_a sig_m sig_xi sig_sigma

% Compute volatilities and welfare gains under different response coefficients

vNU1=[-2:0.5:2];
vNU2=[-2:0.5:2];
for mi=1:1:length(vNU1)
    data1=vNU1(mi);
    for mm=1:1:length(vNU2)
        data2=vNU2(mm);
save info_multi data1 data2
dynare bench_mul.mod
totutility_mean=oo_.mean(41);
totutility_ss=oo_.steady_state(41);
var_ch=oo_.var(1,1);var_h=oo_.var(2,2);
ch_ss=oo_.steady_state(1);
vol_y(mi,mm)=sqrt(oo_.var(26,26))/oo_.steady_state(26);var_y(mi,mm)=oo_.var(26,26);var_pi(mi,mm)=oo_.var(8,8);
vol_pi(mi,mm)=sqrt(oo_.var(8,8))/oo_.steady_state(8);var_gdp(mi,mm)=oo_.var(48,48); var_loan(mi,mm)=oo_.var(57,57);
vol_pr(mi,mm)=sqrt(oo_.var(38,38))/oo_.steady_state(38);var_pr(mi,mm)=oo_.var(38,38);
vol_lls(mi,mm)=sqrt(oo_.var(44,44))/oo_.steady_state(44);var_lls(mi,mm)=oo_.var(44,44);
vol_l(mi,mm)=sqrt(oo_.var(57,57))/oo_.steady_state(57);
vol_yys(mi,mm)=sqrt(oo_.var(45,45))/oo_.steady_state(45);var_yys(mi,mm)=oo_.var(45,45);
W_gain(mi,mm)=(exp((totutility_mean-totutility_ss)*(1-BETA))-1)*10000;
end
end
save Volatility_MP 

%% generate figure 6
clear all
load Volatility_MP 
figure;
subplot(2,2,1)
x=vNU1;y=vNU2;
[xx,yy]=meshgrid(x,y);
surf(xx,yy,vol_y);
zlabel('总产出波动率','FontSize',18,'interpreter','latex')
xlabel('$\psi_\gamma$','FontSize',18,'interpreter','latex')
ylabel('$\psi_\tau$','FontSize',18,'interpreter','latex')

subplot(2,2,2)
surf(xx,yy,vol_pi);
zlabel('通胀波动率','FontSize',18,'interpreter','latex')
xlabel('$\psi_\gamma$','FontSize',18,'interpreter','latex')
ylabel('$\psi_\tau$','FontSize',18,'interpreter','latex')


subplot(2,2,3)
surf(xx,yy,vol_pr);
zlabel('违约率波动率','FontSize',18,'interpreter','latex')
xlabel('$\psi_\gamma$','FontSize',18,'interpreter','latex')
ylabel('$\psi_\tau$','FontSize',18,'interpreter','latex')


subplot(2,2,4)
surf(xx,yy,W_gain);
zlabel('福利收益（%%）','FontSize',18,'interpreter','latex')
xlabel('$\psi_\gamma$','FontSize',18,'interpreter','latex')
ylabel('$\psi_\tau$','FontSize',18,'interpreter','latex')