% Basic DSGE Model for fiscal rules
% With New Keynesian Features
% + Sluggish prices (forward)
% + Inflation smoothing (forward & backward)

close all;
%format long

%----------------------------------------------------------------
% 1. Defining variables
%----------------------------------------------------------------

var 	y
	c
	r
	pi
	w
	l
	a b d q g tc tw tk tr k z i u_q u_b u_l;

varexo	e_a e_q e_r e_g e_tc e_tw e_tk e_l e_b;

parameters	beta alpha sigmaC sigmaL
			rho_r phi_ry phi_rpi theta_p gamma_p
			rho_a  rho_q rho_b rho_l rho_g phi_gb phi_gy rho_tc phi_tcb phi_tcy phi_tctw phi_tctk rho_tw phi_twb phi_twy phi_twtc phi_twtk rho_tk phi_tkb phi_tky phi_tktc phi_tktw phi_tr psi PI Q R TC TW TR CY GY WY BY DY KY delta TK Z phi_rb
			;

%----------------------------------------------------------------
% 2. Calibration
%----------------------------------------------------------------
@#include "calibration_NK_fiscalrules.mod"


%----------------------------------------------------------------
% 3. Model
%----------------------------------------------------------------
model(linear);
    %% Household
	% Euler
	sigmaC*(c(+1)-c)=r+q+u_b(+1)-u_b-pi(+1)+TC/(1+TC)*(tc-tc(+1));
	% Hours supply
	w = sigmaL*l + sigmaC*c + TW/(1-TW)*tw + TC/(1+TC)*tc + u_l;
	k = (1-delta)*k(-1) + delta*i;
	l = -w + z + k(-1);
	1/beta*(r+q-pi(+1)) = (Z*(1-TK)*z(+1))/(1-delta+(1-TK)*Z) - (Z*TK*tk(+1))/(1-delta+(1-TK)*Z);

    % Intermediary firms
	% Production function
	y = a + alpha*k(-1) + (1-alpha)*l;
	y = CY*c+GY*g+delta*KY*i;
	% Price dynamics NKPC backward
	pi = gamma_p/(1+beta*gamma_p)*pi(-1) + beta/(1+beta*gamma_p)*pi(+1) + ((1-theta_p)*(1-theta_p*beta))/((1+beta*gamma_p)*theta_p)*(alpha*z+(1-alpha)*w-a);
	
	% monetary policy
	r = rho_r*r(-1) + (1-rho_r)*(phi_rpi*(pi(-1)) +  phi_ry*(y(-1))) - e_r;
	
	% Debt policy rule
	BY*(b-b(-1)) = DY*d;
	DY*d = 1/beta*(GY*g-TC*CY*(c+tc)-TW*WY*(tw+w+l)-Z*KY*TK*(tk+z+k(-1))+TR*tr)+BY*((1/beta-1)*b(-1)+r+q-pi(+1));
	q = psi*(DY*d(+1)) + u_q;
	% Tax policy rules
	tc = rho_tc*tc(-1) + (1-rho_tc)*(phi_tcb*b(-1)*BY + phi_tcy*y) + e_tc;
	tw = rho_tw*tw(-1) + (1-rho_tw)*(phi_twb*b(-1)*BY + phi_twy*y) + e_tw;
	tk = rho_tw*tk(-1) + (1-rho_tw)*(phi_tkb*b(-1)*BY + phi_tky*y) + e_tk;
	% Expenditure policy rule
	g = rho_g*g(-1) - (1-rho_g)*(phi_gb*b(-1)*BY + phi_gy*y) - e_g;
	% Stabilization rule transfers
	tr = -phi_tr*(b(-1)*BY);

    % Exogenous shocks
	a = rho_a*a(-1) + e_a;
	u_q = rho_q*u_q(-1) + e_q;
	u_b = rho_b*u_b(-1) - e_b; 
	u_l = rho_l*u_l(-1) + e_l;

end;

%----------------------------------------------------------------
% 4. Computation
%----------------------------------------------------------------
check;

initval;
y = 0;
k = 0;
c = u_b;
l = 0;
w = 0;
r = 0;
z = 0;
b = 0;
pi = 0;
a = 0;
d = 0;
q = 0;
g = 0;
tc = 0;
tw = 0;
tk = 0;
tr = 0;
i = 0;
u_q = 0;
u_b = -20;
u_l = 0;
e_a = 0;
e_q = 0;
e_r = 0;
e_g = 0;
e_tc = 0;
e_tw = 0;
e_tk = 0;
e_b = 0;
e_l = 0;
end;

shocks;
var e_a;  stderr 1;
var e_q;  stderr 1;
var e_r;  stderr 1;
var e_g;  stderr 1;
var e_tc;  stderr 1;
var e_tw;  stderr 1;
var e_tk;  stderr 1;
var e_b;  stderr 1;
var e_l;  stderr 1;
end;

steady;

stoch_simul(order=1, irf=40, irf_shocks = (e_b, e_r, e_g, e_tc, e_tw, e_tk));