Merge pull request #1155 from JohannesPfeifer/prefilter_tests

Enable erroneously disabled unit tests

doc/dynare.texi

@@ -2458,7 +2458,7 @@ in @ref{stoch_simul}, if the @code{periods} option is specified. Note that this
 in @ref{forecast} as the initial point at which the forecasts are computed
 
 @item
-in @ref{conditional_forecast} as the initial point at which the conditional forecasts are computed
+in @ref{conditional_forecast} for a calibrated model as the initial point at which the conditional forecasts are computed
 
 @item 
 in @ref{Ramsey} policy, where it also specifies the values of the endogenous states at 
@@ -4480,7 +4480,9 @@ varobs
 @descriptionhead
 
 This block specifies @emph{linear} trends for observed variables as
-functions of model parameters.
+functions of model parameters. In case the @code{loglinear}-option is used,
+this corresponds to a linear trend in the logged observables, i.e. an exponential
+trend in the level of the observables.
 
 Each line inside of the block should be of the form:
 
@@ -4489,7 +4491,7 @@ Each line inside of the block should be of the form:
 @end example
 
 In most cases, variables shouldn't be centered when
-@code{observation_trends} is used.
+@code{observation_trends} is used. 
 
 @examplehead
 
@@ -6661,7 +6663,7 @@ Use @code{plot_conditional_forecast} to graph the results.
 
 @item parameter_set = @code{calibration} | @code{prior_mode} | @code{prior_mean} | @code{posterior_mode} | @code{posterior_mean} | @code{posterior_median}
 Specify the parameter set to use for the forecasting. No default
-value, mandatory option.
+value, mandatory option. Note that in case of estimated models, @code{conditional_forecast} does not support the @code{prefilter}-option.
 
 @item controlled_varexo = (@var{VARIABLE_NAME}@dots{})
 Specify the exogenous variables to use as control variables. No
@@ -6770,6 +6772,8 @@ variables between the first and last specified period. If an intermediate period
 is not specified, a value of 0 is assumed. That is, if you specify only 
 values for periods 1 and 3, the values for period 2 will be 0. Currently, it is not 
 possible to have uncontrolled intermediate periods.
+In case of the presence of @code{observation_trends}, the specified controlled path for 
+these variables needs to include the trend component.
 
 
 @end deffn

license.txt

@@ -98,6 +98,13 @@ Copyright: 1995-2007 Kurt Hornik
            2008-2009 Dynare Team
 License: GPL-3+
 
+Files: matlab/missing/stats/corr.m 
+Copyright: 1993-1996 Kurt Hornik
+           1996-2015 John W. Eaton
+           2013-2015 Julien Bect
+           2016 Dynare Team
+License: GPL-3+
+
 Files: matlab/missing/corrcoef/corrcoef.m matlab/missing/corrcoef/sumskipnan.m
  matlab/missing/corrcoef/flag_implicit_skip_nan.m matlab/missing/corrcoef/tcdf.m
 Copyright: 2000-2005,2008,2009,2011 by Alois Schloegl <alois.schloegl@gmail.com>

matlab/evaluate_smoother.m

@@ -1,4 +1,4 @@
-function oo_=evaluate_smoother(parameters,var_list)
+function [oo_, Smoothed_variables_declaration_order_deviation_form]=evaluate_smoother(parameters,var_list)
 % Evaluate the smoother at parameters.
 %
 % INPUTS
@@ -11,12 +11,17 @@ function oo_=evaluate_smoother(parameters,var_list)
 %    o oo       [structure]  results:
 %                              - SmoothedVariables
 %                              - SmoothedShocks
-%                              - SmoothedVariables
+%                              - FilteredVariablesShockDecomposition
-%                              - SmoothedVariables
+%                              - UpdatedVariables
-%                              - SmoothedVariables
+%                              - FilteredVariables
-%                              - SmoothedVariables
+%                              - SmoothedMeasurementErrors
-%                              - SmoothedVariables
+%                              - FilteredVariablesKStepAhead
-%                              - SmoothedVariables
+%                              - FilteredVariablesKStepAheadVariances
+%    o Smoothed_variables_declaration_order_deviation_form
+%                           Smoothed variables from the Kalman smoother in
+%                           order of declaration of variables (M_.endo_names)
+%                           in deviations from their respective mean, i.e.
+%                           without trend and constant part (used for shock_decomposition)
 %
 % SPECIAL REQUIREMENTS
 %    None
@@ -89,4 +94,8 @@ end
 
 [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp,Trend] = ...
     DsgeSmoother(parameters,dataset_.nobs,transpose(dataset_.data),dataset_info.missing.aindex,dataset_info.missing.state);
-[oo_]=write_smoother_results(M_,oo_,options_,bayestopt_,dataset_,dataset_info,atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,P,PK,decomp,Trend);
+[oo_]=write_smoother_results(M_,oo_,options_,bayestopt_,dataset_,dataset_info,atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,P,PK,decomp,Trend);
+
+if nargout==2
+   Smoothed_variables_declaration_order_deviation_form=atT(oo_.dr.inv_order_var(bayestopt_.smoother_var_list),:);
+end

matlab/imcforecast.m

@@ -68,6 +68,9 @@ else
 end
 
 if estimated_model
+    if options_.prefilter
+        error('imcforecast:: Conditional forecasting does not support the prefiltering option')
+    end
     if ischar(options_cond_fcst.parameter_set)
         switch options_cond_fcst.parameter_set
           case 'posterior_mode'
@@ -107,18 +110,37 @@ if estimated_model
     data_index = dataset_info.missing.aindex;
     gend = dataset_.nobs;
     missing_value = dataset_info.missing.state;
-    [atT,innov,measurement_error,filtered_state_vector,ys,trend_coeff] = DsgeSmoother(xparam,gend,data,data_index,missing_value);
+    [atT,innov,measurement_error,filtered_state_vector,ys,trend_coeff,aK,T,R,P,PK,decomp,trend_addition] = DsgeSmoother(xparam,gend,data,data_index,missing_value);
-    trend = repmat(ys,1,options_cond_fcst.periods+1);
+    %get constant part
-    for i=1:M_.endo_nbr
+    if options_.noconstant
-        j = strmatch(deblank(M_.endo_names(i,:)),options_.varobs,'exact');
+        constant = zeros(size(ys,1),options_cond_fcst.periods+1);
-        if ~isempty(j)
+    else
-            trend(i,:) = trend(i,:)+trend_coeff(j)*(gend+(0:options_cond_fcst.periods));
+        if options_.loglinear
+            constant = repmat(log(ys),1,options_cond_fcst.periods+1);
+        else
+            constant = repmat(ys,1,options_cond_fcst.periods+1);
         end
     end
-    trend = trend(oo_.dr.order_var,:);
+    %get trend part (which also takes care of prefiltering); needs to
+    %include the last period
+    if bayestopt_.with_trend == 1
+        [trend_addition] =compute_trend_coefficients(M_,options_,size(bayestopt_.smoother_mf,1),gend+options_cond_fcst.periods);
+        trend_addition = trend_addition(:,gend:end);
+    else
+        trend_addition=zeros(size(bayestopt_.smoother_mf,1),1+options_cond_fcst.periods);
+    end
+    % add trend to constant
+    for obs_iter=1:length(options_.varobs)
+        j = strmatch(options_.varobs{obs_iter},M_.endo_names,'exact');
+        constant(j,:) = constant(j,:)+trend_addition(obs_iter,:);        
+    end
+    trend = constant(oo_.dr.order_var,:);
     InitState(:,1) = atT(:,end);
 else
     graph_title='Calibration';
+    if ~isfield(oo_.dr,'kstate')
+        error('You need to call stoch_simul before conditional_forecast')
+    end
 end
 
 if isempty(options_.qz_criterium)
@@ -178,7 +200,10 @@ end
 % number of periods with constrained values
 cL = size(constrained_paths,2);
 
-constrained_paths = bsxfun(@minus,constrained_paths,trend(idx,1:cL));
+%transform constrained periods into deviations from steady state; note that
+%trend includes last actual data point and therefore we need to start in
+%period 2
+constrained_paths = bsxfun(@minus,constrained_paths,trend(idx,2:1+cL));
 
 FORCS1_shocks = zeros(n1,cL,options_cond_fcst.replic);
 
@@ -202,17 +227,15 @@ forecasts.controlled_variables = constrained_vars;
 forecasts.instruments = options_cond_fcst.controlled_varexo;
 
 for i = 1:EndoSize
-    eval(['forecasts.cond.Mean.' deblank(M_.endo_names(oo_.dr.order_var(i),:)) ' = mFORCS1(i,:)'';']);
+    forecasts.cond.Mean.(deblank(M_.endo_names(oo_.dr.order_var(i),:)))= mFORCS1(i,:)';
     tmp = sort(squeeze(FORCS1(i,:,:))');
-    eval(['forecasts.cond.ci.' deblank(M_.endo_names(oo_.dr.order_var(i),:)) ...
+    forecasts.cond.ci.(deblank(M_.endo_names(oo_.dr.order_var(i),:))) = [tmp(t1,:)' ,tmp(t2,:)' ]';
-          ' = [tmp(t1,:)'' ,tmp(t2,:)'' ]'';']);
 end
 
 for i = 1:n1
-    eval(['forecasts.controlled_exo_variables.Mean.' deblank(options_cond_fcst.controlled_varexo(i,:)) ' = mFORCS1_shocks(i,:)'';']);
+    forecasts.controlled_exo_variables.Mean.(deblank(options_cond_fcst.controlled_varexo(i,:))) = mFORCS1_shocks(i,:)';
     tmp = sort(squeeze(FORCS1_shocks(i,:,:))');
-    eval(['forecasts.controlled_exo_variables.ci.' deblank(options_cond_fcst.controlled_varexo(i,:))  ...
+    forecasts.controlled_exo_variables.ci.(deblank(options_cond_fcst.controlled_varexo(i,:))) = [tmp(t1,:)' ,tmp(t2,:)' ]';
-          ' = [tmp(t1,:)'' ,tmp(t2,:)'' ]'';']);
 end
 
 clear FORCS1 mFORCS1_shocks;
@@ -231,10 +254,9 @@ end
 mFORCS2 = mean(FORCS2,3);
 
 for i = 1:EndoSize
-    eval(['forecasts.uncond.Mean.' deblank(M_.endo_names(oo_.dr.order_var(i),:)) ' = mFORCS2(i,:)'';']);
+    forecasts.uncond.Mean.(deblank(M_.endo_names(oo_.dr.order_var(i),:)))= mFORCS2(i,:)';
     tmp = sort(squeeze(FORCS2(i,:,:))');
-    eval(['forecasts.uncond.ci.' deblank(M_.endo_names(oo_.dr.order_var(i),:)) ...
+    forecasts.uncond.ci.(deblank(M_.endo_names(oo_.dr.order_var(i),:))) = [tmp(t1,:)' ,tmp(t2,:)' ]';
-          ' = [tmp(t1,:)'' ,tmp(t2,:)'' ]'';']);
 end
 forecasts.graph.title=graph_title;
 forecasts.graph.fname=M_.fname;

matlab/metropolis_hastings_initialization.m

@@ -124,7 +124,7 @@ if ~options_.load_mh_file && ~options_.mh_recover
             trial = 1;
             while validate == 0 && trial <= 10
                 candidate = rand_multivariate_normal( transpose(xparam1), d * options_.mh_init_scale, npar);
-                if all(candidate(:) > mh_bounds.lb) && all(candidate(:) < mh_bounds.ub) 
+                if all(candidate(:) >= mh_bounds.lb) && all(candidate(:) <= mh_bounds.ub) 
                     ix2(j,:) = candidate;
                     ilogpo2(j) = - feval(TargetFun,ix2(j,:)',dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,mh_bounds,oo_);
                     if ~isfinite(ilogpo2(j)) % if returned log-density is
@@ -157,6 +157,7 @@ if ~options_.load_mh_file && ~options_.mh_recover
             end
             if trial > 10 && ~validate
                 disp(['Estimation::mcmc: I''m unable to find a starting value for block ' int2str(j)])
+                fclose(fidlog);
                 return
             end
         end
@@ -166,7 +167,7 @@ if ~options_.load_mh_file && ~options_.mh_recover
     else% Case 2: one chain (we start from the posterior mode)
         fprintf(fidlog,['  Initial values of the parameters:\n']);
         candidate = transpose(xparam1(:));%
-        if all(candidate(:) > mh_bounds.lb) && all(candidate(:) < mh_bounds.ub) 
+        if all(candidate(:) >= mh_bounds.lb) && all(candidate(:) <= mh_bounds.ub) 
             ix2 = candidate;
             ilogpo2 = - feval(TargetFun,ix2',dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,mh_bounds,oo_);
             disp('Estimation::mcmc: Initialization at the posterior mode.')
@@ -179,6 +180,7 @@ if ~options_.load_mh_file && ~options_.mh_recover
         else
             disp('Estimation::mcmc: Initialization failed...')
             disp('Estimation::mcmc: The posterior mode lies outside the prior bounds.')
+            fclose(fidlog);
             return
         end
         fprintf(fidlog,' \n');

matlab/missing/stats/corr.m

@@ -0,0 +1,189 @@
+function retval = corr(x, y) % --*-- Unitary tests --*--
+%@info:
+%! @deftypefn  {Function File} {} corr (@var{x})
+%! @deftypefnx {Function File} {} corr (@var{x}, @var{y})
+%! Compute matrix of correlation coefficients.
+%! @anchor{corr}
+%!
+%! If each row of @var{x} and @var{y} is an observation and each column is
+%! a variable, then the @w{(@var{i}, @var{j})-th} entry of
+%! @code{corr (@var{x}, @var{y})} is the correlation between the
+%! @var{i}-th variable in @var{x} and the @var{j}-th variable in @var{y}.
+%! @tex
+%! $$
+%! {\rm corr}(x,y) = {{\rm cov}(x,y) \over {\rm std}(x) {\rm std}(y)}
+%! $$
+%! @end tex
+%! @ifnottex
+%!
+%! @example
+%! corr (x,y) = cov (x,y) / (std (x) * std (y))
+%! @end example
+%!
+%! @end ifnottex
+%! If called with one argument, compute @code{corr (@var{x}, @var{x})},
+%! the correlation between the columns of @var{x}.
+%! @end deftypefn
+%@eod:
+% 
+% Notes:    - the original Octave code has been rewritten to avoid calling cov, since
+%               there is a long-standing incompatiblity between Matlab's cov and Octave's cov
+%               (see https://savannah.gnu.org/bugs/?40751)
+%           - For compatibility with Matlab, the correlation of a constant
+%               is defined as NaN, not 1
+%
+% Adapted for Matlab (R) from GNU Octave 4.0.1
+% Original files: statistics\base\corr.m, statistics\base\cov.m, and packages\stk-2.3.4\misc\mole\corr\corr.m
+% Original authors: Kurt Hornik <hornik@wu-wien.ac.at> and Julien Bect  <julien.bect@supelec.fr>
+
+% Copyright (C) 1993-1996 Kurt Hornik
+% Copyright (C) 1996-2015 John W. Eaton
+% Copyright (C) 2013-2015 Julien Bect
+% Copyright (C) 2016 Dynare Team
+%
+% This file is part of Dynare.
+%
+% Dynare is free software: you can redistribute it and/or modify
+% it under the terms of the GNU General Public License as published by
+% the Free Software Foundation, either version 3 of the License, or
+% (at your option) any later version.
+%
+% Dynare is distributed in the hope that it will be useful,
+% but WITHOUT ANY WARRANTY; without even the implied warranty of
+% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+% GNU General Public License for more details.
+%
+% You should have received a copy of the GNU General Public License
+% along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
+
+
+if (nargin < 1 || nargin > 2)
+    error('corr needs to be called with 1 or 2 input arguments')
+end
+
+% Input validation
+if (nargin==1 && ~(isnumeric (x) || islogical (x))) || ...
+    (nargin==2 && ~(isnumeric (x) || islogical (x) || isnumeric (y) || islogical (y)))
+    error ('corr: X and Y must be numeric matrices or vectors');
+end
+
+if (nargin==1 && ~ismatrix(x)) || (nargin==2 && (~ismatrix(y) || ~ismatrix(x)))
+error ('corr: X and Y must be 2-D matrices or vectors');
+end
+
+if (nargin == 2)
+    if size(y,1)~=size(x,1)
+        error('corr: X and Y must have the same length')
+    end
+end
+
+% Special case, correlation with scalar is NaN in Matlab
+if isscalar(x)
+    if nargin==1
+        retval = NaN;
+    else
+        retval = NaN(size(y));    
+    end
+    return;
+end
+
+if nargin==2 && isscalar(y)
+    retval = NaN(size(x'));    
+    return;
+end
+
+n = size (x, 1);
+x = x - repmat(mean(x),n,1); %demean
+sx = std(x,1); %standard deviation
+sx(sx==0)=NaN; %take care of constant vectors
+
+if (nargin == 2)
+    y = y - repmat (mean (y), n, 1);
+    sy = std (y, 1);
+    sy (sy == 0) = nan;
+else
+    y = x;
+    sy = sx;
+end
+
+c = x'*y;
+s = sx'*sy;
+retval = c./(n * s);
+
+
+end
+
+
+%@test:1
+%$test
+%$ x = rand (10);
+%$ cc1 = corr (x);
+%$ cc2 = corr (x, x);
+%$ t(1)=dassert(size(cc1),[10, 10]); 
+%$ t(2)=dassert(size (cc2),[10, 10]);
+%$ t(3)=dassert(cc1, cc2, sqrt (eps));
+%@eof:1
+
+%@test:2
+%$ x = [1:3]';
+%$ y = [3:-1:1]';
+%$ t(4)=dassert(corr (x, y), -1, 5*eps);
+%$ t(5)=dassert(corr (x, flipud (y)), 1, 5*eps);
+%$ t(6)=dassert(corr ([x, y]), [1 -1; -1 1], 5*eps);
+%@eof:2
+
+%@test:3
+%$t(7)=dassert (corr (5), NaN);
+%$t(8)=dassert (corr([1 2 3],5),NaN(3,1));
+%$t(9)=dassert (corr(5,[1 2 3]),NaN(1,3));
+%@eof:3
+
+%@test:4
+%$ Test input validation
+%$ try
+%$     corr ()
+%$     t(10) = false;
+%$ catch
+%$     t(10) = true;
+%$ end
+%@eof:4
+
+%@test:5
+%$ Test input validation
+%$ try
+%$     corr (1, 2, 3)
+%$     t(11) = false;
+%$ catch
+%$     t(11) = true;
+%$ end
+%@eof:5
+
+%@test:6
+%$ Test input validation
+%$ try
+%$     corr ([1; 2], ['A', 'B'])
+%$     t(12) = false;
+%$ catch
+%$     t(12) = true;
+%$ end
+%@eof:6
+
+%@test:7
+%$ Test input validation
+%$ try
+%$     corr (ones (2,2,2))
+%$     t(13) = false;
+%$ catch
+%$     t(13) = true;
+%$ end
+%@eof:7
+
+%@test:8
+%$ Test input validation
+%$ try
+%$     corr (ones (2,2), ones (2,2,2))
+%$     t(14) = false;
+%$ catch
+%$     t(14) = true;
+%$ end
+%@eof:8

matlab/plot_icforecast.m

@@ -59,13 +59,9 @@ end
 function build_figure(name,cci1,cci2,mm1,mm2,options_,graphoptions)
 hh = dyn_figure(options_,'Name',['Conditional forecast (' graphoptions.title ,'): ' name '.']);
 H = length(mm1);
-h1 = area(1:H,cci1(2,1:H));
+h1 = area(1:H,cci1(2,1:H),'BaseValue',min([min(cci1(1,:)),min(cci2(1,:))]),'FaceColor',[.9 .9 .9]);
-set(h1,'BaseValue',min([min(cci1(1,:)),min(cci2(1,:))]))
-set(h1,'FaceColor',[.9 .9 .9])
 hold on
-h2 = area(1:H,cci1(1,1:H));
+h2 = area(1:H,cci1(1,1:H),'BaseValue',min([min(cci1(1,:)),min(cci2(1,:))]),'FaceColor',[1 1 1]);
-set(h2,'BaseValue',min([min(cci1(1,:)),min(cci2(1,:))]))
-set(h2,'FaceColor',[1 1 1])
 plot(1:H,mm1,'-k','linewidth',3)
 plot(1:H,mm2,'--k','linewidth',3)
 plot(1:H,cci2(1,:),'--k','linewidth',1)

matlab/shock_decomposition.m

@@ -19,7 +19,7 @@ function oo_ = shock_decomposition(M_,oo_,options_,varlist)
 % SPECIAL REQUIREMENTS
 %    none
 
-% Copyright (C) 2009-2013 Dynare Team
+% Copyright (C) 2009-2016 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -62,7 +62,7 @@ if isempty(parameter_set)
     end
 end
 
-oo = evaluate_smoother(parameter_set,varlist);
+[oo,Smoothed_Variables_deviation_from_mean] = evaluate_smoother(parameter_set,varlist);
 
 % reduced form
 dr = oo.dr;
@@ -77,16 +77,15 @@ A = dr.ghx;
 B = dr.ghu;
 
 % initialization
-gend = eval(['size(oo.SmoothedShocks.' M_.exo_names(1,:),',1)']);
+gend = size(oo.SmoothedShocks.(M_.exo_names(1,:)),1);
 epsilon=NaN(nshocks,gend);
 for i=1:nshocks
-    epsilon(i,:) = eval(['oo.SmoothedShocks.' M_.exo_names(i,:)]);
+    epsilon(i,:) = oo.SmoothedShocks.(M_.exo_names(i,:));
 end
 
 z = zeros(endo_nbr,nshocks+2,gend);
-for i=1:endo_nbr
+
-    z(i,end,:) = eval(['oo.SmoothedVariables.' M_.endo_names(i,:)]);
+z(:,end,:) = Smoothed_Variables_deviation_from_mean;
-end
 
 maximum_lag = M_.maximum_lag;
 lead_lag_incidence = M_.lead_lag_incidence;

tests/Makefile.am

@@ -275,7 +275,7 @@ MODFILES = \
 	optimizers/fs2000_w.mod \
 	differentiate_forward_vars/RBC_differentiate_forward.mod \
 	TeX/fs2000_corr_ME.mod \
-	prior_posterior_function/fs2000_prior_posterior_function.mod
+	prior_posterior_function/fs2000_prior_posterior_function.mod \
 	observation_trends_and_prefiltering/MCMC/Trend_loglin_no_prefilt_first_obs_MC.mod \
 	observation_trends_and_prefiltering/MCMC/Trend_loglin_prefilt_first_obs_MC.mod \
 	observation_trends_and_prefiltering/MCMC/Trend_loglinear_no_prefilter_MC.mod \
@@ -299,8 +299,7 @@ MODFILES = \
 	observation_trends_and_prefiltering/calib_smoother/Tr_no_prefilter_loglin_calib_smoother.mod \
 	observation_trends_and_prefiltering/calib_smoother/Tr_no_prefil_f_obs_loglin_cal_smoother.mod \
 	observation_trends_and_prefiltering/calib_smoother/Tr_prefilter_loglin_calib_smoother.mod \
-	observation_trends_and_prefiltering/calib_smoother/Tr_prefil_f_obs_loglin_cal_smoother.mod \
+	observation_trends_and_prefiltering/calib_smoother/Tr_prefil_f_obs_loglin_cal_smoother.mod
-	reporting/example1.mod
 
 XFAIL_MODFILES = ramst_xfail.mod \
 	estim_param_in_shock_value_xfail.mod \
@@ -524,7 +523,8 @@ EXTRA_DIST = \
 	observation_trends_and_prefiltering/Trend_exp_model_calib_no_prefilter_common.inc \
 	observation_trends_and_prefiltering/Trend_model_calib_prefilter_common.inc \
 	observation_trends_and_prefiltering/Trend_model_calib_no_prefilter_common.inc \
-	observation_trends_and_prefiltering/Trend_load_data_common.inc 
+	observation_trends_and_prefiltering/Trend_load_data_common.inc \
+	observation_trends_and_prefiltering/Trend_no_prefilter_conditional_forecast.inc \
 	optimal_policy/Ramsey/oo_ramsey_policy_initval.mat \
 	optimizers/optimizer_function_wrapper.m \
 	optimizers/fs2000.common.inc \

tests/conditional_forecasts/2/fs2000_est.mod

@@ -107,7 +107,7 @@ varobs gp_obs gy_obs;
 options_.solve_tolf = 1e-12;
 
 // Metropolis replications are too few, this is only for testing purpose
-estimation(order=1,datafile=fsdat_simul,nobs=192,loglinear,mh_replic=10000,mh_nblocks=1,mh_jscale=0.8);
+estimation(order=1,datafile=fsdat_simul,nobs=192,loglinear,mh_replic=2000,mh_nblocks=1,mh_jscale=0.8);
 
 conditional_forecast_paths;
 var gy_obs;

tests/observation_trends_and_prefiltering/MCMC/Trend_loglin_no_prefilt_first_obs_MC.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_exp_model_no_prefilter_common.inc" 
 
-estimation(order=1,datafile='../Exp_AR1_trend_data_with_constant',mh_replic=2000,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='Exp_AR1_trend_data_with_constant',mh_replic=2000,
     mode_compute=4,first_obs=1000,loglinear,smoother,forecast=100,prefilter=0,
     mcmc_jumping_covariance='MCMC_jump_covar',
     filtered_vars, filter_step_ahead = [1,2,4],
     mh_nblocks=1,mh_jscale=0.3) P_obs Y_obs junk2;
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end

tests/observation_trends_and_prefiltering/MCMC/Trend_loglin_prefilt_first_obs_MC.mod

@@ -1,20 +1,23 @@
 @#include "../Trend_exp_model_prefilter_common.inc" 
         
-estimation(order=1,datafile='../Exp_AR1_trend_data_with_constant',mh_replic=2000,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='Exp_AR1_trend_data_with_constant',mh_replic=2000,
     mode_compute=4,first_obs=1000,loglinear,smoother,forecast=100,prefilter=1,
     mcmc_jumping_covariance='MCMC_jump_covar_prefilter',
     filtered_vars, filter_step_ahead = [1,2,4],
     mh_nblocks=1,mh_jscale=1e-4) P_obs Y_obs junk2;
     
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params_prefilter');
+loaded_par=load('orig_params_prefilter');
 
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end
-loaded_par_full=load('../orig_params');
+loaded_par_full=load('orig_params');
 y_forecast_100_periods=loaded_par_full.orig_params(strmatch('const_y',loaded_par_full.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par_full.orig_params(strmatch('g_y',loaded_par_full.param_names,'exact'));
 p_forecast_100_periods=loaded_par_full.orig_params(strmatch('const_p',loaded_par_full.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par_full.orig_params(strmatch('g_p',loaded_par_full.param_names,'exact'));
 

tests/observation_trends_and_prefiltering/MCMC/Trend_loglinear_no_prefilter_MC.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_exp_model_no_prefilter_common.inc" 
 
-estimation(order=1,datafile='../Exp_AR1_trend_data_with_constant',mh_replic=2000,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='Exp_AR1_trend_data_with_constant',mh_replic=2000,
     mode_compute=4,first_obs=1,loglinear,diffuse_filter,smoother,forecast=100,prefilter=0,
     mcmc_jumping_covariance='MCMC_jump_covar',
     filtered_vars, filter_step_ahead = [1,2,4],
     mh_nblocks=1,mh_jscale=0.3) P_obs Y_obs junk2;
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end

tests/observation_trends_and_prefiltering/MCMC/Trend_loglinear_prefilter_MC.mod

@@ -1,20 +1,23 @@
 @#include "../Trend_exp_model_prefilter_common.inc" 
 
-estimation(order=1,datafile='../Exp_AR1_trend_data_with_constant',mh_replic=2000,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='Exp_AR1_trend_data_with_constant',mh_replic=2000,
     mode_compute=4,first_obs=1,loglinear,smoother,forecast=100,prefilter=1,
     mcmc_jumping_covariance='MCMC_jump_covar_prefilter',
     filtered_vars, filter_step_ahead = [1,2,4],
     mh_nblocks=1,mh_jscale=1e-4) P_obs Y_obs junk2;
     
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params_prefilter');
+loaded_par=load('orig_params_prefilter');
 
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end
-loaded_par_full=load('../orig_params');
+loaded_par_full=load('orig_params');
 y_forecast_100_periods=loaded_par_full.orig_params(strmatch('const_y',loaded_par_full.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par_full.orig_params(strmatch('g_y',loaded_par_full.param_names,'exact'));
 p_forecast_100_periods=loaded_par_full.orig_params(strmatch('const_p',loaded_par_full.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par_full.orig_params(strmatch('g_p',loaded_par_full.param_names,'exact'));
 

tests/observation_trends_and_prefiltering/MCMC/Trend_no_prefilter_MC.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_model_no_prefilter_common.inc"
 
-estimation(order=1,datafile='../AR1_trend_data_with_constant',mh_replic=2000,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='AR1_trend_data_with_constant',mh_replic=2000,
             mode_compute=4,first_obs=1,smoother,mh_nblocks=1,mh_jscale=0.3,
             filtered_vars, filter_step_ahead = [1,2,4],
             mcmc_jumping_covariance='MCMC_jump_covar',forecast=100,prefilter=0) P_obs Y_obs junk2;
             
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end

tests/observation_trends_and_prefiltering/MCMC/Trend_no_prefilter_first_obs_MC.mod

@@ -1,15 +1,18 @@
 @#include "../Trend_model_no_prefilter_common.inc"
 
-estimation(order=1,datafile='../AR1_trend_data_with_constant',
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='AR1_trend_data_with_constant',
         mh_replic=2000,mode_compute=4,first_obs=1000,smoother,forecast=100,prefilter=0,
         mcmc_jumping_covariance='MCMC_jump_covar',
         filtered_vars, filter_step_ahead = [1,2,4],
         mh_nblocks=1,mh_jscale=0.3) P_obs Y_obs junk2;
 
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end

tests/observation_trends_and_prefiltering/MCMC/Trend_prefilter_MC.mod

@@ -1,20 +1,23 @@
 @#include "../Trend_model_prefilter_common.inc"
         
-estimation(order=1,datafile='../AR1_trend_data_with_constant',mh_replic=2000,mode_compute=4,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='AR1_trend_data_with_constant',mh_replic=2000,mode_compute=4,
     first_obs=1,smoother,prefilter=1,
     mh_nblocks=1,mh_jscale=1e-4,
     filtered_vars, filter_step_ahead = [1,2,4],
     mcmc_jumping_covariance='MCMC_jump_covar_prefilter',forecast=100) P_obs Y_obs junk2;
 
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params_prefilter');
+loaded_par=load('orig_params_prefilter');
 
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end
-loaded_par_full=load('../orig_params');
+loaded_par_full=load('orig_params');
 y_forecast_100_periods=loaded_par_full.orig_params(strmatch('const_y',loaded_par_full.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par_full.orig_params(strmatch('g_y',loaded_par_full.param_names,'exact'));
 p_forecast_100_periods=loaded_par_full.orig_params(strmatch('const_p',loaded_par_full.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par_full.orig_params(strmatch('g_p',loaded_par_full.param_names,'exact'));
 

tests/observation_trends_and_prefiltering/MCMC/Trend_prefilter_first_obs_MC.mod

@@ -1,20 +1,23 @@
 @#include "../Trend_model_prefilter_common.inc"
 
-estimation(order=1,datafile='../AR1_trend_data_with_constant',mh_replic=2000,mode_compute=4,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='AR1_trend_data_with_constant',mh_replic=2000,mode_compute=4,
         first_obs=1000,smoother,prefilter=1,
         mh_nblocks=1,mh_jscale=1e-4,
         filtered_vars, filter_step_ahead = [1,2,4],
         mcmc_jumping_covariance='MCMC_jump_covar_prefilter',forecast=100) P_obs Y_obs junk2;
 
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params_prefilter');
+loaded_par=load('orig_params_prefilter');
 
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end
-loaded_par_full=load('../orig_params');
+loaded_par_full=load('orig_params');
 y_forecast_100_periods=loaded_par_full.orig_params(strmatch('const_y',loaded_par_full.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par_full.orig_params(strmatch('g_y',loaded_par_full.param_names,'exact'));
 p_forecast_100_periods=loaded_par_full.orig_params(strmatch('const_p',loaded_par_full.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par_full.orig_params(strmatch('g_p',loaded_par_full.param_names,'exact'));
 

tests/observation_trends_and_prefiltering/ML/Trend_loglinear_no_prefilter.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_exp_model_no_prefilter_common.inc" 
 
-estimation(order=1,datafile='../Exp_AR1_trend_data_with_constant',mh_replic=0,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='Exp_AR1_trend_data_with_constant',mh_replic=0,
     mode_compute=4,first_obs=1,
     filtered_vars, filter_step_ahead = [1,2,4],        
     loglinear,smoother,forecast=100,prefilter=0) P_obs Y_obs junk2;
 
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end
@@ -17,3 +20,6 @@ y_forecast_100_periods=loaded_par.orig_params(strmatch('const_y',M_.param_names,
 p_forecast_100_periods=loaded_par.orig_params(strmatch('const_p',M_.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par.orig_params(strmatch('g_p',M_.param_names,'exact'));
 
 @#include "../Trend_diagnostics_ML_common.inc" 
+
+@#include "../Trend_no_prefilter_conditional_forecast.inc"
+

tests/observation_trends_and_prefiltering/ML/Trend_loglinear_no_prefilter_first_obs.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_exp_model_no_prefilter_common.inc" 
 
-estimation(order=1,datafile='../Exp_AR1_trend_data_with_constant',mh_replic=0,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='Exp_AR1_trend_data_with_constant',mh_replic=0,
     mode_compute=4,first_obs=1000,
     filtered_vars, filter_step_ahead = [1,2,4],        
     loglinear,smoother,forecast=100,prefilter=0) P_obs Y_obs junk2;
 
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end
@@ -17,3 +20,5 @@ y_forecast_100_periods=loaded_par.orig_params(strmatch('const_y',M_.param_names,
 p_forecast_100_periods=loaded_par.orig_params(strmatch('const_p',M_.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par.orig_params(strmatch('g_p',M_.param_names,'exact'));
 
 @#include "../Trend_diagnostics_ML_common.inc" 
+
+@#include "../Trend_no_prefilter_conditional_forecast.inc"

tests/observation_trends_and_prefiltering/ML/Trend_loglinear_prefilter.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_exp_model_prefilter_common.inc" 
 
-estimation(order=1,datafile='../Exp_AR1_trend_data_with_constant',mh_replic=0,mode_compute=4,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='Exp_AR1_trend_data_with_constant',mh_replic=0,mode_compute=4,
         first_obs=1,smoother,loglinear,
         filtered_vars, filter_step_ahead = [1,2,4],        
         forecast=100,prefilter=1) P_obs Y_obs junk2;
 
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 
 if max(abs((M_.params-loaded_par.orig_params([1:4,7:8]))./loaded_par.orig_params([1:4,7:8])))>0.03
     error('Parameter estimates do not match')
@@ -18,3 +21,4 @@ y_forecast_100_periods=loaded_par.orig_params(strmatch('const_y',loaded_par.para
 p_forecast_100_periods=loaded_par.orig_params(strmatch('const_p',loaded_par.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par.orig_params(strmatch('g_p',loaded_par.param_names,'exact'));
 
 @#include "../Trend_diagnostics_ML_common.inc" 
+

tests/observation_trends_and_prefiltering/ML/Trend_loglinear_prefilter_first_obs.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_exp_model_prefilter_common.inc" 
 
-estimation(order=1,datafile='../Exp_AR1_trend_data_with_constant',mh_replic=0,mode_compute=4,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='Exp_AR1_trend_data_with_constant',mh_replic=0,mode_compute=4,
         first_obs=1000,smoother,loglinear,
         filtered_vars, filter_step_ahead = [1,2,4],        
         forecast=100,prefilter=1) P_obs  Y_obs junk2;
 
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 
 if max(abs((M_.params-loaded_par.orig_params([1:4,7:8]))./loaded_par.orig_params([1:4,7:8])))>0.03
     error('Parameter estimates do not match')

tests/observation_trends_and_prefiltering/ML/Trend_no_prefilter.mod

@@ -1,13 +1,16 @@
 @#include "../Trend_model_no_prefilter_common.inc"
 
-estimation(order=1,datafile='../AR1_trend_data_with_constant',mh_replic=0,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='AR1_trend_data_with_constant',mh_replic=0,
         mode_compute=4,first_obs=1,
         filtered_vars, filter_step_ahead = [1,2,4],        
         diffuse_filter,smoother,forecast=100,prefilter=0) P_obs Y_obs junk2;
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameter estimates do not match')
 end
@@ -16,3 +19,5 @@ y_forecast_100_periods=loaded_par.orig_params(strmatch('const_y',M_.param_names,
 p_forecast_100_periods=loaded_par.orig_params(strmatch('const_p',M_.param_names,'exact'))+(options_.first_obs+options_.nobs-1+options_.forecast)*loaded_par.orig_params(strmatch('g_p',M_.param_names,'exact'))
 
 @#include "../Trend_diagnostics_ML_common.inc" 
+
+@#include "../Trend_no_prefilter_conditional_forecast.inc"

tests/observation_trends_and_prefiltering/ML/Trend_no_prefilter_first_obs.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_model_no_prefilter_common.inc"
 
-estimation(order=1,datafile='../AR1_trend_data_with_constant',mh_replic=0,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='AR1_trend_data_with_constant',mh_replic=0,
         mode_compute=4,first_obs=1000,
         filtered_vars, filter_step_ahead = [1,2,4],        
         smoother,forecast=100,prefilter=0) P_obs Y_obs junk2;
 
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs(M_.params-loaded_par.orig_params)./loaded_par.orig_params)>0.03
     error('Parameter estimates do not match')
 end

tests/observation_trends_and_prefiltering/ML/Trend_prefilter.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_model_prefilter_common.inc"
 
-estimation(order=1,datafile='../AR1_trend_data_with_constant',mh_replic=0,mode_compute=4,
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='AR1_trend_data_with_constant',mh_replic=0,mode_compute=4,
         first_obs=1,
         filtered_vars, filter_step_ahead = [1,2,4],        
         smoother,forecast=100,prefilter=1) P_obs Y_obs junk2;
 
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
         
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 
 if max(abs((M_.params-loaded_par.orig_params([1:4,7:8]))./loaded_par.orig_params([1:4,7:8])))>0.03
     error('Parameter estimates do not match')

tests/observation_trends_and_prefiltering/ML/Trend_prefilter_first_obs.mod

@@ -1,15 +1,18 @@
 @#include "../Trend_model_prefilter_common.inc"
 
-estimation(order=1,datafile='../AR1_trend_data_with_constant',
+addpath('..');
+generate_trend_stationary_AR1;
+
+estimation(order=1,datafile='AR1_trend_data_with_constant',
     mh_replic=0,mode_compute=4,
     filtered_vars, filter_step_ahead = [1,2,4],        
     first_obs=1000,diffuse_filter,smoother,forecast=100,prefilter=1) P_obs Y_obs junk2;
 
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 
 if max(abs((M_.params-loaded_par.orig_params([1:4,7:8]))./loaded_par.orig_params([1:4,7:8])))>0.03
     error('Parameter estimates do not match')

tests/observation_trends_and_prefiltering/Trend_no_prefilter_conditional_forecast.inc

@@ -0,0 +1,39 @@
+conditional_forecast_paths;
+var P_obs;
+periods  1  2  3:5;
+values   1.01 0.98 1;
+var Y_obs;
+periods  1  2  3:5;
+values   3 3.04 2.98;
+end;
+
+conditional_forecast(periods=100,parameter_set=posterior_mode,replic=1000, controlled_varexo=(e_y,e_p));
+
+plot_conditional_forecast(periods=100) Y_obs P_obs;
+
+%compare unconditional forecasts
+cond_forecast=load('conditional_forecasts.mat');
+if max(abs(cond_forecast.forecasts.uncond.Mean.P_obs(2:end)-oo_.forecast.Mean.P_obs))>1e-8
+    error('Unconditional Forecasts do not match')
+end
+        
+%compare conditional forecasts; histval here sets initval condition for capital different from steady state
+initial_condition_states=NaN(size(oo_.dr.ys));
+initial_condition_states(strmatch('Y_obs',M_.endo_names,'exact')) = oo_.SmoothedVariables.Y_obs(end)-oo_.Smoother.Trend.Y_obs(end);
+initial_condition_states(strmatch('P_obs',M_.endo_names,'exact')) = oo_.SmoothedVariables.P_obs(end)-oo_.Smoother.Trend.P_obs(end);
+initial_condition_states(strmatch('junk1',M_.endo_names,'exact')) = oo_.SmoothedVariables.junk1(end);
+initial_condition_states(strmatch('junk2',M_.endo_names,'exact')) = oo_.SmoothedVariables.junk2(end)-oo_.Smoother.Trend.junk2(end);
+shock_matrix = zeros(options_cond_fcst_.periods ,M_.exo_nbr); %create shock matrix with found controlled shocks
+shock_matrix(1:5,strmatch('e_y',M_.exo_names,'exact')) = cond_forecast.forecasts.controlled_exo_variables.Mean.e_y; %set controlled shocks to their values
+shock_matrix(1:5,strmatch('e_p',M_.exo_names,'exact')) = cond_forecast.forecasts.controlled_exo_variables.Mean.e_p; %set controlled shocks to their values
+
+y_simult = simult_(initial_condition_states,oo_.dr,shock_matrix,1);
+
+if max(abs(y_simult(strmatch('Y_obs',M_.endo_names,'exact'),:)'+(options_.first_obs-1+options_.nobs:options_.first_obs-1+options_.nobs+options_.forecast)'*g_y-cond_forecast.forecasts.cond.Mean.Y_obs))>1e-8
+    error('Conditional Forecasts do not match')
+end
+
+if max(abs(y_simult(strmatch('P_obs',M_.endo_names,'exact'),:)'+(options_.first_obs-1+options_.nobs:options_.first_obs-1+options_.nobs+options_.forecast)'*g_p-cond_forecast.forecasts.cond.Mean.P_obs))>1e-8
+    error('Conditional Forecasts do not match')
+end
+

tests/observation_trends_and_prefiltering/calib_smoother/Tr_no_prefil_f_obs_loglin_cal_smoother.mod

@@ -1,13 +1,16 @@
 @#include "../Trend_exp_model_calib_no_prefilter_common.inc"
 options_.filter_decomposition=1;
 
-calib_smoother(datafile='../Exp_AR1_trend_data_with_constant',prefilter=0,loglinear,first_obs=1000,
+addpath('..');
+generate_trend_stationary_AR1;
+
+calib_smoother(datafile='Exp_AR1_trend_data_with_constant',prefilter=0,loglinear,first_obs=1000,
 //         filter_decomposition,
         filtered_vars, filter_step_ahead = [1,2,4]) P_obs Y_obs junk2;
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameters do not match')
 end

tests/observation_trends_and_prefiltering/calib_smoother/Tr_no_prefilt_first_obs_cal_smooth.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_model_calib_no_prefilter_common.inc"
 options_.filter_decomposition=1;
 
-calib_smoother(datafile='../AR1_trend_data_with_constant',prefilter=0,first_obs=1000,
+addpath('..');
+generate_trend_stationary_AR1;
+
+calib_smoother(datafile='AR1_trend_data_with_constant',prefilter=0,first_obs=1000,
 //         filter_decomposition,
         filtered_vars, filter_step_ahead = [1,2,4]) P_obs Y_obs;
 
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameters do not match')
 end

tests/observation_trends_and_prefiltering/calib_smoother/Tr_no_prefilter_calib_smoother.mod

@@ -1,12 +1,16 @@
 @#include "../Trend_model_calib_no_prefilter_common.inc"
 options_.filter_decomposition=1;
-calib_smoother(datafile='../AR1_trend_data_with_constant',prefilter=0,
+
+addpath('..');
+generate_trend_stationary_AR1;
+
+calib_smoother(datafile='AR1_trend_data_with_constant',prefilter=0,
 //         filter_decomposition,
         filtered_vars, filter_step_ahead = [1,2,4]) P_obs Y_obs;
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameters do not match')
 end

tests/observation_trends_and_prefiltering/calib_smoother/Tr_no_prefilter_loglin_calib_smoother.mod

@@ -1,14 +1,17 @@
 @#include "../Trend_exp_model_calib_no_prefilter_common.inc"
 options_.filter_decomposition=1;
 
-calib_smoother(datafile='../Exp_AR1_trend_data_with_constant',prefilter=0,loglinear,
+addpath('..');
+generate_trend_stationary_AR1;
+
+calib_smoother(datafile='Exp_AR1_trend_data_with_constant',prefilter=0,loglinear,
 //         filter_decomposition,
         filtered_vars, filter_step_ahead = [1,2,4]) P_obs Y_obs junk2;
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params)./loaded_par.orig_params))>0.03
     error('Parameters do not match')
 end

tests/observation_trends_and_prefiltering/calib_smoother/Tr_prefil_f_obs_loglin_cal_smoother.mod

@@ -1,13 +1,16 @@
 @#include "../Trend_exp_model_calib_prefilter_common.inc"
 options_.filter_decomposition=1;
 
-calib_smoother(datafile='../Exp_AR1_trend_data_with_constant',prefilter=1,loglinear,first_obs=1000,
+addpath('..');
+generate_trend_stationary_AR1;
+
+calib_smoother(datafile='Exp_AR1_trend_data_with_constant',prefilter=1,loglinear,first_obs=1000,
 //         filter_decomposition,
         filtered_vars, filter_step_ahead = [1,2,4]) P_obs Y_obs junk2;
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params([1:4,7:8]))./loaded_par.orig_params([1:4,7:8])))>0.03
     error('Parameters do not match')
 end

tests/observation_trends_and_prefiltering/calib_smoother/Tr_prefilt_first_obs_cal_smooth.mod

@@ -1,13 +1,16 @@
 @#include "../Trend_model_calib_prefilter_common.inc"
 options_.filter_decomposition=1;
 
-calib_smoother(datafile='../AR1_trend_data_with_constant',prefilter=1,first_obs=1000,
+addpath('..');
+generate_trend_stationary_AR1;
+
+calib_smoother(datafile='AR1_trend_data_with_constant',prefilter=1,first_obs=1000,
 //         filter_decomposition,
         filtered_vars, filter_step_ahead = [1,2,4]) P_obs Y_obs;
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params([1:4,7:8]))./loaded_par.orig_params([1:4,7:8])))>0.03
     error('Parameters do not match')
 end

tests/observation_trends_and_prefiltering/calib_smoother/Tr_prefilter_calib_smoother.mod

@@ -1,13 +1,16 @@
 @#include "../Trend_model_calib_prefilter_common.inc"
 options_.filter_decomposition=1;
 
-calib_smoother(datafile='../AR1_trend_data_with_constant',prefilter=1,
+addpath('..');
+generate_trend_stationary_AR1;
+
+calib_smoother(datafile='AR1_trend_data_with_constant',prefilter=1,
 //         filter_decomposition,
         filtered_vars, filter_step_ahead = [1,2,4]) P_obs Y_obs;
-load('../AR1_trend_data_with_constant');
+load('AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params([1:4,7:8]))./loaded_par.orig_params([1:4,7:8])))>0.03
     error('Parameters do not match')
 end

tests/observation_trends_and_prefiltering/calib_smoother/Tr_prefilter_loglin_calib_smoother.mod

@@ -1,13 +1,16 @@
 @#include "../Trend_exp_model_calib_prefilter_common.inc"
 options_.filter_decomposition=1;
 
-calib_smoother(datafile='../Exp_AR1_trend_data_with_constant',prefilter=1,loglinear,
+addpath('..');
+generate_trend_stationary_AR1;
+
+calib_smoother(datafile='Exp_AR1_trend_data_with_constant',prefilter=1,loglinear,
 //         filter_decomposition,
         filtered_vars, filter_step_ahead = [1,2,4]) P_obs Y_obs junk2;
-load('../Exp_AR1_trend_data_with_constant');
+load('Exp_AR1_trend_data_with_constant');
 @#include "../Trend_load_data_common.inc" 
 
-loaded_par=load('../orig_params');
+loaded_par=load('orig_params');
 if max(abs((M_.params-loaded_par.orig_params([1:4,7:8]))./loaded_par.orig_params([1:4,7:8])))>0.03
     error('Parameters do not match')
 end

tests/observation_trends_and_prefiltering/generate_trend_stationary_AR1.m

@@ -0,0 +1,85 @@
+function generate_trend_stationary_AR1
+n_periods=10000;
+rho_y=0.5;
+rho_p=0.5;
+g_y=0.0001;
+g_p=-0.0001;
+const_y=2;
+const_p=2;
+sigma_y=0.001;
+sigma_p=0.001;
+
+orig_params=[rho_y rho_p  g_y g_p sigma_y sigma_p]';
+param_names=char('rho_y','rho_p','g_y','g_p','sigma_y','sigma_p');
+
+save orig_params_prefilter orig_params param_names
+
+orig_params=[rho_y rho_p  g_y g_p const_y const_p sigma_y sigma_p]';
+param_names=char('rho_y','rho_p','g_y','g_p','const_y','const_p','sigma_y','sigma_p');
+
+save orig_params orig_params param_names
+
+jumping_covariance=diag([1e-8; 1e-8; 1e-16; 1e-16; 1e-8; 1e-8; 1e-12; 1e-12;])^-1;
+save MCMC_jump_covar jumping_covariance
+jumping_covariance=diag([1e-8; 1e-8; 1e-16; 1e-16; 1e-12; 1e-12;])^-1;
+save MCMC_jump_covar_prefilter jumping_covariance
+
+%% data without constant
+log_P=zeros(1,n_periods);
+log_Y=zeros(1,n_periods);
+junk2_orig=zeros(1,n_periods);
+for ii=2:n_periods
+    log_P(ii)=rho_p*log_P(ii-1)+sigma_p*randn;
+    log_Y(ii)=rho_y*log_Y(ii-1)+sigma_y*randn;
+    junk2_orig(ii)=0.9*junk2_orig(ii-1)+randn;
+end
+%add trend
+log_P=log_P+g_p*(1:n_periods);
+log_Y=log_Y+g_y*(1:n_periods);
+
+Y_obs=exp(log_Y);
+P_obs=exp(log_P);
+junk2=exp(junk2_orig);
+save Exp_AR1_trend_data_no_constant Y_obs P_obs junk2
+% 
+% [b_p,~,~,~,stats_p] = regress(log(P_obs(2:end))',[ones(n_periods-1,1) (2:n_periods)' log(P_obs(1:end-1)')]);
+% [b_y,~,~,~,stats_y] = regress(log(Y_obs(2:end))',[ones(n_periods-1,1) (2:n_periods)' log(Y_obs(1:end-1)')]);
+
+Y_obs=log_Y;
+P_obs=log_P;
+junk2=junk2_orig;  
+save AR1_trend_data_no_constant Y_obs P_obs junk2
+
+% [b_p,~,~,~,stats_p] = regress((P_obs(2:end))',[ones(n_periods-1,1) (2:n_periods)' (P_obs(1:end-1)')]);
+% [b_y,~,~,~,stats_y] = regress((Y_obs(2:end))',[ones(n_periods-1,1) (2:n_periods)' (Y_obs(1:end-1)')]);
+
+%% data with constant
+log_P=zeros(1,n_periods);
+log_Y=zeros(1,n_periods);
+log_P(1,1)=const_p;
+log_Y(1,1)=const_y;
+for ii=2:n_periods
+    log_P(ii)=(1-rho_p)*const_p+rho_p*log_P(ii-1)+sigma_p*randn;
+    log_Y(ii)=(1-rho_y)*const_y+rho_y*log_Y(ii-1)+sigma_y*randn;
+end
+%add trend
+log_P=log_P+g_p*(1:n_periods);
+log_Y=log_Y+g_y*(1:n_periods);
+
+Y_obs=exp(log_Y);
+P_obs=exp(log_P);
+junk2=exp(junk2_orig);
+save Exp_AR1_trend_data_with_constant Y_obs P_obs junk2
+
+% [b,bint,r,rint,stats] = regress(log(P_obs(2:end))',[ones(n_periods-1,1) (2:n_periods)' log(P_obs(1:end-1)')]);
+% [b,bint,r,rint,stats] = regress(log(Y_obs(2:end))',[ones(n_periods-1,1) (2:n_periods)' log(Y_obs(1:end-1)')]);
+
+Y_obs=log_Y;
+P_obs=log_P;
+junk2=junk2_orig;  
+save AR1_trend_data_with_constant Y_obs P_obs junk2
+
+% [b_p,~,~,~,stats_p] = regress((P_obs(2:end))',[ones(n_periods-1,1) (2:n_periods)' (P_obs(1:end-1)')]);
+% [b_y,~,~,~,stats_y] = regress((Y_obs(2:end))',[ones(n_periods-1,1) (2:n_periods)' (Y_obs(1:end-1)')]);
+
+

tests/optimizers/fs2000_2.mod

@@ -1,3 +1,8 @@
 @#include "fs2000.common.inc"
 
-estimation(mode_compute=2,order=1, datafile='../fs2000/fsdat_simul', nobs=192, mh_replic=0);
+estimation(mode_compute=2,order=1, datafile='../fs2000/fsdat_simul', nobs=192, mh_replic=0,
+optim=(
+'MaxIter',5000,
+'TolFun',1e-4,
+'TolX',1e-4)
+);