1function oo_ = initial_condition_decomposition(M_,oo_,options_,varlist,bayestopt_,estim_params_)
2% function oo_ = initial_condition_decomposition(M_,oo_,options_,varlist,bayestopt_,estim_params_)
3% Computes initial condition contribution to a simulated trajectory. The field set is
4% oo_.initval_decomposition. It is a n_var by n_var+2 by nperiods array. The
5% first n_var columns store the respective endogenous initval contribution, column n+1
6% stores the role of the shocks, while column n+2 stores the
7% value of the smoothed variables.  Variables are stored
8% in the order of declaration, i.e. M_.endo_names.
9%
10% INPUTS
11%    M_:            [structure]                Definition of the model
12%    oo_:           [structure]                Storage of results
13%    options_:      [structure]                Options
14%    varlist:       [cell of char array]       List of variables
15%    bayestopt_:    [structure]                Description of the priors
16%    estim_params_: [structure]                Estimated parameters
17%
18% OUTPUTS
19%    oo_:           [structure]                Storage of results
20%
21% SPECIAL REQUIREMENTS
22%    none
23
24% Copyright (C) 2017-2018 Dynare Team
25%
26% This file is part of Dynare.
27%
28% Dynare is free software: you can redistribute it and/or modify
29% it under the terms of the GNU General Public License as published by
30% the Free Software Foundation, either version 3 of the License, or
31% (at your option) any later version.
32%
33% Dynare is distributed in the hope that it will be useful,
34% but WITHOUT ANY WARRANTY; without even the implied warranty of
35% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
36% GNU General Public License for more details.
37%
38% You should have received a copy of the GNU General Public License
39% along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
40
41options_.plot_shock_decomp.colormap = options_.initial_condition_decomp.colormap;
42options_.plot_shock_decomp.nodisplay = options_.initial_condition_decomp.nodisplay;
43options_.plot_shock_decomp.graph_format = options_.initial_condition_decomp.graph_format;
44options_.plot_shock_decomp.fig_name = options_.initial_condition_decomp.fig_name;
45options_.plot_shock_decomp.detail_plot = options_.initial_condition_decomp.detail_plot;
46options_.plot_shock_decomp.init2shocks = options_.initial_condition_decomp.init2shocks;
47options_.plot_shock_decomp.steadystate = options_.initial_condition_decomp.steadystate;
48options_.plot_shock_decomp.write_xls = options_.initial_condition_decomp.write_xls;
49options_.plot_shock_decomp.type = options_.initial_condition_decomp.type;
50options_.plot_shock_decomp.plot_init_date = options_.initial_condition_decomp.plot_init_date;
51options_.plot_shock_decomp.plot_end_date = options_.initial_condition_decomp.plot_end_date;
52options_.plot_shock_decomp.diff = options_.initial_condition_decomp.diff;
53options_.plot_shock_decomp.flip = options_.initial_condition_decomp.flip;
54options_.plot_shock_decomp.max_nrows = options_.initial_condition_decomp.max_nrows;
55
56if isfield(options_.initial_condition_decomp,'init2shocks') % private trap for uimenu calls
57    init2shocks=options_.initial_condition_decomp.init2shocks;
58else
59    init2shocks=[];
60end
61% indices of endogenous variables
62if isempty(varlist)
63    varlist = M_.endo_names(1:M_.orig_endo_nbr);
64end
65
66if ~isequal(varlist,0)
67    [i_var, nvar, index_uniques] = varlist_indices(varlist, M_.endo_names);
68    varlist = varlist(index_uniques);
69end
70
71% number of variables
72endo_nbr = M_.endo_nbr;
73
74% parameter set
75parameter_set = options_.parameter_set;
76if isempty(parameter_set)
77    if isfield(oo_,'posterior_mean')
78        parameter_set = 'posterior_mean';
79    elseif isfield(oo_,'mle_mode')
80        parameter_set = 'mle_mode';
81    elseif isfield(oo_,'posterior')
82        parameter_set = 'posterior_mode';
83    else
84        error(['shock_decomposition: option parameter_set is not specified ' ...
85               'and posterior mode is not available'])
86    end
87end
88
89if ~isfield(oo_,'initval_decomposition') || isequal(varlist,0)
90    if isfield(oo_,'shock_decomposition_info') && isfield(oo_.shock_decomposition_info,'i_var')
91        if isfield (oo_,'realtime_conditional_shock_decomposition') ...
92                || isfield (oo_,'realtime_forecast_shock_decomposition') ...
93                || isfield (oo_,'realtime_shock_decomposition') ...
94                || isfield (oo_,'conditional_shock_decomposition') ...
95                || isfield (oo_,'shock_decomposition')
96            error('initval_decomposition::squeezed shock decompositions are already stored in oo_')
97        end
98    end
99    with_epilogue = options_.initial_condition_decomp.with_epilogue;
100    options_.selected_variables_only = 0; %make sure all variables are stored
101    options_.plot_priors=0;
102    [oo,M,~,~,Smoothed_Variables_deviation_from_mean] = evaluate_smoother(parameter_set,varlist,M_,oo_,options_,bayestopt_,estim_params_);
103
104    % reduced form
105    dr = oo.dr;
106
107    % data reordering
108    order_var = dr.order_var;
109    inv_order_var = dr.inv_order_var;
110
111
112    % coefficients
113    A = dr.ghx;
114    B = dr.ghu;
115
116    % initialization
117    gend = length(oo.SmoothedShocks.(M_.exo_names{1})); %+options_.forecast;
118    z = zeros(endo_nbr,endo_nbr+2,gend);
119    z(:,end,:) = Smoothed_Variables_deviation_from_mean;
120
121    for i=1:endo_nbr
122        z(i,i,1) = Smoothed_Variables_deviation_from_mean(i,1);
123    end
124
125    maximum_lag = M_.maximum_lag;
126
127    k2 = dr.kstate(find(dr.kstate(:,2) <= maximum_lag+1),[1 2]);
128    i_state = order_var(k2(:,1))+(min(i,maximum_lag)+1-k2(:,2))*M_.endo_nbr;
129    for i=1:gend
130        if i > 1 && i <= maximum_lag+1
131            lags = min(i-1,maximum_lag):-1:1;
132        end
133
134        if i > 1
135            tempx = permute(z(:,1:endo_nbr,lags),[1 3 2]);
136            m = min(i-1,maximum_lag);
137            tempx = [reshape(tempx,endo_nbr*m,endo_nbr); zeros(endo_nbr*(maximum_lag-i+1),endo_nbr)];
138            z(:,1:endo_nbr,i) = A(inv_order_var,:)*tempx(i_state,:);
139            lags = lags+1;
140        end
141        z(:,endo_nbr+1,i) = z(:,endo_nbr+2,i) - sum(z(:,1:endo_nbr,i),2);
142
143    end
144
145    if with_epilogue
146        [z, epilogue_steady_state] = epilogue_shock_decomposition(z, M_, oo_);
147        if ~isfield(oo_,'shock_decomposition_info') || ~isfield(oo_.shock_decomposition_info,'epilogue_steady_state')
148            oo_.shock_decomposition_info.epilogue_steady_state = epilogue_steady_state;
149        end
150    end
151    oo_.initval_decomposition = z;
152end
153
154% when varlist==0, we only store results in oo_ and do not make any plot
155if ~isequal(varlist,0)
156
157    % if ~options_.no_graph.shock_decomposition
158    oo=oo_;
159    oo.shock_decomposition = oo_.initval_decomposition;
160    if ~isempty(init2shocks)
161        init2shocks = M_.init2shocks.(init2shocks);
162        n=size(init2shocks,1);
163        for i=1:n
164            j=strmatch(init2shocks{i}{1},M_.endo_names,'exact');
165            oo.shock_decomposition(:,end-1,:)=oo.shock_decomposition(:,j,:)+oo.shock_decomposition(:,end-1,:);
166            oo.shock_decomposition(:,j,:)=0;
167        end
168    end
169    M_.exo_names = M_.endo_names;
170    M_.exo_nbr = M_.endo_nbr;
171    options_.plot_shock_decomp.realtime=0;
172    options_.plot_shock_decomp.screen_shocks=1;
173    options_.plot_shock_decomp.use_shock_groups = '';
174    options_.plot_shock_decomp.init_cond_decomp = 1; % private flag to plotting utilities
175
176    plot_shock_decomposition(M_,oo,options_,varlist);
177end
178