2014-01-03 12:09:39 +01:00
function imcforecast ( constrained_paths, constrained_vars, options_cond_fcst)
2009-10-22 15:03:30 +02:00
% Computes conditional forecasts.
%
% INPUTS
% o constrained_paths [double] m*p array, where m is the number of constrained endogenous variables and p is the number of constrained periods.
% o constrained_vars [char] m*x array holding the names of the controlled endogenous variables.
% o options_cond_fcst [structure] containing the options. The fields are:
% + replic [integer] scalar, number of monte carlo simulations.
% + parameter_set [char] values of the estimated parameters:
% "posterior_mode",
% "posterior_mean",
% "posterior_median",
% "prior_mode" or
% "prior mean".
% [double] np*1 array, values of the estimated parameters.
% + controlled_varexo [char] m*x array, list of controlled exogenous variables.
% + conf_sig [double] scalar in [0,1], probability mass covered by the confidence bands.
%
% OUTPUTS
% None.
%
% SPECIAL REQUIREMENTS
% This routine has to be called after an estimation statement or an estimated_params block.
%
% REMARKS
% [1] Results are stored in a structure which is saved in a mat file called conditional_forecasts.mat.
% [2] Use the function plot_icforecast to plot the results.
2008-08-01 20:53:30 +02:00
2014-01-03 12:09:39 +01:00
% Copyright (C) 2006-2014 Dynare Team
2008-08-01 20:53:30 +02:00
%
% 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/>.
2009-10-22 15:03:30 +02:00
global options_ oo_ M_ bayestopt_
2006-12-15 12:37:24 +01:00
2011-07-30 20:19:15 +02:00
if ~ isfield ( options_cond_fcst , ' parameter_set' ) || isempty ( options_cond_fcst . parameter_set )
2009-10-22 15:03:30 +02:00
options_cond_fcst . parameter_set = ' posterior_mode' ;
2006-12-15 12:37:24 +01:00
end
2009-12-16 18:17:34 +01:00
2009-10-22 17:56:32 +02:00
if ~ isfield ( options_cond_fcst , ' replic' ) || isempty ( options_cond_fcst . replic )
2009-10-22 15:03:30 +02:00
options_cond_fcst . replic = 5000 ;
end
2009-10-22 17:56:32 +02:00
if ~ isfield ( options_cond_fcst , ' periods' ) || isempty ( options_cond_fcst . periods )
options_cond_fcst . periods = 40 ;
2009-10-22 15:03:30 +02:00
end
2009-10-22 17:56:32 +02:00
if ~ isfield ( options_cond_fcst , ' conf_sig' ) || isempty ( options_cond_fcst . conf_sig )
2009-10-22 15:03:30 +02:00
options_cond_fcst . conf_sig = . 8 ;
2006-12-15 12:37:24 +01:00
end
2011-07-30 20:19:15 +02:00
if isequal ( options_cond_fcst . parameter_set , ' calibration' )
estimated_model = 0 ;
else
estimated_model = 1 ;
2009-10-22 15:03:30 +02:00
end
2011-07-30 20:19:15 +02:00
if estimated_model
if ischar ( options_cond_fcst . parameter_set )
switch options_cond_fcst . parameter_set
case ' posterior_mode'
xparam = get_posterior_parameters ( ' mode' ) ;
2013-05-19 23:41:59 +02:00
graph_title = ' Posterior Mode' ;
2011-07-30 20:19:15 +02:00
case ' posterior_mean'
xparam = get_posterior_parameters ( ' mean' ) ;
2013-05-19 23:41:59 +02:00
graph_title = ' Posterior Mean' ;
2011-07-30 20:19:15 +02:00
case ' posterior_median'
xparam = get_posterior_parameters ( ' median' ) ;
2013-05-19 23:41:59 +02:00
graph_title = ' Posterior Median' ;
2011-07-30 20:19:15 +02:00
case ' prior_mode'
xparam = bayestopt_ . p5 ( : ) ;
2013-05-19 23:41:59 +02:00
graph_title = ' Prior Mode' ;
2011-07-30 20:19:15 +02:00
case ' prior_mean'
xparam = bayestopt_ . p1 ;
2013-05-19 23:41:59 +02:00
graph_title = ' Prior Mean' ;
2011-07-30 20:19:15 +02:00
otherwise
disp ( ' imcforecast:: If the input argument is a string, then it has to be equal to:' )
disp ( ' ' ' calibration' ' , ' )
disp ( ' ' ' posterior_mode' ' , ' )
disp ( ' ' ' posterior_mean' ' , ' )
disp ( ' ' ' posterior_median' ' , ' )
disp ( ' ' ' prior_mode' ' or' )
disp ( ' ' ' prior_mean' ' .' )
error ( ' imcforecast:: Wrong argument type!' )
2009-10-22 15:03:30 +02:00
end
else
2011-07-30 20:19:15 +02:00
xparam = options_cond_fcst . parameter_set ;
if length ( xparam ) ~= length ( M_ . params )
error ( ' imcforecast:: The dimension of the vector of parameters doesn' ' t match the number of estimated parameters!' )
end
2009-10-22 15:03:30 +02:00
end
2011-07-30 20:19:15 +02:00
set_parameters ( xparam ) ;
2014-06-26 10:55:54 +02:00
[ dataset_ , dataset_info ] = makedataset ( options_ ) ;
data = transpose ( dataset_ . data ) ;
data_index = dataset_info . missing . aindex ;
gend = dataset_ . nobs ;
missing_value = dataset_info . missing . state ;
2011-11-01 16:44:31 +01:00
[ atT , innov , measurement_error , filtered_state_vector , ys , trend_coeff ] = DsgeSmoother ( xparam , gend , data , data_index , missing_value ) ;
2011-07-30 20:19:15 +02:00
trend = repmat ( ys , 1 , options_cond_fcst . periods + 1 ) ;
for i = 1 : M_ . endo_nbr
j = strmatch ( deblank ( M_ . endo_names ( i , : ) ) , options_ . varobs , ' exact' ) ;
if ~ isempty ( j )
trend ( i , : ) = trend ( i , : ) + trend_coeff ( j ) * ( gend + ( 0 : options_cond_fcst . periods ) ) ;
end
2009-10-06 16:41:56 +02:00
end
2011-07-30 20:19:15 +02:00
trend = trend ( oo_ . dr . order_var , : ) ;
InitState ( : , 1 ) = atT ( : , end ) ;
else
InitState ( : , 1 ) = zeros ( M_ . endo_nbr , 1 ) ;
trend = repmat ( oo_ . steady_state ( oo_ . dr . order_var ) , 1 , options_cond_fcst . periods + 1 ) ;
2013-05-26 08:30:47 +02:00
graph_title = ' Calibration' ;
2009-10-06 16:41:56 +02:00
end
2011-07-30 20:19:15 +02:00
if isempty ( options_ . qz_criterium )
options_ . qz_criterium = 1 + 1e-6 ;
end
2014-07-18 18:36:02 +02:00
2011-09-19 12:19:04 +02:00
[ T , R , ys , info , M_ , options_ , oo_ ] = dynare_resolve ( M_ , options_ , oo_ ) ;
2006-12-15 12:37:24 +01:00
2014-07-18 18:36:02 +02:00
if ~ isdiagonal ( M_ . Sigma_e )
warning ( sprintf ( ' The innovations are correlated (the covariance matrix has non zero off diagonal elements), the results of the conditional forecasts will\ndepend on the ordering of the innovations (as declared after varexo) because a Cholesky decomposition is used to factorize the covariance matrix.\n\n=> It is preferable to declare the correlations in the model block (explicitly imposing the identification restrictions), unless you are satisfied\nwith the implicit identification restrictions implied by the Cholesky decomposition.' ) )
end
2014-03-09 17:05:59 +01:00
sQ = chol ( M_ . Sigma_e , ' lower' ) ;
2006-12-15 12:37:24 +01:00
NumberOfStates = length ( InitState ) ;
2009-10-22 15:03:30 +02:00
FORCS1 = zeros ( NumberOfStates , options_cond_fcst . periods + 1 , options_cond_fcst . replic ) ;
2006-12-15 12:37:24 +01:00
2011-07-30 20:19:15 +02:00
FORCS1 ( : , 1 , : ) = repmat ( InitState , 1 , options_cond_fcst . replic ) ;
2006-12-15 12:37:24 +01:00
2009-10-02 13:47:56 +02:00
EndoSize = M_ . endo_nbr ;
ExoSize = M_ . exo_nbr ;
2006-12-15 12:37:24 +01:00
2009-10-22 15:03:30 +02:00
n1 = size ( constrained_vars , 1 ) ;
n2 = size ( options_cond_fcst . controlled_varexo , 1 ) ;
2006-12-15 12:37:24 +01:00
if n1 ~= n2
2009-10-22 15:03:30 +02:00
error ( [ ' imcforecast:: The number of constrained variables doesn' ' t match the number of controlled shocks' ] )
2006-12-15 12:37:24 +01:00
end
idx = [ ] ;
jdx = [ ] ;
for i = 1 : n1
2011-07-30 20:19:15 +02:00
idx = [ idx ; oo_ . dr . inv_order_var ( constrained_vars ( i , : ) ) ] ;
2009-10-22 15:03:30 +02:00
jdx = [ jdx ; strmatch ( deblank ( options_cond_fcst . controlled_varexo ( i , : ) ) , M_ . exo_names , ' exact' ) ] ;
2006-12-15 12:37:24 +01:00
end
mv = zeros ( n1 , NumberOfStates ) ;
mu = zeros ( ExoSize , n2 ) ;
for i = 1 : n1
mv ( i , idx ( i ) ) = 1 ;
mu ( jdx ( i ) , i ) = 1 ;
end
2012-12-10 21:14:00 +01:00
% number of periods with constrained values
cL = size ( constrained_paths , 2 ) ;
2006-12-15 12:37:24 +01:00
2010-04-01 15:27:42 +02:00
constrained_paths = bsxfun ( @ minus , constrained_paths , trend ( idx , 1 : cL ) ) ;
2011-07-30 20:19:15 +02:00
%randn('state',0);
2006-12-15 12:37:24 +01:00
2009-10-22 15:03:30 +02:00
for b = 1 : options_cond_fcst . replic
shocks = sQ * randn ( ExoSize , options_cond_fcst . periods ) ;
shocks ( jdx , : ) = zeros ( length ( jdx ) , options_cond_fcst . periods ) ;
FORCS1 ( : , : , b ) = mcforecast3 ( cL , options_cond_fcst . periods , constrained_paths , shocks , FORCS1 ( : , : , b ) , T , R , mv , mu ) + trend ;
2006-12-15 12:37:24 +01:00
end
mFORCS1 = mean ( FORCS1 , 3 ) ;
2009-10-22 15:03:30 +02:00
tt = ( 1 - options_cond_fcst . conf_sig ) / 2 ;
t1 = round ( options_cond_fcst . replic * tt ) ;
t2 = round ( options_cond_fcst . replic * ( 1 - tt ) ) ;
2006-12-15 12:37:24 +01:00
2013-10-15 13:31:11 +02:00
forecasts . controlled_variables = constrained_vars ;
2009-10-22 15:03:30 +02:00
forecasts . instruments = options_cond_fcst . controlled_varexo ;
2006-12-15 12:37:24 +01:00
for i = 1 : EndoSize
2013-03-17 22:47:47 +01:00
eval ( [ ' forecasts.cond.Mean.' deblank ( M_ . endo_names ( oo_ . dr . order_var ( i ) , : ) ) ' = mFORCS1(i,:)' ' ;' ] ) ;
2006-12-15 12:37:24 +01:00
tmp = sort ( squeeze ( FORCS1 ( i , : , : ) ) ' ) ;
eval ( [ ' forecasts.cond.ci.' deblank ( M_ . endo_names ( oo_ . dr . order_var ( i ) , : ) ) ...
' = [tmp(t1,:)' ' ,tmp(t2,:)' ' ]' ' ;' ] ) ;
end
clear FORCS1 ;
2009-10-22 15:03:30 +02:00
FORCS2 = zeros ( NumberOfStates , options_cond_fcst . periods + 1 , options_cond_fcst . replic ) ;
for b = 1 : options_cond_fcst . replic
2006-12-15 12:37:24 +01:00
FORCS2 ( : , 1 , b ) = InitState ;
end
2011-10-10 12:17:06 +02:00
%randn('state',0);
2006-12-15 12:37:24 +01:00
2009-10-22 15:03:30 +02:00
for b = 1 : options_cond_fcst . replic
shocks = sQ * randn ( ExoSize , options_cond_fcst . periods ) ;
shocks ( jdx , : ) = zeros ( length ( jdx ) , options_cond_fcst . periods ) ;
FORCS2 ( : , : , b ) = mcforecast3 ( 0 , options_cond_fcst . periods , constrained_paths , shocks , FORCS2 ( : , : , b ) , T , R , mv , mu ) + trend ;
2006-12-15 12:37:24 +01:00
end
mFORCS2 = mean ( FORCS2 , 3 ) ;
for i = 1 : EndoSize
2013-03-17 22:47:47 +01:00
eval ( [ ' forecasts.uncond.Mean.' deblank ( M_ . endo_names ( oo_ . dr . order_var ( i ) , : ) ) ' = mFORCS2(i,:)' ' ;' ] ) ;
2006-12-15 12:37:24 +01:00
tmp = sort ( squeeze ( FORCS2 ( i , : , : ) ) ' ) ;
eval ( [ ' forecasts.uncond.ci.' deblank ( M_ . endo_names ( oo_ . dr . order_var ( i ) , : ) ) ...
' = [tmp(t1,:)' ' ,tmp(t2,:)' ' ]' ' ;' ] ) ;
end
2013-05-19 23:41:59 +02:00
forecasts . graph . title = graph_title ;
forecasts . graph . fname = M_ . fname ;
2006-12-15 12:37:24 +01:00
save ( ' conditional_forecasts.mat' , ' forecasts' ) ;