Var A_e b b_i c c_e c_i c_obs c_p d_p ee_h h_i h_p I I_obs k lambda l_i l_p m_i n pie q_h qh_obs q_k r r_k r_n r_obs rps v w_i w_p x y; 

varexo e_a_e e_mi e_g e_h e_rps e_rn e_p;


parameters beta_p beta_i phi alpha mu r_ss delta thet phi_x gamma iota pi_ss epsilon gamma_p phi_r phi_pie phi_y rho_a rho_mi rho_g rho_eeh rho_rps m_ibar lev nk_ss chi x_ss eeh_bar;   


//% local model parameters: dbar kappa Gbar;

//%Calibrated parameters



x_ss = 1.15;
beta_p = 0.9925;
r_ss = 1/beta_p;
beta_i = 0.97;
delta = 0.025;
lev = 1.02;
alpha = 0.25;
phi = 1.5;
m_ibar = 0.8; 
eeh_bar = 0.12;
phi_x = 1.5;
gamma_p = 0.6;
pi_ss = 1;
%dbar=1.711;


//%estimated parameters (mean):
mu = 0.3;
gamma = 0.94;
iota = 0.691;
phi_r = 0.61;
phi_pie = 1.40;
phi_y = 0.512;

nk_ss = 0.481;
chi = 0.05; 
thet = 0.836;
epsilon=(1-thet)*(1-beta_p*thet)/thet;



//%shocks parameters (mean)

rho_a = 0.5425; 
rho_mi = 0.9; 
rho_g = 0.95;
rho_eeh = 0.958;
rho_rps = 0.8; 


model;

#sp = ((1-alpha)*(1-mu) + x_ss - 1)/x_ss; %income share of patient household
#si = (1-alpha)*mu/x_ss; %income share of impatient household
#llr_n=1/beta_p;
#lr=llr_n;
#lci_Y = (1-beta_i-m_ibar*(beta_p-beta_i))*si/(1-beta_i-m_ibar*(beta_p-beta_i)+eeh_bar*m_ibar*(1-beta_p));
#lqhi_Y = eeh_bar*si/(1-beta_i-m_ibar*(beta_p-beta_i-eeh_bar*(1-beta_p)));
#lr_k=lev/beta_p;
#lk_Y=alpha/(x_ss*(lr_k-(1-delta)));
#lk_n=1/nk_ss;
#lcp_Y  = sp-(1-1/beta_p)*((1-gamma_p)*beta_p*m_ibar*lqhi_Y+lk_Y-lk_Y*nk_ss)/gamma_p;
#lqhp_Y = eeh_bar*lcp_Y/(1-beta_p);
#ll_p=((1-mu)*(1-alpha)/(x_ss*lcp_Y))^(1/(1+phi));
#ll_i=(mu*(1-alpha)/(x_ss*lci_Y))^(1/(1+phi));
#ly=(lk_Y)^(alpha/1-alpha)*(ll_i)^mu*(ll_p)^(1-mu);
#lbi_Y = beta_p*m_ibar*lqhi_Y;
#lb_Y= lk_Y-lk_Y*nk_ss;
#lk=lk_Y*ly;
#lln=nk_ss*lk;
#lc_e=(1-gamma)*lln/gamma;
#lb_i=lbi_Y*ly;
#lc_p=lcp_Y*ly;
#lc_i=lci_Y*ly;
#lb=lb_Y*ly;
#lw_p=1/x_ss*(1-mu)*(1-alpha)*ly/ll_p;
#lw_i=1/x_ss*mu*(1-alpha)*ly/ll_i;
#ld_p=((1-gamma_p)*lb_i+lb)/gamma_p;
#llambda=(beta_p-beta_i)/lc_i;
#lh_i=lc_i*(1-beta_p)/((1-gamma_p)*lc_i*(1-beta_p)+gamma_p*lc_p*(1-beta_i-m_ibar*(beta_p-beta_i)));
#lh_p=(1-lh_i*(1-gamma_p))/gamma_p;
#lq_h=lqhi_Y*ly/lh_i;
#lI=delta*lk_Y*ly;
#lv=lev*lr*lln;
#dbar=(lln-gamma*lv)/(1-gamma);
#kappa=lev*(nk_ss)^chi;
#lc=gamma_p*lc_p+(1-gamma_p)*lc_i;
#Gbar=ly-lc-lI;
#c_ss=log(lc);
#I_ss=log(lI);
#qh_ss=log(lq_h);


//% Patient Households
1=beta_p*exp(c_p)*exp(r_n)/(exp(c_p(+1))*exp(pie(+1)));
exp(w_p)=exp(c_p)*(exp(l_p))^phi;
exp(q_h)/exp(c_p)=exp(ee_h)/exp(h_p)+beta_p*exp(q_h(+1))/exp(c_p(+1));
exp(c_p)+exp(q_h)*(exp(h_p)-exp(h_p(-1)))+ exp(d_p)=exp(w_p)*exp(l_p)+exp(r_n(-1))*exp(d_p(-1))/exp(pie)+(1-1/exp(x))*exp(y);

//% Impatient Households

1/exp(c_i)= beta_i*exp(r_n)/(exp(pie(+1))*exp(c_i(+1)))+exp(lambda)*exp(r_n);
exp(w_i)=exp(c_i)*(exp(l_i))^phi;
exp(q_h)/exp(c_i)=beta_i*exp(q_h(+1))/exp(c_i(+1))+exp(ee_h)/exp(h_i)+exp(lambda)*exp(m_i)*exp(q_h(+1))*exp(pie(+1));
exp(b_i)=exp(m_i)*exp(q_h(+1))*exp(h_i)*exp(pie(+1))/exp(r_n);
exp(c_i)+exp(q_h)*(exp(h_i)-exp(h_i(-1)))+exp(r_n(-1))*exp(b_i(-1))/exp(pie)=exp(w_i)*exp(l_i)+exp(b_i);

//% Entrepreneurs

exp(y)=exp(A_e)*(exp(k(-1)))^alpha*(exp(l_i))^(mu*(1-alpha))*(exp(l_p))^((1-mu)*(1-alpha));

exp(w_i)=mu*(1-alpha)*exp(y)/(exp(x)*exp(l_i));
exp(w_p)=(1-mu)*(1-alpha)*exp(y)/(exp(x)*exp(l_p));
exp(b)=exp(q_k)*exp(k)-exp(n); 

//% Net worth
exp(k)=(1-delta)*exp(k(-1))+(1-phi_x*(exp(I)/exp(I(-1))-1)^2)*exp(I);
1=exp(q_k)*(1-phi_x*(exp(I)/exp(I(-1))-1)^2-2*phi_x*(exp(I)/exp(I(-1))-1)*exp(I)/exp(I(-1)))+exp(pie(+1))/exp(r_n(+1))*exp(q_k(+1))*2*phi_x*(exp(I(+1))/exp(I)-1)*(exp(I(+1))/exp(I))^2;
exp(r_k)=((1/x_ss*alpha*exp(y)/exp(k(-1)))+(1-delta)*exp(q_k))/exp(q_k(-1));
exp(r_k)=kappa*(exp(n(-1))/(exp(k(-1))*exp(q_k(-1))))^(-chi)*exp(r)*exp(rps);
%exp(r_k)=lev*exp(r)*exp(rps);
exp(n)=gamma*exp(v)+(1-gamma)*dbar;
exp(v)=exp(r_k)*exp(q_k(-1))*exp(k(-1))-kappa*(exp(n(-1))/(exp(k(-1))*exp(q_k(-1))))^(-chi)*exp(rps)*exp(r)*(exp(q_k(-1))*exp(k(-1))-exp(n(-1)));
%exp(v)=exp(r_k)*exp(q_k(-1))*exp(k(-1))-lev*exp(rps)*exp(r)*(exp(q_k(-1))*exp(k(-1))-exp(n(-1)));
exp(c_e)=(1-gamma)*exp(n)/gamma;

//%  Inflation
pie-iota*pie(-1)=beta_p*(pie(+1)-iota*pie)-epsilon*(x-log(x_ss))+e_p;

//% interest rate
exp(r(+1))=exp(r_n)/exp(pie(+1));

//% Aggregation

exp(c)=gamma_p*exp(c_p)+(1-gamma_p)*exp(c_i);
%exp(y)=exp(c)+exp(I)+exp(G);
1=gamma_p*exp(h_p)+(1-gamma_p)*exp(h_i);
gamma_p*exp(d_p)=(1-gamma_p)*exp(b_i)+exp(b);

//%  Monetary Policy

%exp(r_n)=(1/beta_p)^(1-phi_r)*((exp(r_n(-1)))^phi_r)*exp(pie)^(phi_pie*(1-phi_r))*(exp(y)/exp(y(-1)))^(phi_y*(1-phi_r))*e_rn;
r_n=phi_r*r_n(-1)+(1-phi_r)*phi_pie*pie+(1-phi_r)*phi_y*(y-y(-1))+(1-phi_r)*log(1/beta_p)- e_rn;


//% shocks

A_e=rho_a*A_e(-1)+ e_a_e;
exp(m_i)=rho_mi*exp(m_i(-1))+(1-rho_mi)*m_ibar + e_mi;
%exp(G)=rho_g*exp(G(-1))+ (1-rho_g)*Gbar+ e_g;
exp(ee_h)=rho_eeh*(exp(ee_h(-1)))+(1-rho_eeh)*eeh_bar + e_h;
exp(rps)=rho_rps*(exp(rps(-1)))+(1-rho_rps)*1+ e_rps;

//%description of variables taken to data

c_obs=c-c_ss;
I_obs=I-I_ss;
qh_obs=q_h-qh_ss;
r_obs=r_n-log(1/beta_p);


end;  %end of model definition

initval;

A_e=     0;
 b=   1.1541;
 b_i=   0.7752;
 c=  -0.4131;
 c_e=  -1.6735;
 c_i=  -1.2967;
 c_p=   -0.0829;
 d_p=   1.9070;
 ee_h=  -2.1203;
 h_i=  -1.2896;
 h_p=   0.3941;
 I=  -1.8789;
 k=   1.8099;
 lambda =  -2.4976;
 l_i =   0.0233;
 l_p =  -0.1233;
 m_i =  -0.2231;
 n =   1.0781;
 pie=      0;
 q_h=   2.2955;
 q_k   =     0;
 r=   0.0075;
 r_k=   0.0273;
 r_n=   0.0075;
 rps=        0;
 v=   1.1054;
 w_i=  -1.2617;
 w_p=  -0.2678;
 x=   0.1398;
 y=   0.3930;
 c_obs = 0;
 I_obs = 0;
 qh_obs = 0;
 r_obs = 0; 
 

end;

check;
resid(1);
%steady;

shocks;


var e_a_e; stderr 0.101;
var e_mi;stderr 0.0041;
%var e_g; stderr 0.06;
var e_h; stderr 0.041;
var e_rps; stderr 0.01;
var e_rn; stderr 0.025;
var e_p; stderr 0.0045;
end;

%stoch_simul(order=1,irf=24);

estimated_params ;

//% Start Values and Priors
mu, 0.3, 0, 1, beta_pdf, 0.35, 0.05;
thet, 0.75, 0, 0.99, beta_pdf, 0.667,0.05;
gamma, 0.97, 0, 0.99, beta_pdf, 0.93, 0.05;
iota, 0.69,0,1, beta_pdf, 0.5, 0.2;
phi_r, 0.65, 0 , inf, beta_pdf, 0.75, 0.1;
phi_pie, 1.35, 0, inf, normal_pdf, 1.5, 0.1;
phi_y, 0.35, 0, inf, normal_pdf, 0, 0.1;
nk_ss, 0.484, 0, inf, beta_pdf, 0.5, 0.15;
chi, 0.03, 0, inf, inv_gamma_pdf, 0.05, 4;
rho_a, 0.95,  , , beta_pdf, 0.8, 0.1;
rho_mi, 0.95, , , beta_pdf, 0.8, 0.1;
rho_eeh, 0.98, , , beta_pdf, 0.8, 0.1;
rho_rps, 0.9, , , beta_pdf, 0.8, 0.1;
rho_g, 0.95, , , beta_pdf, 0.8, 0.1;
stderr e_a_e, 0.01, 0, inf, inv_gamma_pdf, 0.001, 0.01;
stderr e_mi, 0.04, 0, inf, inv_gamma_pdf, 0.001, 0.01;
stderr e_g, 0.011, 0, inf, inv_gamma_pdf, 0.001, 0.01;
stderr e_h, 0.039, 0, inf, inv_gamma_pdf, 0.001, 0.01;
stderr e_rps, 0.007, 0, inf, inv_gamma_pdf, 0.001, 0.01;
stderr e_rn, 0.008, 0, inf, inv_gamma_pdf, 0.001, 0.01;
stderr e_p, 0.039, 0, inf, inv_gamma_pdf, 0.001, 0.01;

end; 


varobs c_obs I_obs pie qh_obs r_obs;

estimation(datafile=datanew,
bayesian_irf,irf=20,
conf_sig=0.95,
smoother,
mh_jscale=0.2,
mode_compute=4,
prior_trunc=1e-100,
mh_replic=5000,
mh_nblocks=1,
lik_init=1)
c_obs I_obs pie qh_obs r_obs;