var c_b c_d h_b h_d N_H N_C p lambda_b mc inflation lambda_d mu C H IH i_hd i_hb B r z q Lambda pi w_h w_c n_dh n_bh n_dc n_bc Y_H Y_C er;
varexo ez eer;

parameters phi,theta,beta_b,beta_d,rho,m,sigma,delta,sig_ez,mi,epsilon_b,epsilon_d,j,lambda_bSS, c_bSS, lambda_dSS, c_dSS, h_bSS, h_dSS, N_HSS, N_CSS, pSS, muSS, CSS, HSS, IHSS, i_hdSS, i_hbSS, BSS, DSS, qSS, piSS, w_hSS, w_cSS, n_dhSS, n_bhSS, n_dcSS, n_bcSS, Y_HSS, Y_CSS;




set_param_value('phi',phi)
set_param_value('sigma',sigma)
set_param_value('theta',theta)
set_param_value('delta',delta)
set_param_value('beta_b',beta_b)
set_param_value('beta_d',beta_d)
set_param_value('rho',rho)
set_param_value('sig_ez',sig_ez)
set_param_value('m',m)
set_param_value('epsilon_d',epsilon_d)
set_param_value('epsilon_b',epsilon_b)
set_param_value('mi',mi)
set_param_value('j',j)
set_param_value('lambda_bSS',lambda_bSS) 
set_param_value('c_bSS',c_bSS) 
set_param_value('lambda_dSS',lambda_dSS) 
set_param_value('c_dSS',c_dSS) 
set_param_value('h_bSS',h_bSS) 
set_param_value('h_dSS',h_dSS) 
set_param_value('N_HSS',N_HSS) 
set_param_value('pSS',pSS) 
set_param_value('muSS',muSS)
set_param_value('CSS',CSS)
set_param_value('HSS',HSS) 
set_param_value('IHSS',IHSS)
set_param_value('i_hdSS',i_hdSS) 
set_param_value('i_hbSS',i_hbSS) 
set_param_value('BSS',BSS) 
set_param_value('qSS',qSS) 
set_param_value('piSS',piSS) 
set_param_value('w_hSS',w_hSS)  
set_param_value('n_dhSS',n_dhSS)
set_param_value('n_bhSS',n_bhSS) 
set_param_value('n_dcSS',n_dcSS) 
set_param_value('n_bcSS',n_bcSS)
set_param_value('Y_HSS',Y_HSS) 

N_HSS=Y_HSS;
N_CSS=CSS;
Y_CSS=CSS;
w_cSS=z;
DSS=BSS;

model;

%---------------------------- BORROWERS ----------------------------------%

% 1) FOC c_b \ C_b
 
lambda_b = (c_b*p)^(-1);
 
% 2) EULER EQUATION BORROWERS (FOC_B) \ mu

mu = lambda_b - beta_b*r(+1)*lambda_b(+1);

% 3) BINDING COLLATERAL CONSTRAINT \ B

B = h_b*m*q;

% 4) FOC h_b \ lambda_b

q*lambda_b = j*h_b^(-1) + m*q*mu + beta_b*lambda_b(+1)*q(+1)*(1-delta);

% 5) FOC n_bc \ n_bc

w_c*lambda_b = (n_bc^(epsilon_b))*(n_bc^(1+epsilon_b) + n_bh^(1+epsilon_b))^(((1+mi)/(1+epsilon_b))-1);

% 6) FOC n_bh \ n_bh

w_h*lambda_b = (n_bh^(epsilon_b))*(n_bc^(1+epsilon_b) + n_bh^(1+epsilon_b))^(((1+mi)/(1+epsilon_b))-1);

% 7) BUDGET CONSTRAINT \ h_b

c_b*p + q*(h_b - (1-delta)*h_b(-1)) = w_h*n_bh + w_c*n_bc - r*B(-1) + B ;

% 8) HOUSING INVESTMENT BORROWERS \ i_hb

i_hb = (h_b - (1-delta)*h_b(-1));

%---------------------------- DEPOSITORS ----------------------------------%

% 9) FOC c_d \ c_d

lambda_d = (c_d*p)^(-1);
 
% 10) EULER EQUATION (FOC_D) \ lambda_d

lambda_d = beta_d*r(+1)*lambda_d(+1);

% 11) FOC h_d \ h_d

q*lambda_d = j*h_d^(-1) + beta_d*lambda_d(+1)*q(+1)*(1-delta);

% 12) FOC n_dc \ n_dc

w_c*lambda_d = (n_dc^(epsilon_d))*(n_dc^(1+epsilon_d) + n_dh^(1+epsilon_d))^(((1+mi)/(1+epsilon_d))-1);

% 13) FOC n_dh \ n_dh

w_h*lambda_d = (n_dh^(epsilon_d))*(n_dc^(1+epsilon_d) + n_dh^(1+epsilon_d))^(((1+mi)/(1+epsilon_d))-1);

% 14) BUDGET CONSTRAINT \ D

c_d*p + q*(h_d - (1-delta)*h_d(-1)) + B = w_h*n_dh + w_c*n_dc + r*B(-1) + pi;

% 15) STOCHASTIC DISCOUNT FACTOR \ Lambda

Lambda = beta_d*lambda_d(+1)*lambda_d^(-1);

% 16) HOUSING INVESTMENT DEPOSITORS \ i_hd

i_hd = (h_d - (1-delta)*h_d(-1));

%---------------------------- FIRMS --------------------------------------%

% 17) PROFITS \ pi

 pi = q*Y_H - w_h*N_H - q*Y_H*(theta/2)*(q*q(-1)^(-1)-1)^(2);

% 18) WAGE IN CONSUMPTION SECTOR \ w_c

w_c = z;

% 19) OPTIMAL PRICE SETTING \ q

((1 - sigma) + sigma*mc) = theta*(inflation - 1)*inflation - theta*((Lambda(+1)/Lambda)*(q(+1)/q)*(Y_H(+1)/Y_H)*(inflation(+1) - 1)*inflation(+1));

%((1 - sigma) + sigma*mc) = theta*(inflation - 1)*inflation - theta*beta_d*((lambda_d(+1)/lambda_d)*(q(+1)/q)*(Y_H(+1)/Y_H)*(inflation(+1) - 1)*inflation(+1));

% 20) SUPPLY H \ Y_H

 Y_H = z*N_H;

% 21) SUPPLY C \ Y_C

 Y_C = z*N_C;

%------------------------- AGGREGATE AND INTERESTING VARIABLES ---------------------------%

% 22) AGGREGATE LABOR IN HOUSING SECTOR \ N_H

N_H = n_bh + n_dh;

% 23) AGGREGATE LABOR IN CONSUMPTION SECTOR \ N_C

N_C = n_bc + n_dc;

% 24) AGGRGATE CONSUMPTION \ C

 C = c_d + c_b;

% 25) AGGREGATE HOUSING \ H

 H = h_d + h_b;

% 26) AGGREGATE HOUSING INVESTMENT \ IH

 IH = i_hb + i_hd;

% 27) MARGINAL COST \ mc

 %mc = w_h/q;
mc = w_h/(q*z);

% 28) INFLATION IN HOUSING PRICES \ inflation

 inflation = q/q(-1);

%------------------------------- MARKET CLEARING -------------------------%

% 29) MARKET CLEARING HOUSING

%Y_C = C;

Y_H = IH + Y_H*theta*4^(-1)*(q*q(-1)^(-1)-1)^(2);

% 30) CLEARING CREDIT MARKET: LOANS = DEPOSITS


%--------------------- EXOGENOUS PROCESS ---------------------------------%

% 32) PRODUCTIVITY LAW OF MOTION \ z

z = (1-rho) + rho*z(-1) + sig_ez*ez;

% 33) INTEREST RATE RULE \ r

%r = r(-1)^(0.5)*(inflation*(IH/(C+IH))+p/p(-1)*(C/(C+IH)))^(phi)*er/beta_d^(0.5);

r = (inflation*(IH/(C+IH))+p/p(-1)*(C/(C+IH)))^(phi)*er/beta_d;

% 34) INTEREST RATE SHOCK \ er

er = (1-rho) + rho*er(-1) + eer;

end;

%------------------------- SS INITIAL VALUES -------------------------%

initval;

r =1/beta_d;
inflation = 1;
Lambda=beta_d;
er = 1;
z=1;
mc=sigma/(sigma-1);
w_c = 1;
p=pSS;
c_b = c_bSS;
h_b = h_bSS;
lambda_b = lambda_bSS; 
mu = muSS ;
q = qSS;
w_h = w_hSS;
n_bc =n_bcSS;
n_bh =n_bhSS;
c_d=c_dSS;
lambda_d =lambda_dSS;
h_d=h_dSS;
H = HSS;
C =CSS;
B = BSS;
n_dc=n_dcSS;
n_dh=n_dhSS;
N_H=N_HSS;
N_C=CSS;
Y_C=CSS;
Y_H = Y_HSS;
i_hb=i_hbSS;
i_hd =i_hdSS;
IH = IHSS;
pi=piSS;


end;

steady(solve_algo=0);
check;

shocks;
var ez; stderr sig_ez;
var eer; stderr sig_ez;
end;

stoch_simul(periods=2100, order=1, drop=400, IRF=20,replic=1000);