Merge branch 'master' into ecb-master

Fixed conflict in preprocessor/ModFile.cc

NEWS

@@ -1,4 +1,48 @@
-Announcement for Dynare 4.5.0 (on 2013-12-16)
+Announcement for Dynare 4.5.1 (on 2017-08-24)
+=============================================
+
+We are pleased to announce the release of Dynare 4.5.1.
+
+This is a bugfix release.
+
+The Windows packages are already available for download at:
+
+ http://www.dynare.org/download/dynare-stable
+
+The Mac and GNU/Linux packages (for Debian and Ubuntu LTS) should follow soon.
+
+This release is compatible with MATLAB versions 7.3 (R2006b) to 9.2 (R2017a)
+and with GNU Octave versions 4.2.
+
+Here is a list of the problems identified in version 4.5.0 and that have been
+fixed in version 4.5.1:
+
+
+ - Fixed out of memory issue with simpsa optimization algorithm.
+
+ - Added missing plots for measurement errors with `generate_trace_plot`
+   command.
+
+ - Posterior moments after MCMC for very big models were not correctly computed
+   and their plotting might crash Dynare.
+
+ - Results of the posterior conditional variance decomposition after MCMC were
+   not correctly computed.
+
+ - Options `use_shock_groups` and `colormap` of the `shock_decomposition`
+   command were not working.
+
+ - Added a clean error message if sensitivity toolbox is used with recursive
+   estimation.
+
+ - Computation of posterior filtered variables was crashing in models with only
+   one variable.
+
+ - Fixed various typos and errors in the reference manual.
+
+
+
+Announcement for Dynare 4.5.0 (on 2017-06-11)
 =============================================
 
 We are pleased to announce the release of Dynare 4.5.0.
@@ -19,11 +63,6 @@ This release is compatible with MATLAB versions ranging from 7.3 (R2006b) to
 Here is the list of major user-visible changes:
 
 
-
-Dynare 4.5
-==========
-
-
  - Ramsey policy
 
    + Added command `ramsey_model` that builds the expanded model with

doc/dynare.texi

@@ -1220,9 +1220,10 @@ or @var{Sigma_e} to name your variable. Not conforming to this rule might yield
 error messages or crashes. Second, to minimize interference with MATLAB or Octave functions that
 may be called by Dynare or user-defined steady state files, it is recommended to avoid using the name
 of MATLAB functions. In particular when working with steady state files, do not use correctly-spelled greek 
-names like @var{alpha}, because there are Matlab functions of the same name. Rather go for @var{alppha} or so.
+names like @var{alpha}, because there are Matlab functions of the same name. Rather go for @var{alppha} or @var{alph}.
 Lastly, please do not name a variable or parameter @var{i}. This may interfere with the imaginary 
-number @var{i} and the index in many loops. Rather, name investment @var{invest}.
+number @var{i} and the index in many loops. Rather, name investment @var{invest}. Using @var{inv} is also not recommended
+as it alread denotes the inverse operator.
 
 Declarations of variables and parameters are made with the following commands:
 
@@ -1984,7 +1985,7 @@ different purposes by allowing the user to attach arbitrary informations to each
 equation and to recover them at runtime. For instance, it is possible to name the
 equations with a @code{name}-tag, using a syntax like:
 @example
-mode;
+model;
     ...
     [name = 'Budget constraint']
     c + k = k^theta*A;
@@ -1996,7 +1997,7 @@ of the model is tagged with a name, the @code{resid} command
 will display the name of the equations (which may be more informative than the
 equation numbers) in addition to the equation number. Several tags for one equation can be separated using a comma. 
 @example
-mode;
+model;
     ...
    [name='Taylor rule',mcp = 'r > -1.94478']
    r = rho*r(-1) + (1-rho)*(gpi*Infl+gy*YGap) + e;
@@ -2362,7 +2363,7 @@ static equilibria, but not necessarily.
 One typical application is to consider an economy at the equilibrium at time 0,
 trigger a shock in first period, and study the trajectory of return to
 the initial equilibrium. To do that, one needs @code{initval} and
-@code{shocks} (@pxref{Shocks on exogenous variables}.
+@code{shocks} (@pxref{Shocks on exogenous variables}).
 
 Another one is to study, how an economy, starting from arbitrary
 initial conditions at time 0 converges toward equilibrium. 
@@ -3233,7 +3234,7 @@ option, @pxref{Model declaration})
 Trust-region algorithm on the entire model.
 
 @item 10
-Levenberg-Marquardt mixed compleproblem (LMMCP) solver
+Levenberg-Marquardt mixed complementarity problem (LMMCP) solver
 (@cite{Kanzow and Petra 2004})
 
 @item 11
@@ -3680,7 +3681,7 @@ to the equation appears currently on the left hand side and where no
 backward looking endogenous variables appear. The block has the form:
 @math{y_{j,t} = f_j(y_t, y_{t+1}, \ldots, y_{t+k})}.
 
-@item SOLVE FORWARD @var{x}
+@item SOLVE BACKWARD @var{x}
 The block contains only equations where endogenous variable attributed
 to the equation does not appear currently on the left hand side and
 where no forward looking endogenous variables appear. The block has
@@ -4386,7 +4387,7 @@ If option @code{periods} is present, sets @code{oo_.skewness},
 MATLAB/Octave vectors of the global workspace with the same name as
 the endogenous variables.
 
-If options @code{irf} is different from zero, sets @code{oo_.irfs}
+If option @code{irf} is different from zero, sets @code{oo_.irfs}
 (see below) and also saves the IRFs in MATLAB/Octave vectors of
 the global workspace (this latter way of accessing the IRFs is
 deprecated and will disappear in a future version).
@@ -4812,7 +4813,7 @@ non-decreasing tuples are stored, @i{i.e.} those for which @math{i_1
 \leq i_2 \leq i_3}. The columns are arranged in the lexicographical
 order of non-decreasing tuples. Also note that for tuples that have
 three distinct indices (@i{i.e.} @math{i_1 \neq i_2} and @math{i_1
-\neq i_3} and @math{i_2 \neq i_3}, since these elements are stored
+\neq i_3} and @math{i_2 \neq i_3}), since these elements are stored
 only once but appears six times in the unfolded @math{G_3} matrix,
 they must be multiplied by 6 when computing the decision
 rules. Similarly, for those tuples that have two equal indices
@@ -4992,7 +4993,7 @@ likelihood estimation. In a Bayesian estimation context, sets a lower bound
 only effective while maximizing the posterior kernel. This lower bound does not
 modify the shape of the prior density, and is only aimed at helping the
 optimizer in identifying the posterior mode (no consequences for the MCMC). For
-some prior densities (namely inverse gamma, gamma, uniform, beta or weibull) it is
+some prior densities (namely inverse gamma, gamma, uniform, beta or Weibull) it is
 possible to shift the support of the prior distributions to the left or the right using
 @ref{prior_3rd_parameter}. In this case the prior density is effectively
 modified (note that the truncated Gaussian density is not implemented in
@@ -5020,7 +5021,7 @@ that @code{inv_gamma_pdf} is equivalent to
 @item @var{PRIOR_3RD_PARAMETER}
 @anchor{prior_3rd_parameter}
 A third parameter of the prior used for generalized beta distribution,
-generalized gamma, generalized weibull and for the uniform distribution. Default: @code{0}
+generalized gamma, generalized Weibull and for the uniform distribution. Default: @code{0}
 
 @item @var{PRIOR_4TH_PARAMETER}
 @anchor{prior_4th_parameter}
@@ -5921,7 +5922,7 @@ distribution of IRFs. The length of the IRFs are controlled by the
 weight of  the DSGE prior  of the VAR model  is calibrated to  the value
 passed (see @cite{Del  Negro and Schorfheide (2004)}). It represents ratio of dummy over actual observations. 
 To assure that the prior is proper, the value must be bigger than @math{(k+n)/T}, 
-where @math{k} is the number of estimated parameters, @math{n} is the number of observables, #
+where @math{k} is the number of estimated parameters, @math{n} is the number of observables,
 and @math{T} is the number of observations. NB:  The previous method
 of   declaring  @code{dsge_prior_weight} as a parameter and then
 calibrating it is now deprecated and will be removed in a future release
@@ -6833,8 +6834,8 @@ observation for which the forecast has been made.
 @defvr {MATLAB/Octave variable} oo_.Smoother.State_uncertainty
 @anchor{oo_.Smoother.State_uncertainty}
 Three-dimensional array set by the @code{estimation} command (if used with the
-@code{smoother}) without Metropolis, 
+@code{smoother} option) without Metropolis, 
-or by the @code{calib_smoother} command, if the @code{o_smoothed_state_uncertainty} option
+or by the @code{calib_smoother} command, if the @code{smoothed_state_uncertainty} option
 has been requested.
 Contains the series of covariance matrices for the state estimate given the full data
 from the Kalman smoother. The @code{M_.endo_nbr} times @code{M_.endo_nbr} times
@@ -7775,7 +7776,7 @@ uncertainty, but ignoring the effect of measurement error on
 observed variables
 
 @item HPDsup
-Lower bound of a 90% HPD interval due to parameter uncertainty, but 
+Upper bound of a 90% HPD forecast interval due to parameter uncertainty, but 
 ignoring the effect of measurement error on 
 observed variables
 
@@ -7785,7 +7786,7 @@ to change the size of the HPD interval} of forecast for observed variables
 due to parameter uncertainty and measurement error 
 
 @item HPDsup_ME
-Lower bound of a 90% HPD interval of forecast for observed variables 
+Upper bound of a 90% HPD interval of forecast for observed variables 
 due to parameter uncertainty and measurement error 
 
 @item Mean
@@ -8385,7 +8386,7 @@ Note that the use of this block requires the use of a constrained optimizer, @i{
 
 @example
 
-osr_param_bounds;
+osr_params_bounds;
 gamma_inf_, 0, 2.5;
 end;
 
@@ -9275,7 +9276,7 @@ named as:
 @code{<namendo>_vs_<namexo>}: for entries of the matrix of the shocks.
 @end itemize
 
-The following files are stored in each directory (we stick with prior sample but similar convetins are used for MC samples):
+The following files are stored in each directory (we stick with prior sample but similar conventions are used for MC samples):
 @itemize
 @item
 @code{<mod_file>_prior_<namendo>_vs_<namexo>.fig}: histogram and CDF plot of the MC sample of the individual entry
@@ -9551,7 +9552,7 @@ In case of estimating a MS-DSGE model,@footnote{An example can be found at @uref
 @table @code
 
 @item parameters = @var{[LIST OF PARAMETERS]}
-This option specifies which parameters are controlled by this Markov Chai
+This option specifies which parameters are controlled by this Markov Chain.
 
 @item number_of_lags = @var{DOUBLE}
 Provides the number of lags that each parameter can take within each regime in this chain. 
@@ -10970,8 +10971,8 @@ This file computes the steady state. It begins with:
 @@#include "modeqs.mod"
 @end example
 Then it initializes parameters (including @code{lab_rat}, excluding
-@math{\alpha}, computes the steady state (using guess values for
+@math{\alpha}), computes the steady state (using guess values for
-endogenous, including @math{\alpha}, then saves values of parameters
+endogenous, including @math{\alpha}), then saves values of parameters
 and endogenous at steady state in a file, using the
 @code{save_params_and_steady_state} command.
 
@@ -11028,7 +11029,7 @@ array (@code{rhos}).
 The advantage of this method is that it uses a shorter syntax, since
 list of values directly given in the loop construct. Note that values
 are given as character strings (the macro-processor does not know
-floating point values. The inconvenient is that you can not reuse an
+floating point values). The inconvenience is that you can not reuse an
 array stored in a MATLAB/Octave variable.
 
 @end table
@@ -11987,7 +11988,7 @@ ans =
 
 @deftypefn{dates} {@var{C} =} gt (@var{A}, @var{B})
 
-Overloads the Matlab/Octave @code{gt} (greater than, @code{>=}) operator. @dates objects @var{A} and @var{B} must have the  same number of elements (say, @code{n}). The returned argument is a @code{n} by @code{1} vector of zeros and ones. The i-th element of @var{C} is equal to @code{1} if and only if the date @code{A(i)} is posterior to the date @code{B(i)}.
+Overloads the Matlab/Octave @code{gt} (greater than, @code{>}) operator. @dates objects @var{A} and @var{B} must have the  same number of elements (say, @code{n}). The returned argument is a @code{n} by @code{1} vector of zeros and ones. The i-th element of @var{C} is equal to @code{1} if and only if the date @code{A(i)} is posterior to the date @code{B(i)}.
 
 @examplehead
 @example
@@ -12136,7 +12137,7 @@ ans =
 
 @deftypefn{dates} {@var{C} =} lt (@var{A}, @var{B})
 
-Overloads the Matlab/Octave @code{lt} (less than, @code{<=}) operator. @dates objects @var{A} and @var{B} must have the  same number of elements (say, @code{n}). The returned argument is a @code{n} by @code{1} vector of zeros and ones. The i-th element of @var{C} is equal to @code{1} if and only if the date @code{A(i)} preceeds the date @code{B(i)}.
+Overloads the Matlab/Octave @code{lt} (less than, @code{<}) operator. @dates objects @var{A} and @var{B} must have the  same number of elements (say, @code{n}). The returned argument is a @code{n} by @code{1} vector of zeros and ones. The i-th element of @var{C} is equal to @code{1} if and only if the date @code{A(i)} preceeds the date @code{B(i)}.
 
 @examplehead
 @example
@@ -14010,7 +14011,7 @@ command. Default: @code{`!'}
 @end defmethod
 
 @anchor{addGraph}
-@defmethod Report addGraph axisShape, data, graphDirName, graphName, graphSize, height, showGrid, showLegend, legendAt, showLegendBox, legendLocation, legendOrientation, legendFontSize, miscTikzAxisOptions, miscTikzPictureOptions, seriesToUse, shade, shadeColor, shadeOpacity, tickFontSize, title, titleFontSize, titleFormat, width, writeCSV, xlabel, ylabel, xAxisTight, xrange, xTicks, xTickLabels, xTickLabelAnchor, xTickLabelRotation, yAxisTight, yTickLabelFixed, yTickLabelPrecision, yTickLabelScaled, yTickLabelZeroFill, yrange, showZeroline, zeroLineColor
+@defmethod Report addGraph axisShape, data, graphDirName, graphName, graphSize, height, showGrid, showLegend, legendAt, showLegendBox, legendLocation, legendOrientation, legendFontSize, miscTikzAxisOptions, miscTikzPictureOptions, seriesToUse, shade, shadeColor, shadeOpacity, tickFontSize, title, titleFontSize, titleFormat, width, writeCSV, xlabel, ylabel, xAxisTight, xrange, xTicks, xTickLabels, xTickLabelAnchor, xTickLabelRotation, yAxisTight, yTickLabelFixed, yTickLabelPrecision, yTickLabelScaled, yTickLabelZeroFill, yrange, showZeroLine, zeroLineColor
 Adds a @code{Graph} to a @code{Section}.
 @optionshead
 @table @code

matlab/discretionary_policy_1.m

@@ -68,6 +68,15 @@ if options_.steadystate_flag
 end
 [U,Uy,W] = feval([M_.fname,'_objective_static'],zeros(endo_nbr,1),[], M_.params);
 if any(any(Uy~=0))
+    non_zero_derivs=find(any(Uy~=0));
+    for ii=1:length(non_zero_derivs)
+        non_zero_deriv_names{ii,1}=deblank(M_.endo_names(non_zero_derivs(ii),:));
+    end
+    disp_string=[non_zero_deriv_names{1,:}];
+    for ii=2:size(non_zero_deriv_names,1)
+        disp_string=[disp_string,', ',non_zero_deriv_names{ii,:}];
+    end
+    fprintf('\nThe derivative of the objective function w.r.t. to variable(s) %s is not 0\n',disp_string)
     error(['discretionary_policy: the objective function must have zero ' ...
            'first order derivatives'])
 end

matlab/dsge_simulated_theoretical_conditional_variance_decomposition.m

@@ -109,10 +109,10 @@ for file = 1:NumberOfDrawsFiles
     for linee = 1:NumberOfDraws
         linea = linea+1;
         if isdrsaved
-            set_parameters(pdraws{linee,1});% Needed to update the covariance matrix of the state innovations.
+            M_=set_parameters_locally(M_,pdraws{linee,1});% Needed to update the covariance matrix of the state innovations.
             dr = pdraws{linee,2};
         else
-            set_parameters(pdraws{linee,1});
+            M_=set_parameters_locally(M_,pdraws{linee,1});
             [dr,info,M_,options_,oo_] = resol(0,M_,options_,oo_);
         end
         if first_call

matlab/dsge_simulated_theoretical_correlation.m

@@ -106,7 +106,7 @@ for file = 1:NumberOfDrawsFiles
         if isdrsaved
             dr = pdraws{linee,2};
         else
-            set_parameters(pdraws{linee,1});
+            M_=set_parameters_locally(M_,pdraws{linee,1});
             [dr,info,M_,options_,oo_] = resol(0,M_,options_,oo_);
         end
         tmp = th_autocovariances(dr,ivar,M_,options_,nodecomposition);

matlab/dsge_simulated_theoretical_covariance.m

@@ -105,7 +105,7 @@ for file = 1:NumberOfDrawsFiles
         if isdrsaved
             dr = pdraws{linee,2};
         else
-            set_parameters(pdraws{linee,1});
+            M_=set_parameters_locally(M_,pdraws{linee,1});
             [dr,info,M_,options_,oo_] = resol(0,M_,options_,oo_);
         end
         tmp = th_autocovariances(dr,ivar,M_,options_,nodecomposition);

matlab/dsge_simulated_theoretical_variance_decomposition.m

@@ -113,7 +113,7 @@ for file = 1:NumberOfDrawsFiles
         if isdrsaved
             dr = pdraws{linee,2};
         else
-            set_parameters(pdraws{linee,1});
+            M_=set_parameters_locally(M_,pdraws{linee,1});
             [dr,info,M_,options_,oo_] = resol(0,M_,options_,oo_);
         end
         if file==1 && linee==1

matlab/generate_trace_plots.m

@@ -50,6 +50,11 @@ for ii=1:size(estim_params_.var_exo,1)
     trace_plot(options_,M_,estim_params_,'StructuralShock',chain_number,parameter_name)
 end
 
+for ii=1:size(estim_params_.var_endo,1)
+    parameter_name=deblank(M_.endo_names(estim_params_.var_endo(ii,1),:));
+    trace_plot(options_,M_,estim_params_,'MeasurementError',chain_number,parameter_name)
+end
+
 for ii=1:size(estim_params_.corrn,1)
     parameter_name_1=deblank(M_.endo_names(estim_params_.corrn(ii,1),:));
     parameter_name_2=deblank(M_.endo_names(estim_params_.corrn(ii,2),:));

matlab/optimization/simpsa.m

@@ -142,8 +142,9 @@ DEFAULT_OPTIONS = simpsaset('TEMP_START',[],...  % starting temperature (if none
 OPTIONS = simpsaset(DEFAULT_OPTIONS,OPTIONS);
 
 % store options in OUTPUT
-
+if nargout>3
-OUTPUT.OPTIONS = OPTIONS;
+    OUTPUT.OPTIONS = OPTIONS;
+end
 
 % initialize simplex
 % ------------------
@@ -175,13 +176,13 @@ end
 
 % initialize OUTPUT structure
 % ---------------------------
-
+if nargout>3
-OUTPUT.TEMPERATURE = zeros(OPTIONS.MAX_ITER_TOTAL,1);
+    OUTPUT.TEMPERATURE = zeros(OPTIONS.MAX_ITER_TOTAL,1);
-OUTPUT.SIMPLEX = zeros(NDIM+1,NDIM,OPTIONS.MAX_ITER_TOTAL);
+    OUTPUT.SIMPLEX = zeros(NDIM+1,NDIM,OPTIONS.MAX_ITER_TOTAL);
-OUTPUT.SIMPLEX_BEST = zeros(OPTIONS.MAX_ITER_TOTAL,NDIM);
+    OUTPUT.SIMPLEX_BEST = zeros(OPTIONS.MAX_ITER_TOTAL,NDIM);
-OUTPUT.COSTS = zeros(OPTIONS.MAX_ITER_TOTAL,NDIM+1);
+    OUTPUT.COSTS = zeros(OPTIONS.MAX_ITER_TOTAL,NDIM+1);
-OUTPUT.COST_BEST = zeros(OPTIONS.MAX_ITER_TOTAL,1);
+    OUTPUT.COST_BEST = zeros(OPTIONS.MAX_ITER_TOTAL,1);
-
+end
 % initialize iteration data
 % -------------------------
 
@@ -304,17 +305,19 @@ while 1
         Y = help(:,2);
         P = help(:,3:end);
 
-        % store temperature at current iteration
+        if nargout>3
-        OUTPUT.TEMPERATURE(nITERATIONS) = TEMP;
+            % store temperature at current iteration
+            OUTPUT.TEMPERATURE(nITERATIONS) = TEMP;
 
-        % store information about simplex at the current iteration
+            % store information about simplex at the current iteration
-        OUTPUT.SIMPLEX(:,:,nITERATIONS) = P;
+            OUTPUT.SIMPLEX(:,:,nITERATIONS) = P;
-        OUTPUT.SIMPLEX_BEST(nITERATIONS,:) = PBEST;
+            OUTPUT.SIMPLEX_BEST(nITERATIONS,:) = PBEST;
-
-        % store cost function value of best vertex in current iteration
-        OUTPUT.COSTS(nITERATIONS,:) = Y;
-        OUTPUT.COST_BEST(nITERATIONS) = YBEST;
 
+            % store cost function value of best vertex in current iteration
+            OUTPUT.COSTS(nITERATIONS,:) = Y;
+            OUTPUT.COST_BEST(nITERATIONS) = YBEST;
+        end
+        
         if strcmp(OPTIONS.DISPLAY,'iter')
             disp(sprintf('%5.0f      %5.0f       %12.6g     %15.6g      %12.6g       %s',nITERATIONS,nFUN_EVALS,Y(1),YBEST,TEMP,ALGOSTEP));
         end
@@ -452,21 +455,23 @@ end
 X = transpose(PBEST);
 FVAL = YBEST;
 
-% store number of function evaluations
+if nargout>3
-OUTPUT.nFUN_EVALS = nFUN_EVALS;
+    % store number of function evaluations
-
+    OUTPUT.nFUN_EVALS = nFUN_EVALS;
-% store number of iterations
+    
-OUTPUT.nITERATIONS = nITERATIONS;
+    % store number of iterations
-
+    OUTPUT.nITERATIONS = nITERATIONS;
-% trim OUTPUT data structure
+    
-OUTPUT.TEMPERATURE(nITERATIONS+1:end) = [];
+    % trim OUTPUT data structure
-OUTPUT.SIMPLEX(:,:,nITERATIONS+1:end) = [];
+    OUTPUT.TEMPERATURE(nITERATIONS+1:end) = [];
-OUTPUT.SIMPLEX_BEST(nITERATIONS+1:end,:) = [];
+    OUTPUT.SIMPLEX(:,:,nITERATIONS+1:end) = [];
-OUTPUT.COSTS(nITERATIONS+1:end,:) = [];
+    OUTPUT.SIMPLEX_BEST(nITERATIONS+1:end,:) = [];
-OUTPUT.COST_BEST(nITERATIONS+1:end) = [];
+    OUTPUT.COSTS(nITERATIONS+1:end,:) = [];
-
+    OUTPUT.COST_BEST(nITERATIONS+1:end) = [];
-% store the amount of time needed in OUTPUT data structure
+    
-OUTPUT.TIME = toc;
+    % store the amount of time needed in OUTPUT data structure
+    OUTPUT.TIME = toc;
+end
 
 return
 

matlab/osr1.m

@@ -126,7 +126,7 @@ elseif isequal(options_.osr.opt_algo,11)
 else
 
     if ~isempty(M_.osr.param_bounds) && ~(ismember(options_.osr.opt_algo,[1,2,5,9]) || ischar(options_.osr.opt_algo))
-        error('OSR: OSR with bounds on parameters requires a constrained optimizer, i.e. 1,2,5, or 9.')
+        error('OSR: OSR with bounds on parameters requires a constrained optimizer, i.e. opt_algo= 1,2,5, or 9.')
     end
     %%do actual optimization
     [p, f, exitflag] = dynare_minimize_objective(str2func('osr_obj'),t0,options_.osr.opt_algo,options_,M_.osr.param_bounds,cellstr(M_.param_names(i_params,:)),[],[], i_params,...

matlab/plot_shock_decomposition.m

@@ -30,7 +30,6 @@ function [z, steady_state] = plot_shock_decomposition(M_,oo_,options_,varlist)
 
 options_.nodisplay = options_.plot_shock_decomp.nodisplay;
 options_.graph_format = options_.plot_shock_decomp.graph_format;
-options_.use_shock_groups = options_.plot_shock_decomp.use_shock_groups;
 
 % indices of endogenous variables
 if size(varlist,1) == 0

matlab/prior_posterior_statistics_core.m

@@ -260,23 +260,23 @@ for b=fpar:B
             else
                 constant_part=repmat(SteadyState(dr.order_var)',[length(options_.filter_step_ahead),1,gend+max(options_.filter_step_ahead)]);
             end
-
+            stock_filter_step_ahead(:,dr.order_var,:,irun(4)) = aK(options_.filter_step_ahead,1:endo_nbr,:) + constant_part;
-            stock_filter_step_ahead(:,dr.order_var,:,irun(4)) = aK(options_.filter_step_ahead,1:endo_nbr,:) + ...
-                constant_part;
-
             %now add trend to observables
             for ii=1:length(options_.filter_step_ahead)
                 if options_.prefilter
-                    stock_filter_step_ahead(ii,IdObs,:,irun(4)) = squeeze(stock_filter_step_ahead(ii,IdObs,:,irun(4)))...
+                    zdim = size(stock_filter_step_ahead(ii,IdObs,:,irun(4)));
-                        +repmat(mean_correction(:,1),1,gend+max(options_.filter_step_ahead))... %constant correction
+                    squeezed = reshape(stock_filter_step_ahead(ii,IdObs,:,irun(4)), [zdim(2:end) 1]);
+                    stock_filter_step_ahead(ii,IdObs,:,irun(4)) = squeezed ...
+                        +repmat(mean_correction(:,1),1,gend+max(options_.filter_step_ahead)) ... %constant correction
                         +[trend_addition repmat(trend_addition(:,end),1,max(options_.filter_step_ahead))+trend_coeff*[1:max(options_.filter_step_ahead)]]; %trend
                 else
-                    stock_filter_step_ahead(ii,IdObs,:,irun(4)) = squeeze(stock_filter_step_ahead(ii,IdObs,:,irun(4)))...
+                    zdim = size(stock_filter_step_ahead(ii,IdObs,:,irun(4)));
+                    squeezed = reshape(stock_filter_step_ahead(ii,IdObs,:,irun(4)), [zdim(2:end) 1]);
+                    stock_filter_step_ahead(ii,IdObs,:,irun(4)) = squeezed ...
                         +[trend_addition repmat(trend_addition(:,end),1,max(options_.filter_step_ahead))+trend_coeff*[1:max(options_.filter_step_ahead)]]; %trend
                 end
             end
         end
-
         if horizon
             yyyy = alphahat(iendo,i_last_obs);
             yf = forcst2a(yyyy,dr,zeros(horizon,exo_nbr));

matlab/set_parameters.m

@@ -31,53 +31,6 @@ function set_parameters(xparam1)
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
 
-global estim_params_ M_
+global M_
 
-nvx = estim_params_.nvx;
+M_ = set_parameters_locally(M_, xparam1);
-ncx = estim_params_.ncx;
-nvn = estim_params_.nvn;
-ncn = estim_params_.ncn;
-np = estim_params_.np;
-Sigma_e = M_.Sigma_e;
-Correlation_matrix = M_.Correlation_matrix;
-offset = 0;
-
-% setting shocks variance on the diagonal of Covariance matrix; used later
-% for updating covariances
-if nvx
-    var_exo = estim_params_.var_exo;
-    for i=1:nvx
-        k = var_exo(i,1);
-        Sigma_e(k,k) = xparam1(i)^2;
-    end
-end
-% and update offset
-offset = offset + nvx + nvn;
-
-% correlations amonx shocks (ncx)
-if ncx
-    corrx = estim_params_.corrx;
-    for i=1:ncx
-        k1 = corrx(i,1);
-        k2 = corrx(i,2);
-        Correlation_matrix(k1,k2) = xparam1(i+offset);
-        Correlation_matrix(k2,k1) = Correlation_matrix(k1,k2);
-    end
-end
-%build covariance matrix from correlation matrix and variances already on
-%diagonal
-Sigma_e = diag(sqrt(diag(Sigma_e)))*Correlation_matrix*diag(sqrt(diag(Sigma_e)));
-if isfield(estim_params_,'calibrated_covariances')
-    Sigma_e(estim_params_.calibrated_covariances.position)=estim_params_.calibrated_covariances.cov_value;
-end
-
-% and update offset
-offset = offset + ncx + ncn;
-
-% structural parameters
-if np
-    M_.params(estim_params_.param_vals(:,1)) = xparam1(offset+1:end);
-end
-
-M_.Sigma_e = Sigma_e;
-M_.Correlation_matrix=Correlation_matrix;

matlab/set_parameters_locally.m

@@ -0,0 +1,85 @@
+function M_=set_parameters_locally(M_,xparam1)
+
+% function M_out=set_parameters(M_,xparam1)
+% Sets parameters value (except measurement errors)
+% This is called for computations such as IRF and forecast
+% when measurement errors aren't taken into account; in contrast to 
+% set_parameters.m, the global M_-structure is not altered
+%
+% INPUTS
+%    xparam1:   vector of parameters to be estimated (initial values)
+%    M_:        Dynare model-structure
+%
+% OUTPUTS
+%    M_:        Dynare model-structure
+%
+% SPECIAL REQUIREMENTS
+%    none
+
+% Copyright (C) 2017 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/>.
+
+global estim_params_
+
+nvx = estim_params_.nvx;
+ncx = estim_params_.ncx;
+nvn = estim_params_.nvn;
+ncn = estim_params_.ncn;
+np = estim_params_.np;
+Sigma_e = M_.Sigma_e;
+Correlation_matrix = M_.Correlation_matrix;
+offset = 0;
+
+% setting shocks variance on the diagonal of Covariance matrix; used later
+% for updating covariances
+if nvx
+    var_exo = estim_params_.var_exo;
+    for i=1:nvx
+        k = var_exo(i,1);
+        Sigma_e(k,k) = xparam1(i)^2;
+    end
+end
+% and update offset
+offset = offset + nvx + nvn;
+
+% correlations amonx shocks (ncx)
+if ncx
+    corrx = estim_params_.corrx;
+    for i=1:ncx
+        k1 = corrx(i,1);
+        k2 = corrx(i,2);
+        Correlation_matrix(k1,k2) = xparam1(i+offset);
+        Correlation_matrix(k2,k1) = Correlation_matrix(k1,k2);
+    end
+end
+%build covariance matrix from correlation matrix and variances already on
+%diagonal
+Sigma_e = diag(sqrt(diag(Sigma_e)))*Correlation_matrix*diag(sqrt(diag(Sigma_e)));
+if isfield(estim_params_,'calibrated_covariances')
+    Sigma_e(estim_params_.calibrated_covariances.position)=estim_params_.calibrated_covariances.cov_value;
+end
+
+% and update offset
+offset = offset + ncx + ncn;
+
+% structural parameters
+if np
+    M_.params(estim_params_.param_vals(:,1)) = xparam1(offset+1:end);
+end
+
+M_.Sigma_e = Sigma_e;
+M_.Correlation_matrix=Correlation_matrix;

preprocessor/DynamicModel.cc

@@ -44,6 +44,10 @@ DynamicModel::DynamicModel(SymbolTable &symbol_table_arg,
   max_endo_lag(0), max_endo_lead(0),
   max_exo_lag(0), max_exo_lead(0),
   max_exo_det_lag(0), max_exo_det_lead(0),
+  max_lag_orig(0), max_lead_orig(0),
+  max_endo_lag_orig(0), max_endo_lead_orig(0),
+  max_exo_lag_orig(0), max_exo_lead_orig(0),
+  max_exo_det_lag_orig(0), max_exo_det_lead_orig(0),
   dynJacobianColsNbr(0),
   global_temporary_terms(true)
 {
@@ -2594,7 +2598,15 @@ DynamicModel::writeOutput(ostream &output, const string &basename, bool block_de
       outstruct = "oo_.";
     }
 
-  output << modstruct << "lead_lag_incidence = [";
+  output << modstruct << "max_endo_lag_orig = " << max_endo_lag_orig << ";" << endl
+         << modstruct << "max_endo_lead_orig = " << max_endo_lead_orig << ";" << endl
+         << modstruct << "max_exo_lag_orig = " << max_exo_lag_orig << ";" << endl
+         << modstruct << "max_exo_lead_orig = " << max_exo_lead_orig << ";" << endl
+         << modstruct << "max_exo_det_lag_orig = " << max_exo_det_lag_orig << ";" << endl
+         << modstruct << "max_exo_det_lead_orig = " << max_exo_det_lead_orig << ";" << endl
+         << modstruct << "max_lag_orig = " << max_lag_orig << ";" << endl
+         << modstruct << "max_lead_orig = " << max_lead_orig << ";" << endl
+         << modstruct << "lead_lag_incidence = [";
   // Loop on endogenous variables
   int nstatic = 0,
     nfwrd   = 0,
@@ -3830,6 +3842,8 @@ DynamicModel::cloneDynamic(DynamicModel &dynamic_model) const
   for (size_t i = 0; i < static_only_equations.size(); i++)
     dynamic_model.addStaticOnlyEquation(static_only_equations[i]->cloneDynamic(dynamic_model),
                                         static_only_equations_lineno[i]);
+
+  dynamic_model.setLeadsLagsOrig();
 }
 
 void
@@ -4008,6 +4022,55 @@ DynamicModel::findUnusedExogenous()
   return unusedExo;
 }
 
+void
+DynamicModel::setLeadsLagsOrig()
+{
+  set<pair<int, int> > dynvars;
+
+  for (int i = 0; i < (int) equations.size(); i++)
+    {
+      equations[i]->collectDynamicVariables(eEndogenous, dynvars);
+      equations[i]->collectDynamicVariables(eExogenous, dynvars);
+      equations[i]->collectDynamicVariables(eExogenousDet, dynvars);
+    }
+
+    for (set<pair<int, int> >::const_iterator it = dynvars.begin();
+         it != dynvars.end(); it++)
+    {
+      int lag = it->second;
+      SymbolType type = symbol_table.getType(it->first);
+
+      if (max_lead_orig < lag)
+        max_lead_orig= lag;
+      else if (-max_lag_orig > lag)
+        max_lag_orig = -lag;
+
+      switch (type)
+        {
+        case eEndogenous:
+          if (max_endo_lead_orig < lag)
+            max_endo_lead_orig = lag;
+          else if (-max_endo_lag_orig > lag)
+            max_endo_lag_orig = -lag;
+          break;
+        case eExogenous:
+          if (max_exo_lead_orig < lag)
+            max_exo_lead_orig = lag;
+          else if (-max_exo_lag_orig > lag)
+            max_exo_lag_orig = -lag;
+          break;
+        case eExogenousDet:
+          if (max_exo_det_lead_orig < lag)
+            max_exo_det_lead_orig = lag;
+          else if (-max_exo_det_lag_orig > lag)
+            max_exo_det_lag_orig = -lag;
+          break;
+        default:
+          break;
+        }
+    }
+}
+
 void
 DynamicModel::computeDerivIDs()
 {

preprocessor/DynamicModel.hh

@@ -61,6 +61,14 @@ private:
   //! Maximum lag and lead over deterministic exogenous variables (positive values)
   /*! Set by computeDerivIDs() */
   int max_exo_det_lag, max_exo_det_lead;
+  //! Maximum lag and lead over all types of variables (positive values) of original model
+  int max_lag_orig, max_lead_orig;
+  //! Maximum lag and lead over endogenous variables (positive values) of original model
+  int max_endo_lag_orig, max_endo_lead_orig;
+  //! Maximum lag and lead over exogenous variables (positive values) of original model
+  int max_exo_lag_orig, max_exo_lead_orig;
+  //! Maximum lag and lead over deterministic exogenous variables (positive values) of original model
+  int max_exo_det_lag_orig, max_exo_det_lead_orig;
 
   //! Cross reference information
   map<int, ExprNode::EquationInfo> xrefs;
@@ -287,6 +295,9 @@ public:
   //! Find exogenous variables not used in model
   set<int> findUnusedExogenous();
 
+  //! Set the max leads/lags of the original model
+  void setLeadsLagsOrig();
+
   //! Copies a dynamic model (only the equations)
   /*! It assumes that the dynamic model given in argument has just been allocated */
   void cloneDynamic(DynamicModel &dynamic_model) const;

preprocessor/DynareBison.yy

@@ -2688,8 +2688,8 @@ plot_shock_decomposition_options_list : plot_shock_decomposition_option COMMA pl
                                       | plot_shock_decomposition_option
                                       ;
 
-plot_shock_decomposition_option : o_psd_use_shock_groups
+plot_shock_decomposition_option : o_use_shock_groups
-                                | o_psd_colormap
+                                | o_colormap
                                 | o_psd_nodisplay
                                 | o_psd_graph_format
                                 | o_psd_detail_plot
@@ -3514,14 +3514,10 @@ o_lmmcp : LMMCP {driver.option_num("lmmcp.status", "1"); };
 o_occbin : OCCBIN {driver.option_num("occbin", "1"); };
 o_function : FUNCTION EQUAL filename { driver.option_str("function", $3); };
 o_sampling_draws : SAMPLING_DRAWS EQUAL INT_NUMBER { driver.option_num("sampling_draws",$3); };
-o_use_shock_groups : USE_SHOCK_GROUPS { driver.option_str("use_shock_groups","default"); }
+o_use_shock_groups : USE_SHOCK_GROUPS { driver.option_str("plot_shock_decomp.use_shock_groups","default"); }
-                   | USE_SHOCK_GROUPS EQUAL symbol { driver.option_str("use_shock_groups", $3); }
+                   | USE_SHOCK_GROUPS EQUAL symbol { driver.option_str("plot_shock_decomp.use_shock_groups", $3); }
                    ;
-o_psd_use_shock_groups : USE_SHOCK_GROUPS { driver.option_str("plot_shock_decomp.use_shock_groups","default"); }
-                       | USE_SHOCK_GROUPS EQUAL symbol { driver.option_str("plot_shock_decomp.use_shock_groups", $3); }
-                       ;
 o_colormap : COLORMAP EQUAL symbol { driver.option_num("plot_shock_decomp.colormap",$3); };
-o_psd_colormap : COLORMAP EQUAL symbol { driver.option_num("plot_shock_decomp.colormap",$3); };
 
 range : symbol ':' symbol
         {

scripts/dynare.el

@@ -28,8 +28,6 @@
 ;;   add this to your .emacs or site-init.el file:
 ;;
 ;;   (require 'dynare)
-;;   (autoload 'dynare-mode "dynare" "Enter dynare mode." t)
-;;   (setq auto-mode-alist (cons '("\\.mod\\'" . dynare-mode) auto-mode-alist))
 
 ;;; Commentary:
 ;;
@@ -107,7 +105,8 @@ For detail, see `comment-dwim'."
     ("(\\(+\\|-\\)[1-9])" . font-lock-constant-face)
     ))
 
-;; define the major mode
+;;; define the major mode
+;;;###autoload
 (define-derived-mode dynare-mode fundamental-mode
   "dynare mode"
   "dynare is a mode for editing mod files used by dynare."
@@ -149,5 +148,9 @@ For detail, see `comment-dwim'."
   (setq dynare-functions-regexp nil)
   )
 
+;;; mode trigger
+;;;###autoload
+(add-to-list 'auto-mode-alist '("\\.mod$" . dynare-mode))
+
 (provide 'dynare)
 ;;; dynare.el ends here

tests/Makefile.am

@@ -302,6 +302,7 @@ MODFILES = \
 	gradient/fs2000_numgrad_5.mod \
 	filter_step_ahead/fs2000_filter_step_ahead_bayesian.mod \
 	filter_step_ahead/fs2000_filter_step_ahead_ML.mod \
+	filter_step_ahead/trend_cycle_decomposition.mod \
 	loglinear/example4_exp.mod \
 	loglinear/example4_loglinear.mod \
 	loglinear/example4_loglinear_lagged_exogenous.mod \
@@ -762,6 +763,7 @@ EXTRA_DIST = \
 	external_function/extFunNoDerivs.m \
 	external_function/extFunWithFirstAndSecondDerivs.m \
 	expectations/expectation_ss_old_steadystate.m \
+	filter_step_ahead/trend_cycle_decomposition_data.m \
 	steady_state/walsh1_old_ss_steadystate.m \
 	data/test.xlsx \
 	gsa/morris/nk_est_data.m \

tests/filter_step_ahead/trend_cycle_decomposition.mod

@@ -0,0 +1,34 @@
+var y yp z mu;
+
+varexo ez eyp emu;
+
+parameters alpha;
+
+alpha = .889;
+
+model(linear);
+  y  = yp + z;
+  yp = mu + yp(-1) + eyp;
+  mu = mu(-1) + emu;
+  z  = alpha*z(-1) + ez;
+end;
+
+initval;
+y=8.655680;
+z=0;
+yp=8.655680;
+mu=0;
+end;
+
+steady(nocheck);
+
+estimated_params;
+stderr emu     , inv_gamma_pdf,  0.002  , inf;
+stderr eyp     , inv_gamma_pdf,  0.002  , inf;
+stderr ez        , inv_gamma_pdf,   0.06 , inf;
+alpha, normal_pdf, 0.9, 0.1;
+end;
+
+varobs y; 
+
+estimation(datafile=trend_cycle_decomposition_data,nobs=82, mh_replic=2000, mode_compute=4, mh_nblocks=1, mh_jscale=0.3, filtered_vars, smoother, diffuse_filter) yp z; 

tests/filter_step_ahead/trend_cycle_decomposition_data.m

@@ -0,0 +1,164 @@
+y=[8.655680
+8.656123
+8.654401
+8.641229
+8.630567
+8.620637
+8.605288
+8.608660
+8.620357
+8.619683
+8.630202
+8.636380
+8.614226
+8.605106
+8.537263
+8.555079
+8.591242
+8.626651
+8.645062
+8.671851
+8.704031
+8.723127
+8.737705
+8.750433
+8.758621
+8.772876
+8.790050
+8.794088
+8.798091
+8.816624
+8.836301
+8.852974
+8.872307
+8.887314
+8.902466
+8.924601
+8.944691
+8.960662
+8.973784
+8.984161
+9.000123
+9.019129
+9.035000
+9.054575
+9.075913
+9.094302
+9.114171
+9.135554
+9.157538
+9.176307
+9.193220
+9.220434
+9.246620
+9.257001
+9.243476
+9.202275
+9.160007
+9.147657
+9.157229
+9.171294
+9.189745
+9.198273
+9.201042
+9.214450
+9.225220
+9.233702
+9.247974
+9.262110
+9.272097
+9.275883
+9.277719
+9.278357
+9.282812
+9.289376
+9.290788
+9.292530
+9.295274
+9.297219
+9.297555
+9.290289
+9.271941
+9.255655];
+ipc=[23.5993149
+24.9374227
+15.5668177
+12.822998
+8.720567
+5.2893706
+3.4198277
+3.274742
+3.664711
+3.0744047
+2.0693272
+1.704983
+1.3786582
+1.3437517
+2.6167173
+4.7126735
+6.3341827
+7.3020707
+6.5659735
+4.3717219
+3.4512364
+4.7505574
+5.6671666
+5.6241466
+5.5933745
+4.7565382
+3.7032831
+3.6624086
+3.8250059
+3.2723449
+3.3568247
+3.7965188
+3.5401286
+2.7513284
+2.5482882
+2.2870368
+2.0836978
+2.6555146
+3.3263235
+3.3905589
+3.2762151
+2.8684465
+2.1706005
+2.3531933
+2.5136208
+2.3339758
+2.0782414
+1.8397957
+1.7220838
+2.3651145
+3.4129652
+3.8212567
+3.9356373
+3.6723844
+3.11287
+2.9099557
+3.1280929
+2.5837809
+1.5138311
+1.191679
+1.286452
+1.8904928
+2.4851907
+2.7408867
+2.2572768
+1.3839678
+0.813925
+0.88364
+1.2723562
+2.0426676
+2.4161978
+1.5967846
+0.9763814
+1.6120089
+2.0962214
+1.6427865
+1.5391928
+2.0537506
+3.1496297
+4.34142
+4.3532191
+2.2026746];