% Model with Unemployment

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

var Y K N C I G Ni X M U Ipsilon Q S K_p LAMBDA R Pi Lambda R_k Q_k P_w A W Ji Mu E W_o B Y_n Zeta;
varexo Epsilon_z Epsilon_i Epsilon_b Epsilon_w Epsilon_p Epsilon_r Epsilon_g Epsilon_eta;

parameters alpha y_c y_i y_g y_ni y_x sigma rho n_s u_s xi beta eta_ni eta_k gamma_z delta nu_a nu_w 
            nu_lambda ji lambda gamma x_s mu fi_a fi_s fi_x fi_b fi_ji s_s gamma_b gamma_o gamma_f i_b
            i_o i_f rho_s rho_pi rho_y k_s a_s w_s b_s y_s c_s g_s i_s kappa r_k q_k sigma_m ipsilon p_w eta e
            psi_ni nu h_s tau tau_1 tau_2 Fi gamma_p dseta fi_eta lambda_p
            tau_p epsilon_p h;

%------------------------------------------------------------------------
% 2. Calibration (Prior values)
%------------------------------------------------------------------------

beta    = 0.99;
delta   = 0.025;
alpha   = 0.33;
y_g     = 0.2;
xi      = 10;
sigma   = 0.5;
sigma_m = 0.5;      %chequear
rho     = 0.895;
s_s     = 0.95;
p_w     = 1;        %chequear 
h       = 0.5;      %chequear

eta_k     = 4;
psi_ni    = 0.5;
eta       = 0.5;
lambda    = 0.75;
lambda_p  = 0.66;
gamma     = 0.5;
epsilon_p = 1.15;
rho_s     = 0.75;
gamma_z   = 1.025;

%Consuption and saving
1 = (beta/gamma_z)*(1 - delta + r_k);

%Capital/employment ratio
r_k = alpha*((k_s/n_s)^(-(1 - alpha)));

%Marginal product
a_s = (1 - alpha)*((k_s/n_s)^alpha);

%Investment
q_k = 1;

%Rates
x_s = 1 - rho;

%Flows
n_s = (s_s/x_s)/(1 + (s_s/x_s));

%Unemployment
u_s = 1 - n_s;

%Matching
s_s = sigma_m*(u_s^(sigma-1))*(ipsilon^(1 - sigma));

%Hiring
1 = (p_w*a_s)/(kappa*x_s) - (w_s/(kappa*x_s)) + beta*(x_s/2) + beta*rho;

%Wage
w_s = ji*(a_s + beta*(kappa/2)*(x_s^2) + beta*kappa*s_s*x_s) + (1 - ji)*b_s;

ji = (eta/(eta + (1 - eta)*mu/e));

mu = 1/(1 - lambda*beta);

e = 1/(1 - rho*lambda*beta);

%Resource constraint
1 = c_s/y_s + g_s/y_s + i_s/y_s + (kappa*(x_s^2)*n_s/y_s);

1/y_s = ((k_s)^(-alpha))*(n_s^(alpha-1));

i_s = ((gamma_z - (1 - delta))*((k_s/n_s)^(1 - alpha)))*y_s;

%From the complete lolinear model
y_c = c_s/y_s;

y_i = i_s/y_s;

y_g = g_s/y_s;

y_ni = r_k*(k_s/y_s);

y_x = (kappa/2)*((x_s^2)*n_s/y_s);

xi = (1 - delta)/gamma_z;

eta_ni = (1 - psi_ni)/psi_ni;

nu = (kappa*x_s)^(-1);

nu_a = nu*p_w*a_s;

nu_w = nu*w_s;

nu_lambda = beta*(1 + beta)/2;

fi_a = ji*p_w*a_s*(w_s^(-1));

fi_x = ji*beta*kappa*(x_s^2)*(w_s^(-1));

fi_b = (1 - ji)*b_s*(w_s^(-1));

fi_s = (1 - ji)*beta*h_s*(w_s^(-1));

fi_ji = ji*((1 - ji)^(-1))*kappa*x_s*(w_s^(-1));

fi_eta = fi_ji*(1 - ji)*((1 - eta)^(-1));

gamma_b = (1 + tau_2)*(Fi^(-1));

gamma_o = dseta*(Fi^(-1));

gamma_f = (tau*(lambda^(-1)) - tau_1)*(Fi^(-1));

Fi = (1 + tau_2) + dseta + (tau*(lambda^(-1)) - tau_1);

dseta = (1 - lambda)*(1 - tau)*(lambda^(-1));

tau_1 = ((nu_w*mu*fi_x) + fi_ji*(1 - ji)*(x_s*beta*lambda)*(nu_w*mu)*mu*(rho*beta) + fi_s*Gamma)*(1 - tau);

tau_2 = -(nu_w*mu)*fi_ji*(1-ji)*(x_s*beta*lambda)*mu*(1 - tau);

Gamma = (1 - eta*x_s*beta*lambda*mu)*(eta^(-1))*mu*nu_w;

i_b = gamma_p*(Fi_p^(-1));

i_o = (dseta_p/tau_p)*(Fi_p^(-1));

i_f = beta*(Fi_p^(-1));

Fi_p = 1 + beta*gamma_p;

dseta_p = (1 - lambda_p)*(1 - lambda_p*beta)*(lambda_p^(-1));

tau_p = 1 + (epsilon_p - 1)*xi;

%------------------------------------------------------------------------
% 3. Model
%------------------------------------------------------------------------

model(linear);

%Technology
Y = alpha*K(-1) + (1 - alpha)*N;

%Resource constraint
Y = y_c*C + y_i*I + y_g*G + y_ni*Ni + y_x*(2*X + N(-1));

%Matching
M = sigma*U + (1 - sigma)*Ipsilon;

%Employment dynamics
N = N(-1) + (1 - rho)*N;

%Transition probabilities
Q = M - Ipsilon;

S = M - U;

%Unemployment
U = -(n_s/u_s)*N(-1);

%Effective capital
K(-1) + Epsilon_z = Ni + K_p(-1);    %K(-1)?

%Physical capital dynamics
K_p = xi*(K_p(-1) - Epsilon_z) + (1 - xi)*(I + Epsilon_i);

%Aggregate vacancies
X = Q + Ipsilon - N(-1);

%Consumption-saving
0 = EXPECTATION(-1)(LAMBDA(+1)) + (R - EXPECTATION(-1)(Pi(+1))) - EXPECTATION(-1)(Epsilon_z(+1));

%Marginal utility
(1 - h/gamma_z)*(1 - beta*(h/gamma_z))*Lambda = (h/gamma_z)*(C(-1) - Epsilon_z) - (1 + beta*((h/gamma_z)^2))*C + beta*(h/gamma_z)*EXPECTATION(-1)(C(+1) + Epsilon_z(+1)) + (1 - (h/gamma_z))*(Epsilon_b - beta*(h/gamma_z)*EXPECTATION(-1)(Epsilon_b(+1)));

LAMBDA = Lambda/Lambda(-1);

%Capital utilization
Ni = eta_ni*R_k;

%Investment
I = (1/(1 + beta))*(I(-1) - Epsilon_z) + (((1/(eta_k*(gamma_z)^2))/(1 + beta))*(Q_k + Epsilon_i)) + ((beta/(1 + beta))*EXPECTATION(-1)(I(+1) + Epsilon_z(+1)));

%Capital renting
P_w + Y - K(-1) = R_k

%Tobinīs q 
Q_k = beta/gamma_z*(1 - delta)*EXPECTATION(-1)(Q_k(+1)) + (1 - (beta/gamma_z)*(1 - delta))*EXPECTATION(-1)(R_k(+1)) - (R - EXPECTATION(-1)(Pi(+1)));

%Aggregate hiring rate
X = nu_a*(P_w + A) - nu_w*W + (nu_lambda*EXPECTATION(-1)(LAMBDA(+1))) + (beta*(EXPECTATION(-1)(X(+1))));

%Marginal product of labor
A = Y - N;

%Weight in Nash bargaining
Ji = -(1 - ji)*(Mu - E);

E = (rho*lambda*beta)*EXPECTATION(-1)(LAMBDA(+1) - Pi(+1) + gamma*Pi + E(+1) - Epsilon_z(+1));

Mu = (x_s*lambda*beta)*EXPECTATION(-1)(X(+1)) - ((x_s*lambda*beta)*(nu_w*mu)*mu*EXPECTATION(-1)(W + gamma*Pi - Pi(+1) - Epsilon_z(+1) - W(+1))) + ((lambda*beta)*EXPECTATION(-1)(Mu(+1) + LAMBDA(+1) + gamma*Pi - Pi(+1) - Epsilon_z(+1)));

%Spillover-free target wage
W_o = fi_a*(P_w + A) + (fi_s + fi_x)*EXPECTATION(-1)(X(+1)) + fi_s*EXPECTATION(-1)(S(+1)) + fi_b*B + (fi_s + (fi_x/2))*EXPECTATION(-1)(LAMBDA(+1)) + fi_ji*(Ji - (rho - s_s)*beta*Ji(+1)) + Epsilon_w;

Epsilon_w = fi_eta*(1 - (rho - s_s)*beta*(rho^eta))*Epsilon_eta;

%Aggregate wage
W = gamma_b*(W(-1) - Pi + (gamma*Pi(-1) - Epsilon_z)) + gamma_o*W_o + gamma_f*EXPECTATION(-1)(W(+1) + Pi(+1) - gamma*Pi + Epsilon_z(+1));

%Phillips curve
Pi = i_b*Pi(-1) + i_o*(P_w + Epsilon_p) + i_f*EXPECTATION(-1)(Pi(+1));

%Taylor rule
R = rho_s*R(-1) + (1 - rho_s)*(rho_pi*Pi + rho_y*(Y - Y_n)) + Epsilon_r;

%Goverment spending
G = Y + ((1 - y_g)/y_g)*Epsilon_g;

%Market tightness
Zeta = Ipsilon - U;

%Benefits
B = K_p;

%Flex Economy

Y_n = 0;

end;

%---------------------------------------------------------------------
% 4. Shocks and Simulation
%---------------------------------------------------------------------

check;

shocks;
var Epsilon_z;      staderr 0*1;
var Epsilon_i;      staderr 0*1;
var Epsilon_w;      staderr 0*1;
var Epsilon_eta;    staderr 0*1;
var Epsilon_p;      staderr 0*1;
var Epsilon_r;      staderr 0*1;
var Epsilon_g;      staderr 0*1;

end;

stoch_simul(order=1,periods=1000);