//* ******  Gali Monacelli (2016) * ******//

/*  Here there are two monetary regimes: e = 0, pi_H = 0 */


%----------------------------------------------
% Define one of two cases:
% - Inflation targeting : pi_H = 0 for all t        corresponds to CU = 0
% - Currency Union: e = 0 for all t
@#define CU = 0
%----------------------------------------------

%----------------------------------------------
% Preamble
%----------------------------------------------


var 
y           ${y}$       
c           
n
a
z
i
w
y_n         % natural output
c_n
n_n
omega       % omega = w - p , real wage
omega_n
y_gap
c_gap
omega_gap
p
p_H
pi
pi_H
pi_w
s
e
s_n
s_gap
i_star
z_1_star
z_2_star
tau
r
i_ann
r_ann
;

varexo eps_z, eps_z_1, eps_z_2, eps_a, eps_tau;

parameters alpha beta phi upsilon rho_a rho_z rho_z_1 rho_z_2 epsilon_p epsilon_w theta_p theta_w rho_tau
// 16 parameters (without the lambdas)
// rho is a constant and thus omitted
;
%----------------------------------------------
% Parameterization: Follows GM(2016)
%----------------------------------------------
alpha = 0.26;
beta = 0.99;
phi = 2.2;
upsilon = 0.3;
rho_a = 0.9;
rho_z = 0.9;
rho_z_1 = 0.9;
rho_z_2 = 0.9;
epsilon_p = 3.8;
epsilon_w = 4.3;
theta_p = 0.8;
theta_w = 0.8;
rho_tau = 0.9;      % = 1 for permanent cut
 

%----------------------------------------------
% Model
%----------------------------------------------

model(linear);

//Composite parameters
#lambda_p = ((1-theta_p)*(1-beta*theta_p))/theta_p*(1-alpha)/(1-alpha + alpha*epsilon_p);
#lambda_w = ((1-theta_w)*(1-beta*theta_w))/(theta_w*(1+epsilon_w*phi));

@#if CU == 0
    pi_H = 0; //inflation targeting
    @#else
      e = 0;  //currency union
@#endif

//AD block
y = (1-upsilon)*c + upsilon*(2-upsilon)*s + upsilon*z_1_star;
c=(1-upsilon)*s+z-z_2_star;
c=c(+1)-(i-pi(+1))+(1-rho_z)*z;
s=e-p_H;
n = (1/(1-alpha))*(y-a);

//AS block
pi_H = beta*pi_H(+1) + (lambda_p*alpha/(1-alpha))*y_gap + lambda_p*omega_gap + lambda_p*upsilon*s_gap + lambda_p*tau;
pi_H = p_H - p_H(-1);
pi = p - p(-1);
p = p_H + upsilon*s;
pi_w = beta*pi_w(+1) + ((lambda_w*phi)/(1-alpha))*y_gap + lambda_w*c_gap - lambda_w*omega_gap;
pi_w = w - w(-1);
omega = w - p;


//gaps
c_gap = c - c_n;
y_gap = y - y_n;
omega_gap = omega - omega_n;
s_gap = s - s_n;


//natural values
n_n = (upsilon/(1+phi))*(z_1_star + z_2_star - z) - (1/(1+phi))*tau;
y_n = a + (1-alpha)*n_n;
s_n = a - z + z_2_star - tau - (alpha + phi)*n_n;
c_n = z +(1-upsilon)*s_n -z_2_star;
omega_n = a - alpha*n_n - tau - upsilon*s_n;


//exogenous processes
a = rho_a*a(-1) + eps_a;
z = rho_z*z(-1) + eps_z;
z_1_star = rho_z_1*z_1_star(-1) + eps_z_1;
z_2_star = rho_z_2*z_2_star(-1) + eps_z_2;
i_star =(1-rho_z_2)*z_2_star;       
tau = rho_tau*tau(-1) - eps_tau;        


r = i - pi_H(+1); //real interest rate

//annualized values
i_ann = 4*i;
r_ann = 4*r;


end;

steady;
check;

%----------------------------------------------
% Define shock variances
%----------------------------------------------
//one-off payroll tax cut
shocks;
var eps_tau = 1;

end;


stoch_simul(order=1,periods=1000, irf = 25, nocorr, nomoments)n i_ann r_ann s;