// Output Gap and Unemployment Model

var Y_US LGDP_US LGDP_US_BAR G_US E1_Y_US GROWTH_US GROWTH4_US GROWTH4_US_BAR;
varexo RES_Y_US RES_LGDP_US_BAR RES_G_US;

parameters beta_us1 tau_us growth_us_ss;

alpha_us1    = 0.5;
alpha_us2    = 0.3;
beta_us1     = 0.75;
tau_us       = 0.1;
growth_us_ss = 2.5;



model;

Y_US = beta_us1*Y_US(-1) + RES_Y_US;
Y_US = LGDP_US - LGDP_US_BAR;
LGDP_US_BAR = LGDP_US_BAR(-1) + G_US/4 + RES_LGDP_US_BAR;
G_US = tau_us*growth_us_ss + (1-tau_us)*G_US(-1) + RES_G_US;

GROWTH_US = 4*(LGDP_US - LGDP_US(-1));
GROWTH4_US = LGDP_US - LGDP_US(-4);
GROWTH4_US_BAR = LGDP_US_BAR - LGDP_US_BAR(-4);

E1_Y_US = Y_US(+1);


end;


shocks;

var RES_Y_US;       stderr 0.30;	  
var RES_LGDP_US_BAR; stderr 0.20;  
var RES_G_US;  stderr 0.05;

end;

unit_root_vars LGDP_US LGDP_US_BAR; 

steady;

check;


estimated_params;

beta_us1,   gamma_pdf,   			0.75, 0.10;
growth_us_ss, normal_pdf,     		2.5, 0.10; 

stderr RES_Y_US, inv_gamma_pdf,       0.25 , inf; 
stderr RES_LGDP_US_BAR, inv_gamma_pdf, 0.10 , inf; 
stderr RES_G_US, inv_gamma_pdf,  0.10  ,inf;


end;

varobs  LGDP_US ;

observation_trends;
//LGDP (growth_ss/4);
//LCPI (pietar_ss/4);
LGDP_US (growth_us_ss/4);
//LCPI_US (pietar_us_ss/4);
end;

//% Here is what you need to produce the TeX output:

//% Declaration of two global variables: 
//global M_.exo_names_tex M_.endo_names_tex;

//% Declaration of the tex names of the exogenous variables 
//% (note that the variables must be declared in the same 
//% order than in the global M_.exo_names_):
M_.exo_names_tex = 'RES\_Y\_US';
M_.exo_names_tex = strvcat(M_.exo_names_tex, 'RES\_LGDP\_US\_BAR');
M_.exo_names_tex = strvcat(M_.exo_names_tex, 'RES\_G\_US');

//% Declaration of the tex names of the endogenous variables 
//% (note that the variables must be declared in the same 
//% order than in the global M_.endo_names_):
M_.endo_names_tex = 'Y\_US';
M_.endo_names_tex = strvcat(M_.endo_names_tex,'LGDP\_US');
M_.endo_names_tex = strvcat(M_.endo_names_tex,'LGDP\_US\_BAR');
M_.endo_names_tex = strvcat(M_.endo_names_tex,'G\_US');
M_.endo_names_tex = strvcat(M_.endo_names_tex,'E1\_Y\_US');
M_.endo_names_tex = strvcat(M_.endo_names_tex,'GROWTH\_US');
M_.endo_names_tex = strvcat(M_.endo_names_tex,'GROWTH4\_US');
M_.endo_names_tex = strvcat(M_.endo_names_tex,'GROWTH4\_US\_BAR');



//
M_.param_names_tex = '\beta_{us1}';
M_.param_names_tex = strvcat(M_.param_names_tex,'\tau_{us}');
M_.param_names_tex = strvcat(M_.param_names_tex,'growth\_us\_ss');


//% Declaration of the tex names of the deep parameters to  
//% be estimated (note that the variables must be declared 
//% in the same order than in the global estim_params_.names):
estim_params_.tex = '\beta_{us1}';
estim_params_.tex = strvcat(estim_params_.tex,'growth\_us_ss');

estim_params_.tex = strvcat(estim_params_.tex,'RES\_Y\_US');
estim_params_.tex = strvcat(estim_params_.tex,'RES\_LGDP\_US\_BAR');
estim_params_.tex = strvcat(estim_params_.tex,'RES\_G\_US');



options_.kalman_algo = 5;
estimation(datafile=data,nobs=56,mode_check,mh_replic=0,mh_nblocks=1,filtered_vars,filter_step_ahead=[1:12]) Y_US;
//estimation(datafile=data,nobs=56,mode_compute=0,mode_file=soe13_mode,mh_replic=0,diffuse_d=11);

//optim_weights;
//YF 1;
//PIE4 1;
//DRS 0.5;
//end;

stoch_simul(order=1) Y_US GROWTH_US GROWTH4_US GROWTH4_US_BAR;

//to make more dynamic for report
data;
nobs = 56;

gl_report
//canada_report_extra__US

%% set firstdate
firstdate = '1994Q4'; 

%% set in-sample ahead forecast steps
fcast_steps = [1 4 8 12];

%% calculate the forecasts
max_step = max(fcast_steps);
[fy,fx,Fy,Fx,Udecomp,ndiffuse] = calc_fcast_all_4a;
//[fy,fx,Fy,Fx,Udecomp,ndiffuse] = calc_fcast_all(max_step);
%% save the forecasts
save([M_.dname '_forecasts.mat'], 'fy', 'fx', 'Fy', 'Fx', 'Udecomp', 'ndiffuse');

disp(sprintf('\n'));disp('Reporting forecast errors');
//report_fcast_errors(fy, Fy, ndiffuse, fcast_steps, firstdate);
report_smoothed_errors_4(fx, Fx, ndiffuse, strvcat('Y_US', 'GROWTH_US', 'GROWTH4_US', 'GROWTH4_US_BAR'), fcast_steps, firstdate);

%% order shocks as they should appear in the decompositions (don't put
%% highly contributing shocks close to each other, since they would get
%% similar colors)
exo_ord = [1 2 3];

disp('Reporting dynamic forecasts 1 period ahead');
report_fcast_dynamic_4a(fx, Fx, Udecomp, ndiffuse, 1, 1, exo_ord, firstdate,var_list_);
disp('Reporting dynamic forecasts 4 periods ahead');
report_fcast_dynamic_4a(fx, Fx, Udecomp, ndiffuse, 1, 4, exo_ord, firstdate,var_list_);
disp('Reporting dynamic forecasts 8 periods ahead');
report_fcast_dynamic_4a(fx, Fx, Udecomp, ndiffuse, 1, 8, exo_ord, firstdate,var_list_);
disp('Reporting dynamic forecasts 12 periods ahead');
report_fcast_dynamic_4a(fx, Fx, Udecomp, ndiffuse, 1, 12, exo_ord, firstdate,var_list_);