% The RBC Model with informal sector

% Informal goods: consumption only
% Formal goods: consumption and investment
% GHH Preference

close all;

%----------------------------------------------------------------------
%1. Defining Variables
%----------------------------------------------------------------------

var  tau kappa lY lY_f lY_i lC_f lC_i ld lK lKg lH lh_f lh_i llambda lx r A_f A_i lG p  lC ltb lca ;

varexo ef ei eg ekappa etau;

parameters  alpha pi etta beta theta gamma sigma r_star R rho_f rho_i sigma_f sigma_i sigma_g sigma_kappa d_bar phi delta  G_ss rho_G omega rho_tau sigma_tau rho_kappa tau_ss kappa_ss;

%----------------------------------------------------------------------
% 2. Calibration
% Taken from Ihrig & Moe 2004 and Mitra 2013
%----------------------------------------------------------------------

alpha = 0.35;  %formal good share of utility
beta = 0.96;  %discount rate

%elasticity of formal output wrt laborgov
theta = 0.43;
pi    = 0.25;
etta = 0.75;

%elasticity of informal output wrt labor
gamma = .495;

sigma = 1.5;
r_star = 0.04;
R = 0.00742;

%data from J. Fernandez's model
rho_f = 0.90;
rho_i = 0.64;
sigma_f = 0.10;
sigma_i = 0.10;
sigma_g = 0.10;
sigma_kappa = 0.10;

phi = 0.028 ;  %coefficient of capital adjustment
delta = 0.1; %capital depreciation

omega = 1.4556789;

G_ss = 0.20;
rho_G = 0.8;
kappa_ss = 0.25;
rho_kappa = .8;

rho_tau = .8;
sigma_tau = .10;
tau_ss = .2;

d_bar = 0.0074;
 
%----------------------------------------------------------------------
%3. Model
%----------------------------------------------------------------------

model;

%Total Consumption
lC = alpha*lC_f+(1-alpha)*lC_i;

%Total Production
exp(lY)= exp(lY_f)+p*exp(lY_i);

%Discount 
exp(llambda) = exp(llambda(+1))*beta*(1+r);

%FOC wrt formal consumption
alpha*(exp(lC_i)/exp(lC_f))^(1-alpha)*(exp(lC_f)^(alpha)*exp(lC_i)^(1-alpha)-((exp(lh_f)+exp(lh_i))^omega)/omega)^(-sigma) = exp(llambda);

%FOC wrt informal consumption
(1-alpha)*(exp(lC_f)/exp(lC_i))^(alpha)*(exp(lC_f)^(alpha)*exp(lC_i)^(1-alpha)-((exp(lh_f) + exp(lh_i))^omega)/omega)^(-sigma) = p*exp(llambda);

%FOC wrt formal labor
(exp(lh_f) + exp(lh_i))^(omega-1)*(exp(lC_f)^(alpha)*exp(lC_i)^(1-alpha)-((exp(lh_f) + exp(lh_i))^omega)/omega)^(-sigma) = (1-(tau))*exp(llambda)*theta*exp(lY_f)/exp(lh_f);

%FOC wrt Informal labor
(exp(lh_f) + exp(lh_i))^(omega-1)*(exp(lC_f)^(alpha)*exp(lC_i)^(1-alpha)-((exp(lh_f) + exp(lh_i))^omega)/omega)^(-sigma) = (1-(tau)*(kappa))*exp(llambda)*p*gamma*exp(lY_i)/exp(lh_i);

%FOC wrt capital
exp(llambda)*(1+phi*(exp(lK)-exp(lK(-1)))) = exp(llambda(+1))*beta*((1-(tau))*(1-theta)*exp(lY_f(+1))/exp(lK)+1-delta+phi*exp((lK(+1))-exp(lK)));

%FOC wrt gov capital
%exp(llambda)*(1+phi*(exp(lKg)-exp(lKg(-1)))) = exp(llambda(+1))*beta*((1-(tau))*(1-theta)*exp(lY_f(+1))/exp(lKg)+1-delta+phi*exp((lKg(+1))-exp(lKg)));

%Interest rate
r = r_star + R*(exp(exp(ld)-d_bar)- 1);

%Budget Constraint
exp(ld) = (1+(r(-1)))*exp(ld(-1))-exp(lY_f)-p*exp(lY_i)+exp(lC_f)+p*exp(lC_i) + exp(lx)+(phi/2)*(exp(lK)-exp(lK(-1)))^2+(1-etta)*exp(lG);

%Formal Production
(lY_f) = A_f+ theta*lh_f +(1-theta-pi)*lK(-1)+(1-pi)*(lKg(-1));
%exp(lY_f) = exp(exp(lA_f))*exp(lh_f)^(theta)*(exp(lK(-1)+exp(lKg)))^(1-theta);

%Informal Production
lY_i = A_i+ gamma*lh_i;

%Investment
exp(lx) = exp(lK) - (1-delta)*exp(lK(-1));

%Government Investment
etta*exp(lG) = exp(lKg) - (1-delta)*exp(lKg(-1));

%Government Budget Constraint
tau*kappa*p*exp(lY_i) + tau*exp(lY_f) = exp(lG);

%Total Labor
exp(lh_f) + exp(lh_i) = exp(lH);

%Informal Sector Identity
lY_i = lC_i ;

%Trade Balance
exp(ltb) = exp(lY_f) - exp(lC_f) - (1-etta)*exp(lG) - exp(lx) - (phi/2)*(exp(lK))^2+phi*exp(lK)*exp(lK(-1))-(phi/2)*exp(lK(-1))^2;

%Current Account
exp(lca) = -(exp(ld)-exp(ld(-1)));

%Stochastic Process for Shocks
A_f = rho_f*(A_f(-1))+ ef;
A_i = rho_i*(A_i(-1)) + ei;

%G-G_ss = rho_G*(G(-1)-G_ss)+eg;

kappa-kappa_ss = rho_kappa*(kappa(-1)-kappa_ss)+ekappa;

tau - tau_ss = rho_tau*(tau(-1)-tau_ss) - etau;

end;


%------------------------------------------------------------------
% 4. Computation
%------------------------------------------------------------------

G_ss = .1;

p_ss = 1.5;
h_ss = 1.1;
%d_bar = 0.007442;
%d_ss = d_bar;

initval;

tau    		= 0.2;
kappa  		 =0.25;
%U      =		 -1.28132;
lY      	=	log(2.8796);
lY_f    	=	 log(1.4447);
lY_i    	=	 log(0.513449);
lC_f    	=	 log(0.772642);
lC_i    	=	 log(0.513449);
ld      	=	 log(.210071);
lK      	=	 log(4.65022);
lH      	=	 log(0.451675);
lh_f    	=	 log(0.191571);
lh_i    	=	 log(0.260104);
llambda 	=	 log(1.1621);
lx      	=	 log(0.465022);
r      	=	 0.0416667;
%A_f    	=	 0;
A_i    	=	 0;
ltb     	=	log(0.00875297);

lG      	=	 log(0.396557);
p      	=	 2.79464;
lC      	=	 log(0.592401);

end;
G_ss = .1;

p_ss = 1.5;
h_ss = 1.1;
d_bar = 0.007442;
d_ss = d_bar

shocks;

%var ef = (sigma_f)^2;
%var ei = (sigma_i)^2;
%var eg = (sigma_g)^2;
var ekappa  = (sigma_kappa)^2;
var etau = (sigma_tau)^2;

%var kappa;
%periods 10:19;
%values .5;
end;

steady;

check;

stoch_simul(hp_filter = 100, order = 1, irf = 100) lY lY_f lY_i lC lC_f lC_i lH lh_f lh_i ltb lca ld r lK lx p lG tau kappa;

%forecast (periods = 100);











