2016-06-01 15:51:13 +02:00
function [fval,info,exit_flag,DLIK,Hess,ys,trend_coeff,Model,DynareOptions,BayesInfo,DynareResults] = non_linear_dsge_likelihood ( xparam1,DynareDataset,DatasetInfo,DynareOptions,Model,EstimatedParameters,BayesInfo,BoundsInfo,DynareResults)
2011-12-26 17:44:41 +01:00
2019-12-14 17:23:47 +01:00
% Evaluates the posterior kernel of a dsge model using a non linear filter.
%
% INPUTS
% - xparam1 [double] n× 1 vector, estimated parameters.
% - DynareDataset [struct] Matlab's structure containing the dataset (initialized by dynare, aka dataset_).
% - DatasetInfo [struct] Matlab's structure describing the dataset (initialized by dynare, aka dataset_info).
% - DynareOptions [struct] Matlab's structure describing the options (initialized by dynare, aka options_).
% - Model [struct] Matlab's structure describing the Model (initialized by dynare, aka M_).
% - EstimatedParameters [struct] Matlab's structure describing the estimated_parameters (initialized by dynare, aka estim_params_).
% - BayesInfo [struct] Matlab's structure describing the priors (initialized by dynare,aka bayesopt_).
% - BoundsInfo [struct] Matlab's structure specifying the bounds on the paramater values (initialized by dynare,aka bayesopt_).
% - DynareResults [struct] Matlab's structure gathering the results (initialized by dynare,aka oo_).
%
% OUTPUTS
% - fval [double] scalar, value of the likelihood or posterior kernel.
% - info [integer] 4× 1 vector, informations resolution of the model and evaluation of the likelihood.
% - exit_flag [integer] scalar, equal to 1 (no issues when evaluating the likelihood) or 0 (not able to evaluate the likelihood).
% - DLIK [double] Empty array.
% - Hess [double] Empty array.
% - ys [double] Empty array.
% - trend_coeff [double] Empty array.
% - Model [struct] Updated Model structure described in INPUTS section.
% - DynareOptions [struct] Updated DynareOptions structure described in INPUTS section.
% - BayesInfo [struct] See INPUTS section.
% - DynareResults [struct] Updated DynareResults structure described in INPUTS section.
2011-12-26 17:44:41 +01:00
2022-04-13 13:15:19 +02:00
% Copyright © 2010-2022 Dynare Team
2011-12-26 17:44:41 +01:00
%
% 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
2021-06-09 17:33:48 +02:00
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
2011-12-26 17:44:41 +01:00
2015-10-09 14:18:31 +02:00
% Initialization of the returned arguments.
fval = [ ] ;
ys = [ ] ;
trend_coeff = [ ] ;
exit_flag = 1 ;
2016-06-01 15:51:13 +02:00
DLIK = [ ] ;
Hess = [ ] ;
2015-10-09 14:18:31 +02:00
2017-02-08 12:35:08 +01:00
% Ensure that xparam1 is a column vector.
2021-01-27 15:41:24 +01:00
% (Don't do the transformation if xparam1 is empty, otherwise it would become a
% 0× 1 matrix, which create issues with older MATLABs when comparing with [] in
% check_bounds_and_definiteness_estimation)
if ~ isempty ( xparam1 )
xparam1 = xparam1 ( : ) ;
end
2017-02-08 12:35:08 +01:00
2015-10-09 14:18:31 +02:00
% Issue an error if loglinear option is used.
if DynareOptions . loglinear
error ( ' non_linear_dsge_likelihood: It is not possible to use a non linear filter with the option loglinear!' )
end
%------------------------------------------------------------------------------
% 1. Get the structural parameters & define penalties
%------------------------------------------------------------------------------
2020-07-17 13:18:39 +02:00
Model = set_all_parameters ( xparam1 , EstimatedParameters , Model ) ;
[ fval , info , exit_flag , Q , H ] = check_bounds_and_definiteness_estimation ( xparam1 , Model , EstimatedParameters , BoundsInfo ) ;
2020-06-22 18:45:13 +02:00
if info ( 1 )
2015-10-09 14:18:31 +02:00
return
end
%------------------------------------------------------------------------------
% 2. call model setup & reduction program
%------------------------------------------------------------------------------
% Linearize the model around the deterministic sdteadystate and extract the matrices of the state equation (T and R).
2021-06-22 14:28:46 +02:00
[ dr , info , Model , DynareResults ] = resol ( 0 , Model , DynareOptions , DynareResults ) ;
2015-10-09 14:18:31 +02:00
2016-06-01 18:22:51 +02:00
if info ( 1 )
2016-07-18 11:11:43 +02:00
if info ( 1 ) == 3 || info ( 1 ) == 4 || info ( 1 ) == 5 || info ( 1 ) == 6 || info ( 1 ) == 19 || ...
2017-05-16 15:10:20 +02:00
info ( 1 ) == 20 || info ( 1 ) == 21 || info ( 1 ) == 23 || info ( 1 ) == 26 || ...
info ( 1 ) == 81 || info ( 1 ) == 84 || info ( 1 ) == 85
2016-06-01 18:22:51 +02:00
%meaningful second entry of output that can be used
fval = Inf ;
info ( 4 ) = info ( 2 ) ;
exit_flag = 0 ;
return
else
fval = Inf ;
info ( 4 ) = 0.1 ;
exit_flag = 0 ;
return
end
2015-10-09 14:18:31 +02:00
end
% Define a vector of indices for the observed variables. Is this really usefull?...
BayesInfo . mf = BayesInfo . mf1 ;
% Get needed informations for kalman filter routines.
start = DynareOptions . presample + 1 ;
2019-12-14 17:24:34 +01:00
Y = transpose ( DynareDataset . data ) ;
2015-10-09 14:18:31 +02:00
%------------------------------------------------------------------------------
% 3. Initial condition of the Kalman filter
%------------------------------------------------------------------------------
2019-12-14 17:30:30 +01:00
mf0 = BayesInfo . mf0 ;
mf1 = BayesInfo . mf1 ;
2019-12-20 11:26:35 +01:00
restrict_variables_idx = dr . restrict_var_list ;
state_variables_idx = restrict_variables_idx ( mf0 ) ;
2019-12-14 17:30:30 +01:00
number_of_state_variables = length ( mf0 ) ;
2015-10-09 14:18:31 +02:00
ReducedForm . steadystate = dr . ys ( dr . order_var ( restrict_variables_idx ) ) ;
ReducedForm . constant = ReducedForm . steadystate + . 5 * dr . ghs2 ( restrict_variables_idx ) ;
ReducedForm . state_variables_steady_state = dr . ys ( dr . order_var ( state_variables_idx ) ) ;
ReducedForm . Q = Q ;
ReducedForm . H = H ;
ReducedForm . mf0 = mf0 ;
ReducedForm . mf1 = mf1 ;
2022-08-30 14:16:14 +02:00
if DynareOptions . order > 3
2019-12-20 11:26:35 +01:00
ReducedForm . use_k_order_solver = true ;
ReducedForm . dr = dr ;
2022-02-04 11:31:09 +01:00
ReducedForm . udr = folded_to_unfolded_dr ( dr , Model , DynareOptions ) ;
2022-08-30 14:16:14 +02:00
if pruning
error ( ' Pruning is not available for orders > 3' ) ;
end
2019-12-20 11:26:35 +01:00
else
ReducedForm . use_k_order_solver = false ;
ReducedForm . ghx = dr . ghx ( restrict_variables_idx , : ) ;
ReducedForm . ghu = dr . ghu ( restrict_variables_idx , : ) ;
ReducedForm . ghxx = dr . ghxx ( restrict_variables_idx , : ) ;
ReducedForm . ghuu = dr . ghuu ( restrict_variables_idx , : ) ;
ReducedForm . ghxu = dr . ghxu ( restrict_variables_idx , : ) ;
2022-08-30 14:16:14 +02:00
if DynareOptions . order == 3
ReducedForm . ghxxx = dr . ghxxx ( restrict_variables_idx , : ) ;
ReducedForm . ghuuu = dr . ghuuu ( restrict_variables_idx , : ) ;
ReducedForm . ghxxu = dr . ghxxu ( restrict_variables_idx , : ) ;
ReducedForm . ghxuu = dr . ghxuu ( restrict_variables_idx , : ) ;
ReducedForm . ghxss = dr . ghxss ( restrict_variables_idx , : ) ;
ReducedForm . ghuss = dr . ghuss ( restrict_variables_idx , : ) ;
end
2019-12-20 11:26:35 +01:00
end
2015-10-09 14:18:31 +02:00
% Set initial condition.
switch DynareOptions . particle . initialization
case 1 % Initial state vector covariance is the ergodic variance associated to the first order Taylor-approximation of the model.
StateVectorMean = ReducedForm . constant ( mf0 ) ;
2020-10-08 12:15:26 +02:00
[ A , B ] = kalman_transition_matrix ( dr , dr . restrict_var_list , dr . restrict_columns , Model . exo_nbr ) ;
StateVectorVariance = lyapunov_symm ( A , B * Q * B ' , DynareOptions . lyapunov_fixed_point_tol , ...
2019-12-20 11:26:35 +01:00
DynareOptions . qz_criterium , DynareOptions . lyapunov_complex_threshold , [ ] , DynareOptions . debug ) ;
2020-10-08 12:15:26 +02:00
StateVectorVariance = StateVectorVariance ( mf0 , mf0 ) ;
2015-10-09 14:18:31 +02:00
case 2 % Initial state vector covariance is a monte-carlo based estimate of the ergodic variance (consistent with a k-order Taylor-approximation of the model).
StateVectorMean = ReducedForm . constant ( mf0 ) ;
old_DynareOptionsperiods = DynareOptions . periods ;
DynareOptions . periods = 5000 ;
2021-01-21 19:26:53 +01:00
old_DynareOptionspruning = DynareOptions . pruning ;
DynareOptions . pruning = DynareOptions . particle . pruning ;
2015-10-09 14:18:31 +02:00
y_ = simult ( DynareResults . steady_state , dr , Model , DynareOptions , DynareResults ) ;
2020-10-08 12:15:26 +02:00
y_ = y_ ( dr . order_var ( state_variables_idx ) , 2001 : 5000 ) ; %state_variables_idx is in dr-order while simult_ is in declaration order
2021-05-27 21:04:35 +02:00
if any ( any ( isnan ( y_ ) ) ) || any ( any ( isinf ( y_ ) ) ) && ~ DynareOptions . pruning
fval = Inf ;
info ( 1 ) = 202 ;
info ( 4 ) = 0.1 ;
exit_flag = 0 ;
return ;
end
StateVectorVariance = cov ( y_ ' ) ;
2015-10-09 14:18:31 +02:00
DynareOptions . periods = old_DynareOptionsperiods ;
2021-01-21 19:26:53 +01:00
DynareOptions . pruning = old_DynareOptionspruning ;
2015-10-09 14:18:31 +02:00
clear ( ' old_DynareOptionsperiods' , ' y_' ) ;
2017-04-02 18:37:45 +02:00
case 3 % Initial state vector covariance is a diagonal matrix (to be used
% if model has stochastic trends).
2015-10-09 14:18:31 +02:00
StateVectorMean = ReducedForm . constant ( mf0 ) ;
StateVectorVariance = DynareOptions . particle . initial_state_prior_std * eye ( number_of_state_variables ) ;
otherwise
error ( ' Unknown initialization option!' )
end
ReducedForm . StateVectorMean = StateVectorMean ;
ReducedForm . StateVectorVariance = StateVectorVariance ;
2020-12-08 22:49:41 +01:00
[ ~ , flag ] = chol ( ReducedForm . StateVectorVariance ) ; %reduced_rank_cholesky(ReducedForm.StateVectorVariance)';
if flag
fval = Inf ;
info ( 1 ) = 201 ;
info ( 4 ) = 0.1 ;
exit_flag = 0 ;
return ;
end
2015-10-09 14:18:31 +02:00
%------------------------------------------------------------------------------
% 4. Likelihood evaluation
%------------------------------------------------------------------------------
DynareOptions . warning_for_steadystate = 0 ;
[ s1 , s2 ] = get_dynare_random_generator_state ( ) ;
2019-12-20 11:26:35 +01:00
LIK = feval ( DynareOptions . particle . algorithm , ReducedForm , Y , start , DynareOptions . particle , DynareOptions . threads , DynareOptions , Model ) ;
2015-10-09 14:18:31 +02:00
set_dynare_random_generator_state ( s1 , s2 ) ;
if imag ( LIK )
2020-12-08 22:25:47 +01:00
fval = Inf ;
2016-06-01 18:22:51 +02:00
info ( 1 ) = 46 ;
info ( 4 ) = 0.1 ;
exit_flag = 0 ;
2020-12-08 22:25:47 +01:00
return
2015-10-09 14:18:31 +02:00
elseif isnan ( LIK )
2020-12-08 22:25:47 +01:00
fval = Inf ;
2016-06-01 18:22:51 +02:00
info ( 1 ) = 45 ;
info ( 4 ) = 0.1 ;
exit_flag = 0 ;
2020-12-08 22:25:47 +01:00
return
2015-10-09 14:18:31 +02:00
else
likelihood = LIK ;
end
DynareOptions . warning_for_steadystate = 1 ;
% ------------------------------------------------------------------------------
% Adds prior if necessary
% ------------------------------------------------------------------------------
lnprior = priordens ( xparam1 ( : ) , BayesInfo . pshape , BayesInfo . p6 , BayesInfo . p7 , BayesInfo . p3 , BayesInfo . p4 ) ;
2019-12-14 17:24:34 +01:00
fval = ( likelihood - lnprior ) ;
2015-10-09 14:18:31 +02:00
if isnan ( fval )
2016-06-01 18:22:51 +02:00
fval = Inf ;
info ( 1 ) = 47 ;
info ( 4 ) = 0.1 ;
2015-10-09 14:18:31 +02:00
exit_flag = 0 ;
return
end
2019-12-14 17:24:34 +01:00
if ~ isreal ( fval )
2016-06-01 18:22:51 +02:00
fval = Inf ;
info ( 1 ) = 48 ;
info ( 4 ) = 0.1 ;
2015-10-09 14:18:31 +02:00
exit_flag = 0 ;
return
end
2016-06-14 17:28:52 +02:00
2019-12-14 17:24:34 +01:00
if isinf ( LIK )
2016-06-14 17:28:52 +02:00
fval = Inf ;
info ( 1 ) = 50 ;
info ( 4 ) = 0.1 ;
exit_flag = 0 ;
return
end