%% Created by Daniel de Lima, 2015

%% Parameters
parameters
betta                   % discount rate
sig                     % intertemporal elasticity of 
varphi                  % inverse Frisch elasticity of labour 
alfa                    % capital share
eta_i                   % inverse elasticity of net investiment 
chi                     % labour utility weight
delta                   % exogenous rate of depreciation
I_ss                    % steady state investment
lambda                  % fraction of capital that can be 
omega                   % proportional transfer to the entering 
theta                   % survival rate of bankers
G_ss
tau
kappa
RkmR_ss
rho_ksi
sigma_ksi
sigma
;

%% Variables
var 
Y                       % output
K                       % capital
L                       % labour
I                       % investment
C                       % consumption
CK
CW
Q                       % market value of capital
varrho                  % marginal utility of consumption
varrhoW
Lambda                  % stochastic discount rate
LambdaW
Rk                      % capital return
R                       % risk-free return rate on deposits
In                      % net investment
N                       % bank's net wealth
Ne
Nn
nu                      % value of banks' capital
eta                     % value of banks' net wealth
phi                     % leverage
W
G
x
z
psi
ksi
B
BW
BK
premium
T
inT
;


%% Exogenous Variables
varexo 
e_ksi
;


%% Calibration
betta=0.990;
sig=1.000;
varphi=3;
alfa=0.3;
eta_i=1.728;
chi=1;
delta=0.025;
G_ss = 0.5508;
I_ss = 0.5460;            
lambda = 0.375;
omega = 0.002;	
theta = 0.97;
tau = 0.001;
kappa = 10;
RkmR_ss = 0.0026961;
rho_ksi = 0.66;
sigma_ksi = 0.05;
sigma = 0.05;


%% Declaring external functions


%% Model
model;

exp(varrho)  =   (exp(CK))^(-sig);        

exp(varrhoW) = (exp(CW) - (chi/(1+varphi))*(exp(L)^(1+varphi)))^(-sig);

exp(Lambda)  =   exp(varrho)/exp(varrho(-1));

exp(LambdaW)  =   exp(varrhoW)/exp(varrhoW(-1));

betta*exp(R)*exp(Lambda(+1))  =  1;

betta*exp(R)*exp(LambdaW(+1))  =  1;

exp(CW) = exp(W)*exp(L) - exp(T) + exp(BW(-1))*exp(R(-1)) - exp(BW);

exp(C) = sigma*exp(CK) + (1-sigma)*exp(CW);

exp(B) = sigma*exp(BK) + (1-sigma)*exp(BW);

chi*exp(L)^varphi    =   exp(W);           

exp(Rk)     =   (alfa*exp(Y)/exp(K(-1))+exp(ksi)*(exp(Q)-delta))/exp(Q(-1));

exp(Y)     =   (exp(ksi)^alfa)*(exp(K(-1))^alfa)*((1-sigma)^(1-alfa))*(exp(L)^(1-alfa));

exp(Q)  =   1+eta_i/2*((In+I_ss)/(In(-1)+I_ss)-1)^2+eta_i*((In+I_ss)/(In(-1)+I_ss)-1)*(In+I_ss)/(In(-1)+I_ss)-betta*exp(Lambda(+1))*eta_i*((In(+1)+I_ss)/(In+I_ss)-1)*((In(+1)+I_ss)/(In+I_ss))^2;

In  =   exp(I)-delta*exp(ksi)*exp(K(-1));

exp(K)  =   exp(ksi)*exp(K(-1))+In;

exp(Y)   =   exp(C)+exp(G)+exp(I)+eta_i/2*((In+I_ss)/(In(-1)+I_ss)-1)^2*(In+I_ss)+tau*psi*exp(K);

(1-sigma)*exp(W) = (1-alfa)*exp(Y)/exp(L);

%%% bloco Financial Acelerator

exp(nu)     =   (1-theta)*betta*exp(Lambda(+1))*(exp(Rk(+1))-exp(R))+betta*exp(Lambda(+1))*theta*exp(x(+1))*exp(nu(+1));

exp(eta)    =   (1-theta)+betta*exp(Lambda(+1))*theta*exp(z(+1))*exp(eta(+1));

exp(phi)    =   1/(1-psi)*exp(eta)/(lambda-exp(nu));

exp(z)      =   (exp(Rk)-exp(R(-1)))*(1-psi(-1))*exp(phi(-1))+exp(R(-1));

exp(x)      =   (exp(phi)*(1-psi)/(exp(phi(-1))*(1-psi(-1))))*exp(z);

exp(Q)*exp(K)     =   exp(phi)*exp(N);

exp(N)      =   exp(Ne)+exp(Nn);

exp(Ne)     =   theta*exp(z)*exp(N(-1));

exp(Nn)    =   omega*(1-psi(-1))*exp(Q)*exp(ksi)*exp(K(-1));

exp(Q)*exp(K) = exp(N) + sigma*exp(B);

exp(premium) = exp(Rk(+1))/exp(R);

%%% fim bloco FA

%% Governo

exp(T) + (exp(Rk) - exp(R(-1)))*exp(Q)*psi(-1)*exp(K(-1)) = exp(G) + tau*psi*exp(Q)*exp(K);

exp(G) = G_ss;

inT = tau*psi*exp(Q)*exp(K);

psi = 0;

%psi = kappa*(Rk(+1)-R-RkmR_ss);

%%% Shocks

ksi=   rho_ksi*ksi(-1)-e_ksi;

end;

%% Initial values
initval;
Y = 1.01309;
K = 3.08382;
L = 0.176927;
I = log(I_ss);
C = 0.505167;
CK = 1.3066;
CW = 0.4382;
Q = 0;
varrho = -1.3066;
varrhoW = -0.0418;
Lambda = 0;
LambdaW = 0;
Rk = 0.0127464;
R = 0.0100503;
In = 0;
N = 1.61388;
Ne = 1.60514;
Nn = -3.13078;
nu = -5.43715;
eta = 0.477444;
phi = 1.46995;
W = 0.530782;
G = log(G_ss);
x = 0.0217232;
z = 0.0217232;
psi = 0;
ksi = 0;
B = 5.81827;
BW = 1.9575;
BK = 8.7938;
premium = 0.00269;
T = log(G_ss);
inT = 0;
end;


%% Shocks
shocks;
var e_ksi=sigma_ksi^2;
end;

check;

steady;

resid(1);

%% Simulation
stoch_simul(order=1, periods=2000, irf=60);

A = [Y_e_ksi K_e_ksi L_e_ksi I_e_ksi C_e_ksi CK_e_ksi CW_e_ksi Q_e_ksi varrho_e_ksi Lambda_e_ksi Rk_e_ksi R_e_ksi In_e_ksi N_e_ksi Ne_e_ksi Nn_e_ksi nu_e_ksi eta_e_ksi phi_e_ksi W_e_ksi G_e_ksi x_e_ksi z_e_ksi psi_e_ksi ksi_e_ksi B_e_ksi premium_e_ksi T_e_ksi inT_e_ksi];
xlswrite('irfs.xlsx',A,'twoHTM');
xlswrite('irfs.xlsx',M_.endo_names_long,'names');