stepan
/
dynare
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Remove global variables from prior_posterior_statistics.m and PosteriorIRF.m

kalman-mex
Johannes Pfeifer 2023-09-14 20:17:01 +02:00
parent 782a2e8d69
commit 5231fc04c1
11 changed files with 259 additions and 296 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

112
matlab/PosteriorIRF.m
View File

@ -1,10 +1,18 @@
function PosteriorIRF(type,dispString) function oo_=PosteriorIRF(type,options_,estim_params_,oo_,M_,bayestopt_,dataset_,dataset_info,dispString)
% PosteriorIRF(type,options_,estim_params_,oo_,M_,bayestopt_,dataset_,dataset_info,dispString)
% Builds posterior IRFs after the MH algorithm. % Builds posterior IRFs after the MH algorithm.
% %
% INPUTS % INPUTS
% o type [char] string specifying the joint density of the % o type [char] string specifying the joint density of the
% deep parameters ('prior','posterior'). % deep parameters ('prior','posterior').
% o dispString [char] string to display in the console. % o options_ [structure] storing the options
% o estim_params_ [structure] storing information about estimated parameters
% o oo_ [structure] storing the results
% o M_ [structure] storing the model information
% o bayestopt_ [structure] storing information about priors
% o dataset_ [structure] storing the dataset
% o dataset_info [structure] Various information about the dataset
% o dispString [char] string to display in the console.
% %
% OUTPUTS % OUTPUTS
% None (oo_ and plots). % None (oo_ and plots).
@ -34,9 +42,6 @@ function PosteriorIRF(type,dispString)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>. % along with Dynare. If not, see <https://www.gnu.org/licenses/>.
global options_ estim_params_ oo_ M_ bayestopt_ dataset_ dataset_info
% Set the number of periods % Set the number of periods
if isempty(options_.irf) || ~options_.irf if isempty(options_.irf) || ~options_.irf
options_.irf = 40; options_.irf = 40;
@ -62,13 +67,7 @@ end
irf_shocks_indx = getIrfShocksIndx(M_, options_); irf_shocks_indx = getIrfShocksIndx(M_, options_);
% Set various parameters & Check or create directories % Set various parameters & Check or create directories
nvx = estim_params_.nvx; npar = estim_params_.nvx+estim_params_.nvn+estim_params_.ncx+estim_params_.ncn+estim_params_.np ;
nvn = estim_params_.nvn;
ncx = estim_params_.ncx;
ncn = estim_params_.ncn;
np = estim_params_.np ;
npar = nvx+nvn+ncx+ncn+np;
offset = npar-np; clear('nvx','nvn','ncx','ncn','np');
nvobs = dataset_.vobs; nvobs = dataset_.vobs;
gend = dataset_.nobs; gend = dataset_.nobs;
@ -118,7 +117,8 @@ elseif strcmpi(type,'gsa')
end end
x=[lpmat0(istable,:) lpmat(istable,:)]; x=[lpmat0(istable,:) lpmat(istable,:)];
clear lpmat istable clear lpmat istable
B=size(x,1); options_.B = B; B=size(x,1);
options_.B = B;
else% type = 'prior' else% type = 'prior'
B = options_.prior_draws; B = options_.prior_draws;
options_.B = B; options_.B = B;
@ -130,23 +130,7 @@ irun2 = 0;
NumberOfIRFfiles_dsge = 1; NumberOfIRFfiles_dsge = 1;
NumberOfIRFfiles_dsgevar = 1; NumberOfIRFfiles_dsgevar = 1;
ifil2 = 1; ifil2 = 1;
% Create arrays
if B <= MAX_nruns
stock_param = zeros(B, npar);
else
stock_param = zeros(MAX_nruns, npar);
end
if B >= MAX_nirfs_dsge
stock_irf_dsge = zeros(options_.irf,nvar,M_.exo_nbr,MAX_nirfs_dsge);
else
stock_irf_dsge = zeros(options_.irf,nvar,M_.exo_nbr,B);
end
if MAX_nirfs_dsgevar if MAX_nirfs_dsgevar
if B >= MAX_nirfs_dsgevar
stock_irf_bvardsge = zeros(options_.irf,nvobs,M_.exo_nbr,MAX_nirfs_dsgevar);
else
stock_irf_bvardsge = zeros(options_.irf,nvobs,M_.exo_nbr,B);
end
NumberOfLags = options_.dsge_varlag; NumberOfLags = options_.dsge_varlag;
NumberOfLagsTimesNvobs = NumberOfLags*nvobs; NumberOfLagsTimesNvobs = NumberOfLags*nvobs;
if options_.noconstant if options_.noconstant
@ -154,10 +138,9 @@ if MAX_nirfs_dsgevar
else else
NumberOfParametersPerEquation = NumberOfLagsTimesNvobs+1; NumberOfParametersPerEquation = NumberOfLagsTimesNvobs+1;
end end
Companion_matrix = diag(ones(nvobs*(NumberOfLags-1),1),-nvobs);
end end
% First block of code executed in parallel. The function devoted to do it is PosteriorIRF_core1.m % First block of code executed in parallel. The function devoted to do it is PosteriorIRF_core1.m
% function. % function.
b = 0; b = 0;
@ -183,7 +166,6 @@ if ~strcmpi(type,'prior')
localVars.x=x; localVars.x=x;
end end
b=0;
if options_.dsge_var if options_.dsge_var
localVars.gend = gend; localVars.gend = gend;
localVars.nvobs = nvobs; localVars.nvobs = nvobs;
@ -202,6 +184,16 @@ localVars.NumberOfIRFfiles_dsgevar=NumberOfIRFfiles_dsgevar;
localVars.ifil2=ifil2; localVars.ifil2=ifil2;
localVars.MhDirectoryName=MhDirectoryName; localVars.MhDirectoryName=MhDirectoryName;
%store main structures
localVars.options_=options_;
localVars.estim_params_= estim_params_;
localVars.M_= M_;
localVars.oo_= oo_;
localVars.bayestopt_= bayestopt_;
localVars.dataset_= dataset_;
localVars.dataset_info= dataset_info;
% Like sequential execution! % Like sequential execution!
if isnumeric(options_.parallel) if isnumeric(options_.parallel)
[fout] = PosteriorIRF_core1(localVars,1,B,0); [fout] = PosteriorIRF_core1(localVars,1,B,0);
@ -225,14 +217,6 @@ else
localVars.NumberOfIRFfiles_dsgevar=NumberOfIRFfiles_dsgevar; localVars.NumberOfIRFfiles_dsgevar=NumberOfIRFfiles_dsgevar;
localVars.ifil2=ifil2; localVars.ifil2=ifil2;
globalVars = struct('M_',M_, ...
'options_', options_, ...
'bayestopt_', bayestopt_, ...
'estim_params_', estim_params_, ...
'oo_', oo_, ...
'dataset_',dataset_, ...
'dataset_info',dataset_info);
% which files have to be copied to run remotely % which files have to be copied to run remotely
NamFileInput(1,:) = {'',[M_.fname '.static.m']}; NamFileInput(1,:) = {'',[M_.fname '.static.m']};
NamFileInput(2,:) = {'',[M_.fname '.dynamic.m']}; NamFileInput(2,:) = {'',[M_.fname '.dynamic.m']};
@ -246,7 +230,7 @@ else
NamFileInput(length(NamFileInput)+1,:)={'',[M_.fname '.steadystate.m']}; NamFileInput(length(NamFileInput)+1,:)={'',[M_.fname '.steadystate.m']};
end end
end end
[fout] = masterParallel(options_.parallel, 1, B,NamFileInput,'PosteriorIRF_core1', localVars, globalVars, options_.parallel_info); [fout] = masterParallel(options_.parallel, 1, B,NamFileInput,'PosteriorIRF_core1', localVars, [], options_.parallel_info);
nosaddle=0; nosaddle=0;
for j=1:length(fout) for j=1:length(fout)
nosaddle = nosaddle + fout(j).nosaddle; nosaddle = nosaddle + fout(j).nosaddle;
@ -260,9 +244,9 @@ if nosaddle
disp(['PosteriorIRF :: Percentage of discarded posterior draws = ' num2str(nosaddle/(B+nosaddle))]) disp(['PosteriorIRF :: Percentage of discarded posterior draws = ' num2str(nosaddle/(B+nosaddle))])
end end
ReshapeMatFiles('irf_dsge',type) ReshapeMatFiles(M_.fname,M_.dname,M_.exo_nbr,M_.endo_nbr,options_,'irf_dsge',type)
if MAX_nirfs_dsgevar if MAX_nirfs_dsgevar
ReshapeMatFiles('irf_bvardsge') ReshapeMatFiles(M_.fname,M_.dname,M_.exo_nbr,M_.endo_nbr,options_,'irf_bvardsge')
end end
if strcmpi(type,'gsa') if strcmpi(type,'gsa')
@ -293,21 +277,21 @@ tit = M_.exo_names;
kdx = 0; kdx = 0;
for file = 1:NumberOfIRFfiles_dsge for file = 1:NumberOfIRFfiles_dsge
load([MhDirectoryName filesep M_.fname '_IRF_DSGEs' int2str(file) '.mat']); temp=load([MhDirectoryName filesep M_.fname '_IRF_DSGEs' int2str(file) '.mat']);
for i = irf_shocks_indx for i = irf_shocks_indx
for j = 1:nvar for j = 1:nvar
for k = 1:size(STOCK_IRF_DSGE,1) for k = 1:size(temp.STOCK_IRF_DSGE,1)
kk = k+kdx; kk = k+kdx;
[MeanIRF(kk,j,i),MedianIRF(kk,j,i),VarIRF(kk,j,i),HPDIRF(kk,:,j,i),... [MeanIRF(kk,j,i),MedianIRF(kk,j,i),VarIRF(kk,j,i),HPDIRF(kk,:,j,i),...
DistribIRF(kk,:,j,i)] = posterior_moments(squeeze(STOCK_IRF_DSGE(k,j,i,:)),0,options_.mh_conf_sig); DistribIRF(kk,:,j,i)] = posterior_moments(squeeze(temp.STOCK_IRF_DSGE(k,j,i,:)),0,options_.mh_conf_sig);
end end
end end
end end
kdx = kdx + size(STOCK_IRF_DSGE,1); kdx = kdx + size(temp.STOCK_IRF_DSGE,1);
end end
clear STOCK_IRF_DSGE; clear temp;
for i = irf_shocks_indx for i = irf_shocks_indx
for j = 1:nvar for j = 1:nvar
@ -332,20 +316,20 @@ if MAX_nirfs_dsgevar
tit = M_.exo_names; tit = M_.exo_names;
kdx = 0; kdx = 0;
for file = 1:NumberOfIRFfiles_dsgevar for file = 1:NumberOfIRFfiles_dsgevar
load([MhDirectoryName filesep M_.fname '_IRF_BVARDSGEs' int2str(file) '.mat']); temp=load([MhDirectoryName filesep M_.fname '_IRF_BVARDSGEs' int2str(file) '.mat']);
for i = irf_shocks_indx for i = irf_shocks_indx
for j = 1:nvar for j = 1:nvar
for k = 1:size(STOCK_IRF_BVARDSGE,1) for k = 1:size(temp.STOCK_IRF_BVARDSGE,1)
kk = k+kdx; kk = k+kdx;
[MeanIRFdsgevar(kk,j,i),MedianIRFdsgevar(kk,j,i),VarIRFdsgevar(kk,j,i),... [MeanIRFdsgevar(kk,j,i),MedianIRFdsgevar(kk,j,i),VarIRFdsgevar(kk,j,i),...
HPDIRFdsgevar(kk,:,j,i),DistribIRFdsgevar(kk,:,j,i)] = ... HPDIRFdsgevar(kk,:,j,i),DistribIRFdsgevar(kk,:,j,i)] = ...
posterior_moments(squeeze(STOCK_IRF_BVARDSGE(k,j,i,:)),0,options_.mh_conf_sig); posterior_moments(squeeze(temp.STOCK_IRF_BVARDSGE(k,j,i,:)),0,options_.mh_conf_sig);
end end
end end
end end
kdx = kdx + size(STOCK_IRF_BVARDSGE,1); kdx = kdx + size(temp.STOCK_IRF_BVARDSGE,1);
end end
clear STOCK_IRF_BVARDSGE; clear temp;
for i = irf_shocks_indx for i = irf_shocks_indx
for j = 1:nvar for j = 1:nvar
name = sprintf('%s_%s', M_.endo_names{IndxVariables(j)}, tit{i}); name = sprintf('%s_%s', M_.endo_names{IndxVariables(j)}, tit{i});
@ -358,25 +342,24 @@ if MAX_nirfs_dsgevar
end end
end end
end end
%
% Finally I build the plots.
%
%% Finally I build the plots.
% Second block of code executed in parallel, with the exception of file % Second block of code executed in parallel, with the exception of file
% .tex generation always run in sequentially. This portion of code is execute in parallel by % .tex generation always run in sequentially. This portion of code is executed in parallel by
% PosteriorIRF_core2.m function. % PosteriorIRF_core2.m function.
if ~options_.nograph && ~options_.no_graph.posterior if ~options_.nograph && ~options_.no_graph.posterior
% Save the local variables. % Save the local variables.
localVars=[]; localVars=[];
localVars.dname=M_.dname;
localVars.fname=M_.fname;
localVars.options_=options_;
Check=options_.TeX; Check=options_.TeX;
if (Check) if (Check)
localVars.varlist_TeX=varlist_TeX; localVars.varlist_TeX=varlist_TeX;
end end
localVars.nvar=nvar; localVars.nvar=nvar;
localVars.MeanIRF=MeanIRF; localVars.MeanIRF=MeanIRF;
localVars.tit=tit; localVars.tit=tit;
@ -445,9 +428,7 @@ if ~options_.nograph && ~options_.no_graph.posterior
if isRemoteOctave if isRemoteOctave
[fout] = PosteriorIRF_core2(localVars,1,M_.exo_nbr,0); [fout] = PosteriorIRF_core2(localVars,1,M_.exo_nbr,0);
else else
globalVars = struct('M_',M_, ... globalVars = [];
'options_', options_);
[fout] = masterParallel(options_.parallel, 1, M_.exo_nbr,NamFileInput,'PosteriorIRF_core2', localVars, globalVars, options_.parallel_info); [fout] = masterParallel(options_.parallel, 1, M_.exo_nbr,NamFileInput,'PosteriorIRF_core2', localVars, globalVars, options_.parallel_info);
end end
end end
@ -455,7 +436,6 @@ if ~options_.nograph && ~options_.no_graph.posterior
[fout] = PosteriorIRF_core2(localVars,1,M_.exo_nbr,0); [fout] = PosteriorIRF_core2(localVars,1,M_.exo_nbr,0);
end end
% END parallel code! % END parallel code!
end end
fprintf('%s: Posterior IRFs, done!\n',dispString); fprintf('%s: Posterior IRFs, done!\n',dispString);

17
matlab/PosteriorIRF_core1.m
View File

@ -1,4 +1,5 @@
function myoutput=PosteriorIRF_core1(myinputs,fpar,B,whoiam, ThisMatlab) function myoutput=PosteriorIRF_core1(myinputs,fpar,B,whoiam, ThisMatlab)
%myoutput=PosteriorIRF_core1(myinputs,fpar,B,whoiam, ThisMatlab)
% Generates and stores Posterior IRFs % Generates and stores Posterior IRFs
% PARALLEL CONTEXT % PARALLEL CONTEXT
% This function perfoms in parallel execution a portion of the PosteriorIRF.m code. % This function perfoms in parallel execution a portion of the PosteriorIRF.m code.
@ -40,9 +41,6 @@ function myoutput=PosteriorIRF_core1(myinputs,fpar,B,whoiam, ThisMatlab)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>. % along with Dynare. If not, see <https://www.gnu.org/licenses/>.
global options_ estim_params_ oo_ M_ bayestopt_ dataset_ dataset_info
if nargin<4 if nargin<4
whoiam=0; whoiam=0;
end end
@ -50,6 +48,14 @@ end
% Reshape 'myinputs' for local computation. % Reshape 'myinputs' for local computation.
% In order to avoid confusion in the name space, the instruction struct2local(myinputs) is replaced by: % In order to avoid confusion in the name space, the instruction struct2local(myinputs) is replaced by:
options_= myinputs.options_;
estim_params_= myinputs.estim_params_;
M_= myinputs.M_;
oo_= myinputs.oo_;
bayestopt_= myinputs.bayestopt_;
dataset_= myinputs.dataset_;
dataset_info= myinputs.dataset_info;
IRUN = myinputs.IRUN; IRUN = myinputs.IRUN;
irun =myinputs.irun; irun =myinputs.irun;
irun2=myinputs.irun2; irun2=myinputs.irun2;
@ -78,13 +84,10 @@ if options_.dsge_var
bounds = prior_bounds(bayestopt_,options_.prior_trunc); bounds = prior_bounds(bayestopt_,options_.prior_trunc);
end end
if whoiam if whoiam
Parallel=myinputs.Parallel; Parallel=myinputs.Parallel;
end end
% MhDirectoryName = myinputs.MhDirectoryName;
if strcmpi(type,'posterior') if strcmpi(type,'posterior')
MhDirectoryName = CheckPath('metropolis',M_.dname); MhDirectoryName = CheckPath('metropolis',M_.dname);
elseif strcmpi(type,'gsa') elseif strcmpi(type,'gsa')
@ -115,12 +118,10 @@ else
h = dyn_waitbar(prct0,'Bayesian (prior) IRFs...'); h = dyn_waitbar(prct0,'Bayesian (prior) IRFs...');
end end
OutputFileName_bvardsge = {}; OutputFileName_bvardsge = {};
OutputFileName_dsge = {}; OutputFileName_dsge = {};
OutputFileName_param = {}; OutputFileName_param = {};
fpar = fpar-1; fpar = fpar-1;
fpar0=fpar; fpar0=fpar;
nosaddle=0; nosaddle=0;

17
matlab/PosteriorIRF_core2.m
View File

@ -1,5 +1,5 @@
function myoutput=PosteriorIRF_core2(myinputs,fpar,npar,whoiam,ThisMatlab) function myoutput=PosteriorIRF_core2(myinputs,fpar,npar,whoiam,ThisMatlab)
% function myoutput=PosteriorIRF_core2(myinputs,fpar,npar,whoiam, ThisMatlab) % myoutput=PosteriorIRF_core2(myinputs,fpar,npar,whoiam, ThisMatlab)
% Generates the Posterior IRFs plot from the IRFs generated in % Generates the Posterior IRFs plot from the IRFs generated in
% PosteriorIRF_core1 % PosteriorIRF_core1
% %
@ -47,8 +47,6 @@ function myoutput=PosteriorIRF_core2(myinputs,fpar,npar,whoiam,ThisMatlab)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>. % along with Dynare. If not, see <https://www.gnu.org/licenses/>.
global options_ M_
if nargin<4 if nargin<4
whoiam=0; whoiam=0;
end end
@ -56,6 +54,7 @@ end
% Reshape 'myinputs' for local computation. % Reshape 'myinputs' for local computation.
% In order to avoid confusion in the name space, the instruction struct2local(myinputs) is replaced by: % In order to avoid confusion in the name space, the instruction struct2local(myinputs) is replaced by:
options_=myinputs.options_;
Check=options_.TeX; Check=options_.TeX;
if (Check) if (Check)
varlist_TeX=myinputs.varlist_TeX; varlist_TeX=myinputs.varlist_TeX;
@ -67,6 +66,9 @@ tit=myinputs.tit;
nn=myinputs.nn; nn=myinputs.nn;
MAX_nirfs_dsgevar=myinputs.MAX_nirfs_dsgevar; MAX_nirfs_dsgevar=myinputs.MAX_nirfs_dsgevar;
HPDIRF=myinputs.HPDIRF; HPDIRF=myinputs.HPDIRF;
dname=myinputs.dname;
fname=myinputs.fname;
if options_.dsge_var if options_.dsge_var
HPDIRFdsgevar=myinputs.HPDIRFdsgevar; HPDIRFdsgevar=myinputs.HPDIRFdsgevar;
MeanIRFdsgevar=myinputs.MeanIRFdsgevar; MeanIRFdsgevar=myinputs.MeanIRFdsgevar;
@ -82,7 +84,7 @@ end
% To save the figures where the function is computed! % To save the figures where the function is computed!
DirectoryName = CheckPath('Output',M_.dname); DirectoryName = CheckPath('Output',dname);
RemoteFlag = 0; RemoteFlag = 0;
if whoiam if whoiam
@ -117,7 +119,6 @@ for i=fpar:npar
h2 = area(1:options_.irf,HPDIRF(:,1,j,i),'FaceColor',[1 1 1],'BaseValue',min(HPDIRF(:,1,j,i))); %white below HPDIinf and minimum of HPDIinf h2 = area(1:options_.irf,HPDIRF(:,1,j,i),'FaceColor',[1 1 1],'BaseValue',min(HPDIRF(:,1,j,i))); %white below HPDIinf and minimum of HPDIinf
end end
plot(1:options_.irf,MeanIRF(:,j,i),'-k','linewidth',3) plot(1:options_.irf,MeanIRF(:,j,i),'-k','linewidth',3)
% plot([1 options_.irf],[0 0],'-r','linewidth',0.5);
box on box on
axis tight axis tight
xlim([1 options_.irf]); xlim([1 options_.irf]);
@ -135,7 +136,6 @@ for i=fpar:npar
plot(1:options_.irf,HPDIRFdsgevar(:,1,j,i),'--k','linewidth',1) plot(1:options_.irf,HPDIRFdsgevar(:,1,j,i),'--k','linewidth',1)
plot(1:options_.irf,HPDIRFdsgevar(:,2,j,i),'--k','linewidth',1) plot(1:options_.irf,HPDIRFdsgevar(:,2,j,i),'--k','linewidth',1)
end end
% plot([1 options_.irf],[0 0],'-r','linewidth',0.5);
box on box on
axis tight axis tight
xlim([1 options_.irf]); xlim([1 options_.irf]);
@ -155,9 +155,9 @@ for i=fpar:npar
if subplotnum == MaxNumberOfPlotPerFigure || (j == nvar && subplotnum> 0) if subplotnum == MaxNumberOfPlotPerFigure || (j == nvar && subplotnum> 0)
figunumber = figunumber+1; figunumber = figunumber+1;
dyn_saveas(hh_fig,[DirectoryName '/' M_.fname '_Bayesian_IRF_' tit{i} '_' int2str(figunumber)],options_.nodisplay,options_.graph_format); dyn_saveas(hh_fig,[DirectoryName '/' fname '_Bayesian_IRF_' tit{i} '_' int2str(figunumber)],options_.nodisplay,options_.graph_format);
if RemoteFlag==1 if RemoteFlag==1
OutputFileName = [OutputFileName; {[DirectoryName,filesep], [M_.fname '_Bayesian_IRF_' deblank(tit(i,:)) '_' int2str(figunumber) '.*']}]; OutputFileName = [OutputFileName; {[DirectoryName,filesep], [fname '_Bayesian_IRF_' deblank(tit(i,:)) '_' int2str(figunumber) '.*']}];
end end
subplotnum = 0; subplotnum = 0;
end end
@ -165,7 +165,6 @@ for i=fpar:npar
if whoiam if whoiam
fprintf('Done! \n'); fprintf('Done! \n');
waitbarString = [ 'Exog. shocks ' int2str(i) '/' int2str(npar) ' done.']; waitbarString = [ 'Exog. shocks ' int2str(i) '/' int2str(npar) ' done.'];
% fMessageStatus((i-fpar+1)/(npar-fpar+1),whoiam,waitbarString, waitbarTitle, Parallel(ThisMatlab));
dyn_waitbar((i-fpar+1)/(npar-fpar+1),[],waitbarString); dyn_waitbar((i-fpar+1)/(npar-fpar+1),[],waitbarString);
end end
end % loop over exo_var end % loop over exo_var

83
matlab/ReshapeMatFiles.m
View File

@ -1,10 +1,15 @@
function ReshapeMatFiles(type, type2) function ReshapeMatFiles(fname, dname, exo_nbr, endo_nbr, options_, type, type2)
% function ReshapeMatFiles(type, type2) % function ReshapeMatFiles(fname, dname, exo_nbr, endo_nbr, options_, type, type2)
% Reshapes and sorts (along the mcmc simulations) the mat files generated by DYNARE. % Reshapes and sorts (along the mcmc simulations) the mat files generated by DYNARE.
% 4D-arrays are splitted along the first dimension. % 4D-arrays are splitted along the first dimension.
% 3D-arrays are splitted along the second dimension. % 3D-arrays are splitted along the second dimension.
% %
% INPUTS: % INPUTS:
% fname: [string] filename
% dname: [string] directory name
% exo_nbr: [integer] number of exogenous variables
% endo_nbr: [integer] number of endogenous variables
% options_: [struct] options structure
% type: statistics type in the repertory: % type: statistics type in the repertory:
% dgse % dgse
% irf_bvardsge % irf_bvardsge
@ -25,7 +30,7 @@ function ReshapeMatFiles(type, type2)
% SPECIAL REQUIREMENTS % SPECIAL REQUIREMENTS
% none % none
% Copyright © 2003-2017 Dynare Team % Copyright © 2003-2023 Dynare Team
% %
% This file is part of Dynare. % This file is part of Dynare.
% %
@ -42,43 +47,41 @@ function ReshapeMatFiles(type, type2)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>. % along with Dynare. If not, see <https://www.gnu.org/licenses/>.
global M_ options_ if nargin==6
MhDirectoryName = [ CheckPath('metropolis',dname) filesep ];
if nargin==1
MhDirectoryName = [ CheckPath('metropolis',M_.dname) filesep ];
else else
if strcmpi(type2,'posterior') if strcmpi(type2,'posterior')
MhDirectoryName = [CheckPath('metropolis',M_.dname) filesep ]; MhDirectoryName = [CheckPath('metropolis',dname) filesep ];
elseif strcmpi(type2,'gsa') elseif strcmpi(type2,'gsa')
if options_.opt_gsa.morris==1 if options_.opt_gsa.morris==1
MhDirectoryName = [CheckPath('gsa/screen',M_.dname) filesep ]; MhDirectoryName = [CheckPath('gsa/screen',dname) filesep ];
elseif options_.opt_gsa.morris==2 elseif options_.opt_gsa.morris==2
MhDirectoryName = [CheckPath('gsa/identif',M_.dname) filesep ]; MhDirectoryName = [CheckPath('gsa/identif',dname) filesep ];
elseif options_.opt_gsa.pprior elseif options_.opt_gsa.pprior
MhDirectoryName = [CheckPath(['gsa' filesep 'prior'],M_.dname) filesep ]; MhDirectoryName = [CheckPath(['gsa' filesep 'prior'],dname) filesep ];
else else
MhDirectoryName = [CheckPath(['gsa' filesep 'mc'],M_.dname) filesep ]; MhDirectoryName = [CheckPath(['gsa' filesep 'mc'],dname) filesep ];
end end
else else
MhDirectoryName = [CheckPath('prior',M_.dname) filesep ]; MhDirectoryName = [CheckPath('prior',dname) filesep ];
end end
end end
switch type switch type
case 'irf_dsge' case 'irf_dsge'
CAPtype = 'IRF_DSGE'; CAPtype = 'IRF_DSGE';
TYPEsize = [ options_.irf , size(options_.varlist,1) , M_.exo_nbr ]; TYPEsize = [ options_.irf , size(options_.varlist,1) , exo_nbr ];
TYPEarray = 4; TYPEarray = 4;
case 'irf_bvardsge' case 'irf_bvardsge'
CAPtype = 'IRF_BVARDSGE'; CAPtype = 'IRF_BVARDSGE';
TYPEsize = [ options_.irf , length(options_.varobs) , M_.exo_nbr ]; TYPEsize = [ options_.irf , length(options_.varobs) , exo_nbr ];
TYPEarray = 4; TYPEarray = 4;
case 'smooth' case 'smooth'
CAPtype = 'SMOOTH'; CAPtype = 'SMOOTH';
TYPEsize = [ M_.endo_nbr , options_.nobs ]; TYPEsize = [ endo_nbr , options_.nobs ];
TYPEarray = 3; TYPEarray = 3;
case 'filter' case 'filter'
CAPtype = 'FILTER'; CAPtype = 'FILTER';
TYPEsize = [ M_.endo_nbr , options_.nobs + 1 ];% TO BE CHECKED! TYPEsize = [ endo_nbr , options_.nobs + 1 ];% TO BE CHECKED!
TYPEarray = 3; TYPEarray = 3;
case 'error' case 'error'
CAPtype = 'ERROR'; CAPtype = 'ERROR';
@ -86,22 +89,22 @@ switch type
TYPEarray = 3; TYPEarray = 3;
case 'innov' case 'innov'
CAPtype = 'INNOV'; CAPtype = 'INNOV';
TYPEsize = [ M_.exo_nbr , options_.nobs ]; TYPEsize = [ exo_nbr , options_.nobs ];
TYPEarray = 3; TYPEarray = 3;
case 'forcst' case 'forcst'
CAPtype = 'FORCST'; CAPtype = 'FORCST';
TYPEsize = [ M_.endo_nbr , options_.forecast ]; TYPEsize = [ endo_nbr , options_.forecast ];
TYPEarray = 3; TYPEarray = 3;
case 'forcst1' case 'forcst1'
CAPtype = 'FORCST1'; CAPtype = 'FORCST1';
TYPEsize = [ M_.endo_nbr , options_.forecast ]; TYPEsize = [ endo_nbr , options_.forecast ];
TYPEarray = 3; TYPEarray = 3;
otherwise otherwise
disp('ReshapeMatFiles :: Unknown argument!') disp('ReshapeMatFiles :: Unknown argument!')
return return
end end
TYPEfiles = dir([MhDirectoryName M_.fname '_' type '*.mat']); TYPEfiles = dir([MhDirectoryName fname '_' type '*.mat']);
NumberOfTYPEfiles = length(TYPEfiles); NumberOfTYPEfiles = length(TYPEfiles);
B = options_.B; B = options_.B;
@ -116,36 +119,29 @@ switch TYPEarray
for f1=1:NumberOfTYPEfiles-foffset for f1=1:NumberOfTYPEfiles-foffset
eval(['STOCK_' CAPtype ' = zeros(NumberOfPeriodsPerTYPEfiles,TYPEsize(2),TYPEsize(3),B);']) eval(['STOCK_' CAPtype ' = zeros(NumberOfPeriodsPerTYPEfiles,TYPEsize(2),TYPEsize(3),B);'])
for f2 = 1:NumberOfTYPEfiles for f2 = 1:NumberOfTYPEfiles
load([MhDirectoryName M_.fname '_' type int2str(f2) '.mat']); load([MhDirectoryName fname '_' type int2str(f2) '.mat']);
eval(['STOCK_' CAPtype '(:,:,1:+size(stock_' type ',3),idx+1:idx+size(stock_' type ',4))=stock_' ... eval(['STOCK_' CAPtype '(:,:,1:+size(stock_' type ',3),idx+1:idx+size(stock_' type ',4))=stock_' ...
type '(jdx+1:jdx+NumberOfPeriodsPerTYPEfiles,:,:,:);']) type '(jdx+1:jdx+NumberOfPeriodsPerTYPEfiles,:,:,:);'])
eval(['idx = idx + size(stock_' type ',4);']) eval(['idx = idx + size(stock_' type ',4);'])
end end
%eval(['STOCK_' CAPtype ' = sort(STOCK_' CAPtype ',4);']) save([MhDirectoryName fname '_' CAPtype 's' int2str(f1) '.mat'],['STOCK_' CAPtype]);
save([MhDirectoryName M_.fname '_' CAPtype 's' int2str(f1) '.mat'],['STOCK_' CAPtype]);
jdx = jdx + NumberOfPeriodsPerTYPEfiles; jdx = jdx + NumberOfPeriodsPerTYPEfiles;
idx = 0; idx = 0;
end end
if reste if reste
eval(['STOCK_' CAPtype ' = zeros(reste,TYPEsize(2),TYPEsize(3),B);']) eval(['STOCK_' CAPtype ' = zeros(reste,TYPEsize(2),TYPEsize(3),B);'])
for f2 = 1:NumberOfTYPEfiles for f2 = 1:NumberOfTYPEfiles
load([MhDirectoryName M_.fname '_' type int2str(f2) '.mat']); load([MhDirectoryName fname '_' type int2str(f2) '.mat']);
eval(['STOCK_' CAPtype '(:,:,:,idx+1:idx+size(stock_' type ',4))=stock_' type '(jdx+1:jdx+reste,:,:,:);']) eval(['STOCK_' CAPtype '(:,:,:,idx+1:idx+size(stock_' type ',4))=stock_' type '(jdx+1:jdx+reste,:,:,:);'])
eval(['idx = idx + size(stock_' type ',4);']) eval(['idx = idx + size(stock_' type ',4);'])
end end
%eval(['STOCK_' CAPtype ' = sort(STOCK_' CAPtype ',4);']) save([MhDirectoryName fname '_' CAPtype 's' int2str(NumberOfTYPEfiles-foffset+1) '.mat'],['STOCK_' CAPtype]);
save([MhDirectoryName M_.fname '_' CAPtype 's' int2str(NumberOfTYPEfiles-foffset+1) '.mat'],['STOCK_' CAPtype]);
end end
else else
load([MhDirectoryName M_.fname '_' type '1.mat']); load([MhDirectoryName fname '_' type '1.mat']);
%eval(['STOCK_' CAPtype ' = sort(stock_' type ',4);'])
eval(['STOCK_' CAPtype ' = stock_' type ';']) eval(['STOCK_' CAPtype ' = stock_' type ';'])
save([MhDirectoryName M_.fname '_' CAPtype 's' int2str(1) '.mat'],['STOCK_' CAPtype ]); save([MhDirectoryName fname '_' CAPtype 's' int2str(1) '.mat'],['STOCK_' CAPtype ]);
end end
% Original file format may be useful in some cases...
% for file = 1:NumberOfTYPEfiles
% delete([MhDirectoryName M_.fname '_' type int2str(file) '.mat'])
% end
case 3 case 3
if NumberOfTYPEfiles>1 if NumberOfTYPEfiles>1
NumberOfPeriodsPerTYPEfiles = ceil( TYPEsize(2)/NumberOfTYPEfiles ); NumberOfPeriodsPerTYPEfiles = ceil( TYPEsize(2)/NumberOfTYPEfiles );
@ -155,31 +151,24 @@ switch TYPEarray
for f1=1:NumberOfTYPEfiles-1 for f1=1:NumberOfTYPEfiles-1
eval(['STOCK_' CAPtype ' = zeros(TYPEsize(1),NumberOfPeriodsPerTYPEfiles,B);']) eval(['STOCK_' CAPtype ' = zeros(TYPEsize(1),NumberOfPeriodsPerTYPEfiles,B);'])
for f2 = 1:NumberOfTYPEfiles for f2 = 1:NumberOfTYPEfiles
load([MhDirectoryName M_.fname '_' type int2str(f2) '.mat']); load([MhDirectoryName fname '_' type int2str(f2) '.mat']);
eval(['STOCK_' CAPtype '(:,:,idx+1:idx+size(stock_ ' type ',3))=stock_' type '(:,jdx+1:jdx+NumberOfPeriodsPerTYPEfiles,:);']) eval(['STOCK_' CAPtype '(:,:,idx+1:idx+size(stock_ ' type ',3))=stock_' type '(:,jdx+1:jdx+NumberOfPeriodsPerTYPEfiles,:);'])
eval(['idx = idx + size(stock_' type ',3);']) eval(['idx = idx + size(stock_' type ',3);'])
end end
%eval(['STOCK_' CAPtype ' = sort(STOCK_' CAPtype ',3);']) save([MhDirectoryName fname '_' CAPtype 's' int2str(f1) '.mat'],['STOCK_' CAPtype]);
save([MhDirectoryName M_.fname '_' CAPtype 's' int2str(f1) '.mat'],['STOCK_' CAPtype]);
jdx = jdx + NumberOfPeriodsPerTYPEfiles; jdx = jdx + NumberOfPeriodsPerTYPEfiles;
idx = 0; idx = 0;
end end
eval(['STOCK_' CAPtype ' = zeros(TYPEsize(1),reste,B);']) eval(['STOCK_' CAPtype ' = zeros(TYPEsize(1),reste,B);'])
for f2 = 1:NumberOfTYPEfiles for f2 = 1:NumberOfTYPEfiles
load([MhDirectoryName M_.fname '_' type int2str(f2) '.mat']); load([MhDirectoryName fname '_' type int2str(f2) '.mat']);
eval(['STOCK_' CAPtype '(:,:,idx+1:idx+size(stock_' type ',3))=stock_' type '(:,jdx+1:jdx+reste,:);']) eval(['STOCK_' CAPtype '(:,:,idx+1:idx+size(stock_' type ',3))=stock_' type '(:,jdx+1:jdx+reste,:);'])
eval(['idx = idx + size(stock_' type ',3);']) eval(['idx = idx + size(stock_' type ',3);'])
end end
%eval(['STOCK_' CAPtype ' = sort(STOCK_' CAPtype ',3);']) save([MhDirectoryName fname '_' CAPtype 's' int2str(NumberOfTYPEfiles) '.mat'],['STOCK_' CAPtype]);
save([MhDirectoryName M_.fname '_' CAPtype 's' int2str(NumberOfTYPEfiles) '.mat'],['STOCK_' CAPtype]);
else else
load([MhDirectoryName M_.fname '_' type '1.mat']); load([MhDirectoryName fname '_' type '1.mat']);
%eval(['STOCK_' CAPtype ' = sort(stock_' type ',3);'])
eval(['STOCK_' CAPtype ' = stock_' type ';']) eval(['STOCK_' CAPtype ' = stock_' type ';'])
save([MhDirectoryName M_.fname '_' CAPtype 's' int2str(1) '.mat'],['STOCK_' CAPtype ]); save([MhDirectoryName fname '_' CAPtype 's' int2str(1) '.mat'],['STOCK_' CAPtype ]);
end end
% Original file format may be useful in some cases...
% for file = 1:NumberOfTYPEfiles
% delete([MhDirectoryName M_.fname '_' type int2str(file) '.mat'])
% end
end end

112
matlab/dynare_estimation_1.m
View File

@ -225,70 +225,70 @@ if ~isequal(options_.mode_compute,0) && ~options_.mh_posterior_mode_estimation
current_optimizer = optimizer_vec{optim_iter}; current_optimizer = optimizer_vec{optim_iter};
[xparam1, fval, ~, hh, options_, Scale, new_rat_hess_info] = dynare_minimize_objective(objective_function,xparam1,current_optimizer,options_,[bounds.lb bounds.ub],bayestopt_.name,bayestopt_,hh,dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,bounds,oo_); [xparam1, fval, ~, hh, options_, Scale, new_rat_hess_info] = dynare_minimize_objective(objective_function,xparam1,current_optimizer,options_,[bounds.lb bounds.ub],bayestopt_.name,bayestopt_,hh,dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,bounds,oo_);
fprintf('\nFinal value of minus the log posterior (or likelihood):%f \n', fval); fprintf('\nFinal value of minus the log posterior (or likelihood):%f \n', fval);
if isnumeric(current_optimizer) if isnumeric(current_optimizer)
if current_optimizer==5 if current_optimizer==5
newratflag = new_rat_hess_info.newratflag; newratflag = new_rat_hess_info.newratflag;
new_rat_hess_info = new_rat_hess_info.new_rat_hess_info; new_rat_hess_info = new_rat_hess_info.new_rat_hess_info;
elseif current_optimizer==6 %save scaling factor elseif current_optimizer==6 %save scaling factor
save([M_.dname filesep 'Output' filesep M_.fname '_optimal_mh_scale_parameter.mat'],'Scale'); save([M_.dname filesep 'Output' filesep M_.fname '_optimal_mh_scale_parameter.mat'],'Scale');
options_.mh_jscale = Scale; options_.mh_jscale = Scale;
bayestopt_.jscale(:) = options_.mh_jscale; bayestopt_.jscale(:) = options_.mh_jscale;
end end
end end
if ~isnumeric(current_optimizer) || ~isequal(current_optimizer,6) %always already computes covariance matrix if ~isnumeric(current_optimizer) || ~isequal(current_optimizer,6) %always already computes covariance matrix
if options_.cova_compute == 1 %user did not request covariance not to be computed if options_.cova_compute == 1 %user did not request covariance not to be computed
if options_.analytic_derivation && strcmp(func2str(objective_function),'dsge_likelihood') if options_.analytic_derivation && strcmp(func2str(objective_function),'dsge_likelihood')
ana_deriv_old = options_.analytic_derivation; ana_deriv_old = options_.analytic_derivation;
options_.analytic_derivation = 2; options_.analytic_derivation = 2;
[~,~,~,~,hh] = feval(objective_function,xparam1, ... [~,~,~,~,hh] = feval(objective_function,xparam1, ...
dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,bounds,oo_); dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,bounds,oo_);
options_.analytic_derivation = ana_deriv_old; options_.analytic_derivation = ana_deriv_old;
elseif ~isnumeric(current_optimizer) || ~(isequal(current_optimizer,5) && newratflag~=1 && strcmp(func2str(objective_function),'dsge_likelihood')) elseif ~isnumeric(current_optimizer) || ~(isequal(current_optimizer,5) && newratflag~=1 && strcmp(func2str(objective_function),'dsge_likelihood'))
% enter here if i) not mode_compute_5, ii) if mode_compute_5 and newratflag==1; % enter here if i) not mode_compute_5, ii) if mode_compute_5 and newratflag==1;
% with flag==0 or 2 and dsge_likelihood, we force to use % with flag==0 or 2 and dsge_likelihood, we force to use
% the hessian from outer product gradient of optimizer 5 below % the hessian from outer product gradient of optimizer 5 below
if options_.hessian.use_penalized_objective if options_.hessian.use_penalized_objective
penalized_objective_function = str2func('penalty_objective_function'); penalized_objective_function = str2func('penalty_objective_function');
hh = hessian(penalized_objective_function, xparam1, options_.gstep, objective_function, fval, dataset_, dataset_info, options_, M_, estim_params_, bayestopt_, bounds,oo_); hh = hessian(penalized_objective_function, xparam1, options_.gstep, objective_function, fval, dataset_, dataset_info, options_, M_, estim_params_, bayestopt_, bounds,oo_);
else else
hh = hessian(objective_function, xparam1, options_.gstep, dataset_, dataset_info, options_, M_, estim_params_, bayestopt_, bounds,oo_); hh = hessian(objective_function, xparam1, options_.gstep, dataset_, dataset_info, options_, M_, estim_params_, bayestopt_, bounds,oo_);
end
hh = reshape(hh, nx, nx);
elseif isnumeric(current_optimizer) && isequal(current_optimizer,5)
% other numerical hessian options available with optimizer
% 5 and dsge_likelihood
%
% if newratflag == 0
% compute outer product gradient of optimizer 5
%
% if newratflag == 2
% compute 'mixed' outer product gradient of optimizer 5
% with diagonal elements computed with numerical second order derivatives
%
% uses univariate filters, so to get max # of available
% densities for outer product gradient
kalman_algo0 = options_.kalman_algo;
compute_hessian = 1;
if ~((options_.kalman_algo == 2) || (options_.kalman_algo == 4))
options_.kalman_algo=2;
if options_.lik_init == 3
options_.kalman_algo=4;
end
elseif newratflag==0 % hh already contains outer product gradient with univariate filter
compute_hessian = 0;
end
if compute_hessian
crit = options_.newrat.tolerance.f;
newratflag = newratflag>0;
hh = reshape(mr_hessian(xparam1,objective_function,fval,newratflag,crit,new_rat_hess_info,[bounds.lb bounds.ub],bayestopt_.p2,0,dataset_, dataset_info, options_,M_,estim_params_,bayestopt_,bounds,oo_), nx, nx);
end
options_.kalman_algo = kalman_algo0;
end end
hh = reshape(hh, nx, nx);
elseif isnumeric(current_optimizer) && isequal(current_optimizer,5)
% other numerical hessian options available with optimizer
% 5 and dsge_likelihood
%
% if newratflag == 0
% compute outer product gradient of optimizer 5
%
% if newratflag == 2
% compute 'mixed' outer product gradient of optimizer 5
% with diagonal elements computed with numerical second order derivatives
%
% uses univariate filters, so to get max # of available
% densities for outer product gradient
kalman_algo0 = options_.kalman_algo;
compute_hessian = 1;
if ~((options_.kalman_algo == 2) || (options_.kalman_algo == 4))
options_.kalman_algo=2;
if options_.lik_init == 3
options_.kalman_algo=4;
end
elseif newratflag==0 % hh already contains outer product gradient with univariate filter
compute_hessian = 0;
end
if compute_hessian
crit = options_.newrat.tolerance.f;
newratflag = newratflag>0;
hh = reshape(mr_hessian(xparam1,objective_function,fval,newratflag,crit,new_rat_hess_info,[bounds.lb bounds.ub],bayestopt_.p2,0,dataset_, dataset_info, options_,M_,estim_params_,bayestopt_,bounds,oo_), nx, nx);
end
options_.kalman_algo = kalman_algo0;
end end
end end
parameter_names = bayestopt_.name; end
parameter_names = bayestopt_.name;
end end
if options_.cova_compute || current_optimizer==5 || current_optimizer==6 if options_.cova_compute || current_optimizer==5 || current_optimizer==6
save([M_.dname filesep 'Output' filesep M_.fname '_mode.mat'],'xparam1','hh','parameter_names','fval'); save([M_.dname filesep 'Output' filesep M_.fname '_mode.mat'],'xparam1','hh','parameter_names','fval');
@ -472,7 +472,7 @@ if (any(bayestopt_.pshape >0 ) && options_.mh_replic) || ...
if error_flag if error_flag
error('%s: I cannot compute the posterior IRFs!',dispString) error('%s: I cannot compute the posterior IRFs!',dispString)
end end
PosteriorIRF('posterior',dispString); oo_=PosteriorIRF('posterior',options_,estim_params_,oo_,M_,bayestopt_,dataset_,dataset_info,dispString);
end end
if options_.moments_varendo if options_.moments_varendo
if error_flag if error_flag
@ -504,7 +504,7 @@ if (any(bayestopt_.pshape >0 ) && options_.mh_replic) || ...
error('%s: I cannot compute the posterior statistics!',dispString) error('%s: I cannot compute the posterior statistics!',dispString)
end end
if options_.order==1 && ~options_.particle.status if options_.order==1 && ~options_.particle.status
prior_posterior_statistics('posterior',dataset_,dataset_info,dispString); %get smoothed and filtered objects and forecasts oo_=prior_posterior_statistics('posterior',dataset_,dataset_info,M_,oo_,options_,estim_params_,bayestopt_,dispString); %get smoothed and filtered objects and forecasts
else else
error('%s: Particle Smoothers are not yet implemented.',dispString) error('%s: Particle Smoothers are not yet implemented.',dispString)
end end

4
matlab/dynare_sensitivity.m
View File

@ -410,9 +410,9 @@ if options_gsa.rmse
end end
end end
prior_posterior_statistics('gsa',dataset_, dataset_info,'gsa::mcmc'); oo_=prior_posterior_statistics('gsa',dataset_, dataset_info,M_,oo_,options_,estim_params_,bayestopt_,'gsa::mcmc');
if options_.bayesian_irf if options_.bayesian_irf
PosteriorIRF('gsa','gsa::mcmc'); oo_=PosteriorIRF('gsa',options_,estim_params_,oo_,M_,bayestopt_,dataset_,dataset_info,'gsa::mcmc');
end end
options_gsa.load_rmse=0; options_gsa.load_rmse=0;
% else % else

32
matlab/pm3.m
View File

@ -1,4 +1,4 @@
function pm3(n1,n2,ifil,B,tit1,tit2,tit3,tit_tex,names1,names2,name3,DirectoryName,var_type,dispString) function oo_=pm3(M_,options_,oo_,n1,n2,ifil,B,tit1,tit2,tit_tex,names1,names2,name3,DirectoryName,var_type,dispString)
% Computes, stores and plots the posterior moment statistics. % Computes, stores and plots the posterior moment statistics.
% %
@ -8,8 +8,7 @@ function pm3(n1,n2,ifil,B,tit1,tit2,tit3,tit_tex,names1,names2,name3,DirectoryNa
% ifil [scalar] number of moment files to load % ifil [scalar] number of moment files to load
% B [scalar] number of subdraws % B [scalar] number of subdraws
% tit1 [string] Figure title % tit1 [string] Figure title
% tit2 [string] not used % tit2 [string] Save name for figure
% tit3 [string] Save name for figure
% tit_tex [cell array] TeX-Names for Variables % tit_tex [cell array] TeX-Names for Variables
% names1 [cell array] Names of all variables in the moment matrix from % names1 [cell array] Names of all variables in the moment matrix from
% which names2 is selected % which names2 is selected
@ -20,6 +19,10 @@ function pm3(n1,n2,ifil,B,tit1,tit2,tit3,tit_tex,names1,names2,name3,DirectoryNa
% var_type [string] suffix of the filename from which to load moment % var_type [string] suffix of the filename from which to load moment
% matrix % matrix
% dispString [string] string to be displayes in the command window % dispString [string] string to be displayes in the command window
%
% OUTPUTS
% oo_ [structure] storing the results
% PARALLEL CONTEXT % PARALLEL CONTEXT
% See also the comment in posterior_sampler.m funtion. % See also the comment in posterior_sampler.m funtion.
@ -42,8 +45,6 @@ function pm3(n1,n2,ifil,B,tit1,tit2,tit3,tit_tex,names1,names2,name3,DirectoryNa
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>. % along with Dynare. If not, see <https://www.gnu.org/licenses/>.
global options_ M_ oo_
nn = 3; nn = 3;
MaxNumberOfPlotsPerFigure = nn^2; % must be square MaxNumberOfPlotsPerFigure = nn^2; % must be square
varlist = names2; varlist = names2;
@ -284,9 +285,7 @@ if strcmp(var_type,'_trend_coeff') || max(max(abs(Mean(:,:))))<=10^(-6) || all(a
fprintf(['%s: ' tit1 ', done!\n'],dispString); fprintf(['%s: ' tit1 ', done!\n'],dispString);
return %not do plots return %not do plots
end end
%%
%% Finally I build the plots. %% Finally I build the plots.
%%
if ~options_.nograph && ~options_.no_graph.posterior if ~options_.nograph && ~options_.no_graph.posterior
% Block of code executed in parallel, with the exception of file % Block of code executed in parallel, with the exception of file
@ -297,8 +296,6 @@ if ~options_.nograph && ~options_.no_graph.posterior
% %
% %%% The file .TeX! are not saved in parallel. % %%% The file .TeX! are not saved in parallel.
% Store the variable mandatory for local/remote parallel computing. % Store the variable mandatory for local/remote parallel computing.
localVars=[]; localVars=[];
@ -313,8 +310,14 @@ if ~options_.nograph && ~options_.no_graph.posterior
end end
localVars.MaxNumberOfPlotsPerFigure=MaxNumberOfPlotsPerFigure; localVars.MaxNumberOfPlotsPerFigure=MaxNumberOfPlotsPerFigure;
localVars.name3=name3; localVars.name3=name3;
localVars.tit3=tit3; localVars.tit2=tit2;
localVars.Mean=Mean; localVars.Mean=Mean;
localVars.TeX=options_.TeX;
localVars.nodisplay=options_.nodisplay;
localVars.graph_format=options_.graph_format;
localVars.dname=M_.dname;
localVars.fname=M_.fname;
% Like sequential execution! % Like sequential execution!
nvar0=nvar; nvar0=nvar;
@ -332,16 +335,13 @@ if ~options_.nograph && ~options_.no_graph.posterior
if isRemoteOctave if isRemoteOctave
fout = pm3_core(localVars,1,nvar,0); fout = pm3_core(localVars,1,nvar,0);
else else
globalVars = struct('M_',M_, ... globalVars = [];
'options_', options_, ...
'oo_', oo_);
[fout, nvar0, totCPU] = masterParallel(options_.parallel, 1, nvar, [],'pm3_core', localVars,globalVars, options_.parallel_info); [fout, nvar0, totCPU] = masterParallel(options_.parallel, 1, nvar, [],'pm3_core', localVars,globalVars, options_.parallel_info);
end end
end end
else else
% For the time being in Octave enviroment the pm3.m is executed only in % For the time being in Octave enviroment the pm3.m is executed only in
% serial modality, to avoid problem with the plots. % serial modality, to avoid problem with the plots.
fout = pm3_core(localVars,1,nvar,0); fout = pm3_core(localVars,1,nvar,0);
end end
@ -365,8 +365,8 @@ if ~options_.nograph && ~options_.no_graph.posterior
if subplotnum == MaxNumberOfPlotsPerFigure || i == nvar if subplotnum == MaxNumberOfPlotsPerFigure || i == nvar
fprintf(fidTeX,'\\begin{figure}[H]\n'); fprintf(fidTeX,'\\begin{figure}[H]\n');
fprintf(fidTeX,'\\centering \n'); fprintf(fidTeX,'\\centering \n');
fprintf(fidTeX,['\\includegraphics[width=%2.2f\\textwidth]{%s/Output/%s_' name3 '_%s}\n'],options_.figures.textwidth*min(subplotnum/nn,1),M_.dname,M_.fname, tit3{i}); fprintf(fidTeX,['\\includegraphics[width=%2.2f\\textwidth]{%s/Output/%s_' name3 '_%s}\n'],options_.figures.textwidth*min(subplotnum/nn,1),M_.dname,M_.fname, tit2{i});
fprintf(fidTeX,'\\label{Fig:%s:%s}\n',name3,tit3{i}); fprintf(fidTeX,'\\label{Fig:%s:%s}\n',name3,tit2{i});
fprintf(fidTeX,'\\caption{%s}\n',tit1); fprintf(fidTeX,'\\caption{%s}\n',tit1);
fprintf(fidTeX,'\\end{figure}\n'); fprintf(fidTeX,'\\end{figure}\n');
fprintf(fidTeX,' \n'); fprintf(fidTeX,' \n');

40
matlab/pm3_core.m
View File

@ -1,13 +1,23 @@
function myoutput=pm3_core(myinputs,fpar,nvar,whoiam, ThisMatlab) function myoutput=pm3_core(myinputs,fpar,nvar,whoiam, ThisMatlab)
% myoutput=pm3_core(myinputs,fpar,nvar,whoiam, ThisMatlab)
% PARALLEL CONTEXT % PARALLEL CONTEXT
% Core functionality for pm3.m function, which can be parallelized. % Core functionality for pm3.m function, which can be parallelized.
% INPUTS % INPUTS
% See the comment in posterior_sampler_core.m funtion. % o myimput [struc] The mandatory variables for local/remote
% parallel computing obtained from prior_posterior_statistics.m
% function.
% o fpar and nvar [integer] first variable and number of variables
% o whoiam [integer] In concurrent programming a modality to refer to the different threads running in parallel is needed.
% The integer whoaim is the integer that
% allows us to distinguish between them. Then it is the index number of this CPU among all CPUs in the
% cluster.
% o ThisMatlab [integer] Allows us to distinguish between the
% 'main' Matlab, the slave Matlab worker, local Matlab, remote Matlab,
% ... Then it is the index number of this slave machine in the cluster.
%
% OUTPUTS % OUTPUTS
% o myoutput [struc] % o myoutput [struct] Contains file names
% %
% %
% SPECIAL REQUIREMENTS. % SPECIAL REQUIREMENTS.
@ -45,16 +55,18 @@ varlist=myinputs.varlist;
MaxNumberOfPlotsPerFigure=myinputs.MaxNumberOfPlotsPerFigure; MaxNumberOfPlotsPerFigure=myinputs.MaxNumberOfPlotsPerFigure;
name3=myinputs.name3; name3=myinputs.name3;
tit3=myinputs.tit3; tit2=myinputs.tit2;
Mean=myinputs.Mean; Mean=myinputs.Mean;
options_.TeX=myinputs.TeX;
options_.nodisplay=myinputs.nodisplay;
options_.graph_format=myinputs.graph_format;
fname=myinputs.fname;
dname=myinputs.dname;
if whoiam if whoiam
Parallel=myinputs.Parallel; Parallel=myinputs.Parallel;
end end
global options_ M_ oo_
if options_.TeX if options_.TeX
varlist_TeX=myinputs.varlist_TeX; varlist_TeX=myinputs.varlist_TeX;
end end
@ -64,8 +76,6 @@ if whoiam
h = dyn_waitbar(prct0,'Parallel plots pm3 ...'); h = dyn_waitbar(prct0,'Parallel plots pm3 ...');
end end
figunumber = 0; figunumber = 0;
subplotnum = 0; subplotnum = 0;
hh_fig = dyn_figure(options_.nodisplay,'Name',[tit1 ' ' int2str(figunumber+1)]); hh_fig = dyn_figure(options_.nodisplay,'Name',[tit1 ' ' int2str(figunumber+1)]);
@ -110,14 +120,14 @@ for i=fpar:nvar
if whoiam if whoiam
if Parallel(ThisMatlab).Local==0 if Parallel(ThisMatlab).Local==0
DirectoryName = CheckPath('Output',M_.dname); DirectoryName = CheckPath('Output',dname);
end end
end end
if subplotnum == MaxNumberOfPlotsPerFigure || i == nvar if subplotnum == MaxNumberOfPlotsPerFigure || i == nvar
dyn_saveas(hh_fig,[M_.dname '/Output/' M_.fname '_' name3 '_' tit3{i}],options_.nodisplay,options_.graph_format); dyn_saveas(hh_fig,[dname '/Output/' fname '_' name3 '_' tit2{i}],options_.nodisplay,options_.graph_format);
if RemoteFlag==1 if RemoteFlag==1
OutputFileName = [OutputFileName; {[M_.dname, filesep, 'Output',filesep], [M_.fname '_' name3 '_' deblank(tit3(i,:)) '.*']}]; OutputFileName = [OutputFileName; {[dname, filesep, 'Output',filesep], [fname '_' name3 '_' deblank(tit2(i,:)) '.*']}];
end end
subplotnum = 0; subplotnum = 0;
figunumber = figunumber+1; figunumber = figunumber+1;
@ -127,12 +137,8 @@ for i=fpar:nvar
end end
if whoiam if whoiam
% waitbarString = [ 'Variable ' int2str(i) '/' int2str(nvar) ' done.'];
% fMessageStatus((i-fpar+1)/(nvar-fpar+1),whoiam,waitbarString, waitbarTitle, Parallel(ThisMatlab));
dyn_waitbar((i-fpar+1)/(nvar-fpar+1),h); dyn_waitbar((i-fpar+1)/(nvar-fpar+1),h);
end end
end end
if whoiam if whoiam

129
matlab/prior_posterior_statistics.m
View File

@ -1,17 +1,19 @@
function prior_posterior_statistics(type,dataset,dataset_info,dispString) function oo_=prior_posterior_statistics(type,dataset_,dataset_info,M_,oo_,options_,estim_params_,bayestopt_,dispString)
% function prior_posterior_statistics(type,dataset,dataset_info,dispString) % oo_=prior_posterior_statistics(type,dataset_,dataset_info,M_,oo_,options_,estim_params_,bayestopt_,dispString))
% Computes Monte Carlo filter smoother and forecasts % Computes Monte Carlo filter smoother and forecasts
% %
% INPUTS % INPUTS
% type: posterior % type [string] posterior, prior, or gsa
% prior % o dataset_ [structure] storing the dataset
% gsa % o dataset_info [structure] Various information about the dataset
% dataset: data structure % o M_ [structure] storing the model information
% dataset_info: dataset structure % o oo_ [structure] storing the results
% dispString: string to display in the command window % o options_ [structure] storing the options
% % o estim_params_ [structure] storing information about estimated parameters
% o bayestopt_ [structure] storing information about priors
% dispString: [string] display info in the command window
% OUTPUTS % OUTPUTS
% none % oo_: [structure] storing the results
% %
% SPECIAL REQUIREMENTS % SPECIAL REQUIREMENTS
% none % none
@ -37,36 +39,23 @@ function prior_posterior_statistics(type,dataset,dataset_info,dispString)
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>. % along with Dynare. If not, see <https://www.gnu.org/licenses/>.
if nargin < 4 if nargin < 9
dispString = 'prior_posterior_statistics'; dispString = 'prior_posterior_statistics';
end end
global options_ estim_params_ oo_ M_ bayestopt_
localVars=[]; localVars=[];
Y = transpose(dataset.data); Y = transpose(dataset_.data);
gend = dataset.nobs; gend = dataset_.nobs;
data_index = dataset_info.missing.aindex;
missing_value = dataset_info.missing.state;
mean_varobs = dataset_info.descriptive.mean;
nvx = estim_params_.nvx;
nvn = estim_params_.nvn; nvn = estim_params_.nvn;
ncx = estim_params_.ncx; npar = estim_params_.nvx+nvn+estim_params_.ncx+estim_params_.ncn+estim_params_.np;
ncn = estim_params_.ncn;
np = estim_params_.np ;
npar = nvx+nvn+ncx+ncn+np;
naK = length(options_.filter_step_ahead); naK = length(options_.filter_step_ahead);
MaxNumberOfBytes=options_.MaxNumberOfBytes; MaxNumberOfBytes=options_.MaxNumberOfBytes;
endo_nbr=M_.endo_nbr; endo_nbr=M_.endo_nbr;
exo_nbr=M_.exo_nbr; exo_nbr=M_.exo_nbr;
meas_err_nbr=length(M_.Correlation_matrix_ME); meas_err_nbr=length(M_.Correlation_matrix_ME);
iendo = 1:endo_nbr;
horizon = options_.forecast; horizon = options_.forecast;
IdObs = bayestopt_.mfys;
if horizon if horizon
i_last_obs = gend+(1-M_.maximum_endo_lag:0); i_last_obs = gend+(1-M_.maximum_endo_lag:0);
end end
@ -130,13 +119,6 @@ varlist = options_.varlist;
if isempty(varlist) if isempty(varlist)
varlist = sort(M_.endo_names(1:M_.orig_endo_nbr)); varlist = sort(M_.endo_names(1:M_.orig_endo_nbr));
end end
nvar = length(varlist);
SelecVariables = [];
for i=1:nvar
if ~isempty(strmatch(varlist{i}, M_.endo_names, 'exact'))
SelecVariables = [SelecVariables; strmatch(varlist{i}, M_.endo_names, 'exact')];
end
end
n_variables_to_fill=13; n_variables_to_fill=13;
@ -171,17 +153,17 @@ localVars.filter_covariance=filter_covariance;
localVars.smoothed_state_uncertainty=smoothed_state_uncertainty; localVars.smoothed_state_uncertainty=smoothed_state_uncertainty;
localVars.gend=gend; localVars.gend=gend;
localVars.Y=Y; localVars.Y=Y;
localVars.data_index=data_index; localVars.data_index=dataset_info.missing.aindex;
localVars.missing_value=missing_value; localVars.missing_value=dataset_info.missing.state;
localVars.varobs=options_.varobs; localVars.varobs=options_.varobs;
localVars.mean_varobs=mean_varobs; localVars.mean_varobs=dataset_info.descriptive.mean;
localVars.irun=irun; localVars.irun=irun;
localVars.endo_nbr=endo_nbr; localVars.endo_nbr=endo_nbr;
localVars.nvn=nvn; localVars.nvn=nvn;
localVars.naK=naK; localVars.naK=naK;
localVars.horizon=horizon; localVars.horizon=horizon;
localVars.iendo=iendo; localVars.iendo=1:endo_nbr;
localVars.IdObs=IdObs; localVars.IdObs=bayestopt_.mfys;
if horizon if horizon
localVars.i_last_obs=i_last_obs; localVars.i_last_obs=i_last_obs;
localVars.MAX_nforc1=MAX_nforc1; localVars.MAX_nforc1=MAX_nforc1;
@ -212,6 +194,13 @@ localVars.MAX_momentsno = MAX_momentsno;
localVars.ifil=ifil; localVars.ifil=ifil;
localVars.DirectoryName = DirectoryName; localVars.DirectoryName = DirectoryName;
localVars.M_=M_;
localVars.oo_=oo_;
localVars.options_=options_;
localVars.estim_params_=estim_params_;
localVars.bayestopt_=bayestopt_;
if strcmpi(type,'posterior') if strcmpi(type,'posterior')
record=load_last_mh_history_file(DirectoryName, M_.fname); record=load_last_mh_history_file(DirectoryName, M_.fname);
[nblck, npar] = size(record.LastParameters); [nblck, npar] = size(record.LastParameters);
@ -290,11 +279,7 @@ else
end end
end end
localVars.ifil = ifil; localVars.ifil = ifil;
globalVars = struct('M_',M_, ... globalVars = [];
'options_', options_, ...
'bayestopt_', bayestopt_, ...
'estim_params_', estim_params_, ...
'oo_', oo_);
% which files have to be copied to run remotely % which files have to be copied to run remotely
NamFileInput(1,:) = {'',[M_.fname '.static.m']}; NamFileInput(1,:) = {'',[M_.fname '.static.m']};
NamFileInput(2,:) = {'',[M_.fname '.dynamic.m']}; NamFileInput(2,:) = {'',[M_.fname '.dynamic.m']};
@ -330,28 +315,28 @@ if ~isnumeric(options_.parallel)
end end
if options_.smoother if options_.smoother
pm3(endo_nbr,gend,ifil(1),B,'Smoothed variables',... oo_=pm3(M_,options_,oo_,endo_nbr,gend,ifil(1),B,'Smoothed variables',...
'',varlist, M_.endo_names_tex,M_.endo_names,... varlist, M_.endo_names_tex,M_.endo_names,...
varlist,'SmoothedVariables',DirectoryName,'_smooth',dispString); varlist,'SmoothedVariables',DirectoryName,'_smooth',dispString);
pm3(exo_nbr,gend,ifil(2),B,'Smoothed shocks',... oo_=pm3(M_,options_,oo_,exo_nbr,gend,ifil(2),B,'Smoothed shocks',...
'',M_.exo_names,M_.exo_names_tex,M_.exo_names,... M_.exo_names,M_.exo_names_tex,M_.exo_names,...
M_.exo_names,'SmoothedShocks',DirectoryName,'_inno',dispString); M_.exo_names,'SmoothedShocks',DirectoryName,'_inno',dispString);
pm3(endo_nbr,1,ifil(9),B,'Trend_coefficients',... oo_=pm3(M_,options_,oo_,endo_nbr,1,ifil(9),B,'Trend_coefficients',...
'',varlist,M_.endo_names_tex,M_.endo_names,... varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'TrendCoeff',DirectoryName,'_trend_coeff',dispString); varlist,'TrendCoeff',DirectoryName,'_trend_coeff',dispString);
pm3(endo_nbr,gend,ifil(10),B,'Smoothed constant',... oo_=pm3(M_,options_,oo_,endo_nbr,gend,ifil(10),B,'Smoothed constant',...
'',varlist,M_.endo_names_tex,M_.endo_names,... varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'Constant',DirectoryName,'_smoothed_constant',dispString); varlist,'Constant',DirectoryName,'_smoothed_constant',dispString);
pm3(endo_nbr,gend,ifil(11),B,'Smoothed trend',... oo_=pm3(M_,options_,oo_,endo_nbr,gend,ifil(11),B,'Smoothed trend',...
'',varlist,M_.endo_names_tex,M_.endo_names,... varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'Trend',DirectoryName,'_smoothed_trend',dispString); varlist,'Trend',DirectoryName,'_smoothed_trend',dispString);
pm3(endo_nbr,gend,ifil(1),B,'Updated Variables',... oo_=pm3(M_,options_,oo_,endo_nbr,gend,ifil(1),B,'Updated Variables',...
'',varlist,M_.endo_names_tex,M_.endo_names,... varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'UpdatedVariables',DirectoryName, ... varlist,'UpdatedVariables',DirectoryName, ...
'_update',dispString); '_update',dispString);
if smoothed_state_uncertainty if smoothed_state_uncertainty
pm3(endo_nbr,endo_nbr,ifil(13),B,'State Uncertainty',... oo_=pm3(M_,options_,oo_,endo_nbr,endo_nbr,ifil(13),B,'State Uncertainty',...
'',varlist,M_.endo_names_tex,M_.endo_names,... varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'StateUncertainty',DirectoryName,'_state_uncert',dispString); varlist,'StateUncertainty',DirectoryName,'_state_uncert',dispString);
end end
@ -360,24 +345,24 @@ if options_.smoother
meas_error_names{obs_iter,1}=['SE_EOBS_' M_.endo_names{strmatch(options_.varobs{obs_iter},M_.endo_names,'exact')}]; meas_error_names{obs_iter,1}=['SE_EOBS_' M_.endo_names{strmatch(options_.varobs{obs_iter},M_.endo_names,'exact')}];
texnames{obs_iter,1}=['\sigma^{ME}_' M_.endo_names_tex{strmatch(options_.varobs{obs_iter},M_.endo_names,'exact')}]; texnames{obs_iter,1}=['\sigma^{ME}_' M_.endo_names_tex{strmatch(options_.varobs{obs_iter},M_.endo_names,'exact')}];
end end
pm3(meas_err_nbr,gend,ifil(3),B,'Smoothed measurement errors',... oo_=pm3(M_,options_,oo_,meas_err_nbr,gend,ifil(3),B,'Smoothed measurement errors',...
'',meas_error_names,texnames,meas_error_names,... meas_error_names,texnames,meas_error_names,...
meas_error_names,'SmoothedMeasurementErrors',DirectoryName,'_error',dispString) meas_error_names,'SmoothedMeasurementErrors',DirectoryName,'_error',dispString);
end end
end end
if options_.filtered_vars if options_.filtered_vars
pm3(endo_nbr,gend,ifil(4),B,'One step ahead forecast (filtered variables)',... oo_=pm3(M_,options_,oo_,endo_nbr,gend,ifil(4),B,'One step ahead forecast (filtered variables)',...
'',varlist,M_.endo_names_tex,M_.endo_names,... varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'FilteredVariables',DirectoryName,'_filter_step_ahead',dispString); varlist,'FilteredVariables',DirectoryName,'_filter_step_ahead',dispString);
end end
if options_.forecast if options_.forecast
pm3(endo_nbr,horizon,ifil(6),B,'Forecasted variables (mean)',... oo_=pm3(M_,options_,oo_,endo_nbr,horizon,ifil(6),B,'Forecasted variables (mean)',...
'',varlist,M_.endo_names_tex,M_.endo_names,... varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'MeanForecast',DirectoryName,'_forc_mean',dispString); varlist,'MeanForecast',DirectoryName,'_forc_mean',dispString);
pm3(endo_nbr,horizon,ifil(7),B,'Forecasted variables (point)',... oo_=pm3(M_,options_,oo_,endo_nbr,horizon,ifil(7),B,'Forecasted variables (point)',...
'',varlist,M_.endo_names_tex,M_.endo_names,... varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'PointForecast',DirectoryName,'_forc_point',dispString); varlist,'PointForecast',DirectoryName,'_forc_point',dispString);
if ~isequal(M_.H,0) && ~isempty(intersect(options_.varobs,varlist)) if ~isequal(M_.H,0) && ~isempty(intersect(options_.varobs,varlist))
texnames = cell(length(options_.varobs), 1); texnames = cell(length(options_.varobs), 1);
@ -387,15 +372,15 @@ if options_.forecast
texnames{obs_iter}=M_.endo_names_tex{strmatch(options_.varobs{obs_iter},M_.endo_names,'exact')}; texnames{obs_iter}=M_.endo_names_tex{strmatch(options_.varobs{obs_iter},M_.endo_names,'exact')};
end end
varlist_forecast_ME=intersect(options_.varobs,varlist); varlist_forecast_ME=intersect(options_.varobs,varlist);
pm3(meas_err_nbr,horizon,ifil(12),B,'Forecasted variables (point) with ME',... oo_=pm3(M_,options_,oo_,meas_err_nbr,horizon,ifil(12),B,'Forecasted variables (point) with ME',...
'',varlist_forecast_ME,texnames,obs_names,... varlist_forecast_ME,texnames,obs_names,...
varlist_forecast_ME,'PointForecastME',DirectoryName,'_forc_point_ME',dispString) varlist_forecast_ME,'PointForecastME',DirectoryName,'_forc_point_ME',dispString);
end end
end end
if options_.filter_covariance if options_.filter_covariance
pm3(endo_nbr,endo_nbr,ifil(8),B,'Filtered covariances',... oo_=pm3(M_,options_,oo_,endo_nbr,endo_nbr,ifil(8),B,'Filtered covariances',...
'',varlist,M_.endo_names_tex,M_.endo_names,... varlist,M_.endo_names_tex,M_.endo_names,...
varlist,'FilterCovariance',DirectoryName,'_filter_covar',dispString); varlist,'FilterCovariance',DirectoryName,'_filter_covar',dispString);
end end

7
matlab/prior_posterior_statistics_core.m
View File

@ -47,14 +47,17 @@ function myoutput=prior_posterior_statistics_core(myinputs,fpar,B,whoiam, ThisMa
% You should have received a copy of the GNU General Public License % You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <https://www.gnu.org/licenses/>. % along with Dynare. If not, see <https://www.gnu.org/licenses/>.
global options_ oo_ M_ bayestopt_ estim_params_
if nargin<4 if nargin<4
whoiam=0; whoiam=0;
end end
% Reshape 'myinputs' for local computation. % Reshape 'myinputs' for local computation.
% In order to avoid confusion in the name space, the instruction struct2local(myinputs) is replaced by: % In order to avoid confusion in the name space, the instruction struct2local(myinputs) is replaced by:
M_=myinputs.M_;
oo_=myinputs.oo_;
options_=myinputs.options_;
estim_params_=myinputs.estim_params_;
bayestopt_=myinputs.bayestopt_;
type=myinputs.type; type=myinputs.type;
run_smoother=myinputs.run_smoother; run_smoother=myinputs.run_smoother;

2
tests/estimation/fs2000.mod
View File

@ -132,4 +132,4 @@ oo_.MarginalDensity.LaplaceApproximation = Laplace; %reset correct Laplace
%test prior sampling %test prior sampling
options_.prior_draws=100; options_.prior_draws=100;
options_.no_graph.posterior=0; options_.no_graph.posterior=0;
prior_posterior_statistics('prior',dataset_,dataset_info); %get smoothed and filtered objects and forecasts oo_=prior_posterior_statistics('prior',dataset_,dataset_info,M_,oo_,options_,estim_params_,bayestopt_); %get smoothed and filtered objects and forecasts