diff --git a/matlab/conditional_variance_decomposition.m b/matlab/conditional_variance_decomposition.m index 047200a382f0335f4895bd40d0b6391324a30098..0a0de48451c6daf2d97302358349ed7bd7f1d977 100644 --- a/matlab/conditional_variance_decomposition.m +++ b/matlab/conditional_variance_decomposition.m @@ -35,6 +35,7 @@ function PackedConditionalVarianceDecomposition = conditional_variance_decomposi ConditionalVariance = zeros(StateSpaceModel.number_of_state_equations,StateSpaceModel.number_of_state_equations); ConditionalVariance = repmat(ConditionalVariance,[1 1 length(Steps) StateSpaceModel.number_of_state_innovations]); BB = StateSpaceModel.impulse_matrix*transpose(StateSpaceModel.impulse_matrix); + for h = 1:length(Steps) for t = 0:Steps(h) for i=1:StateSpaceModel.number_of_state_innovations @@ -44,6 +45,7 @@ function PackedConditionalVarianceDecomposition = conditional_variance_decomposi end end end + ConditionalVariance = ConditionalVariance(SubsetOfVariables,SubsetOfVariables,:,:); NumberOfVariables = length(SubsetOfVariables); PackedConditionalVarianceDecomposition = zeros(NumberOfVariables*(NumberOfVariables+1)/2,length(Steps),StateSpaceModel.number_of_state_innovations); diff --git a/matlab/conditional_variance_decomposition_mc_analysis.m b/matlab/conditional_variance_decomposition_mc_analysis.m index ee0a749be677609f8ec3639f25ad2caf8904fb5a..c0701c6927536835b03415064dba7fcf57ea8816 100644 --- a/matlab/conditional_variance_decomposition_mc_analysis.m +++ b/matlab/conditional_variance_decomposition_mc_analysis.m @@ -45,7 +45,7 @@ function oo_ = conditional_variance_decomposition_mc_analysis(NumberOfSimulation if isfield(temporary_structure,'dsge') eval(['temporary_structure = oo_.' TYPE 'TheoreticalMoments.dsge;']) if isfield(temporary_structure,'ConditionalVarianceDecomposition') - eval(['temporary_structure = oo_.' TYPE 'TheoreticalMoments.dsge.VarianceDecomposition.mean;']) + eval(['temporary_structure = oo_.' TYPE 'TheoreticalMoments.dsge.ConditionalVarianceDecomposition.mean;']) if isfield(temporary_structure,name) if sum(Steps-temporary_structure.(name)(1,:)) == 0 % Nothing (new) to do here... diff --git a/matlab/dsge_simulated_theoretical_conditional_variance_decomposition.m b/matlab/dsge_simulated_theoretical_conditional_variance_decomposition.m index 25f8a7e9915f6e6286407b2c198aa4aef13bdd84..2b993e21b13c86cf46f42fb88287b22ad98db140 100644 --- a/matlab/dsge_simulated_theoretical_conditional_variance_decomposition.m +++ b/matlab/dsge_simulated_theoretical_conditional_variance_decomposition.m @@ -40,6 +40,7 @@ if strcmpi(type,'posterior') posterior = 1; elseif strcmpi(type,'prior') DrawsFiles = dir([M_.dname '/prior/draws/' type '_draws*' ]); + CheckPath('prior/moments'); posterior = 0; else disp('dsge_simulated_theoretical_conditional_variance_decomposition:: Unknown type!') @@ -53,7 +54,7 @@ if ~posterior end options_.varlist = options_.prior_analysis_endo_var_list; end -[ivar,vartan, options_] = set_stationary_variables_list(options_, M_); +[ivar,vartan ] = set_stationary_variables_list(options_, M_); if ~posterior if exist('temp','var') options_.varlist = temp; @@ -65,9 +66,6 @@ nvar = length(ivar); nar = options_.ar; options_.ar = 0; - -NumberOfDrawsFiles = length(DrawsFiles); - NumberOfDrawsFiles = rows(DrawsFiles); NumberOfSavedElementsPerSimulation = nvar*(nvar+1)/2*M_.exo_nbr*length(Steps); MaXNumberOfConditionalDecompLines = ceil(options_.MaxNumberOfBytes/NumberOfSavedElementsPerSimulation/8); @@ -81,21 +79,13 @@ else NumberOfConditionalDecompFiles = ceil(SampleSize/MaXNumberOfConditionalDecompLines); end -NumberOfConditionalDecompLines = rows(Conditional_decomposition_array); -ConditionalDecompFileNumber = 1; +NumberOfConditionalDecompLines = size(Conditional_decomposition_array,4); +ConditionalDecompFileNumber = 0; StateSpaceModel.number_of_state_equations = M_.endo_nbr; StateSpaceModel.number_of_state_innovations = M_.exo_nbr; -endo_nbr = M_.endo_nbr; -nstatic = oo_.dr.nstatic; -npred = oo_.dr.npred; -iv = (1:endo_nbr)'; -ic = [ nstatic+(1:npred) endo_nbr+(1:size(oo_.dr.ghx,2)-npred) ]'; -aux = oo_.dr.transition_auxiliary_variables; -k = find(aux(:,2) > npred); -aux(:,2) = aux(:,2) + nstatic; -aux(k,2) = aux(k,2) + oo_.dr.nfwrd; +first_call = 1; linea = 0; for file = 1:NumberOfDrawsFiles @@ -115,11 +105,28 @@ for file = 1:NumberOfDrawsFiles set_parameters(pdraws{linee,1}); [dr,info] = resol(oo_.steady_state,0); end + if first_call + endo_nbr = M_.endo_nbr; + nstatic = dr.nstatic; + npred = dr.npred; + iv = (1:endo_nbr)'; + ic = [ nstatic+(1:npred) endo_nbr+(1:size(dr.ghx,2)-npred) ]'; + aux = dr.transition_auxiliary_variables; + k = find(aux(:,2) > npred); + aux(:,2) = aux(:,2) + nstatic; + aux(k,2) = aux(k,2) + dr.nfwrd; + StateSpaceModel.number_of_state_equations = M_.endo_nbr+rows(aux); + StateSpaceModel.number_of_state_innovations = M_.exo_nbr; + first_call = 0; + clear('endo_nbr','nstatic','npred','k'); + end [StateSpaceModel.transition_matrix,StateSpaceModel.impulse_matrix] = kalman_transition_matrix(dr,iv,ic,aux,M_.exo_nbr); StateSpaceModel.state_innovations_covariance_matrix = M_.Sigma_e; clear('dr'); - Conditional_decomposition_array(:,:,:,linea) = conditional_variance_decomposition(StateSpaceModel, Steps, ivar); + Conditional_decomposition_array(:,:,:,linea) = conditional_variance_decomposition(StateSpaceModel, Steps, ivar); if linea == NumberOfConditionalDecompLines + ConditionalDecompFileNumber = ConditionalDecompFileNumber + 1; + linea = 0; if posterior save([M_.dname '/metropolis/' M_.fname '_PosteriorConditionalVarianceDecomposition' int2str(ConditionalDecompFileNumber) '.mat' ], ... 'Conditional_decomposition_array'); @@ -127,14 +134,10 @@ for file = 1:NumberOfDrawsFiles save([M_.dname '/prior/moments/' M_.fname '_PriorConditionalVarianceDecomposition' int2str(ConditionalDecompFileNumber) '.mat' ], ... 'Conditional_decomposition_array'); end - ConditionalDecompFileNumber = ConditionalDecompFileNumber + 1; - linea = 0; - test = ConditionalDecompFileNumber-NumberOfConditionalDecompFiles; - if ~test% Prepare the last round... - Conditional_decomposition_array = zeros(nvar*(nvar+1)/2,length(Steps),M_.exo_nbr,NumberOfLinesInTheLastConditionalDecompFile); + if (ConditionalDecompFileNumber==NumberOfConditionalDecompFiles-1)% Prepare last round. + Conditional_decomposition_array = zeros(nvar*(nvar+1)/2, length(Steps),M_.exo_nbr,NumberOfLinesInTheLastConditionalDecompFile) ; NumberOfConditionalDecompLines = NumberOfLinesInTheLastConditionalDecompFile; - ConditionalDecompFileNumber = ConditionalDecompFileNumber - 1; - elseif test<0; + elseif ConditionalDecompFileNumber<NumberOfConditionalDecompFiles-1 Conditional_decomposition_array = zeros(nvar*(nvar+1)/2,length(Steps),M_.exo_nbr,MaXNumberOfConditionalDecompLines); else clear('Conditional_decomposition_array'); diff --git a/matlab/dsge_simulated_theoretical_correlation.m b/matlab/dsge_simulated_theoretical_correlation.m index c03eb2d419ecfecc594254ba6b9cd7f8e071b567..30ce0c516419e4a93a7fc4e1c82aa2af1d693eb1 100644 --- a/matlab/dsge_simulated_theoretical_correlation.m +++ b/matlab/dsge_simulated_theoretical_correlation.m @@ -40,6 +40,7 @@ if strcmpi(type,'posterior') posterior = 1; elseif strcmpi(type,'prior') DrawsFiles = dir([M_.dname '/prior/draws/' type '_draws*' ]); + CheckPath('prior/moments'); posterior = 0; else disp('dsge_simulated_theoretical_correlation:: Unknown type!'); diff --git a/matlab/dsge_simulated_theoretical_covariance.m b/matlab/dsge_simulated_theoretical_covariance.m index 79f58a3b9680e740ed16b98a075ac3bb6b932541..d647db5993df51f0741bc5793b62156a4552807c 100644 --- a/matlab/dsge_simulated_theoretical_covariance.m +++ b/matlab/dsge_simulated_theoretical_covariance.m @@ -40,6 +40,7 @@ if strcmpi(type,'posterior') posterior = 1; elseif strcmpi(type,'prior') DrawsFiles = dir([M_.dname '/prior/draws/' type '_draws*' ]); + CheckPath('prior/moments'); posterior = 0; else disp('dsge_simulated_theoretical_covariance:: Unknown type!') diff --git a/matlab/dsge_simulated_theoretical_variance_decomposition.m b/matlab/dsge_simulated_theoretical_variance_decomposition.m index 9a2afaa681277e2ef267a604ab995ad2dcec1656..53a70588edc625b230bbc1b00af0ed5fc2ecbb80 100644 --- a/matlab/dsge_simulated_theoretical_variance_decomposition.m +++ b/matlab/dsge_simulated_theoretical_variance_decomposition.m @@ -41,6 +41,7 @@ if strcmpi(type,'posterior') posterior = 1; elseif strcmpi(type,'prior') DrawsFiles = dir([M_.dname '/prior/draws/' type '_draws*' ]); + CheckPath('prior/moments'); posterior = 0; else disp('dsge_simulated_theoretical_variance_decomposition:: Unknown type!') diff --git a/matlab/prior_sampler.m b/matlab/prior_sampler.m index 3f8c98141bbecfc1a0a4aae0817956031f3bbc3b..407eaefc7402ee217aa8073d160863fa20e7efa5 100644 --- a/matlab/prior_sampler.m +++ b/matlab/prior_sampler.m @@ -34,7 +34,7 @@ function results = prior_sampler(drsave,M_,bayestopt_,options_,oo_) prior_draw(1,bayestopt_); PriorDirectoryName = CheckPath('prior/draws'); work = ~drsave; - iteration = 1; + iteration = 0; loop_indx = 0; file_indx = []; count_bk_indeterminacy = 0; @@ -55,48 +55,39 @@ function results = prior_sampler(drsave,M_,bayestopt_,options_,oo_) if NumberOfSimulations <= NumberOfElementsPerFile TableOfInformations = [ 1 , NumberOfSimulations , 1] ; else - NumberOfFiles = fix(NumberOfSimulations/NumberOfElementsPerFile) ; + NumberOfFiles = ceil(NumberOfSimulations/NumberOfElementsPerFile) ; NumberOfElementsInTheLastFile = NumberOfSimulations - NumberOfElementsPerFile*(NumberOfFiles-1) ; - if ~isint(NumberOfSimulations/NumberOfElementsPerFile) - NumberOfFiles = NumberOfFiles + 1 ; - end - TableOfInformations = NaN(NumberOfFiles,3); + TableOfInformations = NaN(NumberOfFiles,3) ; TableOfInformations(:,1) = transpose(1:NumberOfFiles) ; TableOfInformations(1:NumberOfFiles-1,2) = NumberOfElementsPerFile*ones(NumberOfFiles-1,1) ; TableOfInformations(NumberOfFiles,2) = NumberOfElementsInTheLastFile ; TableOfInformations(1,3) = 1; - TableOfInformations(2:end,3) = cumsum(TableOfInformations(2:end,2))+1; + TableOfInformations(2:end,3) = cumsum(TableOfInformations(1:end-1,2))+1; end pdraws = cell(TableOfInformations(1,2),drsave+1) ; sampled_prior_expectation = zeros(NumberOfParameters,1); sampled_prior_covariance = zeros(NumberOfParameters,NumberOfParameters); - file_line_number = 1; - file_indx_number = 1; - + file_line_number = 0; + file_indx_number = 0; + % Simulations. - while iteration <= NumberOfSimulations + while iteration < NumberOfSimulations loop_indx = loop_indx+1; - if ( (file_line_number-1)==TableOfInformations(file_indx_number,2)) - save([ PriorDirectoryName '/prior_draws' int2str(file_indx_number) '.mat' ],'pdraws'); - file_line_number = 1; - file_indx_number = file_indx_number + 1; - pdraws = cell(TableOfInformations(file_indx_number,2),drsave+1) ; - end params = prior_draw(); set_all_parameters(params); [dr,INFO] = resol(oo_.steady_state,work); - switch INFO(1) + switch INFO(1) case 0 + file_line_number = file_line_number + 1 ; + iteration = iteration + 1; pdraws(file_line_number,1) = {params}; if drsave pdraws(file_line_number,2) = {dr}; end - iteration = iteration+1; - file_line_number = file_line_number+1; [sampled_prior_expectation,sampled_prior_covariance] = ... - recursive_prior_moments(sampled_prior_expectation,sampled_prior_covariance,params,iteration) ; + recursive_prior_moments(sampled_prior_expectation,sampled_prior_covariance,params,iteration); case 1 count_static_undefined = count_static_undefined + 1; case 2 @@ -118,6 +109,14 @@ function results = prior_sampler(drsave,M_,bayestopt_,options_,oo_) otherwise count_unknown_problem = count_unknown_problem + 1 ; end + if ( file_line_number==TableOfInformations(file_indx_number+1,2) ) + file_indx_number = file_indx_number + 1; + save([ PriorDirectoryName '/prior_draws' int2str(file_indx_number) '.mat' ],'pdraws'); + if file_indx_number<NumberOfFiles + pdraws = cell(TableOfInformations(file_indx_number+1,2),drsave+1); + end + file_line_number = 0; + end end % Get informations about BK conditions and other things... @@ -141,8 +140,7 @@ function results = prior_sampler(drsave,M_,bayestopt_,options_,oo_) results.prior.mean = sampled_prior_expectation; results.prior.variance = sampled_prior_covariance; results.prior.mass = 1-results.garbage_share; - - + function [mu,sigma] = recursive_prior_moments(m0,s0,newobs,iter) % Recursive estimation of order one and two moments (expectation and % covariance matrix). newobs should be a row vector. I do not use the