close all;

var c c_h im s y y_star ii i_star inflation inflation_star n w b a v nx 
inflation_h q ex u z;

var c_obs ex_obs y_obs inflation_obs ii_obs i_star_obs;

varexo Taylor_shock output_shock philip_shock world_output_shock 
inflation_star_shock is_shock;

parameters theta phiy phip phiy_star phip_star kappa_a kappa_y lambda 
disc beta epsilon rhoi rhoy alpha phi psi M mu theta_gamma gamma
etha delta rhoy_star rho_inflation lambda_ex etha_star rho_inflation_1
rho_is phii_star phii;


lambda_ex = 0.1;      
etha_star = 2.14;     
etha = 0.9;           
%delta = 0.25;        
delta = 0.22;         
disc = 0.01;          
alpha = 0.58;         
                      
gamma = 0.2572;       
phi = 1.6;            
psi = 0.01;           
rho_inflation_1 = 0.8; 
rhoi = 0.66;           
rhoy = 0.86;           
rhoy_star = 0.45;      
rho_inflation = 0.635; 
rho_is = 0.8;         
beta = exp(-disc);    
                      
phip = 1.50;          
phiy = 0.125;          
phip_star = 1.5;      
phiy_star = 0.5;      
theta = 2;            
epsilon = 6;          
M = epsilon/(epsilon-1);
mu = log(M);
omega = 0.75;        
                     
theta_gamma = gamma*(theta-1);
lambda = ((1-omega)*(1-beta*omega)/omega)*((1-alpha)/(1-alpha+alpha*epsilon));  
kappa_y = lambda*((phi+alpha)/(1-alpha));                                       
kappa_a = lambda*((phi+1)/(1-alpha));                                           
phii_star = 0.50;       
phii = 0.43;            
%dsge_prior_weight = 1.5;

model(linear);

-(theta_gamma + theta)*c + z = -theta_gamma*c(-1) - theta*c(+1) + ii - inflation(+1);  
w = phi*n - theta_gamma*c(-1) + theta*c;                                    
b = (1+ii(-1)-inflation)*b(-1)+y(-1)-c(-1)-delta*s(-1);                     
ex = y_star + lambda_ex*etha_star*s + (1-lambda_ex)*(ex(-1)-y_star(-1));    
ii = i_star - psi*b + (s(+1)-s) + inflation_h(+1) - inflation_star(+1);     
%ii = i_star - psi*b + (s(+1)-s) + inflation(+1) - inflation_star(+1);     
q = (1-delta)*s;                                                            
inflation = inflation_h + delta*(s-s(-1));                                  
c = (1-delta)*c_h + delta*im;                                               
c_h = etha*delta*s + c;                                                     
y = ((1-delta)*etha + lambda_ex*etha_star)*delta*s + (1-delta)*c + delta*y_star + delta*(1-lambda_ex)*(ex(-1)-y_star(-1));               
ii = phii*ii(-1)+ phip*inflation + phiy*y + v;                              
i_star = phii_star*i_star(-1) + phip_star*inflation_star + phiy_star*y_star;  
inflation_h = beta*inflation_h(+1) + kappa_y*y - kappa_a*a - lambda*theta_gamma*c(-1) + lambda*theta*c + lambda*delta*s + u; 
y = a + (1-alpha)*n;                                                         
nx = delta*(ex-im);                                                          
inflation_star = rho_inflation*inflation_star(-1) + inflation_star_shock;    
v = rhoi*v(-1) + Taylor_shock;                                               
a = rhoy*a(-1) + output_shock;                                               
z = rho_is*z(-1) + is_shock;                                                 
y_star = rhoy_star*y_star(-1) + world_output_shock;                          
u = rho_inflation_1*u(-1) + philip_shock;                                    
c_obs = c-c(-1);
ex_obs = ex-ex(-1);
y_obs = y-y(-1);
inflation_obs=inflation;
ii_obs=ii;
i_star_obs=i_star;
end;

write_latex_dynamic_model;

initval;
c=0;
c_h=0;
im=0;
s=0;
y=0;
y_star=0;
ii=0;
i_star=0;
inflation=0;
inflation_star=0;
n=0;
w=0;
b=0;
a=0;
v=0;
nx=0;
inflation_h=0;
q=0;
ex=0;
u=0;
z=0;
c_obs = 0;
ex_obs = 0;
y_obs = 0;
inflation_obs=0;
ii_obs=0;
i_star_obs=0;
end;

shocks; %Volt. Historicas
 var is_shock = (0.10)^2;
 var Taylor_shock = (0.10)^2;
 var output_shock = (0.10)^2;
 var philip_shock = (0.10)^2;
 var world_output_shock = (0.10)^2;
 var inflation_star_shock = (0.10)^2;
end;

%check;
%steady(solve_algo=2, maxit=5000);
%stoch_simul(order=2,periods=500,irf=20) w n c ex y ii inflation i_star;

estimated_params;
   stderr is_shock,0.10,,,inv_gamma_pdf,0.01,inf;
   stderr Taylor_shock,0.10,,,inv_gamma_pdf,0.01,inf;
   stderr output_shock,0.10,,,inv_gamma_pdf,0.01,inf;
   stderr philip_shock,0.10,,,inv_gamma_pdf,0.01,inf;
   stderr world_output_shock,0.10,,,inv_gamma_pdf,0.01,inf;
   stderr inflation_star_shock,0.10,,,inv_gamma_pdf,0.01,inf;
   alpha,0.10,,,beta_pdf,0.5,0.2;
   phi,0.10,,,gamma_pdf,2.5,1.0;
   gamma,0.10,,,beta_pdf,0.5,0.2;
   rhoi,0.10,,,beta_pdf,0.50,0.2;  
   rhoy,0.10,,,beta_pdf,0.50,0.2;  
   rhoy_star,0.10,,,beta_pdf,0.5,0.2; 
   rho_inflation,0.10,,,beta_pdf,0.5,0.2; 
   rho_inflation_1,0.10,,,beta_pdf,0.5,0.2;
   rho_is,0.10,,,beta_pdf,0.5,0.2;
   phii_star,0.10,,,beta_pdf,0.50,0.20;  
   phii,0.10,,,beta_pdf,0.50,0.20;     
   dsge_prior_weight,0.2728,,,gamma_pdf,0.5,0.4;    
end; 

varobs c_obs ex_obs y_obs ii_obs inflation_obs i_star_obs;
estimated_params_init(use_calibration);
end;

rand('seed',122); randn('seed',123);
estimation(datafile=DATABmr,dsge_var,dsge_varlag=1,xls_range=B1:G54,first_obs=2,prefilter=1,mh_nblocks=1,tex,mode_compute=4,mh_replic=10000,forecast=6,mode_check,moments_varendo,bayesian_irf) c_obs ex_obs y_obs ii_obs inflation_obs i_star_obs;