From 58cb6e42b6645592a0939860899556ce60c1d1bb Mon Sep 17 00:00:00 2001 From: adjemian Date: Mon, 6 Mar 2006 08:38:11 +0000 Subject: [PATCH] TeX formatting and various things... git-svn-id: https://www.dynare.org/svn/dynare/dynare_v4@656 ac1d8469-bf42-47a9-8791-bf33cf982152 --- matlab/dynare_estimation.m | 43 +++++++++++--------------------------- matlab/metropolis.m | 41 ++++++++++++++++++------------------ matlab/posterior_moments.m | 1 - 3 files changed, 32 insertions(+), 53 deletions(-) diff --git a/matlab/dynare_estimation.m b/matlab/dynare_estimation.m index 240152fcb..6200645a6 100644 --- a/matlab/dynare_estimation.m +++ b/matlab/dynare_estimation.m @@ -191,8 +191,7 @@ rawdata = read_variables(options_.datafile,options_.varobs,[],options_.xls_sheet k = []; k1 = []; for i=1:n_varobs - k = [k strmatch(deblank(options_.varobs(i,:)),M_.endo_names(dr.order_var,:), ... - 'exact')]; + k = [k strmatch(deblank(options_.varobs(i,:)),M_.endo_names(dr.order_var,:),'exact')]; k1 = [k1 strmatch(deblank(options_.varobs(i,:)),M_.endo_names, 'exact')]; end @@ -257,29 +256,9 @@ if options_.mode_compute > 0 & options_.posterior_mode_estimation disp(sprintf('Objective function at mode: %f',fval)) disp(sprintf('Objective function at mode: %f',DsgeLikelihood(xparam1,gend,data))) elseif options_.mode_compute == 5 - if isfield(options_,'hess') - flag = options_.hess; - else - flag = 1; - end - if ~exist('igg'), % by M. Ratto - hh=[]; - gg=[]; - igg=[]; - end % by M. Ratto - if isfield(options_,'ftol') - crit = options_.ftol; - else - crit = 1.e-7; - end - if isfield(options_,'nit') - nit = options_.nit; - else - nit=1000; - end - %[xparam1, hh, gg, fval] = newrat('DsgeLikelihood',xparam1,[],[],flag,gend,data); - [xparam1, hh, gg, fval, invhess] = newrat('DsgeLikelihood',xparam1,hh,gg,igg,crit,nit,flag,gend,data); - eval(['save ' M_.fname '_mode xparam1 hh gg fval invhess;']); + flag = 0; + [xparam1, hh, gg, fval] = newrat('DsgeLikelihood',xparam1,[],[],flag,gend,data); + eval(['save ' M_.fname '_mode xparam1 hh gg fval;']); end if options_.mode_compute ~= 5 hh = reshape(hessian('DsgeLikelihood',xparam1,gend,data),nx,nx); @@ -483,9 +462,11 @@ elseif ~any(bayestopt_.pshape > 0) & options_.posterior_mode_estimation end +OutputDirectoryName = CheckPath('Output'); + if any(bayestopt_.pshape > 0) & options_.TeX %% Bayesian estimation (posterior mode) Latex output if np - filename = [M_.fname '_Posterior_Mode_1.TeX']; + filename = [OutputDirectoryName '\' M_.fname '_Posterior_Mode_1.TeX']; fidTeX = fopen(filename,'w'); fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n'); fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (parameters)\n'); @@ -502,7 +483,7 @@ if any(bayestopt_.pshape > 0) & options_.TeX %% Bayesian estimation (posterior m ip = nvx+nvn+ncx+ncn+1; for i=1:np fprintf(fidTeX,'$%s$ & %s & %7.3f & %6.4f & %8.4f & %7.4f \\\\ \n',... - M_.param_names_tex(estim_params_.param_vals(i,1),:),...%deblank(estim_params_.tex(i,:)),... + M_.param_names_tex(estim_params_.param_vals(i,1),:),... deblank(pnames(bayestopt_.pshape(ip)+1,:)),... bayestopt_.pmean(ip),... estim_params_.param_vals(i,6),... @@ -520,7 +501,7 @@ if any(bayestopt_.pshape > 0) & options_.TeX %% Bayesian estimation (posterior m fclose(fidTeX); end if nvx - TeXfile = [M_.fname '_Posterior_Mode_2.TeX']; + TeXfile = [OutputDirectoryName '\' M_.fname '_Posterior_Mode_2.TeX']; fidTeX = fopen(TeXfile,'w'); fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n'); fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (standard deviation of structural shocks)\n'); @@ -556,7 +537,7 @@ if any(bayestopt_.pshape > 0) & options_.TeX %% Bayesian estimation (posterior m fclose(fidTeX); end if nvn - TeXfile = [M_.fname '_Posterior_Mode_3.TeX']; + TeXfile = [OutputDirectoryName '\' M_.fname '_Posterior_Mode_3.TeX']; fidTeX = fopen(TeXfile,'w'); fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n'); fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (standard deviation of measurement errors)\n'); @@ -590,7 +571,7 @@ if any(bayestopt_.pshape > 0) & options_.TeX %% Bayesian estimation (posterior m fclose(fidTeX); end if ncx - TeXfile = [M_.fname '_Posterior_Mode_4.TeX']; + TeXfile = [OutputDirectoryName '\' M_.fname '_Posterior_Mode_4.TeX']; fidTeX = fopen(TeXfile,'w'); fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n'); fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (correlation of structural shocks)\n'); @@ -625,7 +606,7 @@ if any(bayestopt_.pshape > 0) & options_.TeX %% Bayesian estimation (posterior m fclose(fidTeX); end if ncn - TeXfile = [M_.fname '_Posterior_Mode_5.TeX']; + TeXfile = [OutputDirectoryName '\' M_.fname '_Posterior_Mode_5.TeX']; fidTeX = fopen(TeXfile,'w'); fprintf(fidTeX,'%% TeX-table generated by dynare_estimation (Dynare).\n'); fprintf(fidTeX,'%% RESULTS FROM POSTERIOR MAXIMIZATION (correlation of measurement errors)\n'); diff --git a/matlab/metropolis.m b/matlab/metropolis.m index 153781934..ac9dd31f0 100644 --- a/matlab/metropolis.m +++ b/matlab/metropolis.m @@ -4,7 +4,7 @@ global M_ oo_ options_ bayestopt_ estim_params_ bayestopt_.penalty = 1e8; -DirectoryName = CheckPath('metropolis'); +MhDirectoryName = CheckPath('metropolis'); nblck = options_.mh_nblck; nruns = ones(nblck,1)*options_.mh_replic; @@ -23,9 +23,9 @@ if options_.load_mh_file == 0 disp('MH: One Chain mode.') end % Delete old mh files... - files = eval(['dir(''' DirectoryName '/' M_.fname '_mh*_blck*.mat'');']); + files = dir([ MhDirectoryName '/' M_.fname '_mh*_blck*.mat']); if length(files) - delete([ DirectoryName '/' M_.fname '_mh*_blck*.mat']); + delete([ MhDirectoryName '/' M_.fname '_mh*_blck*.mat']); disp('MH: Old _mh files succesfully erased!') end % Initial values... @@ -78,9 +78,9 @@ if options_.load_mh_file == 0 fline = ones(nblck,1); NewFile = ones(nblck,1); % Creation of the mh-history file: - file = eval(['dir(''' DirectoryName '/' M_.fname '_mh_history.mat'');']); + file = dir([MhDirectoryName '/' M_.fname '_mh_history.mat']); if length(files) - delete([ DirectoryName '/' M_.fname '_mh_history.mat']); + delete([ MhDirectoryName '/' M_.fname '_mh_history.mat']); disp('MH: Old mh_history file succesfully erased!') end AnticipatedNumberOfFiles = floor(nruns(1)/MAX_nruns); @@ -102,8 +102,8 @@ if options_.load_mh_file == 0 save([DirectoryName '/' M_.fname '_mh_history'],'record'); elseif options_.load_mh_file == 1% Here we consider previous mh files (previous mh did not crash). disp('MH: I''m loading past metropolis-hastings simulations...') - file = eval(['dir(''' DirectoryName '/' M_.fname '_mh_history.mat'');']); - files = eval(['dir(''' DirectoryName '/' M_.fname '_mh*.mat'');']); + file = dir([ MhDirectoryName '/' M_.fname '_mh_history.mat' ]); + files = dir([ MhDirectoryName '/' M_.fname '_mh*.mat']); if ~length(files) disp('MH:: FAILURE! there is no MH file to load here!') return @@ -112,7 +112,7 @@ elseif options_.load_mh_file == 1% Here we consider previous mh files (previous disp('MH:: FAILURE! there is no MH-history file!') return else - load([ DirectoryName '/' M_.fname '_mh_history']) + load([ MhDirectoryName '/' M_.fname '_mh_history']) end past_number_of_blocks = record.Nblck; if past_number_of_blocks ~= nblck @@ -148,7 +148,7 @@ elseif options_.load_mh_file == 1% Here we consider previous mh files (previous record.MhDraws(end,3) = AnticipatedNumberOfLinesInTheLastFile; randn('state',record.Seeds.Normal); rand('state',record.Seeds.Unifor); - save([DirectoryName '/' M_.fname '_mh_history'],'record'); + save([MhDirectoryName '/' M_.fname '_mh_history'],'record'); disp(['MH: ... It''s done. I''ve loaded ' int2str(NumberOfPreviousSimulations) ' simulations.']) disp(' ') elseif options_.load_mh_file == -1% The previous metropolis-hastings @@ -156,13 +156,13 @@ elseif options_.load_mh_file == -1% The previous metropolis-hastings % recover the saved draws... disp('MH: Recover mode!') disp(' ') - file = eval(['dir(''' DirectoryName '/' M_.fname '_mh_history.mat'');']); + file = dir([MhDirectoryName '/' M_.fname '_mh_history.mat']); if ~length(file) disp('MH:: FAILURE! there is no MH-history file!') return else - load([ DirectoryName '/' M_.fname '_mh_history']) - end + load([ MhDirectoryName '/' M_.fname '_mh_history']) + end nblck = record.Nblck; options_.mh_nblck = nblck; if size(record.MhDraws,1) == 1 @@ -202,12 +202,12 @@ elseif options_.load_mh_file == -1% The previous metropolis-hastings ExpectedNumberOfMhFilesPerBlock = sum(record.MhDraws(:,2),1); ExpectedNumberOfMhFiles = ExpectedNumberOfMhFilesPerBlock*nblck; % I count the total number of saved mh files... - AllMhFiles = eval(['dir(''' DirectoryName '/' M_.fname '_mh*_blck*.mat'');']); + AllMhFiles = dir([MhDirectoryName '/' M_.fname '_mh*_blck*.mat']); TotalNumberOfMhFiles = length(AllMhFiles); % I count the number of saved mh files per block NumberOfMhFilesPerBlock = zeros(nblck,1); for i = 1:nblck - BlckMhFiles = eval(['dir(''' DirectoryName '/' M_.fname '_mh*_blck' int2str(i) '.mat'');']); + BlckMhFiles = dir([ MhDirectoryName '/' M_.fname '_mh*_blck' int2str(i) '.mat']); NumberOfMhFilesPerBlock(i) = length(BlckMhFiles); end tmp = NumberOfMhFilesPerBlock(1); b = 1; @@ -229,8 +229,7 @@ elseif options_.load_mh_file == -1% The previous metropolis-hastings % (if there was a complete session before the crash) ? if OldMh ante = sum(record.MhDraws(1:end-1,2),1); - load(['./' DirectoryName '/' M_.fname '_mh' int2str(ante) '_blck' ... - int2str(CrashedBlck) '.mat'],'logpo2'); + load([MhDirectoryName '/' M_.fname '_mh' int2str(ante) '_blck' int2str(CrashedBlck) '.mat'],'logpo2'); if length(logpo2) == MAX_nruns IsTheLastFileOfThePreviousMhFull = 1; else @@ -245,7 +244,7 @@ elseif options_.load_mh_file == -1% The previous metropolis-hastings MhFileNumber = ante; while MhFileExist MhFileNumber = MhFileNumber + 1; - if ~exist(['./' DirectoryName '/' M_.fname '_mh' int2str(MhFileNumber) '_blck' int2str(CrashedBlck) '.mat']) + if ~exist([MhDirectoryName '/' M_.fname '_mh' int2str(MhFileNumber) '_blck' int2str(CrashedBlck) '.mat']) MhFileExist = 0; end end @@ -267,7 +266,7 @@ elseif options_.load_mh_file == -1% The previous metropolis-hastings nruns(CrashedBlck) = nruns(CrashedBlck)-NumberOfSavedDraws; % I initialize with the last saved mh file of the inccomplete % block: - load(['./' DirectoryName '/' M_.fname '_mh' int2str(MhFileNumber-1) '_blck' int2str(CrashedBlck) '.mat']); + load([MhDirectoryName '/' M_.fname '_mh' int2str(MhFileNumber-1) '_blck' int2str(CrashedBlck) '.mat']); ilogpo2(CrashedBlck) = logpo2(end); ix2(CrashedBlck,:) = x2(end,:); NewFile(CrashedBlck) = MhFileNumber; @@ -319,7 +318,7 @@ for b = fblck:nblck prtfrc = j/nruns(b); waitbar(prtfrc,hh,[ '(' int2str(b) '/' int2str(nblck) ') ' sprintf('%f done, acceptation rate %f',prtfrc,isux/j)]); if (irun == InitSizeArray(b)) | (j == nruns(b)) % Now I save the simulations - save([DirectoryName '/' M_.fname '_mh' int2str(NewFile(b)) '_blck' int2str(b)],'x2','logpo2'); + save([MhDirectoryName '/' M_.fname '_mh' int2str(NewFile(b)) '_blck' int2str(b)],'x2','logpo2'); InitSizeArray(b) = min(nruns(b)-j,MAX_nruns); if j == nruns(b) % I record the last draw... record.LastParameters(b,:) = x2(end,:); @@ -342,8 +341,8 @@ for b = fblck:nblck end% End of the loop over the mh-blocks. record.Seeds.Normal = randn('state'); record.Seeds.Unifor = rand('state'); -save([DirectoryName '/' M_.fname '_mh_history'],'record'); +save([MhDirectoryName '/' M_.fname '_mh_history'],'record'); disp(['MH: Number of mh files : ' int2str(NewFile(1)) ' per block.']) disp(['MH: Total number of generated files : ' int2str(NewFile(1)*nblck) '.']) disp(['MH: Total number of iterations : ' int2str((NewFile(1)-1)*MAX_nruns+irun-1) '.']) -disp(' ') +disp(' ') \ No newline at end of file diff --git a/matlab/posterior_moments.m b/matlab/posterior_moments.m index 0284dc672..122744a20 100644 --- a/matlab/posterior_moments.m +++ b/matlab/posterior_moments.m @@ -39,6 +39,5 @@ else bandwidth = 0; % Rule of thumb optimal bandwidth parameter. kernel_function = 'gaussian'; % Gaussian kernel for Fast Fourrier Transform approximaton. optimal_bandwidth = mh_optimal_bandwidth(xx,length(xx),bandwidth,kernel_function); - % [abscissa,f] = kernel_density_estimate(,,,) [density(:,1),density(:,2)] = kernel_density_estimate(xx,number_of_grid_points,optimal_bandwidth,kernel_function); end \ No newline at end of file