Updated version following new revision of dsgelikelihood.m

git-svn-id: https://www.dynare.org/svn/dynare/trunk@3318 ac1d8469-bf42-47a9-8791-bf33cf982152
time-shift
ratto 2010-01-08 10:06:21 +00:00
parent e7adf28fbd
commit 16794e4ec9
1 changed files with 268 additions and 259 deletions

View File

@ -1,4 +1,4 @@
function [fval,llik,cost_flag,ys,trend_coeff,info] = DsgeLikelihood_hh(xparam1,gend,data,data_index,number_of_observations,no_more_missing_observations) function [fval,llik,cost_flag,ys,trend_coeff,info] = DsgeLikelihood(xparam1,gend,data,data_index,number_of_observations,no_more_missing_observations)
% function [fval,cost_flag,ys,trend_coeff,info] = DsgeLikelihood(xparam1,gend,data,data_index,number_of_observations,no_more_missing_observations) % function [fval,cost_flag,ys,trend_coeff,info] = DsgeLikelihood(xparam1,gend,data,data_index,number_of_observations,no_more_missing_observations)
% Evaluates the posterior kernel of a dsge model. % Evaluates the posterior kernel of a dsge model.
% %
@ -42,9 +42,9 @@ global bayestopt_ estim_params_ options_ trend_coeff_ M_ oo_
fval = []; fval = [];
ys = []; ys = [];
trend_coeff = []; trend_coeff = [];
llik = NaN;
cost_flag = 1; cost_flag = 1;
nobs = size(options_.varobs,1); nobs = size(options_.varobs,1);
llik=NaN;
%------------------------------------------------------------------------------ %------------------------------------------------------------------------------
% 1. Get the structural parameters & define penalties % 1. Get the structural parameters & define penalties
%------------------------------------------------------------------------------ %------------------------------------------------------------------------------
@ -128,26 +128,26 @@ M_.H = H;
[T,R,SteadyState,info] = dynare_resolve(bayestopt_.restrict_var_list,... [T,R,SteadyState,info] = dynare_resolve(bayestopt_.restrict_var_list,...
bayestopt_.restrict_columns,... bayestopt_.restrict_columns,...
bayestopt_.restrict_aux); bayestopt_.restrict_aux);
if info(1) == 1 | info(1) == 2 | info(1) == 5 if info(1) == 1 || info(1) == 2 || info(1) == 5
fval = bayestopt_.penalty+1; fval = bayestopt_.penalty+1;
cost_flag = 0; cost_flag = 0;
return return
elseif info(1) == 3 | info(1) == 4 | info(1) == 20 elseif info(1) == 3 || info(1) == 4 || info(1)==6 ||info(1) == 19 || info(1) == 20 || info(1) == 21
fval = bayestopt_.penalty+info(2);%^2; % penalty power raised in DR1.m and resol already. GP July'08 fval = bayestopt_.penalty+info(2);
cost_flag = 0; cost_flag = 0;
return return
end end
bayestopt_.mf = bayestopt_.mf1; bayestopt_.mf = bayestopt_.mf1;
if ~options_.noconstant if options_.noconstant
if options_.loglinear == 1 constant = zeros(nobs,1);
else
if options_.loglinear
constant = log(SteadyState(bayestopt_.mfys)); constant = log(SteadyState(bayestopt_.mfys));
else else
constant = SteadyState(bayestopt_.mfys); constant = SteadyState(bayestopt_.mfys);
end end
else
constant = zeros(nobs,1);
end end
if bayestopt_.with_trend == 1 if bayestopt_.with_trend
trend_coeff = zeros(nobs,1); trend_coeff = zeros(nobs,1);
t = options_.trend_coeffs; t = options_.trend_coeffs;
for i=1:length(t) for i=1:length(t)
@ -242,16 +242,6 @@ elseif options_.lik_init == 3 % Diffuse Kalman filter
end end
end end
if kalman_algo == 2 if kalman_algo == 2
no_correlation_flag = 1;
if length(H)==1
H = zeros(nobs,1);
else
if all(all(abs(H-diag(diag(H)))<1e-14))% ie, the covariance matrix is diagonal...
H = diag(H);
else
no_correlation_flag = 1;
end
end
end end
kalman_tol = options_.kalman_tol; kalman_tol = options_.kalman_tol;
riccati_tol = options_.riccati_tol; riccati_tol = options_.riccati_tol;
@ -273,6 +263,16 @@ if (kalman_algo==1)% Multivariate Kalman Filter
end end
end end
if (kalman_algo==2)% Univariate Kalman Filter if (kalman_algo==2)% Univariate Kalman Filter
no_correlation_flag = 1;
if length(H)==1 & H == 0
H = zeros(nobs,1);
else
if all(all(abs(H-diag(diag(H)))<1e-14))% ie, the covariance matrix is diagonal...
H = diag(H);
else
no_correlation_flag = 0;
end
end
if no_correlation_flag if no_correlation_flag
[LIK, lik] = univariate_kalman_filter(T,R,Q,H,Pstar,Y,start,mf,kalman_tol,riccati_tol,data_index,number_of_observations,no_more_missing_observations); [LIK, lik] = univariate_kalman_filter(T,R,Q,H,Pstar,Y,start,mf,kalman_tol,riccati_tol,data_index,number_of_observations,no_more_missing_observations);
else else
@ -281,29 +281,38 @@ if (kalman_algo==2)% Univariate Kalman Filter
end end
if (kalman_algo==3)% Multivariate Diffuse Kalman Filter if (kalman_algo==3)% Multivariate Diffuse Kalman Filter
if no_missing_data_flag if no_missing_data_flag
data1 = data - trend; [LIK, lik] = diffuse_kalman_filter(ST,R1,Q,H,Pinf,Pstar,Y,start,Z,kalman_tol, ...
if any(any(H ~= 0)) riccati_tol);
[LIK, lik] = DiffuseLikelihoodH1_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,start);
else else
[LIK, lik] = DiffuseLikelihood1_Z(ST,Z,R1,Q,Pinf,Pstar,data1,start); [LIK, lik] = missing_observations_diffuse_kalman_filter(ST,R1,Q,H,Pinf, ...
Pstar,Y,start,Z,kalman_tol,riccati_tol,...
data_index,number_of_observations,...
no_more_missing_observations);
end end
if isinf(LIK) if isinf(LIK)
kalman_algo = 4; kalman_algo = 4;
end end
else
error(['The diffuse filter is not yet implemented for models with missing observations'])
end
end end
if (kalman_algo==4)% Univariate Diffuse Kalman Filter if (kalman_algo==4)% Univariate Diffuse Kalman Filter
data1 = data - trend; no_correlation_flag = 1;
if any(any(H ~= 0)) if length(H)==1 & H == 0
if ~estim_params_.ncn H = zeros(nobs,1);
[LIK, lik] = DiffuseLikelihoodH3_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,trend,start);
else else
[LIK, lik] = DiffuseLikelihoodH3corr_Z(ST,Z,R1,Q,H,Pinf,Pstar,data1,trend,start); if all(all(abs(H-diag(diag(H)))<1e-14))% ie, the covariance matrix is diagonal...
H = diag(H);
else
no_correlation_flag = 0;
end end
end
if no_correlation_flag
[LIK, lik] = univariate_diffuse_kalman_filter(ST,R1,Q,H,Pinf,Pstar,Y, ...
start,Z,kalman_tol,riccati_tol,data_index,...
number_of_observations,no_more_missing_observations);
else else
[LIK, lik] = DiffuseLikelihood3_Z(ST,Z,R1,Q,Pinf,Pstar,data1,start); [LIK, lik] = univariate_diffuse_kalman_filter_corr(ST,R1,Q,H,Pinf,Pstar, ...
Y,start,Z,kalman_tol,riccati_tol,...
data_index,number_of_observations,...
no_more_missing_observations);
end end
end end
if imag(LIK) ~= 0 if imag(LIK) ~= 0