0001 function dynare_estimation_1(var_list_,dname)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032 global M_ options_ oo_ estim_params_ bayestopt_ dataset_
0033
0034
0035 if options_.order > 1
0036 if options_.particle.status && options_.order==2
0037 disp(' ')
0038 disp('Estimation using a non linear filter!')
0039 disp(' ')
0040 elseif options_.particle.status && options_.order>2
0041 error(['Non linear filter are not implemented with order ' int2str(options_.order) ' approximation of the model!'])
0042 elseif ~options_.particle.status && options_.order==2
0043 disp('If you want to estimate the model with a second order approximation using a non linear filter, set options_.particle.status=1;')
0044 disp('I set order=1!')
0045 options_.order=1;
0046 else
0047 error(['Cannot estimate a model with an order ' int2str(options_.order) ' approximation!'])
0048 end
0049 end
0050
0051 if ~options_.dsge_var
0052 if options_.particle.status
0053 objective_function = str2func('non_linear_dsge_likelihood');
0054 else
0055 objective_function = str2func('dsge_likelihood');
0056 end
0057 else
0058 objective_function = str2func('DsgeVarLikelihood');
0059 end
0060
0061 [dataset_,xparam1, M_, options_, oo_, estim_params_,bayestopt_] = dynare_estimation_init(var_list_, dname, [], M_, options_, oo_, estim_params_, bayestopt_);
0062
0063 data = dataset_.data;
0064 rawdata = dataset_.rawdata;
0065 data_index = dataset_.missing.aindex;
0066 missing_value = dataset_.missing.state;
0067
0068
0069 gend = options_.nobs;
0070
0071 n_varobs = size(options_.varobs,1);
0072
0073 nvx = estim_params_.nvx;
0074 nvn = estim_params_.nvn;
0075 ncx = estim_params_.ncx;
0076 ncn = estim_params_.ncn;
0077 np = estim_params_.np ;
0078 nx = nvx+nvn+ncx+ncn+np;
0079
0080 pnames = [' ';'beta ';'gamm ';'norm ';'invg ';'unif ';'invg2'];
0081
0082 lb = bayestopt_.lb;
0083 ub = bayestopt_.ub;
0084
0085 dr = oo_.dr;
0086
0087
0088 if ~isempty(options_.mode_file) && ~options_.mh_posterior_mode_estimation
0089 load(options_.mode_file);
0090 end
0091
0092 if ~isempty(estim_params_)
0093 set_parameters(xparam1);
0094 end
0095
0096
0097 if options_.dsge_var
0098 if dataset_.missing.state
0099 error('I cannot estimate a DSGE-VAR model with missing observations!')
0100 end
0101 if options_.noconstant
0102 evalin('base',...
0103 ['[mYY,mXY,mYX,mXX,Ydata,Xdata] = ' ...
0104 'var_sample_moments(options_.first_obs,' ...
0105 'options_.first_obs+options_.nobs-1,options_.dsge_varlag,-1,' ...
0106 'options_.datafile, options_.varobs,options_.xls_sheet,options_.xls_range);'])
0107 else
0108 evalin('base',['[mYY,mXY,mYX,mXX,Ydata,Xdata] = ' ...
0109 'var_sample_moments(options_.first_obs,' ...
0110 'options_.first_obs+options_.nobs-1,options_.dsge_varlag,0,' ...
0111 'options_.datafile, options_.varobs,options_.xls_sheet,options_.xls_range);'])
0112 end
0113 end
0114
0115
0116 oo_ = initial_estimation_checks(objective_function,xparam1,dataset_,M_,estim_params_,options_,bayestopt_,oo_);
0117
0118 if isequal(options_.mode_compute,0) && isempty(options_.mode_file) && options_.mh_posterior_mode_estimation==0
0119 if options_.smoother == 1
0120 [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp] = DsgeSmoother(xparam1,gend,data,data_index,missing_value);
0121 oo_.Smoother.SteadyState = ys;
0122 oo_.Smoother.TrendCoeffs = trend_coeff;
0123 if options_.filter_covariance
0124 oo_.Smoother.variance = P;
0125 end
0126 i_endo = bayestopt_.smoother_saved_var_list;
0127 if options_.nk ~= 0
0128 oo_.FilteredVariablesKStepAhead = ...
0129 aK(options_.filter_step_ahead,i_endo,:);
0130 if ~isempty(PK)
0131 oo_.FilteredVariablesKStepAheadVariances = ...
0132 PK(options_.filter_step_ahead,i_endo,i_endo,:);
0133 end
0134 if ~isempty(decomp)
0135 oo_.FilteredVariablesShockDecomposition = ...
0136 decomp(options_.filter_step_ahead,i_endo,:,:);
0137 end
0138 end
0139 for i=bayestopt_.smoother_saved_var_list'
0140 i1 = dr.order_var(bayestopt_.smoother_var_list(i));
0141 eval(['oo_.SmoothedVariables.' deblank(M_.endo_names(i1,:)) ...
0142 ' = atT(i,:)'';']);
0143 if options_.nk > 0
0144 eval(['oo_.FilteredVariables.' deblank(M_.endo_names(i1,:)) ...
0145 ' = squeeze(aK(1,i,:));']);
0146 end
0147 eval(['oo_.UpdatedVariables.' deblank(M_.endo_names(i1,:)) ' = updated_variables(i,:)'';']);
0148 end
0149 for i=1:M_.exo_nbr
0150 eval(['oo_.SmoothedShocks.' deblank(M_.exo_names(i,:)) ' = innov(i,:)'';']);
0151 end
0152 end
0153 return
0154 end
0155
0156 if isequal(options_.mode_compute,6)
0157
0158 delete([M_.fname '_optimal_mh_scale_parameter.mat'])
0159 end
0160
0161
0162
0163 if ~isequal(options_.mode_compute,0) && ~options_.mh_posterior_mode_estimation
0164 switch options_.mode_compute
0165 case 1
0166 optim_options = optimset('display','iter','LargeScale','off', ...
0167 'MaxFunEvals',100000,'TolFun',1e-8,'TolX',1e-6);
0168 if isfield(options_,'optim_opt')
0169 eval(['optim_options = optimset(optim_options,' options_.optim_opt ');']);
0170 end
0171 if options_.analytic_derivation,
0172 optim_options = optimset(optim_options,'GradObj','on');
0173 end
0174 [xparam1,fval,exitflag,output,lamdba,grad,hessian_fmincon] = ...
0175 fmincon(objective_function,xparam1,[],[],[],[],lb,ub,[],optim_options,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0176 case 2
0177 error('ESTIMATION: mode_compute=2 option (Lester Ingber''s Adaptive Simulated Annealing) is no longer available')
0178 case 3
0179 optim_options = optimset('display','iter','MaxFunEvals',100000,'TolFun',1e-8,'TolX',1e-6);
0180 if isfield(options_,'optim_opt')
0181 eval(['optim_options = optimset(optim_options,' options_.optim_opt ');']);
0182 end
0183 if options_.analytic_derivation,
0184 optim_options = optimset(optim_options,'GradObj','on');
0185 end
0186 [xparam1,fval,exitflag] = fminunc(objective_function,xparam1,optim_options,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0187 case 4
0188 H0 = 1e-4*eye(nx);
0189 crit = 1e-7;
0190 nit = 1000;
0191 verbose = 2;
0192 if options_.analytic_derivation,
0193 analytic_grad=1;
0194 else
0195 analytic_grad=[];
0196 end
0197
0198 [fval,xparam1,grad,hessian_csminwel,itct,fcount,retcodehat] = ...
0199 csminwel1(objective_function,xparam1,H0,analytic_grad,crit,nit,options_.gradient_method,options_.gradient_epsilon,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0200 disp(sprintf('Objective function at mode: %f',fval))
0201 case 5
0202 if isfield(options_,'hess')
0203 flag = options_.hess;
0204 else
0205 flag = 1;
0206 end
0207 if ~exist('igg','var'),
0208 hh=[];
0209 gg=[];
0210 igg=[];
0211 end
0212 if isfield(options_,'ftol')
0213 crit = options_.ftol;
0214 else
0215 crit = 1.e-5;
0216 end
0217 if isfield(options_,'nit')
0218 nit = options_.nit;
0219 else
0220 nit=1000;
0221 end
0222 [xparam1,hh,gg,fval,invhess] = newrat(objective_function,xparam1,hh,gg,igg,crit,nit,flag,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0223 parameter_names = bayestopt_.name;
0224 save([M_.fname '_mode.mat'],'xparam1','hh','gg','fval','invhess','parameter_names');
0225 case 6
0226 fval = feval(objective_function,xparam1,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0227 OldMode = fval;
0228 if ~exist('MeanPar','var')
0229 MeanPar = xparam1;
0230 end
0231 if exist('hh','var')
0232 CovJump = inv(hh);
0233 else
0234
0235 varinit = 'prior';
0236 if strcmpi(varinit,'prior')
0237 stdev = bayestopt_.p2;
0238 indx = find(isinf(stdev));
0239 stdev(indx) = ones(length(indx),1)*sqrt(10);
0240 vars = stdev.^2;
0241 CovJump = diag(vars);
0242 elseif strcmpi(varinit,'eye')
0243 vars = ones(length(bayestopt_.p2),1)*0.1;
0244 CovJump = diag(vars);
0245 else
0246 disp('gmhmaxlik :: Error!')
0247 return
0248 end
0249 end
0250 OldPostVar = CovJump;
0251 Scale = options_.mh_jscale;
0252 for i=1:options_.Opt6Iter
0253 if i == 1
0254 if options_.Opt6Iter > 1
0255 flag = '';
0256 else
0257 flag = 'LastCall';
0258 end
0259 [xparam1,PostVar,Scale,PostMean] = ...
0260 gmhmaxlik(objective_function,xparam1,[lb ub],options_.Opt6Numb,Scale,flag,MeanPar,CovJump,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0261 fval = feval(objective_function,xparam1,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0262 options_.mh_jscale = Scale;
0263 mouvement = max(max(abs(PostVar-OldPostVar)));
0264 disp(' ')
0265 disp('========================================================== ')
0266 disp([' Change in the covariance matrix = ' num2str(mouvement) '.'])
0267 disp([' Mode improvement = ' num2str(abs(OldMode-fval))])
0268 disp([' New value of jscale = ' num2str(Scale)])
0269 disp('========================================================== ')
0270 OldMode = fval;
0271 else
0272 OldPostVar = PostVar;
0273 if i<options_.Opt6Iter
0274 flag = '';
0275 else
0276 flag = 'LastCall';
0277 end
0278 [xparam1,PostVar,Scale,PostMean] = ...
0279 gmhmaxlik(objective_function,xparam1,[lb ub],...
0280 options_.Opt6Numb,Scale,flag,PostMean,PostVar,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0281 fval = feval(objective_function,xparam1,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0282 options_.mh_jscale = Scale;
0283 mouvement = max(max(abs(PostVar-OldPostVar)));
0284 fval = dsge_likelihood(xparam1,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0285 disp(['Change in the covariance matrix = ' num2str(mouvement) '.'])
0286 disp(['Mode improvement = ' num2str(abs(OldMode-fval))])
0287 OldMode = fval;
0288 end
0289 hh = inv(PostVar);
0290 save([M_.fname '_mode.mat'],'xparam1','hh');
0291 save([M_.fname '_optimal_mh_scale_parameter.mat'],'Scale');
0292 bayestopt_.jscale = ones(length(xparam1),1)*Scale;
0293 end
0294 case 7
0295
0296 optim_options = optimset('display','iter','MaxFunEvals',1000000,'MaxIter',6000,'TolFun',1e-8,'TolX',1e-6);
0297 if isfield(options_,'optim_opt')
0298 eval(['optim_options = optimset(optim_options,' options_.optim_opt ');']);
0299 end
0300 [xparam1,fval,exitflag] = fminsearch(objective_function,xparam1,optim_options,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0301 case 8
0302
0303 [xparam1,fval,exitflag] = simplex_optimization_routine(objective_function,xparam1,options_.simplex,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0304 case 9
0305 H0 = 1e-4*ones(nx,1);
0306 warning('off','CMAES:NonfinitenessRange');
0307 warning('off','CMAES:InitialSigma');
0308 [x, fval, COUNTEVAL, STOPFLAG, OUT, BESTEVER] = cmaes(func2str(objective_function),xparam1,H0,options_.cmaes,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0309 xparam1=BESTEVER.x;
0310 disp(sprintf('\n Objective function at mode: %f',fval))
0311 case 101
0312 myoptions=soptions;
0313 myoptions(2)=1e-6;
0314 myoptions(3)=1e-6;
0315 myoptions(5)= 1.0;
0316 [xparam1,fval]=solvopt(xparam1,objective_function,[],myoptions,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0317 case 102
0318
0319
0320
0321 LB = xparam1 - 1;
0322 UB = xparam1 + 1;
0323 neps=10;
0324
0325 maxy=0;
0326 epsilon=1.0e-9;
0327 rt_=.10;
0328 t=15.0;
0329 ns=10;
0330 nt=10;
0331 maxevl=100000000;
0332 idisp =1;
0333 npar=length(xparam1);
0334
0335 disp(['size of param',num2str(length(xparam1))])
0336 c=.1*ones(npar,1);
0337
0338
0339 vm=1*ones(npar,1);
0340 disp(['number of parameters= ' num2str(npar) 'max= ' num2str(maxy) 't= ' num2str(t)]);
0341 disp(['rt_= ' num2str(rt_) 'eps= ' num2str(epsilon) 'ns= ' num2str(ns)]);
0342 disp(['nt= ' num2str(nt) 'neps= ' num2str(neps) 'maxevl= ' num2str(maxevl)]);
0343
0344 disp ' ';
0345 disp ' ';
0346 disp(['starting values(x) ' num2str(xparam1')]);
0347 disp(['initial step length(vm) ' num2str(vm')]);
0348 disp(['lower bound(lb)', 'initial conditions', 'upper bound(ub)' ]);
0349 disp([LB xparam1 UB]);
0350 disp(['c vector ' num2str(c')]);
0351
0352 [xparam1, fval, nacc, nfcnev, nobds, ier, t, vm] = sa(objective_function,xparam1,maxy,rt_,epsilon,ns,nt ...
0353 ,neps,maxevl,LB,UB,c,idisp ,t,vm,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0354 case 'prior'
0355 hh = diag(bayestopt_.p2.^2);
0356 otherwise
0357 if ischar(options_.mode_compute)
0358 [xparam1, fval, retcode ] = feval(options_.mode_compute,objective_function,xparam1,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0359 else
0360 error(['dynare_estimation:: mode_compute = ' int2str(options_.mode_compute) ' option is unknown!'])
0361 end
0362 end
0363 if ~isequal(options_.mode_compute,6) && ~isequal(options_.mode_compute,'prior')
0364 if options_.cova_compute == 1
0365 if options_.analytic_derivation && strcmp(func2str(objective_function),'dsge_likelihood'),
0366 ana_deriv = options_.analytic_derivation;
0367 options_.analytic_derivation = 2;
0368 [junk1, junk2, hh] = feval(objective_function,xparam1, ...
0369 dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0370 options_.analytic_derivation = ana_deriv;
0371
0372 else
0373 hh = reshape(hessian(objective_function,xparam1, ...
0374 options_.gstep,dataset_,options_,M_,estim_params_,bayestopt_,oo_),nx,nx);
0375 end
0376 end
0377 end
0378 parameter_names = bayestopt_.name;
0379 if options_.cova_compute
0380 save([M_.fname '_mode.mat'],'xparam1','hh','parameter_names');
0381 else
0382 save([M_.fname '_mode.mat'],'xparam1','parameter_names');
0383 end
0384 end
0385
0386 if options_.cova_compute == 0
0387 hh = [];
0388 end
0389
0390 if ~options_.mh_posterior_mode_estimation && options_.cova_compute
0391 try
0392 chol(hh);
0393 catch
0394 disp(' ')
0395 disp('POSTERIOR KERNEL OPTIMIZATION PROBLEM!')
0396 disp(' (minus) the hessian matrix at the "mode" is not positive definite!')
0397 disp('=> posterior variance of the estimated parameters are not positive.')
0398 disp('You should try to change the initial values of the parameters using')
0399 disp('the estimated_params_init block, or use another optimization routine.')
0400 warning('The results below are most likely wrong!');
0401 end
0402 end
0403
0404 if options_.mode_check == 1 && ~options_.mh_posterior_mode_estimation
0405 mode_check(objective_function,xparam1,hh,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0406 end
0407
0408 oo_.posterior.optimization.mode = xparam1;
0409 oo_.posterior.optimization.variance = [];
0410 if ~options_.mh_posterior_mode_estimation
0411 if options_.cova_compute
0412 invhess = inv(hh);
0413 stdh = sqrt(diag(invhess));
0414 oo_.posterior.optimization.variance = invhess;
0415 end
0416 else
0417 variances = bayestopt_.p2.*bayestopt_.p2;
0418 idInf = isinf(variances);
0419 variances(idInf) = 1;
0420 invhess = options_.mh_posterior_mode_estimation*diag(variances);
0421 xparam1 = bayestopt_.p5;
0422 idNaN = isnan(xparam1);
0423 xparam1(idNaN) = bayestopt_.p1(idNaN);
0424 xparam1 = transpose(xparam1);
0425 end
0426
0427 if ~options_.cova_compute
0428 stdh = NaN(length(xparam1),1);
0429 end
0430
0431 if any(bayestopt_.pshape > 0) && ~options_.mh_posterior_mode_estimation
0432 disp(' ')
0433 disp('RESULTS FROM POSTERIOR MAXIMIZATION')
0434 tstath = zeros(nx,1);
0435 for i = 1:nx
0436 tstath(i) = abs(xparam1(i))/stdh(i);
0437 end
0438
0439 header_width = row_header_width(M_,estim_params_,bayestopt_);
0440
0441 tit1 = sprintf('%-*s %7s %8s %7s %6s %4s %6s\n',header_width-2,' ','prior mean', ...
0442 'mode','s.d.','t-stat','prior','pstdev');
0443 if np
0444 ip = nvx+nvn+ncx+ncn+1;
0445 disp('parameters')
0446 disp(tit1)
0447 for i=1:np
0448 name = bayestopt_.name{ip};
0449 disp(sprintf('%-*s %7.3f %8.4f %7.4f %7.4f %4s %6.4f', ...
0450 header_width,name, ...
0451 bayestopt_.p1(ip),xparam1(ip),stdh(ip),tstath(ip), ...
0452 pnames(bayestopt_.pshape(ip)+1,:), ...
0453 bayestopt_.p2(ip)));
0454 eval(['oo_.posterior_mode.parameters.' name ' = xparam1(ip);']);
0455 eval(['oo_.posterior_std.parameters.' name ' = stdh(ip);']);
0456 ip = ip+1;
0457 end
0458 end
0459 if nvx
0460 ip = 1;
0461 disp('standard deviation of shocks')
0462 disp(tit1)
0463 for i=1:nvx
0464 k = estim_params_.var_exo(i,1);
0465 name = deblank(M_.exo_names(k,:));
0466 disp(sprintf('%-*s %7.3f %8.4f %7.4f %7.4f %4s %6.4f', ...
0467 header_width,name,bayestopt_.p1(ip),xparam1(ip), ...
0468 stdh(ip),tstath(ip),pnames(bayestopt_.pshape(ip)+1,:), ...
0469 bayestopt_.p2(ip)));
0470 M_.Sigma_e(k,k) = xparam1(ip)*xparam1(ip);
0471 eval(['oo_.posterior_mode.shocks_std.' name ' = xparam1(ip);']);
0472 eval(['oo_.posterior_std.shocks_std.' name ' = stdh(ip);']);
0473 ip = ip+1;
0474 end
0475 end
0476 if nvn
0477 disp('standard deviation of measurement errors')
0478 disp(tit1)
0479 ip = nvx+1;
0480 for i=1:nvn
0481 name = deblank(options_.varobs(estim_params_.var_endo(i,1),:));
0482 disp(sprintf('%-*s %7.3f %8.4f %7.4f %7.4f %4s %6.4f', ...
0483 header_width,name,bayestopt_.p1(ip), ...
0484 xparam1(ip),stdh(ip),tstath(ip), ...
0485 pnames(bayestopt_.pshape(ip)+1,:), ...
0486 bayestopt_.p2(ip)));
0487 eval(['oo_.posterior_mode.measurement_errors_std.' name ' = xparam1(ip);']);
0488 eval(['oo_.posterior_std.measurement_errors_std.' name ' = stdh(ip);']);
0489 ip = ip+1;
0490 end
0491 end
0492 if ncx
0493 disp('correlation of shocks')
0494 disp(tit1)
0495 ip = nvx+nvn+1;
0496 for i=1:ncx
0497 k1 = estim_params_.corrx(i,1);
0498 k2 = estim_params_.corrx(i,2);
0499 name = [deblank(M_.exo_names(k1,:)) ',' deblank(M_.exo_names(k2,:))];
0500 NAME = [deblank(M_.exo_names(k1,:)) '_' deblank(M_.exo_names(k2,:))];
0501 disp(sprintf('%-*s %7.3f %8.4f %7.4f %7.4f %4s %6.4f', ...
0502 header_width,name,bayestopt_.p1(ip),xparam1(ip),stdh(ip),tstath(ip), ...
0503 pnames(bayestopt_.pshape(ip)+1,:), bayestopt_.p2(ip)));
0504 M_.Sigma_e(k1,k2) = xparam1(ip)*sqrt(M_.Sigma_e(k1,k1)*M_.Sigma_e(k2,k2));
0505 M_.Sigma_e(k2,k1) = M_.Sigma_e(k1,k2);
0506 eval(['oo_.posterior_mode.shocks_corr.' NAME ' = xparam1(ip);']);
0507 eval(['oo_.posterior_std.shocks_corr.' NAME ' = stdh(ip);']);
0508 ip = ip+1;
0509 end
0510 end
0511 if ncn
0512 disp('correlation of measurement errors')
0513 disp(tit1)
0514 ip = nvx+nvn+ncx+1;
0515 for i=1:ncn
0516 k1 = estim_params_.corrn(i,1);
0517 k2 = estim_params_.corrn(i,2);
0518 name = [deblank(M_.endo_names(k1,:)) ',' deblank(M_.endo_names(k2,:))];
0519 NAME = [deblank(M_.endo_names(k1,:)) '_' deblank(M_.endo_names(k2,:))];
0520 disp(sprintf('%-*s %7.3f %8.4f %7.4f %7.4f %4s %6.4f', ...
0521 header_width,name,bayestopt_.p1(ip),xparam1(ip),stdh(ip),tstath(ip), ...
0522 pnames(bayestopt_.pshape(ip)+1,:), bayestopt_.p2(ip)));
0523 eval(['oo_.posterior_mode.measurement_errors_corr.' NAME ' = xparam1(ip);']);
0524 eval(['oo_.posterior_std.measurement_errors_corr.' NAME ' = stdh(ip);']);
0525 ip = ip+1;
0526 end
0527 end
0528
0529 if options_.cova_compute
0530 estim_params_nbr = size(xparam1,1);
0531 scale_factor = -sum(log10(diag(invhess)));
0532 log_det_invhess = -estim_params_nbr*log(scale_factor)+log(det(scale_factor*invhess));
0533 likelihood = feval(objective_function,xparam1,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0534 oo_.MarginalDensity.LaplaceApproximation = .5*estim_params_nbr*log(2*pi) + .5*log_det_invhess - likelihood;
0535 disp(' ')
0536 disp(sprintf('Log data density [Laplace approximation] is %f.',oo_.MarginalDensity.LaplaceApproximation))
0537 disp(' ')
0538 end
0539 elseif ~any(bayestopt_.pshape > 0) && ~options_.mh_posterior_mode_estimation
0540 disp(' ')
0541 disp('RESULTS FROM MAXIMUM LIKELIHOOD')
0542 tstath = zeros(nx,1);
0543 for i = 1:nx
0544 tstath(i) = abs(xparam1(i))/stdh(i);
0545 end
0546 header_width = row_header_width(M_,estim_params_,bayestopt_);
0547 tit1 = sprintf('%-*s %10s %7s %6s\n',header_width-2,' ','Estimate','s.d.','t-stat');
0548 if np
0549 ip = nvx+nvn+ncx+ncn+1;
0550 disp('parameters')
0551 disp(tit1)
0552 for i=1:np
0553 name = bayestopt_.name{ip};
0554 disp(sprintf('%-*s %8.4f %7.4f %7.4f', ...
0555 header_width,name,xparam1(ip),stdh(ip),tstath(ip)));
0556 eval(['oo_.mle_mode.parameters.' name ' = xparam1(ip);']);
0557 eval(['oo_.mle_std.parameters.' name ' = stdh(ip);']);
0558 ip = ip+1;
0559 end
0560 end
0561 if nvx
0562 ip = 1;
0563 disp('standard deviation of shocks')
0564 disp(tit1)
0565 for i=1:nvx
0566 k = estim_params_.var_exo(i,1);
0567 name = deblank(M_.exo_names(k,:));
0568 disp(sprintf('%-*s %8.4f %7.4f %7.4f',header_width,name,xparam1(ip),stdh(ip),tstath(ip)));
0569 M_.Sigma_e(k,k) = xparam1(ip)*xparam1(ip);
0570 eval(['oo_.mle_mode.shocks_std.' name ' = xparam1(ip);']);
0571 eval(['oo_.mle_std.shocks_std.' name ' = stdh(ip);']);
0572 ip = ip+1;
0573 end
0574 end
0575 if nvn
0576 disp('standard deviation of measurement errors')
0577 disp(tit1)
0578 ip = nvx+1;
0579 for i=1:nvn
0580 name = deblank(options_.varobs(estim_params_.var_endo(i,1),:));
0581 disp(sprintf('%-*s %8.4f %7.4f %7.4f',header_width,name,xparam1(ip),stdh(ip),tstath(ip)))
0582 eval(['oo_.mle_mode.measurement_errors_std.' name ' = xparam1(ip);']);
0583 eval(['oo_.mle_std.measurement_errors_std.' name ' = stdh(ip);']);
0584 ip = ip+1;
0585 end
0586 end
0587 if ncx
0588 disp('correlation of shocks')
0589 disp(tit1)
0590 ip = nvx+nvn+1;
0591 for i=1:ncx
0592 k1 = estim_params_.corrx(i,1);
0593 k2 = estim_params_.corrx(i,2);
0594 name = [deblank(M_.exo_names(k1,:)) ',' deblank(M_.exo_names(k2,:))];
0595 NAME = [deblank(M_.exo_names(k1,:)) '_' deblank(M_.exo_names(k2,:))];
0596 disp(sprintf('%-*s %8.4f %7.4f %7.4f', header_width,name,xparam1(ip),stdh(ip),tstath(ip)));
0597 M_.Sigma_e(k1,k2) = xparam1(ip)*sqrt(M_.Sigma_e(k1,k1)*M_.Sigma_e(k2,k2));
0598 M_.Sigma_e(k2,k1) = M_.Sigma_e(k1,k2);
0599 eval(['oo_.mle_mode.shocks_corr.' NAME ' = xparam1(ip);']);
0600 eval(['oo_.mle_std.shocks_corr.' NAME ' = stdh(ip);']);
0601 ip = ip+1;
0602 end
0603 end
0604 if ncn
0605 disp('correlation of measurement errors')
0606 disp(tit1)
0607 ip = nvx+nvn+ncx+1;
0608 for i=1:ncn
0609 k1 = estim_params_.corrn(i,1);
0610 k2 = estim_params_.corrn(i,2);
0611 name = [deblank(M_.endo_names(k1,:)) ',' deblank(M_.endo_names(k2,:))];
0612 NAME = [deblank(M_.endo_names(k1,:)) '_' deblank(M_.endo_names(k2,:))];
0613 disp(sprintf('%-*s %8.4f %7.4f %7.4f',header_width,name,xparam1(ip),stdh(ip),tstath(ip)));
0614 eval(['oo_.mle_mode.measurement_error_corr.' NAME ' = xparam1(ip);']);
0615 eval(['oo_.mle_std.measurement_error_corr.' NAME ' = stdh(ip);']);
0616 ip = ip+1;
0617 end
0618 end
0619 end
0620
0621
0622 OutputDirectoryName = CheckPath('Output',M_.dname);
0623
0624 if any(bayestopt_.pshape > 0) && options_.TeX
0625 if np
0626 filename = [OutputDirectoryName '/' M_.fname '_Posterior_Mode_1.TeX'];
0627 fidTeX = fopen(filename,'w');
0628 fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n');
0629 fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (parameters)\n');
0630 fprintf(fidTeX,['%% ' datestr(now,0)]);
0631 fprintf(fidTeX,' \n');
0632 fprintf(fidTeX,' \n');
0633 fprintf(fidTeX,'{\\tiny \n');
0634 fprintf(fidTeX,'\\begin{table}\n');
0635 fprintf(fidTeX,'\\centering\n');
0636 fprintf(fidTeX,'\\caption{Results from posterior maximization (parameters)}\n ');
0637 fprintf(fidTeX,'\\label{Table:Posterior:1}\n');
0638 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0639 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0640 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n');
0641 fprintf(fidTeX,'\\hline \\\\ \n');
0642 ip = nvx+nvn+ncx+ncn+1;
0643 for i=1:np
0644 fprintf(fidTeX,'$%s$ & %s & %7.3f & %6.4f & %8.4f & %7.4f \\\\ \n',...
0645 M_.param_names_tex(estim_params_.param_vals(i,1),:),...
0646 deblank(pnames(bayestopt_.pshape(ip)+1,:)),...
0647 bayestopt_.p1(ip),...
0648 bayestopt_.p2(ip),...
0649 xparam1(ip),...
0650 stdh(ip));
0651 ip = ip + 1;
0652 if ~mod(i,50) && i<np,
0653 fprintf(fidTeX,'\\hline \n');
0654 fprintf(fidTeX,'\\multicolumn{6}{c}{(Table continues next page ...)} \\\\ \n')';
0655 fprintf(fidTeX,'\\end{tabular}\n ');
0656 fprintf(fidTeX,'\\end{table}\n');
0657 fprintf(fidTeX,'\\begin{table}\n');
0658 fprintf(fidTeX,'\\centering\n');
0659 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0660 fprintf(fidTeX,'\\multicolumn{6}{c}{( ... Table continued)} \\\\ \n')';
0661 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0662 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n');
0663 fprintf(fidTeX,'\\hline \\\\ \n');
0664 end
0665 end
0666 fprintf(fidTeX,'\\hline\\hline \n');
0667 fprintf(fidTeX,'\\end{tabular}\n ');
0668 fprintf(fidTeX,'\\end{table}\n');
0669 fprintf(fidTeX,'} \n');
0670 fprintf(fidTeX,'%% End of TeX file.\n');
0671 fclose(fidTeX);
0672 end
0673 if nvx
0674 TeXfile = [OutputDirectoryName '/' M_.fname '_Posterior_Mode_2.TeX'];
0675 fidTeX = fopen(TeXfile,'w');
0676 fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n');
0677 fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (standard deviation of structural shocks)\n');
0678 fprintf(fidTeX,['%% ' datestr(now,0)]);
0679 fprintf(fidTeX,' \n');
0680 fprintf(fidTeX,' \n');
0681 fprintf(fidTeX,'{\\tiny \n');
0682 fprintf(fidTeX,'\\begin{table}\n');
0683 fprintf(fidTeX,'\\centering\n');
0684 fprintf(fidTeX,'\\caption{Results from posterior maximization (standard deviation of structural shocks)}\n ');
0685 fprintf(fidTeX,'\\label{Table:Posterior:2}\n');
0686 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0687 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0688 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n');
0689 fprintf(fidTeX,'\\hline \\\\ \n');
0690 ip = 1;
0691 for i=1:nvx
0692 k = estim_params_.var_exo(i,1);
0693 fprintf(fidTeX,[ '$%s$ & %4s & %7.3f & %6.4f & %8.4f & %7.4f \\\\ \n'],...
0694 deblank(M_.exo_names_tex(k,:)),...
0695 deblank(pnames(bayestopt_.pshape(ip)+1,:)),...
0696 bayestopt_.p1(ip),...
0697 bayestopt_.p2(ip),...
0698 xparam1(ip), ...
0699 stdh(ip));
0700 ip = ip+1;
0701 if ~mod(i,50) && i<nvx,
0702 fprintf(fidTeX,'\\hline \n');
0703 fprintf(fidTeX,'\\multicolumn{6}{c}{(Table continues next page ...)} \\\\ \n')';
0704 fprintf(fidTeX,'\\end{tabular}\n ');
0705 fprintf(fidTeX,'\\end{table}\n');
0706 fprintf(fidTeX,'\\begin{table}\n');
0707 fprintf(fidTeX,'\\centering\n');
0708 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0709 fprintf(fidTeX,'\\multicolumn{6}{c}{( ... Table continued)} \\\\ \n')';
0710 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0711 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n');
0712 fprintf(fidTeX,'\\hline \\\\ \n');
0713 end
0714 end
0715 fprintf(fidTeX,'\\hline\\hline \n');
0716 fprintf(fidTeX,'\\end{tabular}\n ');
0717 fprintf(fidTeX,'\\end{table}\n');
0718 fprintf(fidTeX,'} \n');
0719 fprintf(fidTeX,'%% End of TeX file.\n');
0720 fclose(fidTeX);
0721 end
0722 if nvn
0723 TeXfile = [OutputDirectoryName '/' M_.fname '_Posterior_Mode_3.TeX'];
0724 fidTeX = fopen(TeXfile,'w');
0725 fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n');
0726 fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (standard deviation of measurement errors)\n');
0727 fprintf(fidTeX,['%% ' datestr(now,0)]);
0728 fprintf(fidTeX,' \n');
0729 fprintf(fidTeX,' \n');
0730 fprintf(fidTeX,'\\begin{table}\n');
0731 fprintf(fidTeX,'\\centering\n');
0732 fprintf(fidTeX,'\\caption{Results from posterior maximization (standard deviation of measurement errors)}\n ');
0733 fprintf(fidTeX,'\\label{Table:Posterior:3}\n');
0734 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0735 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0736 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n')
0737 fprintf(fidTeX,'\\hline \\\\ \n');
0738 ip = nvx+1;
0739 for i=1:nvn
0740 idx = strmatch(options_.varobs(estim_params_.var_endo(i,1),:),M_.endo_names);
0741 fprintf(fidTeX,'$%s$ & %4s & %7.3f & %6.4f & %8.4f & %7.4f \\\\ \n',...
0742 deblank(M_.endo_names_tex(idx,:)), ...
0743 deblank(pnames(bayestopt_.pshape(ip)+1,:)), ...
0744 bayestopt_.p1(ip), ...
0745 bayestopt_.p2(ip),...
0746 xparam1(ip),...
0747 stdh(ip));
0748 ip = ip+1;
0749 if ~mod(i,50) && i<nvn,
0750 fprintf(fidTeX,'\\hline \n');
0751 fprintf(fidTeX,'\\multicolumn{6}{c}{(Table continues next page ...)} \\\\ \n')';
0752 fprintf(fidTeX,'\\end{tabular}\n ');
0753 fprintf(fidTeX,'\\end{table}\n');
0754 fprintf(fidTeX,'\\begin{table}\n');
0755 fprintf(fidTeX,'\\centering\n');
0756 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0757 fprintf(fidTeX,'\\multicolumn{6}{c}{( ... Table continued)} \\\\ \n')';
0758 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0759 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n');
0760 fprintf(fidTeX,'\\hline \\\\ \n');
0761 end
0762 end
0763 fprintf(fidTeX,'\\hline\\hline \n');
0764 fprintf(fidTeX,'\\end{tabular}\n ');
0765 fprintf(fidTeX,'\\end{table}\n');
0766 fprintf(fidTeX,'%% End of TeX file.\n');
0767 fclose(fidTeX);
0768 end
0769 if ncx
0770 TeXfile = [OutputDirectoryName '/' M_.fname '_Posterior_Mode_4.TeX'];
0771 fidTeX = fopen(TeXfile,'w');
0772 fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n');
0773 fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (correlation of structural shocks)\n');
0774 fprintf(fidTeX,['%% ' datestr(now,0)]);
0775 fprintf(fidTeX,' \n');
0776 fprintf(fidTeX,' \n');
0777 fprintf(fidTeX,'\\begin{table}\n');
0778 fprintf(fidTeX,'\\centering\n');
0779 fprintf(fidTeX,'\\caption{Results from posterior parameters (correlation of structural shocks)}\n ');
0780 fprintf(fidTeX,'\\label{Table:Posterior:4}\n');
0781 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0782 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0783 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n')
0784 fprintf(fidTeX,'\\hline \\\\ \n');
0785 ip = nvx+nvn+1;
0786 for i=1:ncx
0787 k1 = estim_params_.corrx(i,1);
0788 k2 = estim_params_.corrx(i,2);
0789 fprintf(fidTeX,[ '$%s$ & %s & %7.3f & %6.4f & %8.4f & %7.4f \\\\ \n'],...
0790 [deblank(M_.exo_names_tex(k1,:)) ',' deblank(M_.exo_names_tex(k2,:))], ...
0791 deblank(pnames(bayestopt_.pshape(ip)+1,:)), ...
0792 bayestopt_.p1(ip), ...
0793 bayestopt_.p2(ip), ...
0794 xparam1(ip), ...
0795 stdh(ip));
0796 ip = ip+1;
0797 if ~mod(i,50) && i<ncx,
0798 fprintf(fidTeX,'\\hline \n');
0799 fprintf(fidTeX,'\\multicolumn{6}{c}{(Table continues next page ...)} \\\\ \n')';
0800 fprintf(fidTeX,'\\end{tabular}\n ');
0801 fprintf(fidTeX,'\\end{table}\n');
0802 fprintf(fidTeX,'\\begin{table}\n');
0803 fprintf(fidTeX,'\\centering\n');
0804 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0805 fprintf(fidTeX,'\\multicolumn{6}{c}{( ... Table continued)} \\\\ \n')';
0806 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0807 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n');
0808 fprintf(fidTeX,'\\hline \\\\ \n');
0809 end
0810 end
0811 fprintf(fidTeX,'\\hline\\hline \n');
0812 fprintf(fidTeX,'\\end{tabular}\n ');
0813 fprintf(fidTeX,'\\end{table}\n');
0814 fprintf(fidTeX,'%% End of TeX file.\n');
0815 fclose(fidTeX);
0816 end
0817 if ncn
0818 TeXfile = [OutputDirectoryName '/' M_.fname '_Posterior_Mode_5.TeX'];
0819 fidTeX = fopen(TeXfile,'w');
0820 fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n');
0821 fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (correlation of measurement errors)\n');
0822 fprintf(fidTeX,['%% ' datestr(now,0)]);
0823 fprintf(fidTeX,' \n');
0824 fprintf(fidTeX,' \n');
0825 fprintf(fidTeX,'\\begin{table}\n');
0826 fprintf(fidTeX,'\\centering\n');
0827 fprintf(fidTeX,'\\caption{Results from posterior parameters (correlation of measurement errors)}\n ');
0828 fprintf(fidTeX,'\\label{Table:Posterior:5}\n');
0829 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0830 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0831 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n')
0832 fprintf(fidTeX,'\\hline \\\\ \n');
0833 ip = nvx+nvn+ncx+1;
0834 for i=1:ncn
0835 k1 = estim_params_.corrn(i,1);
0836 k2 = estim_params_.corrn(i,2);
0837 fprintf(fidTeX,'$%s$ & %s & %7.3f & %6.4f & %8.4f & %7.4f \\\\ \n',...
0838 [deblank(M_.endo_names_tex(k1,:)) ',' deblank(M_.endo_names_tex(k2,:))], ...
0839 pnames(bayestopt_.pshape(ip)+1,:), ...
0840 bayestopt_.p1(ip), ...
0841 bayestopt_.p2(ip), ...
0842 xparam1(ip), ...
0843 stdh(ip));
0844 ip = ip+1;
0845 if ~mod(i,50) && i<ncn,
0846 fprintf(fidTeX,'\\hline \n');
0847 fprintf(fidTeX,'\\multicolumn{6}{c}{(Table continues next page ...)} \\\\ \n')';
0848 fprintf(fidTeX,'\\end{tabular}\n ');
0849 fprintf(fidTeX,'\\end{table}\n');
0850 fprintf(fidTeX,'\\begin{table}\n');
0851 fprintf(fidTeX,'\\centering\n');
0852 fprintf(fidTeX,'\\begin{tabular}{l|lcccc} \n');
0853 fprintf(fidTeX,'\\multicolumn{6}{c}{( ... Table continued)} \\\\ \n')';
0854 fprintf(fidTeX,'\\hline\\hline \\\\ \n');
0855 fprintf(fidTeX,' & Prior distribution & Prior mean & Prior s.d. & Posterior mode & s.d. \\\\ \n');
0856 fprintf(fidTeX,'\\hline \\\\ \n');
0857 end
0858 end
0859 fprintf(fidTeX,'\\hline\\hline \n');
0860 fprintf(fidTeX,'\\end{tabular}\n ');
0861 fprintf(fidTeX,'\\end{table}\n');
0862 fprintf(fidTeX,'%% End of TeX file.\n');
0863 fclose(fidTeX);
0864 end
0865 end
0866
0867 if np > 0
0868 pindx = estim_params_.param_vals(:,1);
0869 save([M_.fname '_params.mat'],'pindx');
0870 end
0871
0872 if (any(bayestopt_.pshape >0 ) && options_.mh_replic) || ...
0873 (any(bayestopt_.pshape >0 ) && options_.load_mh_file)
0874 bounds = prior_bounds(bayestopt_,options_);
0875 bounds(:,1)=max(bounds(:,1),lb);
0876 bounds(:,2)=min(bounds(:,2),ub);
0877 bayestopt_.lb = bounds(:,1);
0878 bayestopt_.ub = bounds(:,2);
0879 if any(xparam1 < bounds(:,1)) || any(xparam1 > bounds(:,2))
0880 find(xparam1 < bounds(:,1))
0881 find(xparam1 > bounds(:,2))
0882 error('Mode values are outside prior bounds. Reduce prior_trunc.')
0883 end
0884
0885 if options_.mh_replic
0886 if options_.load_mh_file && options_.use_mh_covariance_matrix
0887 invhess = compute_mh_covariance_matrix;
0888 end
0889 ana_deriv = options_.analytic_derivation;
0890 options_.analytic_derivation = 0;
0891 if options_.cova_compute
0892 feval(options_.posterior_sampling_method,objective_function,options_.proposal_distribution,xparam1,invhess,bounds,dataset_,options_,M_,estim_params_,bayestopt_,oo_);
0893 else
0894 error('I Cannot start the MCMC because the hessian of the posterior kernel at the mode was not computed.')
0895 end
0896 options_.analytic_derivation = ana_deriv;
0897 end
0898 if options_.mh_posterior_mode_estimation
0899 CutSample(M_, options_, estim_params_);
0900 return
0901 else
0902 if ~options_.nodiagnostic && options_.mh_replic > 1000 && options_.mh_nblck > 1
0903 McMCDiagnostics(options_, estim_params_, M_);
0904 end
0905
0906 CutSample(M_, options_, estim_params_);
0907
0908 if options_.mh_replic
0909 [marginal,oo_] = marginal_density(M_, options_, estim_params_, oo_);
0910 oo_ = GetPosteriorParametersStatistics(estim_params_, M_, options_, bayestopt_, oo_);
0911 oo_ = PlotPosteriorDistributions(estim_params_, M_, options_, bayestopt_, oo_);
0912 [oo_.posterior.metropolis.mean,oo_.posterior.metropolis.variance] ...
0913 = GetPosteriorMeanVariance(M_,options_.mh_drop);
0914 else
0915 load([M_.fname '_results'],'oo_');
0916 end
0917 metropolis_draw(1);
0918 if options_.bayesian_irf
0919 PosteriorIRF('posterior');
0920 end
0921 if options_.moments_varendo
0922 oo_ = compute_moments_varendo('posterior',options_,M_,oo_,var_list_);
0923 end
0924 if options_.smoother || ~isempty(options_.filter_step_ahead) || options_.forecast
0925 prior_posterior_statistics('posterior',dataset_);
0926 end
0927 xparam = get_posterior_parameters('mean');
0928 set_all_parameters(xparam);
0929 end
0930 end
0931
0932 if (~((any(bayestopt_.pshape > 0) && options_.mh_replic) || (any(bayestopt_.pshape ...
0933 > 0) && options_.load_mh_file)) ...
0934 || ~options_.smoother ) && options_.partial_information == 0
0935
0936 [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp] = DsgeSmoother(xparam1,dataset_.info.ntobs,dataset_.data,dataset_.missing.aindex,dataset_.missing.state);
0937 oo_.Smoother.SteadyState = ys;
0938 oo_.Smoother.TrendCoeffs = trend_coeff;
0939 oo_.Smoother.variance = P;
0940 i_endo = bayestopt_.smoother_saved_var_list;
0941 if options_.nk ~= 0
0942 oo_.FilteredVariablesKStepAhead = aK(options_.filter_step_ahead, ...
0943 i_endo,:);
0944 if isfield(options_,'kalman_algo')
0945 if ~isempty(PK)
0946 oo_.FilteredVariablesKStepAheadVariances = ...
0947 PK(options_.filter_step_ahead,i_endo,i_endo,:);
0948 end
0949 if ~isempty(decomp)
0950 oo_.FilteredVariablesShockDecomposition = ...
0951 decomp(options_.filter_step_ahead,i_endo,:,:);
0952 end
0953 end
0954 end
0955 for i=bayestopt_.smoother_saved_var_list'
0956 i1 = dr.order_var(bayestopt_.smoother_var_list(i));
0957 eval(['oo_.SmoothedVariables.' deblank(M_.endo_names(i1,:)) ' = ' ...
0958 'atT(i,:)'';']);
0959 if options_.nk > 0
0960 eval(['oo_.FilteredVariables.' deblank(M_.endo_names(i1,:)) ...
0961 ' = squeeze(aK(1,i,:));']);
0962 end
0963 eval(['oo_.UpdatedVariables.' deblank(M_.endo_names(i1,:)) ...
0964 ' = updated_variables(i,:)'';']);
0965 end
0966 if ~options_.nograph,
0967 [nbplt,nr,nc,lr,lc,nstar] = pltorg(M_.exo_nbr);
0968 if options_.TeX
0969 fidTeX = fopen([M_.fname '_SmoothedShocks.TeX'],'w');
0970 fprintf(fidTeX,'%% TeX eps-loader file generated by dynare_estimation.m (Dynare).\n');
0971 fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
0972 fprintf(fidTeX,' \n');
0973 end
0974 for plt = 1:nbplt,
0975 hh = dyn_figure(options_,'Name','Smoothed shocks');
0976 NAMES = [];
0977 if options_.TeX, TeXNAMES = []; end
0978 nstar0=min(nstar,M_.exo_nbr-(plt-1)*nstar);
0979 for i=1:nstar0,
0980 k = (plt-1)*nstar+i;
0981 subplot(nr,nc,i);
0982 plot([1 gend],[0 0],'-r','linewidth',.5)
0983 hold on
0984 plot(1:gend,innov(k,:),'-k','linewidth',1)
0985 hold off
0986 name = deblank(M_.exo_names(k,:));
0987 if isempty(NAMES)
0988 NAMES = name;
0989 else
0990 NAMES = char(NAMES,name);
0991 end
0992 if ~isempty(options_.XTick)
0993 set(gca,'XTick',options_.XTick)
0994 set(gca,'XTickLabel',options_.XTickLabel)
0995 end
0996 xlim([1 gend])
0997 if options_.TeX
0998 texname = M_.exo_names_tex(k,:);
0999 if isempty(TeXNAMES)
1000 TeXNAMES = ['$ ' deblank(texname) ' $'];
1001 else
1002 TeXNAMES = char(TeXNAMES,['$ ' deblank(texname) ' $']);
1003 end
1004 end
1005 title(name,'Interpreter','none')
1006 eval(['oo_.SmoothedShocks.' deblank(name) ' = innov(k,:)'';']);
1007 end
1008 dyn_saveas(hh,[M_.fname '_SmoothedShocks' int2str(plt)],options_);
1009 if options_.TeX
1010 fprintf(fidTeX,'\\begin{figure}[H]\n');
1011 for jj = 1:nstar0
1012 fprintf(fidTeX,'\\psfrag{%s}[1][][0.5][0]{%s}\n',deblank(NAMES(jj,:)),deblank(TeXNAMES(jj,:)));
1013 end
1014 fprintf(fidTeX,'\\centering \n');
1015 fprintf(fidTeX,'\\includegraphics[scale=0.5]{%s_SmoothedShocks%s}\n',M_.fname,int2str(plt));
1016 fprintf(fidTeX,'\\caption{Smoothed shocks.}');
1017 fprintf(fidTeX,'\\label{Fig:SmoothedShocks:%s}\n',int2str(plt));
1018 fprintf(fidTeX,'\\end{figure}\n');
1019 fprintf(fidTeX,'\n');
1020 end
1021 end
1022 if options_.TeX
1023 fprintf(fidTeX,'\n');
1024 fprintf(fidTeX,'%% End of TeX file.\n');
1025 fclose(fidTeX);
1026 end
1027 end
1028
1029
1030
1031 if options_.prefilter == 1
1032 yf = atT(bayestopt_.mf,:)+repmat(dataset_.descriptive.mean',1,gend);
1033 elseif options_.loglinear == 1
1034 yf = atT(bayestopt_.mf,:)+repmat(log(ys(bayestopt_.mfys)),1,gend)+...
1035 trend_coeff*[1:gend];
1036 else
1037 yf = atT(bayestopt_.mf,:)+repmat(ys(bayestopt_.mfys),1,gend)+...
1038 trend_coeff*[1:gend];
1039 end
1040 if nvn
1041 number_of_plots_to_draw = 0;
1042 index = [];
1043 for i=1:n_varobs
1044 if max(abs(measurement_error)) > 0.000000001
1045 number_of_plots_to_draw = number_of_plots_to_draw + 1;
1046 index = cat(1,index,i);
1047 end
1048 eval(['oo_.SmoothedMeasurementErrors.' deblank(options_.varobs(i,:)) ...
1049 ' = measurement_error(i,:)'';']);
1050 end
1051 if ~options_.nograph
1052 [nbplt,nr,nc,lr,lc,nstar] = pltorg(number_of_plots_to_draw);
1053 if options_.TeX
1054 fidTeX = fopen([M_.fname '_SmoothedObservationErrors.TeX'],'w');
1055 fprintf(fidTeX,'%% TeX eps-loader file generated by dynare_estimation.m (Dynare).\n');
1056 fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
1057 fprintf(fidTeX,' \n');
1058 end
1059 for plt = 1:nbplt
1060 hh = dyn_figure(options_,'Name','Smoothed observation errors');
1061 NAMES = [];
1062 if options_.TeX, TeXNAMES = []; end
1063 for i=1:min(nstar,number_of_plots_to_draw-(nbplt-1)*nstar)
1064 k = (plt-1)*nstar+i;
1065 subplot(nr,nc,i);
1066 plot([1 gend],[0 0],'-r','linewidth',.5)
1067 hold on
1068 plot(1:gend,measurement_error(index(k),:),'-k','linewidth',1)
1069 hold off
1070 name = deblank(options_.varobs(index(k),:));
1071 if isempty(NAMES)
1072 NAMES = name;
1073 else
1074 NAMES = char(NAMES,name);
1075 end
1076 if ~isempty(options_.XTick)
1077 set(gca,'XTick',options_.XTick)
1078 set(gca,'XTickLabel',options_.XTickLabel)
1079 end
1080 if options_.TeX
1081 idx = strmatch(options_.varobs(index(k),:),M_.endo_names,'exact');
1082 texname = M_.endo_names_tex(idx,:);
1083 if isempty(TeXNAMES)
1084 TeXNAMES = ['$ ' deblank(texname) ' $'];
1085 else
1086 TeXNAMES = char(TeXNAMES,['$ ' deblank(texname) ' $']);
1087 end
1088 end
1089 title(name,'Interpreter','none')
1090 end
1091 dyn_saveas(hh,[M_.fname '_SmoothedObservationErrors' int2str(plt)],options_);
1092 if options_.TeX
1093 fprintf(fidTeX,'\\begin{figure}[H]\n');
1094 for jj = 1:nstar
1095 fprintf(fidTeX,'\\psfrag{%s}[1][][0.5][0]{%s}\n',deblank(NAMES(jj,:)),deblank(TeXNAMES(jj,:)));
1096 end
1097 fprintf(fidTeX,'\\centering \n');
1098 fprintf(fidTeX,'\\includegraphics[scale=0.5]{%s_SmoothedObservationErrors%s}\n',M_.fname,int2str(plt));
1099 fprintf(fidTeX,'\\caption{Smoothed observation errors.}');
1100 fprintf(fidTeX,'\\label{Fig:SmoothedObservationErrors:%s}\n',int2str(plt));
1101 fprintf(fidTeX,'\\end{figure}\n');
1102 fprintf(fidTeX,'\n');
1103 end
1104 end
1105 if options_.TeX
1106 fprintf(fidTeX,'\n');
1107 fprintf(fidTeX,'%% End of TeX file.\n');
1108 fclose(fidTeX);
1109 end
1110 end
1111 end
1112
1113
1114
1115 if ~options_.nograph
1116 [nbplt,nr,nc,lr,lc,nstar] = pltorg(n_varobs);
1117 if options_.TeX
1118 fidTeX = fopen([M_.fname '_HistoricalAndSmoothedVariables.TeX'],'w');
1119 fprintf(fidTeX,'%% TeX eps-loader file generated by dynare_estimation.m (Dynare).\n');
1120 fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
1121 fprintf(fidTeX,' \n');
1122 end
1123 for plt = 1:nbplt,
1124 hh = dyn_figure(options_,'Name','Historical and smoothed variables');
1125 NAMES = [];
1126 if options_.TeX, TeXNAMES = []; end
1127 nstar0=min(nstar,n_varobs-(plt-1)*nstar);
1128 for i=1:nstar0,
1129 k = (plt-1)*nstar+i;
1130 subplot(nr,nc,i);
1131 plot(1:gend,yf(k,:),'--r','linewidth',1)
1132 hold on
1133 plot(1:gend,rawdata(:,k),'--k','linewidth',1)
1134 hold off
1135 name = deblank(options_.varobs(k,:));
1136 if isempty(NAMES)
1137 NAMES = name;
1138 else
1139 NAMES = char(NAMES,name);
1140 end
1141 if ~isempty(options_.XTick)
1142 set(gca,'XTick',options_.XTick)
1143 set(gca,'XTickLabel',options_.XTickLabel)
1144 end
1145 xlim([1 gend])
1146 if options_.TeX
1147 idx = strmatch(options_.varobs(k,:),M_.endo_names,'exact');
1148 texname = M_.endo_names_tex(idx,:);
1149 if isempty(TeXNAMES)
1150 TeXNAMES = ['$ ' deblank(texname) ' $'];
1151 else
1152 TeXNAMES = char(TeXNAMES,['$ ' deblank(texname) ' $']);
1153 end
1154 end
1155 title(name,'Interpreter','none')
1156 end
1157 dyn_saveas(hh,[M_.fname '_HistoricalAndSmoothedVariables' int2str(plt)],options_);
1158 if options_.TeX
1159 fprintf(fidTeX,'\\begin{figure}[H]\n');
1160 for jj = 1:nstar0,
1161 fprintf(fidTeX,'\\psfrag{%s}[1][][0.5][0]{%s}\n',deblank(NAMES(jj,:)),deblank(TeXNAMES(jj,:)));
1162 end
1163 fprintf(fidTeX,'\\centering \n');
1164 fprintf(fidTeX,'\\includegraphics[scale=0.5]{%s_HistoricalAndSmoothedVariables%s}\n',M_.fname,int2str(plt));
1165 fprintf(fidTeX,'\\caption{Historical and smoothed variables.}');
1166 fprintf(fidTeX,'\\label{Fig:HistoricalAndSmoothedVariables:%s}\n',int2str(plt));
1167 fprintf(fidTeX,'\\end{figure}\n');
1168 fprintf(fidTeX,'\n');
1169 end
1170 end
1171 if options_.TeX
1172 fprintf(fidTeX,'\n');
1173 fprintf(fidTeX,'%% End of TeX file.\n');
1174 fclose(fidTeX);
1175 end
1176 end
1177 end
1178
1179 if options_.forecast > 0 && options_.mh_replic == 0 && ~options_.load_mh_file
1180 dyn_forecast(var_list_,'smoother');
1181 end
1182
1183 if np > 0
1184 pindx = estim_params_.param_vals(:,1);
1185 save([M_.fname '_pindx.mat'] ,'pindx');
1186 end
1187