dynare_estimation_init

 function dynare_estimation_init(var_list_, gsa_flag)
 preforms initialization tasks before estimation or
 global sensitivity analysis

 INPUTS
   var_list_:  selected endogenous variables vector
   dname:      alternative directory name
   gsa_flag:   flag for GSA operation (optional)

 OUTPUTS
   data:    data after required transformation
   rawdata:  data as in the data file
   xparam1:    initial value of estimated parameters as returned by
               set_prior()

 SPECIAL REQUIREMENTS
   none

0001 function [dataset_,xparam1, M_, options_, oo_, estim_params_,bayestopt_, fake] = dynare_estimation_init(var_list_, dname, gsa_flag, M_, options_, oo_, estim_params_, bayestopt_)
0002 
0003 % function dynare_estimation_init(var_list_, gsa_flag)
0004 % preforms initialization tasks before estimation or
0005 % global sensitivity analysis
0006 %
0007 % INPUTS
0008 %   var_list_:  selected endogenous variables vector
0009 %   dname:      alternative directory name
0010 %   gsa_flag:   flag for GSA operation (optional)
0011 %
0012 % OUTPUTS
0013 %   data:    data after required transformation
0014 %   rawdata:  data as in the data file
0015 %   xparam1:    initial value of estimated parameters as returned by
0016 %               set_prior()
0017 %
0018 % SPECIAL REQUIREMENTS
0019 %   none
0020 
0021 % Copyright (C) 2003-2011 Dynare Team
0022 %
0023 % This file is part of Dynare.
0024 %
0025 % Dynare is free software: you can redistribute it and/or modify
0026 % it under the terms of the GNU General Public License as published by
0027 % the Free Software Foundation, either version 3 of the License, or
0028 % (at your option) any later version.
0029 %
0030 % Dynare is distributed in the hope that it will be useful,
0031 % but WITHOUT ANY WARRANTY; without even the implied warranty of
0032 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0033 % GNU General Public License for more details.
0034 %
0035 % You should have received a copy of the GNU General Public License
0036 % along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
0037 
0038 if isempty(gsa_flag)
0039     gsa_flag = 0;
0040 else% Decide if a DSGE or DSGE-VAR has to be estimated.
0041     if ~isempty(strmatch('dsge_prior_weight',M_.param_names))
0042         options_.dsge_var = 1;
0043     end
0044     var_list_ = check_list_of_variables(options_, M_, var_list_);
0045     options_.varlist = var_list_;
0046 end
0047 
0048 % Get the indices of the observed variables in M_.endo_names.
0049 options_.lgyidx2varobs = zeros(size(M_.endo_names,1),1);
0050 for i = 1:size(M_.endo_names,1)
0051     tmp = strmatch(deblank(M_.endo_names(i,:)),options_.varobs,'exact');
0052     if ~isempty(tmp)
0053         if length(tmp)>1
0054             disp(' ')
0055             error(['Multiple declarations of ' deblank(M_.endo_names(i,:)) ' as an observed variable is not allowed!'])
0056         end
0057         options_.lgyidx2varobs(i) = tmp;
0058     end
0059 end
0060 
0061 if options_.order>2
0062     error(['I cannot estimate a model with a ' int2str(options_.order) ' order approximation of the model!'])
0063 end
0064 
0065 % Set options_.lik_init equal to 3 if diffuse filter is used or
0066 % kalman_algo refers to a diffuse filter algorithm.
0067 if (options_.diffuse_filter==1) || (options_.kalman_algo > 2)
0068     if options_.lik_init == 2
0069         error(['options diffuse_filter, lik_init and/or kalman_algo have ' ...
0070                'contradictory settings'])
0071     else
0072         options_.lik_init = 3;
0073     end
0074 end
0075 
0076 % If options_.lik_init == 1
0077 %  set by default options_.qz_criterium to 1-1e-6
0078 %  and check options_.qz_criterium < 1-eps if options_.lik_init == 1
0079 % Else set by default options_.qz_criterium to 1+1e-6
0080 if options_.lik_init == 1
0081     if isempty(options_.qz_criterium)
0082         options_.qz_criterium = 1-1e-6;
0083     elseif options_.qz_criterium > 1-eps
0084         error(['estimation: option qz_criterium is too large for estimating ' ...
0085                'a stationary model. If your model contains unit roots, use ' ...
0086                'option diffuse_filter'])
0087     end
0088 elseif isempty(options_.qz_criterium)
0089     options_.qz_criterium = 1+1e-6;
0090 end
0091 
0092 % If the data are prefiltered then there must not be constants in the
0093 % measurement equation of the DSGE model or in the DSGE-VAR model.
0094 if options_.prefilter == 1
0095     options_.noconstant = 1;
0096 end
0097 
0098 % Set options related to filtered variables.
0099 if ~isequal(options_.filtered_vars,0) && isempty(options_.filter_step_ahead)
0100     options_.filter_step_ahead = 1;
0101 end
0102 if ~isequal(options_.filtered_vars,0) && isequal(options_.filter_step_ahead,0)
0103     options_.filter_step_ahead = 1;
0104 end
0105 if ~isequal(options_.filter_step_ahead,0)
0106     options_.nk = max(options_.filter_step_ahead);
0107 end
0108 
0109 % Set the name of the directory where (intermediary) results will be saved.
0110 if isempty(dname)
0111     M_.dname = M_.fname;
0112 else
0113     M_.dname = dname;
0114 end
0115 
0116 % Set the number of observed variables.
0117 n_varobs = size(options_.varobs,1);
0118 
0119 % Set priors over the estimated parameters.
0120 if ~isempty(estim_params_)
0121     [xparam1,estim_params_,bayestopt_,lb,ub,M_] = set_prior(estim_params_,M_,options_);
0122     if any(bayestopt_.pshape > 0)
0123         % Plot prior densities.
0124         if ~options_.nograph && options_.plot_priors
0125             plot_priors(bayestopt_,M_,estim_params_,options_)
0126         end
0127         % Set prior bounds
0128         bounds = prior_bounds(bayestopt_,options_);
0129         bounds(:,1)=max(bounds(:,1),lb);
0130         bounds(:,2)=min(bounds(:,2),ub);
0131     else
0132         % No priors are declared so Dynare will estimate the model by
0133         % maximum likelihood with inequality constraints for the parameters.
0134         options_.mh_replic = 0;% No metropolis.
0135         bounds(:,1) = lb;
0136         bounds(:,2) = ub;
0137     end
0138     % Test if initial values of the estimated parameters are all between
0139     % the prior lower and upper bounds.
0140     if any(xparam1 < bounds(:,1)) || any(xparam1 > bounds(:,2))
0141         find(xparam1 < bounds(:,1))
0142         find(xparam1 > bounds(:,2))
0143         error('Initial parameter values are outside parameter bounds')
0144     end
0145     lb = bounds(:,1);
0146     ub = bounds(:,2);
0147     bayestopt_.lb = lb;
0148     bayestopt_.ub = ub;
0149 else% If estim_params_ is empty...
0150     xparam1 = [];
0151     bayestopt_.lb = [];
0152     bayestopt_.ub = [];
0153     bayestopt_.jscale = [];
0154     bayestopt_.pshape = [];
0155     bayestopt_.p1 = [];
0156     bayestopt_.p2 = [];
0157     bayestopt_.p3 = [];
0158     bayestopt_.p4 = [];
0159     bayestopt_.p5 = [];
0160     bayestopt_.p6 = [];
0161     bayestopt_.p7 = [];
0162     estim_params_.nvx = 0;
0163     estim_params_.nvn = 0;
0164     estim_params_.ncx = 0;
0165     estim_params_.ncn = 0;
0166     estim_params_.np = 0;
0167 end
0168 
0169 % storing prior parameters in results
0170 oo_.prior.mean = bayestopt_.p1;
0171 oo_.prior.variance = diag(bayestopt_.p2.^2);
0172 
0173 % Is there a linear trend in the measurement equation?
0174 if ~isfield(options_,'trend_coeffs') % No!
0175     bayestopt_.with_trend = 0;
0176 else% Yes!
0177     bayestopt_.with_trend = 1;
0178     bayestopt_.trend_coeff = {};
0179     trend_coeffs = options_.trend_coeffs;
0180     nt = length(trend_coeffs);
0181     for i=1:n_varobs
0182         if i > length(trend_coeffs)
0183             bayestopt_.trend_coeff{i} = '0';
0184         else
0185             bayestopt_.trend_coeff{i} = trend_coeffs{i};
0186         end
0187     end
0188 end
0189 
0190 % Set the "size" of penalty.
0191 bayestopt_.penalty = 1e8;
0192 
0193 % Get informations about the variables of the model.
0194 dr = set_state_space(oo_.dr,M_);
0195 oo_.dr = dr;
0196 nstatic = dr.nstatic;          % Number of static variables.
0197 npred = dr.npred;              % Number of predetermined variables.
0198 nspred = dr.nspred;            % Number of predetermined variables in the state equation.
0199 
0200 % Test if observed variables are declared.
0201 if isempty(options_.varobs)
0202     error('VAROBS is missing')
0203 end
0204 
0205 % Setting resticted state space (observed + predetermined variables)
0206 var_obs_index = [];
0207 k1 = [];
0208 for i=1:n_varobs
0209     var_obs_index = [var_obs_index; strmatch(deblank(options_.varobs(i,:)),M_.endo_names(dr.order_var,:),'exact')];
0210     k1 = [k1; strmatch(deblank(options_.varobs(i,:)),M_.endo_names, 'exact')];
0211 end
0212 
0213 % Define union of observed and state variables
0214 k2 = union(var_obs_index,[dr.nstatic+1:dr.nstatic+dr.npred]', 'rows');
0215 % Set restrict_state to postion of observed + state variables in expanded state vector.
0216 oo_.dr.restrict_var_list = k2;
0217 bayestopt_.restrict_var_list = k2;
0218 % set mf0 to positions of state variables in restricted state vector for likelihood computation.
0219 [junk,bayestopt_.mf0] = ismember([dr.nstatic+1:dr.nstatic+dr.npred]',k2);
0220 % Set mf1 to positions of observed variables in restricted state vector for likelihood computation.
0221 [junk,bayestopt_.mf1] = ismember(var_obs_index,k2);
0222 % Set mf2 to positions of observed variables in expanded state vector for filtering and smoothing.
0223 bayestopt_.mf2  = var_obs_index;
0224 bayestopt_.mfys = k1;
0225 
0226 [junk,ic] = intersect(k2,nstatic+(1:npred)');
0227 oo_.dr.restrict_columns = [ic; length(k2)+(1:nspred-npred)'];
0228 
0229 k3 = [];
0230 k3p = [];
0231 if options_.selected_variables_only
0232     for i=1:size(var_list_,1)
0233         k3 = [k3; strmatch(var_list_(i,:),M_.endo_names(dr.order_var,:), ...
0234                            'exact')];
0235         k3p = [k3; strmatch(var_list_(i,:),M_.endo_names, ...
0236                            'exact')];
0237     end
0238 else
0239     k3 = (1:M_.endo_nbr)';
0240     k3p = (1:M_.endo_nbr)';
0241 end
0242 
0243 % Define union of observed and state variables
0244 if options_.block == 1
0245     k1 = k1';
0246     [k2, i_posA, i_posB] = union(k1', M_.state_var', 'rows');
0247     % Set restrict_state to postion of observed + state variables in expanded state vector.
0248     oo_.dr.restrict_var_list  = [k1(i_posA) M_.state_var(sort(i_posB))];
0249     % set mf0 to positions of state variables in restricted state vector for likelihood computation.
0250     [junk,bayestopt_.mf0] = ismember(M_.state_var',oo_.dr.restrict_var_list);
0251     % Set mf1 to positions of observed variables in restricted state vector for likelihood computation.
0252     [junk,bayestopt_.mf1] = ismember(k1,oo_.dr.restrict_var_list);
0253     % Set mf2 to positions of observed variables in expanded state vector for filtering and smoothing.
0254     bayestopt_.mf2  = var_obs_index;
0255     bayestopt_.mfys = k1;
0256     oo_.dr.restrict_columns = [size(i_posA,1)+(1:size(M_.state_var,2))];
0257 
0258     [k2, i_posA, i_posB] = union(k3p, M_.state_var', 'rows');
0259     bayestopt_.smoother_var_list = [k3p(i_posA); M_.state_var(sort(i_posB))'];
0260     [junk,junk,bayestopt_.smoother_saved_var_list] = intersect(k3p,bayestopt_.smoother_var_list(:));
0261     [junk,ic] = intersect(bayestopt_.smoother_var_list,M_.state_var);
0262     bayestopt_.smoother_restrict_columns = ic;
0263     [junk,bayestopt_.smoother_mf] = ismember(k1, ...
0264                                              bayestopt_.smoother_var_list);
0265 else
0266     k2 = union(var_obs_index,[dr.nstatic+1:dr.nstatic+dr.npred]', 'rows');
0267     % Set restrict_state to postion of observed + state variables in expanded state vector.
0268     oo_.dr.restrict_var_list = k2;
0269     % set mf0 to positions of state variables in restricted state vector for likelihood computation.
0270     [junk,bayestopt_.mf0] = ismember([dr.nstatic+1:dr.nstatic+dr.npred]',k2);
0271     % Set mf1 to positions of observed variables in restricted state vector for likelihood computation.
0272     [junk,bayestopt_.mf1] = ismember(var_obs_index,k2);
0273     % Set mf2 to positions of observed variables in expanded state vector for filtering and smoothing.
0274     bayestopt_.mf2  = var_obs_index;
0275     bayestopt_.mfys = k1;
0276     [junk,ic] = intersect(k2,nstatic+(1:npred)');
0277     oo_.dr.restrict_columns = [ic; length(k2)+(1:nspred-npred)'];
0278 
0279     bayestopt_.smoother_var_list = union(k2,k3);
0280     [junk,junk,bayestopt_.smoother_saved_var_list] = intersect(k3,bayestopt_.smoother_var_list(:));
0281     [junk,ic] = intersect(bayestopt_.smoother_var_list,nstatic+(1:npred)');
0282     bayestopt_.smoother_restrict_columns = ic;
0283     [junk,bayestopt_.smoother_mf] = ismember(var_obs_index, ...
0284                                              bayestopt_.smoother_var_list);
0285 end;
0286 
0287 if options_.analytic_derivation,
0288     if ~(exist('sylvester3','file')==2),        
0289         dynareroot = strrep(which('dynare'),'dynare.m','');
0290         addpath([dynareroot 'gensylv'])
0291     end
0292 end
0293 
0294 % Test if the data file is declared.
0295 if isempty(options_.datafile)
0296     if gsa_flag
0297         dataset_ = [];
0298 %         rawdata = [];
0299 %         data_info = [];
0300         return
0301     else
0302         error('datafile option is missing')
0303     end
0304 end
0305 
0306 % If jscale isn't specified for an estimated parameter, use global option options_.jscale, set to 0.2, by default.
0307 k = find(isnan(bayestopt_.jscale));
0308 bayestopt_.jscale(k) = options_.mh_jscale;
0309 
0310 % Load and transform data.
0311 transformation = [];
0312 if options_.loglinear && ~options_.logdata
0313     transformation = @log;
0314 end
0315 xls.sheet = options_.xls_sheet;
0316 xls.range = options_.xls_range;
0317 
0318 if ~isfield(options_,'nobs')
0319     options_.nobs = [];
0320 end
0321 
0322 dataset_ = initialize_dataset(options_.datafile,options_.varobs,options_.first_obs,options_.nobs,transformation,options_.prefilter,xls);
0323 
0324 options_.nobs = dataset_.info.ntobs;
0325 
0326 % setting noconstant option
0327 if options_.diffuse_filter
0328     steadystate_check_flag = 0;
0329 else
0330     steadystate_check_flag = 1;
0331 end
0332 
0333 M = M_;
0334 nvx = estim_params_.nvx;
0335 ncx = estim_params_.ncx;
0336 nvn = estim_params_.nvn;
0337 ncn = estim_params_.ncn;
0338 if isfield(estim_params_,'param_vals')
0339   M.params(estim_params_.param_vals(:,1)) = xparam1(nvx+ncx+nvn+ncn+1:end);
0340 end;
0341 oo_.steady_state = evaluate_steady_state(oo_.steady_state,M,options_,oo_,steadystate_check_flag);
0342 if all(abs(oo_.steady_state(bayestopt_.mfys))<1e-9)
0343     options_.noconstant = 1;
0344 else
0345     options_.noconstant = 0;
0346 end
0347 
0348