Added routines for PAC expectations.
Stéphane Adjemian (Charybdis)
parent
d347f57bab
commit
a135c49469
|
|
@ -0,0 +1,33 @@
|
|||
function equation(pacname)
|
||||
|
||||
% Updates the parameters of a PAC equation.
|
||||
%
|
||||
% INPUTS
|
||||
% - pacname [string] Name of the pac equation.
|
||||
%
|
||||
% OUTPUTS
|
||||
% - none
|
||||
%
|
||||
% SPECIAL REQUIREMENTS
|
||||
% none
|
||||
|
||||
% Copyright (C) 2018 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/>.
|
||||
|
||||
global M_ oo_
|
||||
|
||||
M_ = pac.update.parameters(pacname, M_, oo_);
|
||||
|
|
@ -0,0 +1,113 @@
|
|||
function DynareModel = update_pac_parameters(pacname, DynareModel, DynareOutput)
|
||||
|
||||
% Updates the parameters of a PAC equation.
|
||||
%
|
||||
% INPUTS
|
||||
% - pacname [string] Name of the pac equation.
|
||||
% - DynareModel [struct] M_ global structure (model properties)
|
||||
% - DynareOutput [struct] oo_ global structure (model results)
|
||||
%
|
||||
% OUTPUTS
|
||||
% - none
|
||||
%
|
||||
% SPECIAL REQUIREMENTS
|
||||
% none
|
||||
|
||||
% Copyright (C) 2018 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/>.
|
||||
|
||||
% Check that the first input is a row character array.
|
||||
if ~isrow(pacname)==1 || ~ischar(pacname)
|
||||
error('First input argument must be a row character array!')
|
||||
end
|
||||
|
||||
% Check the name of the PAC model.
|
||||
if ~isfield(DynareModel.pac_expectation, pacname)
|
||||
error('PAC model %s is not defined in the model block!', pacname)
|
||||
end
|
||||
|
||||
% Get PAC model description
|
||||
pacmodel = DynareModel.pac_expectation.(pacname);
|
||||
|
||||
% Get the name of the associated VAR model and test its existence.
|
||||
if ~isfield(DynareModel.var, pacmodel.var_model_name)
|
||||
error('Unknown VAR (%s) in PAC model (%s)!', pacmodel.var_model_name, pacname)
|
||||
end
|
||||
|
||||
varmodel = DynareModel.var.(pacmodel.var_model_name);
|
||||
|
||||
% Check that we have the values of the VAR matrices.
|
||||
if ~isfield(DynareOutput.var, pacmodel.var_model_name)
|
||||
error('VAR model %s has to be estimated first!', pacmodel.var_model_name)
|
||||
end
|
||||
|
||||
varcalib = DynareOutput.var.(pacmodel.var_model_name);
|
||||
|
||||
if ~isfield(varcalib, 'CompanionMatrix') || any(isnan(varcalib.CompanionMatrix(:)))
|
||||
error('VAR model %s has to be estimated first.', pacmode.var_model_name)
|
||||
end
|
||||
|
||||
% Build the vector of PAC parameters (ECM parameter + autoregressive parameters).
|
||||
pacvalues = DynareModel.params(pacmodel.equation_params);
|
||||
|
||||
% Get the indices for the stationary/nonstationary variables in the VAR system.
|
||||
if any(varmodel.nonstationary)
|
||||
idns = find(varmodel.nonstationary);
|
||||
if length(idns)<length(varmodel.eqn)
|
||||
ids = find(~varmodel.nonstationary);
|
||||
else
|
||||
% All the variables in the system are non stationary.
|
||||
ids = [];
|
||||
end
|
||||
else
|
||||
% All the variables in the system are stationary.
|
||||
idns = [];
|
||||
ids = 1:length(varmodel.eqn);
|
||||
end
|
||||
|
||||
% Get the value of the discount factor.
|
||||
beta = DynareModel.params(pacmodel.discount_index);
|
||||
|
||||
% Is growth argument passed to pac_expectation?
|
||||
if isfield(pacmodel, 'growth_index')
|
||||
growth_flag = true;
|
||||
growth_type = pacmodel.growth_type;
|
||||
else
|
||||
growth_flag = false;
|
||||
end
|
||||
|
||||
% Get h0 and h1 vectors (plus the parameter for the growth neutrality correction).
|
||||
if growth_flag
|
||||
[h0, h1, growthneutrality] = hVectors([pacvalues; beta], varcalib.CompanionMatrix, ids, idns);
|
||||
else
|
||||
[h0, h1] = hVectors([pacvalues; beta], varcalib.CompanionMatrix, ids, idns);
|
||||
end
|
||||
|
||||
% Update the parameters related to the stationary components.
|
||||
if length(ids)
|
||||
DynareModel.params(pacmodel.h0_param_indices) = h0;
|
||||
end
|
||||
|
||||
% Update the parameters related to the nonstationary components.
|
||||
if length(idns)
|
||||
DynareModel.params(pacmodel.h1_param_indices) = h1;
|
||||
end
|
||||
|
||||
% Update the parameter related to the growth neutrality correction.
|
||||
if growth_flag
|
||||
DynareModel.params(pacmodel.growth_neutrality_param_index) = growthneutrality;
|
||||
end
|
||||
|
|
@ -65,6 +65,7 @@ p = {'/distributions/' ; ...
|
|||
'/lmmcp/' ; ...
|
||||
'/optimization/' ; ...
|
||||
'/ols/'; ...
|
||||
'/pac-tools/'; ...
|
||||
'/modules/dates/src/' ; ...
|
||||
'/modules/dseries/src/' ; ...
|
||||
'/utilities/doc/' ; ...
|
||||
|
|
|
|||
|
|
@ -0,0 +1,39 @@
|
|||
function alpha = a2alpha(a)
|
||||
|
||||
% Computes the m alpha coefficients from the m a coefficients of the PAC model.
|
||||
%
|
||||
% INPUTS
|
||||
% - a [double] m*1 vector of coefficients.
|
||||
%
|
||||
% OUTPUTS
|
||||
% - alpha [double] m*1 vector of coefficients.
|
||||
%
|
||||
% NOTES
|
||||
%
|
||||
% Given the current estimate of the PAC parameters a_0, a_1, ..., a_{m-1}, the routine does the following:
|
||||
%
|
||||
% \alpha_{m} = a_{m-1}
|
||||
% \alpha_{m-1} = a_{m-2}-a_{m-1}
|
||||
% \alpha_{m-2} = a_{m-3}-a_{m-2}
|
||||
% ...
|
||||
% \alpha_3 = a_2-a_3
|
||||
% \alpha_2 = a_1-a_2
|
||||
% \alpha_1 = a_0-a_1-1
|
||||
%
|
||||
% Note that the last elements of input a are (a_0, a_1, ..., a_{m-1}).
|
||||
|
||||
% Return an error if the input is not a vector
|
||||
if ~isvector(a)
|
||||
error('Input argument has to be a vector of doubles!')
|
||||
end
|
||||
|
||||
% Get the number of PAC parameters (without the discount factor)
|
||||
m = length(a);
|
||||
|
||||
% Initialize the vector of transformed PAC parameters.
|
||||
alpha = zeros(m, 1);
|
||||
|
||||
% Compute the transformed parameters
|
||||
alpha(m) = a(m);
|
||||
alpha(2:m-1) = a(2:m-1)-a(3:m);
|
||||
alpha(1) = a(1)-a(2)-1;
|
||||
|
|
@ -0,0 +1,40 @@
|
|||
function G = buildGmatrix(alpha, beta)
|
||||
|
||||
% Builds the G matrix needed for PAC.
|
||||
%
|
||||
% INPUTS
|
||||
% - alpha [double] m*1 vector of PAC parameters (lag polynomial parameters).
|
||||
% - beta [double] scalar, discount factor.
|
||||
%
|
||||
% OUTPUTS
|
||||
% - G [double] (m+1)*(m+1) matrix.
|
||||
|
||||
if nargin<2
|
||||
error('Two input arguments are required (vector of alpha parameters and beta discount parameter)!')
|
||||
end
|
||||
|
||||
% Return an error if the first input is not a real vector.
|
||||
if ~isnumeric(alpha) || ~isreal(alpha) || ~isvector(alpha)
|
||||
error('First input argument has to be a vector of doubles!')
|
||||
end
|
||||
|
||||
% Return an error if the second input argument is not a discount factor
|
||||
if ~isnumeric(beta) || ~isreal(beta) || ~isscalar(beta) || beta<eps || beta>1-eps
|
||||
error('Second input argument has to be a discount factor!')
|
||||
end
|
||||
|
||||
% Get the number of parameters
|
||||
m = length(alpha);
|
||||
|
||||
% Return an error if params is too small.
|
||||
if m<1
|
||||
error('First input argument has to be a vector with at least one element.')
|
||||
end
|
||||
|
||||
% Initialize the returned G matrix.
|
||||
G = zeros(m);
|
||||
|
||||
% Fill the returned G matrix.
|
||||
G(1:m-1,2:m) = eye(m-1);
|
||||
G(m, :) = -flip(transpose(alpha));
|
||||
G(m, :) = G(m, :).*flip(cumprod(beta*ones(1,m)));
|
||||
|
|
@ -0,0 +1,41 @@
|
|||
function [G, alpha, beta] = buildGmatrixWithAlphaAndBeta(params)
|
||||
|
||||
% Builds the G matrix needed for PAC.
|
||||
%
|
||||
% INPUTS
|
||||
% - params [double] (m+1)*1 vector of PAC parameters.
|
||||
%
|
||||
% OUTPUTS
|
||||
% - G [double] (m+1)*(m+1) matrix.
|
||||
% - alpha [double] m*1 vector of PAC parameters.
|
||||
% - beta [double] scalar, discount factor.
|
||||
|
||||
% Return an error if the input is not a vector.
|
||||
if ~isvector(params) || ~isnumeric(params) || ~isreal(params)
|
||||
error('Input argument has to be a vector of doubles!')
|
||||
end
|
||||
|
||||
% Get the number of parameters
|
||||
m = length(params)-1;
|
||||
|
||||
% Return an error if params is too small.
|
||||
if m<1
|
||||
error('Input argument has to be a vector with at least two elements. The last element is the discount factor.')
|
||||
end
|
||||
|
||||
% Get the transformed PAC parameters and discount factor.
|
||||
alpha = flip(a2alpha(params(1:m)));
|
||||
beta = params(end);
|
||||
|
||||
% Return an error if beta is not a discount factor
|
||||
if beta<eps || beta>1-eps
|
||||
error('beta has to be a discount factor!')
|
||||
end
|
||||
|
||||
% Initialize the returned G matrix.
|
||||
G = zeros(m);
|
||||
|
||||
% Fill the returned G matrix.
|
||||
G(1:m-1,2:m) = eye(m-1);
|
||||
G(m, :) = -transpose(alpha);
|
||||
G(m, :) = G(m, :).*flip(cumprod(beta*ones(1,m)));
|
||||
|
|
@ -0,0 +1,36 @@
|
|||
function [h0, h1, longruncorrectionparameter] = hVectors(params, H, ids, idns)
|
||||
|
||||
% INPUTS
|
||||
% - params [double] (m+1)*1 vector, PAC parameters (lag polynomial coefficients and discount factor).
|
||||
% - H [double] (np*np) matrix, companion matrix of the VAR(p) model.
|
||||
% - ids [integer] n*1 vector, selection of the stationary components of the VAR.
|
||||
% - idns [integer] n*1 vector, selection of the non stationary components of the VAR.
|
||||
%
|
||||
% OUTPUTS
|
||||
% - h0 [double] 1*n vector.
|
||||
% - h1 [double] 1*n vector.
|
||||
|
||||
if nargout>1 && nargin<4
|
||||
error('Wrong number of Inputs/Outputs!')
|
||||
end
|
||||
|
||||
[G, alpha, beta] = buildGmatrixWithAlphaAndBeta(params);
|
||||
|
||||
A_1 = polyval(alpha, 1.0);
|
||||
A_b = polyval(alpha, beta);
|
||||
|
||||
m = length(alpha);
|
||||
n = length(H);
|
||||
|
||||
tmp = eye(n*m)-kron(G, transpose(H));
|
||||
|
||||
h0 = A_1*A_b*((kron(iota(m, m), H))'*(tmp\kron(iota(m, m), iota(n, ids))));
|
||||
|
||||
if nargout>1
|
||||
h1 = A_1*A_b*(kron(iota(m, m)'*inv(eye(m)-G), H')*(tmp\kron(iota(m, m), iota(n, idns))));
|
||||
end
|
||||
|
||||
if nargout>2
|
||||
d = A_1*A_b*(iota(m, m)'*inv((eye(m)-G)*(eye(m)-G))*iota(m, m));
|
||||
longruncorrectionparameter = (1-sum(params(2:end-2))-d);
|
||||
end
|
||||
|
|
@ -0,0 +1,33 @@
|
|||
function i = iota(n, idx)
|
||||
|
||||
% Returns a selection vector.
|
||||
%
|
||||
% INPUTS
|
||||
% - n [integer] scalar, dimension of the returned vector.
|
||||
% - idx [integer] vector or scalar, non zero entries indices in the returned vector.
|
||||
%
|
||||
% OUTPUTS
|
||||
% - i [integer] n*1 vector. All elements are zero except those specified in idx.
|
||||
|
||||
if ~isscalar(n) || ~isvector(idx)
|
||||
error('First input has to be a scalar and second input a vector!')
|
||||
end
|
||||
|
||||
if (n-round(n))>eps(0) || any((idx-round(idx))>eps(0))
|
||||
error('Inputs must have integer values!')
|
||||
end
|
||||
|
||||
if n<1
|
||||
error('First input argument must be a strictly positive integer!')
|
||||
end
|
||||
|
||||
if any(idx>n)
|
||||
error('Elements in second argument cannot be greater than the first argument!')
|
||||
end
|
||||
|
||||
if any(idx<1)
|
||||
error('Elements in the second argument must be positive!')
|
||||
end
|
||||
|
||||
i = zeros(n, 1);
|
||||
i(idx) = 1;
|
||||
Loading…
Reference in New Issue