Added the possibility to have exogenous variables in the optimal part of PAC.
Works with iterative ols and nls.time-shift
Stéphane Adjemian (Charybdis)
parent
ba416f12ad
commit
6997e0a4a6
GPG Key ID:
A6D44CB9C64CE77B
|
|
@ -71,6 +71,9 @@ stack = dbstack;
|
|||
ipnames__ = ipnames_; % The user provided order.
|
||||
if isequal(stack(2).name, 'iterative_ols')
|
||||
ipnames_ = [M_.pac.(pacmodl).equations.(eqtag).ec.params; M_.pac.(pacmodl).equations.(eqtag).ar.params'];
|
||||
if isfield(M_.pac.(pacmodl).equations.(eqtag), 'optim_additive')
|
||||
ipnames_ = [ipnames_; M_.pac.(pacmodl).equations.(eqtag).optim_additive.params(~isnan(M_.pac.(pacmodl).equations.(eqtag).optim_additive.params))'];
|
||||
end
|
||||
if isfield(M_.pac.(pacmodl).equations.(eqtag), 'additive')
|
||||
ipnames_ = [ipnames_; M_.pac.(pacmodl).equations.(eqtag).additive.params(~isnan(M_.pac.(pacmodl).equations.(eqtag).additive.params))'];
|
||||
end
|
||||
|
|
|
|||
|
|
@ -57,6 +57,13 @@ else
|
|||
is_exogenous_variables = false;
|
||||
end
|
||||
|
||||
% Set flag for models with exogenous variables (in the optimizing agents part)
|
||||
if isfield(M_.pac.(pacmodl).equations.(eqtag), 'optim_additive')
|
||||
is_optim_exogenous_variables = length(M_.pac.(pacmodl).equations.(eqtag).optim_additive.vars)>0;
|
||||
else
|
||||
is_optim_exogenous_variables = false;
|
||||
end
|
||||
|
||||
if is_non_optimizing_agents
|
||||
non_optimizing_behaviour = M_.pac.(pacmodl).equations.(eqtag).non_optimizing_behaviour;
|
||||
non_optimizing_behaviour_params = NaN(length(non_optimizing_behaviour.params), 1);
|
||||
|
|
@ -100,14 +107,17 @@ if is_exogenous_variables
|
|||
additive.scaling_factor(residual_jd) = [];
|
||||
additive.estimation = ismember(additive.params, ipnames_);
|
||||
else
|
||||
additive.params = [];
|
||||
additive.vars = [];
|
||||
additive.isendo = [];
|
||||
additive.lags = [];
|
||||
additive.scaling_factor = [];
|
||||
additive.estimation = [];
|
||||
additive = struct('params', [], 'vars', [], 'isendo', [], 'lags', [], 'scaling_factor', [], 'estimation', []);
|
||||
end
|
||||
|
||||
if is_optim_exogenous_variables
|
||||
optim_additive = M_.pac.(pacmodl).equations.(eqtag).optim_additive;
|
||||
optim_additive.estimation = ismember(optim_additive.params, ipnames_);
|
||||
else
|
||||
optim_additive = struct('params', [], 'vars', [], 'isendo', [], 'lags', [], 'scaling_factor', [], 'estimation', []);
|
||||
end
|
||||
|
||||
|
||||
% Build PAC expectation matrix expression.
|
||||
dataForPACExpectation0 = dseries();
|
||||
listofvariables0 = {};
|
||||
|
|
@ -207,6 +217,43 @@ else
|
|||
dataForExogenousVariables_ = dseries();
|
||||
end
|
||||
|
||||
% Build data for exogenous variables (in the optimizing behaviour term).
|
||||
if is_optim_exogenous_variables
|
||||
listofvariables4 = {}; j = 0;
|
||||
dataForOptimExogenousVariables = dseries(); % Estimated parameters
|
||||
dataForOptimExogenousVariables_ = 0; % Calibrated parameters
|
||||
is_any_calibrated_parameter_optim_x = false;
|
||||
is_any_estimated_parameter_optim_x = false;
|
||||
for i=1:length(optim_additive.vars)
|
||||
if optim_additive.isendo(i)
|
||||
variable = M_.endo_names{optim_additive.vars(i)};
|
||||
else
|
||||
variable = M_.exo_names{optim_additive.vars(i)};
|
||||
end
|
||||
if optim_additive.estimation(i)
|
||||
j = j+1;
|
||||
is_any_estimated_parameter_optim_x = true;
|
||||
listofvariables4{j} = variable;
|
||||
dataForOptimExogenousVariables = [dataForOptimExogenousVariables, optim_additive.scaling_factor(i)*data{variable}.lag(optim_additive.lags(i))];
|
||||
else
|
||||
is_any_calibrated_parameter_optim_x = true;
|
||||
tmp = data{variable}.lag(optim_additive.lags(i)).data;
|
||||
if ~isnan(optim_additive.params(i))
|
||||
tmp = M_.params(optim_additive.params(i))*tmp;
|
||||
end
|
||||
tmp = optim_additive.scaling_factor(i)*tmp;
|
||||
dataForOptimExogenousVariables_ = dataForOptimExogenousVariables_+tmp;
|
||||
end
|
||||
end
|
||||
if is_any_calibrated_parameter_optim_x
|
||||
dataForOptimExogenousVariables_ = dseries(dataForOptimExogenousVariables_, data.dates(1), 'exogenous_variables_associated_with_calibrated_parameters');
|
||||
end
|
||||
else
|
||||
listofvariables4 = {};
|
||||
dataForOptimExogenousVariables = dseries();
|
||||
dataForOptimExogenousVariables_ = dseries();
|
||||
end
|
||||
|
||||
% Reorder ec.vars locally if necessary. Second variable must be the
|
||||
% endogenous variable, while the first must be the associated trend.
|
||||
if M_.pac.(pacmodl).equations.(eqtag).ec.isendo(2)
|
||||
|
|
@ -232,6 +279,10 @@ end
|
|||
|
||||
XDATA = DataForOLS{listofvariables3{:}}(range).data;
|
||||
|
||||
if is_optim_exogenous_variables && is_any_estimated_parameter_optim_x
|
||||
XDATA = [XDATA, dataForOptimExogenousVariables{listofvariables4{:}}(range).data];
|
||||
end
|
||||
|
||||
if is_exogenous_variables && is_any_estimated_parameter_x
|
||||
XDATA = [XDATA, dataForExogenousVariables{listofvariables2{:}}(range).data];
|
||||
end
|
||||
|
|
@ -245,30 +296,25 @@ else
|
|||
end
|
||||
|
||||
% Get indices in params0 for EC and AR parameters
|
||||
[~, params_id_1] = setdiff(ipnames_, [share_of_optimizing_agents_index, additive.params]);
|
||||
[~, params_id_1] = setdiff(ipnames_, [share_of_optimizing_agents_index, optim_additive.params, additive.params, ]);
|
||||
|
||||
% Get indices in params0 for other parameters (optimizing agents
|
||||
% share plus parameters related to exogenous variables).
|
||||
% Get indices in params0 for EC and AR parameters plus parameters related to exogenous variables in the optimal part.
|
||||
[~, params_id_5] = setdiff(ipnames_, [share_of_optimizing_agents_index, additive.params]);
|
||||
|
||||
% Get indices in params0 for other parameters (optimizing agents share plus parameters related to exogenous variables).
|
||||
[~, params_id_2] = setdiff(1:length(ipnames_), params_id_1);
|
||||
|
||||
% Get indices in params0 for the parameters associated to the
|
||||
% exogenous variables.
|
||||
% Get indices in params0 for the parameters associated to the exogenous variables.
|
||||
params_id_3 = setdiff(params_id_2, params_id_0);
|
||||
[~, params_id_3_o] = ismember(optim_additive.params, ipnames_);
|
||||
params_id_3_no = setdiff(params_id_3, params_id_3_o);
|
||||
|
||||
% Get indices in params0 for EC/AR parameters and parameters
|
||||
% associated to the exogenous variables (if any).
|
||||
% Get indices in params0 for EC/AR parameters and parameters associated to the exogenous variables (if any).
|
||||
params_id_4 = [params_id_1; params_id_3];
|
||||
|
||||
[~, params_id_1] = setdiff(ipnames_, [share_of_optimizing_agents_index, additive.params]);
|
||||
|
||||
|
||||
% Get values for EC-AR parameters plus the parameters associated to
|
||||
% the exogenous variables (if any).
|
||||
% Get values for EC-AR parameters plus the parameters associated to the exogenous variables (if any).
|
||||
params0_ = params0([params_id_1; params_id_3]);
|
||||
|
||||
% Get values for parameters associated to the exogenous variables.
|
||||
params0__ = params0(params_id_3);
|
||||
|
||||
% Get value of the share of optimizing agents.
|
||||
if estimate_non_optimizing_agents_share
|
||||
share_of_optimizing_agents = params0(params_id_0);
|
||||
|
|
@ -301,11 +347,17 @@ while noconvergence
|
|||
if is_exogenous_variables && is_any_calibrated_parameter_x
|
||||
YDATA = YDATA-dataForExogenousVariables_(range).data;
|
||||
end
|
||||
if is_optim_exogenous_variables && is_any_calibrated_parameter_optim_x
|
||||
YDATA = YDATA-share_of_optimizing_agents*dataForOptimExogenousVariables_(range).data;
|
||||
end
|
||||
% Run OLS to estimate PAC parameters (autoregressive parameters and error correction parameter).
|
||||
params1_ = (XDATA\YDATA);
|
||||
params1_(1:length(params_id_1)) = params1_(1:length(params_id_1))/share_of_optimizing_agents;
|
||||
params1_(1:length(params_id_5)) = params1_(1:length(params_id_5))/share_of_optimizing_agents;
|
||||
% Compute residuals and sum of squareed residuals.
|
||||
r = YDATA-XDATA*(params1_*share_of_optimizing_agents);
|
||||
r = YDATA-XDATA(:,1:length(params_id_5))*params1_(1:length(params_id_5))*share_of_optimizing_agents;
|
||||
if is_optim_exogenous_variables && is_any_estimated_parameter_optim_x
|
||||
r = r-XDATA(:,length(params_id_5)+1:end)*params1_(length(params_id_5)+1:end);
|
||||
end
|
||||
ssr = r'*r;
|
||||
% Update convergence dummy variable and display iteration info.
|
||||
noconvergence = max(abs(params0_-params1_))>1e-6;
|
||||
|
|
@ -329,7 +381,10 @@ while noconvergence
|
|||
end
|
||||
end
|
||||
% Set vector for regressor.
|
||||
ZDATA = XDATA(:,1:length(params_id_1))*params0_(1:length(params_id_1))+PacExpectations-dataForNonOptimizingBehaviour(range).data*non_optimizing_behaviour_params;
|
||||
ZDATA = XDATA(:,1:length(params_id_5))*params0_(1:length(params_id_5))+PacExpectations-dataForNonOptimizingBehaviour(range).data*non_optimizing_behaviour_params;
|
||||
if is_optim_exogenous_variables && is_any_calibrated_parameter_optim_x
|
||||
ZDATA = ZDATA+dataForOptimExogenousVariables_(range).data;
|
||||
end
|
||||
% Update the (estimated) share of optimizing agents by running OLS
|
||||
share_of_optimizing_agents = (ZDATA\YDATA);
|
||||
% Force the share of optimizing agents to be in [0,1].
|
||||
|
|
|
|||
|
|
@ -432,6 +432,7 @@ MODFILES = \
|
|||
pac/trend-component-22/example.mod \
|
||||
pac/trend-component-23/example.mod \
|
||||
pac/trend-component-24/example.mod \
|
||||
pac/trend-component-25/example.mod \
|
||||
write/example1.mod \
|
||||
ecb/backward-models/irf/solow_1.mod \
|
||||
ecb/backward-models/irf/solow_2.mod \
|
||||
|
|
|
|||
|
|
@ -0,0 +1,8 @@
|
|||
#!/bin/sh
|
||||
|
||||
rm -rf example
|
||||
rm -rf +example
|
||||
rm -f example.log
|
||||
rm -f *.mat
|
||||
rm -f *.m
|
||||
rm -f *.dat
|
||||
|
|
@ -0,0 +1,185 @@
|
|||
// --+ options: json=compute, stochastic +--
|
||||
|
||||
var x1 x2 x1bar x2bar z y x u v s;
|
||||
|
||||
varexo ex1 ex2 ex1bar ex2bar ez ey ex eu ev es;
|
||||
|
||||
parameters
|
||||
rho_1 rho_2 rho_3 rho_4
|
||||
a_x1_0 a_x1_1 a_x1_2 a_x1_x2_1 a_x1_x2_2
|
||||
a_x2_0 a_x2_1 a_x2_2 a_x2_x1_1 a_x2_x1_2
|
||||
e_c_m c_z_1 c_z_2 c_z_dx2 c_z_u c_z_dv c_z_s beta
|
||||
lambda;
|
||||
|
||||
rho_1 = .9;
|
||||
rho_2 = -.2;
|
||||
rho_3 = .4;
|
||||
rho_4 = -.3;
|
||||
|
||||
|
||||
a_x1_0 = -.9;
|
||||
a_x1_1 = .4;
|
||||
a_x1_2 = .3;
|
||||
a_x1_x2_1 = .1;
|
||||
a_x1_x2_2 = .2;
|
||||
|
||||
a_x2_0 = -.9;
|
||||
a_x2_1 = .2;
|
||||
a_x2_2 = -.1;
|
||||
a_x2_x1_1 = -.1;
|
||||
a_x2_x1_2 = .2;
|
||||
|
||||
beta = .2;
|
||||
e_c_m = .5;
|
||||
c_z_1 = .2;
|
||||
c_z_2 = -.1;
|
||||
c_z_dx2 = .3;
|
||||
c_z_u = .3;
|
||||
c_z_dv = .4;
|
||||
c_z_s = -.2;
|
||||
|
||||
lambda = 0.5; // Share of optimizing agents.
|
||||
|
||||
trend_component_model(model_name=toto, eqtags=['eq:x1', 'eq:x2', 'eq:x1bar', 'eq:x2bar'], targets=['eq:x1bar', 'eq:x2bar']);
|
||||
|
||||
pac_model(auxiliary_model_name=toto, discount=beta, model_name=pacman);
|
||||
|
||||
model;
|
||||
|
||||
[name='eq:u']
|
||||
s = .3*s(-1) - .1*s(-2) + es;
|
||||
|
||||
[name='eq:diff(v)']
|
||||
diff(v) = .5*diff(v(-1)) + ev;
|
||||
|
||||
[name='eq:u']
|
||||
u = .5*u(-1) - .2*u(-2) + eu;
|
||||
|
||||
[name='eq:y']
|
||||
y = rho_1*y(-1) + rho_2*y(-2) + ey;
|
||||
|
||||
[name='eq:x']
|
||||
x = rho_3*x(-1) + rho_4*x(-2) + ex;
|
||||
|
||||
[name='eq:x1']
|
||||
diff(x1) = a_x1_0*(x1(-1)-x1bar(-1)) + a_x1_1*diff(x1(-1)) + a_x1_2*diff(x1(-2)) + a_x1_x2_1*diff(x2(-1)) + a_x1_x2_2*diff(x2(-2)) + ex1;
|
||||
|
||||
[name='eq:x2']
|
||||
diff(x2) = a_x2_0*(x2(-1)-x2bar(-1)) + a_x2_1*diff(x1(-1)) + a_x2_2*diff(x1(-2)) + a_x2_x1_1*diff(x2(-1)) + a_x2_x1_2*diff(x2(-2)) + ex2;
|
||||
|
||||
[name='eq:x1bar']
|
||||
x1bar = x1bar(-1) + ex1bar;
|
||||
|
||||
[name='eq:x2bar']
|
||||
x2bar = x2bar(-1) + ex2bar;
|
||||
|
||||
[name='zpac']
|
||||
diff(z) = lambda*(e_c_m*(x1(-1)-z(-1)) + c_z_1*diff(z(-1)) + c_z_2*diff(z(-2)) + pac_expectation(pacman) + c_z_s*s + c_z_dv*diff(v) ) + (1-lambda)*( y + x) + c_z_u*u + c_z_dx2*diff(x2) + ez;
|
||||
|
||||
end;
|
||||
|
||||
shocks;
|
||||
var ex1 = 1.0;
|
||||
var ex2 = 1.0;
|
||||
var ex1bar = 1.0;
|
||||
var ex2bar = 1.0;
|
||||
var ez = 1.0;
|
||||
var ey = 0.1;
|
||||
var ex = 0.1;
|
||||
var eu = 0.05;
|
||||
var ev = 0.05;
|
||||
var es = 0.07;
|
||||
end;
|
||||
|
||||
// Initialize the PAC model (build the Companion VAR representation for the auxiliary model).
|
||||
pac.initialize('pacman');
|
||||
|
||||
// Update the parameters of the PAC expectation model (h0 and h1 vectors).
|
||||
pac.update.expectation('pacman');
|
||||
|
||||
// Set initial conditions to zero. Please use more sensible values if any...
|
||||
initialconditions = dseries(zeros(10, M_.endo_nbr+M_.exo_nbr), 2000Q1, vertcat(M_.endo_names,M_.exo_names));
|
||||
|
||||
// Simulate the model for 500 periods
|
||||
TrueData = simul_backward_model(initialconditions, 5000);
|
||||
|
||||
// Define a structure describing the parameters to be estimated (with initial conditions).
|
||||
clear eparams
|
||||
eparams.e_c_m = .9;
|
||||
eparams.c_z_1 = .5;
|
||||
eparams.c_z_2 = .2;
|
||||
eparams.c_z_s = .9;
|
||||
eparams.c_z_dx2 = .1;
|
||||
eparams.lambda = .7;
|
||||
|
||||
// Define the dataset used for estimation
|
||||
edata = TrueData;
|
||||
edata.ez = dseries(NaN(TrueData.nobs, 1), 2000Q1, 'ez');
|
||||
pac.estimate.iterative_ols('zpac', eparams, edata, 2005Q1:2005Q1+4000);
|
||||
|
||||
e_c_m_iterative_ols = M_.params(strmatch('e_c_m', M_.param_names, 'exact'));
|
||||
c_z_1_iterative_ols = M_.params(strmatch('c_z_1', M_.param_names, 'exact'));
|
||||
c_z_2_iterative_ols = M_.params(strmatch('c_z_2', M_.param_names, 'exact'));
|
||||
c_z_s_iterative_ols = M_.params(strmatch('c_z_s', M_.param_names, 'exact'));
|
||||
c_z_dx2_iterative_ols = M_.params(strmatch('c_z_dx2', M_.param_names, 'exact'));
|
||||
lambda_iterative_ols = M_.params(strmatch('lambda', M_.param_names, 'exact'));
|
||||
|
||||
disp(sprintf('Estimate of e_c_m: %f', e_c_m_iterative_ols))
|
||||
disp(sprintf('Estimate of c_z_1: %f', c_z_1_iterative_ols))
|
||||
disp(sprintf('Estimate of c_z_2: %f', c_z_2_iterative_ols))
|
||||
disp(sprintf('Estimate of c_z_s: %f', c_z_s_iterative_ols))
|
||||
disp(sprintf('Estimate of c_z_dx2: %f', c_z_dx2_iterative_ols))
|
||||
disp(sprintf('Estimate of lambda: %f', lambda_iterative_ols))
|
||||
|
||||
skipline(2)
|
||||
|
||||
// Define a structure describing the parameters to be estimated (with initial conditions).
|
||||
clear eparams
|
||||
eparams.e_c_m = .9;
|
||||
eparams.c_z_1 = .5;
|
||||
eparams.c_z_2 = .2;
|
||||
eparams.c_z_s = .9;
|
||||
eparams.c_z_dx2 = .1;
|
||||
eparams.lambda = .7;
|
||||
|
||||
edata = TrueData;
|
||||
edata.ez = dseries(NaN(TrueData.nobs, 1), 2000Q1, 'ez');
|
||||
pac.estimate.nls('zpac', eparams, edata, 2005Q1:2005Q1+4000, 'fmincon');
|
||||
|
||||
e_c_m_nls = M_.params(strmatch('e_c_m', M_.param_names, 'exact'));
|
||||
c_z_1_nls = M_.params(strmatch('c_z_1', M_.param_names, 'exact'));
|
||||
c_z_2_nls = M_.params(strmatch('c_z_2', M_.param_names, 'exact'));
|
||||
c_z_s_nls = M_.params(strmatch('c_z_s', M_.param_names, 'exact'));
|
||||
c_z_dx2_nls = M_.params(strmatch('c_z_dx2', M_.param_names, 'exact'));
|
||||
lambda_nls = M_.params(strmatch('lambda', M_.param_names, 'exact'));
|
||||
|
||||
disp(sprintf('Estimate of e_c_m: %f', e_c_m_nls))
|
||||
disp(sprintf('Estimate of c_z_1: %f', c_z_1_nls))
|
||||
disp(sprintf('Estimate of c_z_2: %f', c_z_2_nls))
|
||||
disp(sprintf('Estimate of c_z_s: %f', c_z_s_nls))
|
||||
disp(sprintf('Estimate of c_z_dx2: %f', c_z_dx2_nls))
|
||||
disp(sprintf('Estimate of lambda: %f', lambda_nls))
|
||||
|
||||
if abs(e_c_m_nls-e_c_m_iterative_ols)>.01
|
||||
error('Iterative OLS and NLS do not provide consistent estimates (e_c_m)')
|
||||
end
|
||||
|
||||
if abs(c_z_1_nls-c_z_1_iterative_ols)>.01
|
||||
error('Iterative OLS and NLS do not provide consistent estimates (c_z_1)')
|
||||
end
|
||||
|
||||
if abs(c_z_2_nls-c_z_2_iterative_ols)>.01
|
||||
error('Iterative OLS and NLS do not provide consistent estimates (c_z_2)')
|
||||
end
|
||||
|
||||
if abs(c_z_s_nls-c_z_s_iterative_ols)>.01
|
||||
error('Iterative OLS and NLS do not provide consistent estimates (c_z_s)')
|
||||
end
|
||||
|
||||
if abs(c_z_dx2_nls-c_z_dx2_iterative_ols)>.01
|
||||
error('Iterative OLS and NLS do not provide consistent estimates (c_z_dx2)')
|
||||
end
|
||||
|
||||
if abs(lambda_nls-lambda_iterative_ols)>.01
|
||||
error('Iterative OLS and NLS do not provide consistent estimates (lambda)')
|
||||
end
|
||||
Loading…
Reference in New Issue