Home > matlab > stoch_simul.m

stoch_simul

PURPOSE ^

Copyright (C) 2001-2011 Dynare Team

SYNOPSIS ^

function info=stoch_simul(var_list)

DESCRIPTION ^

 Copyright (C) 2001-2011 Dynare Team

 This file is part of Dynare.

 Dynare is free software: you can redistribute it and/or modify
 it under the terms of the GNU General Public License as published by
 the Free Software Foundation, either version 3 of the License, or
 (at your option) any later version.

 Dynare is distributed in the hope that it will be useful,
 but WITHOUT ANY WARRANTY; without even the implied warranty of
 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 GNU General Public License for more details.

 You should have received a copy of the GNU General Public License
 along with Dynare.  If not, see <http://www.gnu.org/licenses/>.

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function info=stoch_simul(var_list)
0002 
0003 % Copyright (C) 2001-2011 Dynare Team
0004 %
0005 % This file is part of Dynare.
0006 %
0007 % Dynare is free software: you can redistribute it and/or modify
0008 % it under the terms of the GNU General Public License as published by
0009 % the Free Software Foundation, either version 3 of the License, or
0010 % (at your option) any later version.
0011 %
0012 % Dynare is distributed in the hope that it will be useful,
0013 % but WITHOUT ANY WARRANTY; without even the implied warranty of
0014 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0015 % GNU General Public License for more details.
0016 %
0017 % You should have received a copy of the GNU General Public License
0018 % along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
0019 
0020 global M_ options_ oo_ it_
0021 
0022 test_for_deep_parameters_calibration(M_);
0023 
0024 dr = oo_.dr;
0025 
0026 options_old = options_;
0027 if options_.linear
0028     options_.order = 1;
0029 end
0030 if options_.order == 1
0031     options_.replic = 1;
0032 elseif options_.order == 3
0033     options_.k_order_solver = 1;
0034 end
0035 
0036 if isempty(options_.qz_criterium)
0037     options_.qz_criterium = 1+1e-6;
0038 end
0039 
0040 if options_.partial_information == 1 || options_.ACES_solver == 1
0041     PI_PCL_solver = 1;
0042     if options_.order ~= 1
0043         warning('STOCH_SIMUL: forcing order=1 since you are using partial_information or ACES solver')
0044         options_.order = 1;
0045     end
0046 else
0047     PI_PCL_solver = 0;
0048 end
0049 
0050 TeX = options_.TeX;
0051 
0052 if size(var_list,1) == 0
0053     var_list = M_.endo_names(1:M_.orig_endo_nbr, :);
0054 end
0055 
0056 [i_var,nvar] = varlist_indices(var_list,M_.endo_names);
0057 
0058 iter_ = max(options_.periods,1);
0059 if M_.exo_nbr > 0
0060     oo_.exo_simul= ones(iter_ + M_.maximum_lag + M_.maximum_lead,1) * oo_.exo_steady_state';
0061 end
0062 
0063 check_model;
0064 
0065 oo_.dr=set_state_space(dr,M_);
0066 
0067 if PI_PCL_solver
0068     [oo_.dr, info] = PCL_resol(oo_.steady_state,0);
0069 elseif options_.discretionary_policy
0070     if ~options_.linear
0071         error(['discretionary_policy solves only linear_quadratic ' ...
0072                'problems']);
0073     end
0074     [oo_.dr,ys,info] = discretionary_policy_1(oo_,options_.instruments);
0075 else
0076     [oo_.dr,info,M_,options_,oo_] = resol(0,M_,options_,oo_);
0077 end
0078 
0079 if info(1)
0080     options_ = options_old;
0081     print_info(info, options_.noprint);
0082     return
0083 end
0084 
0085 if ~options_.noprint
0086     disp(' ')
0087     disp('MODEL SUMMARY')
0088     disp(' ')
0089     disp(['  Number of variables:         ' int2str(M_.endo_nbr)])
0090     disp(['  Number of stochastic shocks: ' int2str(M_.exo_nbr)])
0091     if (options_.block)
0092         disp(['  Number of state variables:   ' int2str(oo_.dr.npred+oo_.dr.nboth)])
0093         disp(['  Number of jumpers:           ' int2str(oo_.dr.nfwrd+oo_.dr.nboth)])
0094     else
0095         disp(['  Number of state variables:   ' ...
0096               int2str(length(find(oo_.dr.kstate(:,2) <= M_.maximum_lag+1)))])
0097         disp(['  Number of jumpers:           ' ...
0098               int2str(length(find(oo_.dr.kstate(:,2) == M_.maximum_lag+2)))])
0099     end;
0100     disp(['  Number of static variables:  ' int2str(oo_.dr.nstatic)])
0101     my_title='MATRIX OF COVARIANCE OF EXOGENOUS SHOCKS';
0102     labels = deblank(M_.exo_names);
0103     headers = char('Variables',labels);
0104     lh = size(labels,2)+2;
0105     dyntable(my_title,headers,labels,M_.Sigma_e,lh,10,6);
0106     if options_.partial_information
0107         disp(' ')
0108         disp('SOLUTION UNDER PARTIAL INFORMATION')
0109         disp(' ')
0110 
0111         if isfield(options_,'varobs')&& ~isempty(options_.varobs)
0112             PCL_varobs=options_.varobs;
0113             disp('OBSERVED VARIABLES')
0114         else
0115             PCL_varobs=M_.endo_names;
0116             disp(' VAROBS LIST NOT SPECIFIED')
0117             disp(' ASSUMED OBSERVED VARIABLES')
0118         end
0119         for i=1:size(PCL_varobs,1)
0120             disp(['    ' PCL_varobs(i,:)])
0121         end
0122     end
0123     disp(' ')
0124     if options_.order <= 2 && ~PI_PCL_solver
0125         disp_dr(oo_.dr,options_.order,var_list);
0126     end
0127 end
0128 
0129 if options_.periods > 0 && ~PI_PCL_solver
0130     if options_.periods <= options_.drop
0131         disp(['STOCH_SIMUL error: The horizon of simulation is shorter' ...
0132               ' than the number of observations to be DROPed'])
0133         options_ =options_old;
0134         return
0135     end
0136     if isempty(M_.endo_histval)
0137         y0 = oo_.dr.ys;
0138     else
0139         y0 = M_.endo_histval;
0140     end
0141     oo_.endo_simul = simult(y0,oo_.dr);
0142     dyn2vec;
0143 end
0144 
0145 if options_.nomoments == 0
0146     if PI_PCL_solver
0147         PCL_Part_info_moments (0, PCL_varobs, oo_.dr, i_var);
0148     elseif options_.periods == 0
0149         % There is no code for theoretical moments at 3rd order
0150         if options_.order <= 2
0151             disp_th_moments(oo_.dr,var_list);
0152         end
0153     else
0154         disp_moments(oo_.endo_simul,var_list);
0155     end
0156 end
0157 
0158 
0159 if options_.irf
0160     var_listTeX = M_.endo_names_tex(i_var,:);
0161 
0162     if TeX
0163         fidTeX = fopen([M_.fname '_IRF.TeX'],'w');
0164         fprintf(fidTeX,'%% TeX eps-loader file generated by stoch_simul.m (Dynare).\n');
0165         fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
0166         fprintf(fidTeX,' \n');
0167     end
0168     SS(M_.exo_names_orig_ord,M_.exo_names_orig_ord)=M_.Sigma_e+1e-14*eye(M_.exo_nbr);
0169     cs = transpose(chol(SS));
0170     tit(M_.exo_names_orig_ord,:) = M_.exo_names;
0171     if TeX
0172         titTeX(M_.exo_names_orig_ord,:) = M_.exo_names_tex;
0173     end
0174     irf_shocks_indx = getIrfShocksIndx();
0175     for i=irf_shocks_indx
0176         if SS(i,i) > 1e-13
0177             if PI_PCL_solver
0178                 y=PCL_Part_info_irf (0, PCL_varobs, i_var, M_, oo_.dr, options_.irf, i);
0179             else
0180                 y=irf(oo_.dr,cs(M_.exo_names_orig_ord,i), options_.irf, options_.drop, ...
0181                       options_.replic, options_.order);
0182             end
0183             if options_.relative_irf
0184                 y = 100*y/cs(i,i);
0185             end
0186             irfs   = [];
0187             mylist = [];
0188             if TeX
0189                 mylistTeX = [];
0190             end
0191             for j = 1:nvar
0192                 assignin('base',[deblank(M_.endo_names(i_var(j),:)) '_' deblank(M_.exo_names(i,:))],...
0193                          y(i_var(j),:)');
0194                 eval(['oo_.irfs.' deblank(M_.endo_names(i_var(j),:)) '_' ...
0195                       deblank(M_.exo_names(i,:)) ' = y(i_var(j),:);']);
0196                 if max(y(i_var(j),:)) - min(y(i_var(j),:)) > 1e-10
0197                     irfs  = cat(1,irfs,y(i_var(j),:));
0198                     if isempty(mylist)
0199                         mylist = deblank(var_list(j,:));
0200                     else
0201                         mylist = char(mylist,deblank(var_list(j,:)));
0202                     end
0203                     if TeX
0204                         if isempty(mylistTeX)
0205                             mylistTeX = deblank(var_listTeX(j,:));
0206                         else
0207                             mylistTeX = char(mylistTeX,deblank(var_listTeX(j,:)));
0208                         end
0209                     end
0210                 end
0211             end
0212             if options_.nograph == 0
0213                 number_of_plots_to_draw = size(irfs,1);
0214                 [nbplt,nr,nc,lr,lc,nstar] = pltorg(number_of_plots_to_draw);
0215                 if nbplt == 0
0216                 elseif nbplt == 1
0217                     if options_.relative_irf
0218                         hh = dyn_figure(options_,'Name',['Relative response to' ...
0219                                             ' orthogonalized shock to ' tit(i,:)]);
0220                     else
0221                         hh = dyn_figure(options_,'Name',['Orthogonalized shock to' ...
0222                                             ' ' tit(i,:)]);
0223                     end
0224                     for j = 1:number_of_plots_to_draw
0225                         subplot(nr,nc,j);
0226                         plot(1:options_.irf,transpose(irfs(j,:)),'-k','linewidth',1);
0227                         hold on
0228                         plot([1 options_.irf],[0 0],'-r','linewidth',0.5);
0229                         hold off
0230                         xlim([1 options_.irf]);
0231                         title(deblank(mylist(j,:)),'Interpreter','none');
0232                     end
0233                     dyn_saveas(hh,[M_.fname '_IRF_' deblank(tit(i,:))],options_);
0234                     if TeX
0235                         fprintf(fidTeX,'\\begin{figure}[H]\n');
0236                         for j = 1:number_of_plots_to_draw
0237                             fprintf(fidTeX,['\\psfrag{%s}[1][][0.5][0]{$%s$}\n'],deblank(mylist(j,:)),deblank(mylistTeX(j,:)));
0238                         end
0239                         fprintf(fidTeX,'\\centering \n');
0240                         fprintf(fidTeX,'\\includegraphics[scale=0.5]{%s_IRF_%s}\n',M_.fname,deblank(tit(i,:)));
0241                         fprintf(fidTeX,'\\caption{Impulse response functions (orthogonalized shock to $%s$).}',titTeX(i,:));
0242                         fprintf(fidTeX,'\\label{Fig:IRF:%s}\n',deblank(tit(i,:)));
0243                         fprintf(fidTeX,'\\end{figure}\n');
0244                         fprintf(fidTeX,' \n');
0245                     end
0246                 else
0247                     for fig = 1:nbplt-1
0248                         if options_.relative_irf
0249                             hh = dyn_figure(options_,'Name',['Relative response to orthogonalized shock' ...
0250                                                 ' to ' tit(i,:) ' figure ' int2str(fig)]);
0251                         else
0252                             hh = dyn_figure(options_,'Name',['Orthogonalized shock to ' tit(i,:) ...
0253                                                 ' figure ' int2str(fig)]);
0254                         end
0255                         for plt = 1:nstar
0256                             subplot(nr,nc,plt);
0257                             plot(1:options_.irf,transpose(irfs((fig-1)*nstar+plt,:)),'-k','linewidth',1);
0258                             hold on
0259                             plot([1 options_.irf],[0 0],'-r','linewidth',0.5);
0260                             hold off
0261                             xlim([1 options_.irf]);
0262                             title(deblank(mylist((fig-1)*nstar+plt,:)),'Interpreter','none');
0263                         end
0264                         dyn_saveas(hh,[ M_.fname '_IRF_' deblank(tit(i,:)) int2str(fig)],options_);
0265                         if TeX
0266                             fprintf(fidTeX,'\\begin{figure}[H]\n');
0267                             for j = 1:nstar
0268                                 fprintf(fidTeX,['\\psfrag{%s}[1][][0.5][0]{$%s$}\n'],deblank(mylist((fig-1)*nstar+j,:)),deblank(mylistTeX((fig-1)*nstar+j,:)));
0269                             end
0270                             fprintf(fidTeX,'\\centering \n');
0271                             fprintf(fidTeX,'\\includegraphics[scale=0.5]{%s_IRF_%s%s}\n',M_.fname,deblank(tit(i,:)),int2str(fig));
0272                             if options_.relative_irf
0273                                 fprintf(fidTeX,['\\caption{Relative impulse response' ...
0274                                                 ' functions (orthogonalized shock to $%s$).}'],deblank(titTeX(i,:)));
0275                             else
0276                                 fprintf(fidTeX,['\\caption{Impulse response functions' ...
0277                                                 ' (orthogonalized shock to $%s$).}'],deblank(titTeX(i,:)));
0278                             end
0279                             fprintf(fidTeX,'\\label{Fig:BayesianIRF:%s:%s}\n',deblank(tit(i,:)),int2str(fig));
0280                             fprintf(fidTeX,'\\end{figure}\n');
0281                             fprintf(fidTeX,' \n');
0282                         end
0283                     end
0284                     hh = dyn_figure(options_,'Name',['Orthogonalized shock to ' tit(i,:) ' figure ' int2str(nbplt) '.']);
0285                     m = 0;
0286                     for plt = 1:number_of_plots_to_draw-(nbplt-1)*nstar;
0287                         m = m+1;
0288                         subplot(lr,lc,m);
0289                         plot(1:options_.irf,transpose(irfs((nbplt-1)*nstar+plt,:)),'-k','linewidth',1);
0290                         hold on
0291                         plot([1 options_.irf],[0 0],'-r','linewidth',0.5);
0292                         hold off
0293                         xlim([1 options_.irf]);
0294                         title(deblank(mylist((nbplt-1)*nstar+plt,:)),'Interpreter','none');
0295                     end
0296                     dyn_saveas(hh,[ M_.fname '_IRF_' deblank(tit(i,:)) int2str(nbplt) ],options_);
0297                     if TeX
0298                         fprintf(fidTeX,'\\begin{figure}[H]\n');
0299                         for j = 1:m
0300                             fprintf(fidTeX,['\\psfrag{%s}[1][][0.5][0]{$%s$}\n'],deblank(mylist((nbplt-1)*nstar+j,:)),deblank(mylistTeX((nbplt-1)*nstar+j,:)));
0301                         end
0302                         fprintf(fidTeX,'\\centering \n');
0303                         fprintf(fidTeX,'\\includegraphics[scale=0.5]{%s_IRF_%s%s}\n',M_.fname,deblank(tit(i,:)),int2str(nbplt));
0304                         if options_.relative_irf
0305                             fprintf(fidTeX,['\\caption{Relative impulse response functions' ...
0306                                             ' (orthogonalized shock to $%s$).}'],deblank(titTeX(i,:)));
0307                         else
0308                             fprintf(fidTeX,['\\caption{Impulse response functions' ...
0309                                             ' (orthogonalized shock to $%s$).}'],deblank(titTeX(i,:)));
0310                         end
0311                         fprintf(fidTeX,'\\label{Fig:IRF:%s:%s}\n',deblank(tit(i,:)),int2str(nbplt));
0312                         fprintf(fidTeX,'\\end{figure}\n');
0313                         fprintf(fidTeX,' \n');
0314                     end
0315                 end
0316             end
0317         end
0318     end
0319     if TeX
0320         fprintf(fidTeX,' \n');
0321         fprintf(fidTeX,'%% End Of TeX file. \n');
0322         fclose(fidTeX);
0323     end
0324 end
0325 
0326 if options_.SpectralDensity.trigger == 1
0327     [omega,f] = UnivariateSpectralDensity(oo_.dr,var_list);
0328 end
0329 
0330 
0331 options_ = options_old;
0332 % temporary fix waiting for local options
0333 options_.partial_information = 0;

Generated on Tue 22-May-2012 02:40:23 by m2html © 2005