1function [nvar,vartan,NumberOfDecompFiles] = ... 2 dsge_simulated_theoretical_variance_decomposition(SampleSize,M_,options_,oo_,type) 3% function [nvar,vartan,NumberOfDecompFiles] = ... 4% dsge_simulated_theoretical_variance_decomposition(SampleSize,M_,options_,oo_,type) 5% This function computes the posterior or prior distribution of the variance 6% decomposition of the observed endogenous variables. 7% 8% INPUTS 9% SampleSize [integer] scalar, number of simulations. 10% M_ [structure] Dynare structure describing the model. 11% options_ [structure] Dynare structure defining global options. 12% oo_ [structure] Dynare structure where the results are saved. 13% type [string] 'prior' or 'posterior' 14% 15% 16% OUTPUTS 17% nvar [integer] nvar is the number of stationary variables. 18% vartan [char] array of characters (with nvar rows). 19% CovarFileNumber [integer] scalar, number of prior or posterior data files (for covariance). 20 21% Copyright (C) 2007-2020 Dynare Team 22% 23% This file is part of Dynare. 24% 25% Dynare is free software: you can redistribute it and/or modify 26% it under the terms of the GNU General Public License as published by 27% the Free Software Foundation, either version 3 of the License, or 28% (at your option) any later version. 29% 30% Dynare is distributed in the hope that it will be useful, 31% but WITHOUT ANY WARRANTY; without even the implied warranty of 32% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 33% GNU General Public License for more details. 34% 35% You should have received a copy of the GNU General Public License 36% along with Dynare. If not, see <http://www.gnu.org/licenses/>. 37 38nodecomposition = 0; 39 40% Get informations about the _posterior_draws files. 41if strcmpi(type,'posterior') 42 NumberOfDrawsFiles = length(dir([M_.dname '/metropolis/' M_.fname '_' type '_draws*' ])); 43 posterior = 1; 44elseif strcmpi(type,'prior') 45 NumberOfDrawsFiles = length(dir([M_.dname '/prior/draws/' type '_draws*' ])); 46 CheckPath('prior/moments',M_.dname); 47 posterior = 0; 48else 49 disp('dsge_simulated_theoretical_variance_decomposition:: Unknown type!') 50 error() 51end 52 53%delete old stale files before creating new ones 54if posterior 55 delete_stale_file([M_.dname '/metropolis/' M_.fname '_PosteriorVarianceDecomposition*']); 56 delete_stale_file([M_.dname '/metropolis/' M_.fname '_PosteriorVarianceDecompME*']); 57else 58 delete_stale_file([M_.dname '/prior/moments/' M_.fname '_PriorVarianceDecomposition*']); 59 delete_stale_file([M_.dname '/prior/moments/' M_.fname '_PriorVarianceDecompME*']);end 60 61% Set varlist (vartan) 62if ~posterior 63 if isfield(options_,'varlist') 64 temp = options_.varlist; 65 end 66 options_.varlist = options_.prior_analysis_endo_var_list; 67end 68[ivar,vartan,options_] = get_variables_list(options_,M_); 69if ~posterior 70 if exist('temp','var') 71 options_.varlist = temp; 72 end 73end 74nvar = length(ivar); 75 76% Set the size of the auto-correlation function to zero. 77nar = options_.ar; 78options_.ar = 0; 79 80 81 82nexo = M_.exo_nbr; 83 84NumberOfSavedElementsPerSimulation = nvar*(nexo+1); 85MaXNumberOfDecompLines = ceil(options_.MaxNumberOfBytes/NumberOfSavedElementsPerSimulation/8); 86 87ME_present=0; 88if ~all(diag(M_.H)==0) 89 if (isoctave && octave_ver_less_than('6')) || (~isoctave && matlab_ver_less_than('8.1')) 90 [observable_pos_requested_vars,index_subset,index_observables]=intersect_stable(ivar,options_.varobs_id); 91 else 92 [observable_pos_requested_vars,index_subset,index_observables]=intersect(ivar,options_.varobs_id,'stable'); 93 end 94 if ~isempty(observable_pos_requested_vars) 95 ME_present=1; 96 nobs_ME=length(observable_pos_requested_vars); 97 NumberOfSavedElementsPerSimulation_ME = nobs_ME*(nexo+1); 98 MaXNumberOfDecompLines_ME = ceil(options_.MaxNumberOfBytes/NumberOfSavedElementsPerSimulation_ME/8); 99 end 100end 101 102if SampleSize<=MaXNumberOfDecompLines 103 Decomposition_array = zeros(SampleSize,nvar*nexo); 104 NumberOfDecompFiles = 1; 105else 106 Decomposition_array = zeros(MaXNumberOfDecompLines,nvar*nexo); 107 NumberOfLinesInTheLastDecompFile = mod(SampleSize,MaXNumberOfDecompLines); 108 NumberOfDecompFiles = ceil(SampleSize/MaXNumberOfDecompLines); 109end 110 111NumberOfDecompLines = rows(Decomposition_array); 112DecompFileNumber = 1; 113 114if ME_present 115 if SampleSize<=MaXNumberOfDecompLines_ME 116 Decomposition_array_ME = zeros(SampleSize,nobs_ME*(nexo+1)); 117 NumberOfDecompFiles_ME = 1; 118 else 119 Decomposition_array_ME = zeros(MaXNumberOfDecompLines_ME,nobs_ME*(nexo+1)); 120 NumberOfLinesInTheLastDecompFile_ME = mod(SampleSize,MaXNumberOfDecompLines_ME); 121 NumberOfDecompFiles_ME = ceil(SampleSize/MaXNumberOfDecompLines_ME); 122 end 123 NumberOfDecompLines_ME = rows(Decomposition_array_ME); 124 DecompFileNumber_ME = 1; 125end 126% Compute total variances (covariances are not saved) and variances 127% implied by each structural shock. 128linea = 0; 129linea_ME = 0; 130only_non_stationary_vars=0; 131for file = 1:NumberOfDrawsFiles 132 if posterior 133 temp=load([M_.dname '/metropolis/' M_.fname '_' type '_draws' num2str(file) ]); 134 else 135 temp=load([M_.dname '/prior/draws/' type '_draws' num2str(file) ]); 136 end 137 isdrsaved = columns(temp.pdraws)-1; 138 NumberOfDraws = rows(temp.pdraws); 139 for linee = 1:NumberOfDraws 140 linea = linea+1; 141 linea_ME = linea_ME+1; 142 if isdrsaved 143 M_=set_parameters_locally(M_,temp.pdraws{linee,1});% Needed to update the covariance matrix of the state innovations. 144 dr = temp.pdraws{linee,2}; 145 else 146 M_=set_parameters_locally(M_,temp.pdraws{linee,1}); 147 [dr,info,M_,options_,oo_] = resol(0,M_,options_,oo_); 148 end 149 if file==1 && linee==1 150 [tmp, stationary_vars] = th_autocovariances(dr,ivar,M_,options_,nodecomposition); 151 if isempty(stationary_vars) 152 fprintf('\ndsge_simulated_theoretical_variance_decomposition:: All requested endogenous variables have a unit root and thus infinite variance.\n') 153 fprintf('dsge_simulated_theoretical_variance_decomposition:: No decomposition is performed.\n') 154 only_non_stationary_vars=1; 155 end 156 end 157 if only_non_stationary_vars 158 Decomposition_array(linea,:) = NaN; 159 else 160 tmp = th_autocovariances(dr,ivar,M_,options_,nodecomposition); 161 for i=1:nvar 162 for j=1:nexo 163 Decomposition_array(linea,(i-1)*nexo+j) = tmp{2}(i,j); 164 end 165 end 166 if ME_present 167 M_ = set_measurement_errors(temp.pdraws{linee,1},temp.estim_params_,M_); 168 ME_Variance=diag(M_.H); 169 tmp_ME=NaN(nobs_ME,nexo+1); 170 tmp_ME(:,1:end-1)=tmp{2}(index_subset,:).*repmat(diag(tmp{1}(index_subset,index_subset))./(diag(tmp{1}(index_subset,index_subset))+ME_Variance(index_observables)),1,nexo); 171 tmp_ME(:,end)=1-sum(tmp_ME(:,1:end-1),2); 172 for i=1:nobs_ME 173 for j=1:nexo+1 174 Decomposition_array_ME(linea,(i-1)*(nexo+1)+j) = tmp_ME(i,j); 175 end 176 end 177 end 178 end 179 if linea == NumberOfDecompLines 180 if posterior 181 save([M_.dname '/metropolis/' M_.fname '_PosteriorVarianceDecomposition' int2str(DecompFileNumber) '.mat' ],'Decomposition_array'); 182 else 183 save([M_.dname '/prior/moments/' M_.fname '_PriorVarianceDecomposition' int2str(DecompFileNumber) '.mat' ],'Decomposition_array'); 184 end 185 DecompFileNumber = DecompFileNumber + 1; 186 linea = 0; 187 test = DecompFileNumber-NumberOfDecompFiles; 188 if ~test% Prepare the last round... 189 Decomposition_array = zeros(NumberOfLinesInTheLastDecompFile,nvar*nexo); 190 NumberOfDecompLines = NumberOfLinesInTheLastDecompFile; 191 elseif test<0 192 Decomposition_array = zeros(MaXNumberOfDecompLines,nvar*nexo); 193 else 194 clear('Decomposition_array'); 195 end 196 end 197 if ME_present 198 if linea_ME == NumberOfDecompLines_ME 199 if posterior 200 save([M_.dname '/metropolis/' M_.fname '_PosteriorVarianceDecompME' int2str(DecompFileNumber_ME) '.mat' ],'Decomposition_array_ME'); 201 else 202 save([M_.dname '/prior/moments/' M_.fname '_PriorVarianceDecompME' int2str(DecompFileNumber_ME) '.mat' ],'Decomposition_array_ME'); 203 end 204 DecompFileNumber_ME = DecompFileNumber_ME + 1; 205 linea_ME = 0; 206 test = DecompFileNumber_ME-NumberOfDecompFiles_ME; 207 if ~test% Prepare the last round... 208 Decomposition_array_ME = zeros(NumberOfLinesInTheLastDecompFile_ME,nobs_ME*(nexo+1)); 209 NumberOfDecompLines_ME = NumberOfLinesInTheLastDecompFile_ME; 210 elseif test<0 211 Decomposition_array_ME = zeros(MaXNumberOfDecompLines_ME,nobs_ME*(nexo+1)); 212 else 213 clear('Decomposition_array_ME'); 214 end 215 end 216 end 217 end 218end 219 220options_.ar = nar; 221