Rewrote evaluate routine.
- Can handle more than one equation.
- Can handle identities.
- Forbids dynamic equations.
- Can handle following LHS y, diff(y), diff(diff(y)), log(y), diff(log(y)) and
diff(diff(log(y))), other transformations will result in an error.
- Added integration tests.
Remark 1. In the integration tests I compare the values returned by the
evaluate routine with the values computed with the simulation
routines. Normally the discrepancies should be small, but this is not
the case when the endogenous variable appear under a log on the
LHS. My current conclusion is that this has more to do with the
cumulation of the accuracy errors in the simulation routine (a
sequence of Newton algorithms) rather than with the evaluate routine.
Remark 2. Currently the only allowed nonlinear transformation on the LHS
endogenous variable is the log. It is not difficult to generalise, at
some point I had all the matlab functions allowed by Dynare,
but this would complicate the code for not much gain.
time-shift
Stéphane Adjemia (Scylla)
parent
de458022e5
commit
8740355407
GPG Key ID:
A6D44CB9C64CE77B
|
|
@ -1,4 +1,4 @@
|
|||
function nds = evaluate(ds, eqtag)
|
||||
function ds = evaluate(ds, eqtags, firstperiod, lastperiod)
|
||||
|
||||
% Copyright (C) 2019 Dynare Team
|
||||
%
|
||||
|
|
@ -19,47 +19,136 @@ function nds = evaluate(ds, eqtag)
|
|||
|
||||
global M_
|
||||
|
||||
% Get equation
|
||||
[LHS, RHS] = get_lhs_and_rhs(eqtag, M_, true);
|
||||
|
||||
% Parse equation and return list of parameters, endogenous and exogenous variables.
|
||||
[pnames, enames, xnames] = get_variables_and_parameters_in_equation(LHS, RHS, M_);
|
||||
|
||||
% Load parameter values.
|
||||
commands = sprintf('%s = %s;', pnames{1}, num2str(M_.params(strcmp(pnames{1},M_.param_names)), 16));
|
||||
for i=2:length(pnames)
|
||||
commands = sprintf('%s %s = %s;', commands, pnames{i}, ...
|
||||
num2str(M_.params(strcmp(pnames{i},M_.param_names)), 16));
|
||||
if ischar(eqtags)
|
||||
eqtags = {eqtags};
|
||||
end
|
||||
|
||||
list_of_expression_tokens = {'+', '-', '*', '/', '^', ...
|
||||
'exp(', 'log(', 'sqrt(', 'abs(', 'sign(', ...
|
||||
'sin(', 'cos(', 'tan(', 'asin(', 'acos(', 'atan(', ...
|
||||
'min(', 'max(', ...
|
||||
'normcdf(', 'normpdf(', 'erf(', ...
|
||||
'diff(', 'adl(', ')'};
|
||||
|
||||
if nargin<3
|
||||
range = ds.dates(1):ds.dates(end);
|
||||
elseif nargin<4
|
||||
range = firstperiod:ds.dates(end);
|
||||
else
|
||||
range = firstperiod:lastperiod;
|
||||
end
|
||||
eval(commands)
|
||||
|
||||
|
||||
% Substitute endogenous variable x with ds.x
|
||||
enames = unique(enames);
|
||||
for i=1:length(enames)
|
||||
if ismember(enames{i}, ds.name)
|
||||
RHS = regexprep(RHS, sprintf('\\<(%s)\\>', enames{i}), sprintf('ds.%s', enames{i}));
|
||||
else
|
||||
error('Endogenous variable %s is unknown in dseries objet.', enames{i})
|
||||
for i=1:length(eqtags)
|
||||
% Get equation
|
||||
[LHS, RHS] = get_lhs_and_rhs(eqtags{i}, M_, true);
|
||||
% Parse equation and return list of parameters, endogenous and exogenous variables.
|
||||
[pnames, enames, xnames] = get_variables_and_parameters_in_equation(LHS, RHS, M_);
|
||||
% Load parameter values.
|
||||
if ~isempty(pnames)
|
||||
commands = sprintf('%s = %s;', pnames{1}, num2str(M_.params(strcmp(pnames{1},M_.param_names)), 16));
|
||||
for j=2:length(pnames)
|
||||
commands = sprintf('%s %s = %s;', commands, pnames{i}, num2str(M_.params(strcmp(pnames{i},M_.param_names)), 16));
|
||||
end
|
||||
eval(commands)
|
||||
end
|
||||
end
|
||||
|
||||
% Substitute exogenous variable x with ds.x, except if
|
||||
if ~isfield(M_, 'simulation_exo_names')
|
||||
M_.simulation_exo_names = M_.exo_names;
|
||||
end
|
||||
xnames = unique(xnames);
|
||||
for i=1:length(xnames)
|
||||
if ismember(xnames{i}, M_.simulation_exo_names)
|
||||
if ismember(xnames{i}, ds.name)
|
||||
RHS = regexprep(RHS, sprintf('\\<(%s)\\>', xnames{i}), sprintf('ds.%s', xnames{i}));
|
||||
% Remove repetitions in enames
|
||||
enames = unique(enames);
|
||||
% Test if LHS is an endogenous variable
|
||||
is_lhs_expression = ~ismember(LHS, enames);
|
||||
if is_lhs_expression
|
||||
variable = strsplit(LHS, list_of_expression_tokens);
|
||||
variable(cellfun(@(x) all(isempty(x)), variable)) = [];
|
||||
if length(variable)>1
|
||||
error('It is not possible to have an expression with more than one variable on the LHS (%s).', LHS)
|
||||
else
|
||||
RHS = regexprep(RHS, sprintf('\\<(%s)\\>', xnames{i}), '0');
|
||||
warning('Exogenous variable %s is unknown in dseries objet. Assign zero value.', xnames{i})
|
||||
if isequal(LHS, sprintf('log(%s)', variable{1}))
|
||||
transform = {'exp'};
|
||||
elseif isequal(LHS, sprintf('diff(%s)', variable{1}))
|
||||
transform = {'cumsum'};
|
||||
elseif isequal(LHS, sprintf('diff(log(%s))', variable{1}))
|
||||
transform = {'cumsum', 'exp'};
|
||||
elseif isequal(LHS, sprintf('diff(diff(%s))', variable{1}))
|
||||
transform = {'cumsum', 'cumsum'};
|
||||
elseif isequal(LHS, sprintf('diff(diff(log(%s)))', variable{1}))
|
||||
transform = {'cumsum', 'cumsum', 'exp'};
|
||||
else
|
||||
error('Cannot proceed with provided LHS (%s in %s)', LHS, eqtags{i})
|
||||
end
|
||||
lhs = variable{1};
|
||||
end
|
||||
else
|
||||
RHS = regexprep(RHS, sprintf('(\\ *)(+)(\\ *)%s', xnames{i}), '');
|
||||
lhs = LHS;
|
||||
transform = {};
|
||||
end
|
||||
end
|
||||
|
||||
nds = eval(RHS);
|
||||
% Throw an error if the equation is dynamic.
|
||||
if exactcontains(RHS, lhs)
|
||||
error('RHS cannot contain LHS variable (%s in %s)', lhs, eqtags{i})
|
||||
end
|
||||
% Substitute endogenous variable x with ds.
|
||||
for j=1:length(enames)
|
||||
if ismember(enames{j}, ds.name)
|
||||
RHS = exactstrrep(RHS, enames{j}, sprintf('ds(range).%s', enames{j}));
|
||||
else
|
||||
error('Endogenous variable %s is unknown in dseries objet.', enames{j})
|
||||
end
|
||||
end
|
||||
% Substitute exogenous variable x with ds.x, except if
|
||||
if ~isfield(M_, 'simulation_exo_names')
|
||||
M_.simulation_exo_names = M_.exo_names;
|
||||
end
|
||||
xnames = unique(xnames);
|
||||
for j=1:length(xnames)
|
||||
if ismember(xnames{j}, M_.simulation_exo_names)
|
||||
if ismember(xnames{j}, ds.name)
|
||||
RHS = exactstrrep(RHS, xnames{j}, sprintf('ds(range).%s', xnames{j}));
|
||||
else
|
||||
RHS = exactstrrep(RHS, xnames{j}, '0');
|
||||
warning('Exogenous variable %s is unknown in dseries objet. Assign zero value.', xnames{j})
|
||||
end
|
||||
else
|
||||
RHS = regexprep(RHS, sprintf('(\\ *)(+)(\\ *)%s', xnames{j}), '');
|
||||
end
|
||||
end
|
||||
if isempty(transform)
|
||||
ds{LHS}(range) = eval(RHS);
|
||||
else
|
||||
tmp = eval(RHS);
|
||||
switch length(transform)
|
||||
case 1
|
||||
if isequal(transform{1}, 'cumsum')
|
||||
ds{lhs}(range) = cumsum(tmp)+ds{lhs}(range(1)-1).data;
|
||||
else
|
||||
ds{lhs}(range) = feval(transform{1}, tmp);
|
||||
end
|
||||
case 2
|
||||
if isequal(transform{2}, 'cumsum')
|
||||
% Squared first difference.
|
||||
t2 = zeros(length(range), 1);
|
||||
for t = 1:length(range)
|
||||
t2(t) = 2*ds{lhs}(range(t)-1).data-ds{lhs}(range(t)-2).data+tmp(range(t)).data;
|
||||
end
|
||||
ds{lhs}(range) = dseries(t2, range(1));
|
||||
else
|
||||
t2 = zeros(length(range), 1);
|
||||
for t = 1:length(range)
|
||||
t1 = feval(transform{2}, log(ds{lhs}(range(t)-1))+tmp(range(t)).data );
|
||||
t2(t) = t1.data;
|
||||
end
|
||||
ds{lhs}(range) = dseries(t2, range(1));
|
||||
% $$$ % The commented version below is more efficient but the discrepancy with what is returned by simulating
|
||||
% $$$ % the model is much bigger (see pac/trend-component-28/example4.mod).
|
||||
% $$$ tmp = cumsum(tmp)+log(ds{lhs}(range(1)-1).data);
|
||||
% $$$ ds{lhs}(range) = feval(transform{2}, tmp);
|
||||
end
|
||||
case 3
|
||||
t2 = zeros(length(range), 1);
|
||||
for t = 1:length(range)
|
||||
t2(t) = feval(transform{3}, 2*log(ds{lhs}(range(t)-1).data)-log(ds{lhs}(range(t)-2).data)+tmp(range(t)).data);
|
||||
end
|
||||
ds{lhs}(range) = dseries(t2, range(1));
|
||||
otherwise
|
||||
error('More than 3 unary ops. in LHS not implemented.')
|
||||
end
|
||||
end
|
||||
end
|
||||
|
|
@ -0,0 +1,22 @@
|
|||
function b = exactcontains(str, word)
|
||||
|
||||
% Same as contains but with exact word matching.
|
||||
|
||||
% Copyright (C) 2019 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/>.
|
||||
|
||||
b = ~isempty(regexp(str, sprintf('\\<%s\\>', word)));
|
||||
|
|
@ -444,6 +444,10 @@ ECB_MODFILES = \
|
|||
pac/trend-component-27/example.mod \
|
||||
pac/trend-component-28/example1.mod \
|
||||
pac/trend-component-28/example2.mod \
|
||||
pac/trend-component-28/example3.mod \
|
||||
pac/trend-component-28/example4.mod \
|
||||
pac/trend-component-28/example5.mod \
|
||||
pac/trend-component-28/example6.mod \
|
||||
pac/trend-component-29/example2.mod \
|
||||
pac/trend-component-29/example2.mod \
|
||||
write/example1.mod \
|
||||
|
|
|
|||
|
|
@ -1,15 +1,27 @@
|
|||
#!/bin/sh
|
||||
|
||||
rm -f *.mat
|
||||
rm -f *.m
|
||||
rm -f *.dat
|
||||
rm -f *.log
|
||||
|
||||
rm -rf example1
|
||||
rm -rf +example1
|
||||
rm -f example1.log
|
||||
rm -f *.mat
|
||||
rm -f *.m
|
||||
rm -f *.dat
|
||||
|
||||
rm -rf example2
|
||||
rm -rf +example2
|
||||
rm -f example2.log
|
||||
rm -f *.mat
|
||||
rm -f *.m
|
||||
rm -f *.dat
|
||||
|
||||
rm -rf example3
|
||||
rm -rf +example3
|
||||
|
||||
rm -rf example4
|
||||
rm -rf +example4
|
||||
|
||||
rm -rf example5
|
||||
rm -rf +example5
|
||||
|
||||
rm -rf example6
|
||||
rm -rf +example6
|
||||
|
||||
|
||||
|
||||
|
|
|
|||
|
|
@ -1,6 +1,8 @@
|
|||
// --+ options: json=compute, stochastic +--
|
||||
|
||||
var x1 x2 x1bar x2bar z y x u v s;
|
||||
PLOTS = true;
|
||||
|
||||
var x1 x2 x3 x1bar x2bar z y x u v s;
|
||||
|
||||
varexo ex1
|
||||
ex2
|
||||
|
|
@ -18,7 +20,8 @@ parameters
|
|||
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 cx cy beta
|
||||
lambda;
|
||||
lambda
|
||||
px3;
|
||||
|
||||
rho_1 = .9;
|
||||
rho_2 = -.2;
|
||||
|
|
@ -52,6 +55,8 @@ cy = 1.0;
|
|||
|
||||
lambda = 0.5; // Share of optimizing agents.
|
||||
|
||||
px3 = -.1;
|
||||
|
||||
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);
|
||||
|
|
@ -79,6 +84,9 @@ diff(x1) = a_x1_0*(x1(-1)-x1bar(-1)) + a_x1_1*diff(x1(-1)) + a_x1_2*diff(x1(-2))
|
|||
[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:x3']
|
||||
x3 = px3*x + y ;
|
||||
|
||||
[name='eq:x1bar']
|
||||
x1bar = x1bar(-1) + ex1bar;
|
||||
|
||||
|
|
@ -115,9 +123,26 @@ initialconditions = dseries(init, 2000Q1, vertcat(M_.endo_names,M_.exo_names));
|
|||
|
||||
// Simulate the model for 500 periods
|
||||
TrueData = simul_backward_model(initialconditions, 500);
|
||||
TrueData_ = copy(TrueData);
|
||||
|
||||
TrueData.x1bis = equation.evaluate(TrueData, 'eq:x1');
|
||||
TrueData = equation.evaluate(TrueData, 'eq:x3', 2004Q1);
|
||||
|
||||
if max(abs(TrueData.x1bis.data(5:end)-diff(TrueData.x1.data(4:end))))>1e-8
|
||||
if PLOTS
|
||||
|
||||
figure(1)
|
||||
plot(TrueData_.x.data(12:end), TrueData_.x3.data(12:end)-TrueData.x3.data(12:end), '.k', 'linewidth', 2);
|
||||
|
||||
figure(2)
|
||||
plot(TrueData_.x.data(12:end), '-k', 'linewidth', 2);
|
||||
|
||||
figure(3)
|
||||
plot([TrueData_.x3.data(12:end)-TrueData.x3.data(12:end)], '-k', 'linewidth', 2);
|
||||
|
||||
end
|
||||
|
||||
|
||||
fprintf('Max. abs. error is %s.\n', num2str(max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end))), 16));
|
||||
|
||||
if max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end)))>1e-12
|
||||
error('equation.evaluate() returned wrong values.')
|
||||
end
|
||||
|
|
@ -1,9 +1,11 @@
|
|||
// --+ options: json=compute, stochastic +--
|
||||
|
||||
var x1 x2 x1bar x2bar z y x u v s;
|
||||
PLOTS = true;
|
||||
|
||||
var x1 x2 x3 x1bar x2bar z y x u v s;
|
||||
|
||||
varexo ex1
|
||||
ex2 (used='estimationonly')
|
||||
ex2
|
||||
ex1bar (used='estimationonly')
|
||||
ex2bar (used='estimationonly')
|
||||
ez
|
||||
|
|
@ -18,7 +20,8 @@ parameters
|
|||
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 cx cy beta
|
||||
lambda;
|
||||
lambda
|
||||
px3;
|
||||
|
||||
rho_1 = .9;
|
||||
rho_2 = -.2;
|
||||
|
|
@ -52,6 +55,8 @@ cy = 1.0;
|
|||
|
||||
lambda = 0.5; // Share of optimizing agents.
|
||||
|
||||
px3 = -.1;
|
||||
|
||||
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);
|
||||
|
|
@ -79,6 +84,9 @@ diff(x1) = a_x1_0*(x1(-1)-x1bar(-1)) + a_x1_1*diff(x1(-1)) + a_x1_2*diff(x1(-2))
|
|||
[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:x3']
|
||||
log(x3) = x ;
|
||||
|
||||
[name='eq:x1bar']
|
||||
x1bar = x1bar(-1) + ex1bar;
|
||||
|
||||
|
|
@ -110,16 +118,30 @@ pac.initialize('pacman');
|
|||
pac.update.expectation('pacman');
|
||||
|
||||
// Set initial conditions to zero for non logged variables, and one for logged variables
|
||||
init = zeros(10, M_.endo_nbr+M_.exo_nbr);
|
||||
init = rand(10, M_.endo_nbr+M_.exo_nbr);
|
||||
initialconditions = dseries(init, 2000Q1, vertcat(M_.endo_names,M_.exo_names));
|
||||
|
||||
// Simulate the model for 500 periods
|
||||
TrueData = simul_backward_model(initialconditions, 500);
|
||||
TrueData_ = copy(TrueData);
|
||||
|
||||
TrueData.diff_x2_fit = equation.evaluate(TrueData, 'eq:x2');
|
||||
TrueData = equation.evaluate(TrueData, 'eq:x3', 2004Q1);
|
||||
|
||||
dx2 = diff(TrueData.x2);
|
||||
if PLOTS
|
||||
|
||||
if max(abs((TrueData.diff_x2_fit.data+TrueData.ex2.data-dx2.data)))>1e-8
|
||||
figure(1)
|
||||
plot(TrueData_.x.data(12:end), TrueData_.x3.data(12:end)-TrueData.x3.data(12:end), '.k', 'linewidth', 2);
|
||||
|
||||
figure(2)
|
||||
plot(TrueData_.x.data(12:end), '-k', 'linewidth', 2);
|
||||
|
||||
figure(3)
|
||||
plot([TrueData_.x3.data(12:end)-TrueData.x3.data(12:end)], '-k', 'linewidth', 2);
|
||||
|
||||
end
|
||||
|
||||
fprintf('Max. abs. error is %s.\n', num2str(max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end))), 16));
|
||||
|
||||
if max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end)))>1e-5
|
||||
error('equation.evaluate() returned wrong values.')
|
||||
end
|
||||
|
|
@ -0,0 +1,147 @@
|
|||
// --+ options: json=compute, stochastic +--
|
||||
|
||||
PLOTS = true;
|
||||
|
||||
var x1 x2 x3 x1bar x2bar z y x u v s;
|
||||
|
||||
varexo ex1
|
||||
ex2
|
||||
ex1bar (used='estimationonly')
|
||||
ex2bar (used='estimationonly')
|
||||
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 cx cy beta
|
||||
lambda
|
||||
px3;
|
||||
|
||||
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;
|
||||
cx = 1.0;
|
||||
cy = 1.0;
|
||||
|
||||
|
||||
lambda = 0.5; // Share of optimizing agents.
|
||||
|
||||
px3 = -.1;
|
||||
|
||||
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:x3']
|
||||
diff(x3) = x ;
|
||||
|
||||
[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)*( cy*y + cx*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 for non logged variables, and one for logged variables
|
||||
init = randn(10, M_.endo_nbr+M_.exo_nbr);
|
||||
initialconditions = dseries(init, 2000Q1, vertcat(M_.endo_names,M_.exo_names));
|
||||
|
||||
// Simulate the model for 500 periods
|
||||
TrueData = simul_backward_model(initialconditions, 500);
|
||||
TrueData_ = copy(TrueData);
|
||||
|
||||
TrueData = equation.evaluate(TrueData, 'eq:x3', 2004Q1);
|
||||
|
||||
if PLOTS
|
||||
|
||||
figure(1)
|
||||
plot(TrueData_.x.data(12:end), TrueData_.x3.data(12:end)-TrueData.x3.data(12:end), '.k', 'linewidth', 2);
|
||||
|
||||
figure(2)
|
||||
plot(TrueData_.x.data(12:end), '-k', 'linewidth', 2);
|
||||
|
||||
figure(3)
|
||||
plot([TrueData_.x3.data(12:end)-TrueData.x3.data(12:end)], '-k', 'linewidth', 2);
|
||||
|
||||
end
|
||||
|
||||
fprintf('Max. abs. error is %s.\n', num2str(max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end))), 16));
|
||||
|
||||
if max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end)))>1e-12
|
||||
error('equation.evaluate() returned wrong values.')
|
||||
end
|
||||
|
|
@ -0,0 +1,147 @@
|
|||
// --+ options: json=compute, stochastic +--
|
||||
|
||||
PLOTS = true;
|
||||
|
||||
var x1 x2 x3 x1bar x2bar z y x u v s;
|
||||
|
||||
varexo ex1
|
||||
ex2
|
||||
ex1bar (used='estimationonly')
|
||||
ex2bar (used='estimationonly')
|
||||
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 cx cy beta
|
||||
lambda
|
||||
px3;
|
||||
|
||||
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;
|
||||
cx = 1.0;
|
||||
cy = 1.0;
|
||||
|
||||
|
||||
lambda = 0.5; // Share of optimizing agents.
|
||||
|
||||
px3 = -.1;
|
||||
|
||||
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:x3']
|
||||
diff(log(x3)) = x ;
|
||||
|
||||
[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)*( cy*y + cx*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.01;
|
||||
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 for non logged variables, and one for logged variables
|
||||
init = ones(10, M_.endo_nbr+M_.exo_nbr);
|
||||
initialconditions = dseries(init, 2000Q1, vertcat(M_.endo_names,M_.exo_names));
|
||||
|
||||
// Simulate the model for 500 periods
|
||||
TrueData = simul_backward_model(initialconditions, 5000);
|
||||
TrueData_ = copy(TrueData);
|
||||
|
||||
TrueData = equation.evaluate(TrueData, 'eq:x3', 2004Q1);
|
||||
|
||||
if PLOTS
|
||||
|
||||
figure(1)
|
||||
plot(TrueData_.x.data(12:end), TrueData_.x3.data(12:end)-TrueData.x3.data(12:end), '.k', 'linewidth', 2);
|
||||
|
||||
figure(2)
|
||||
plot(TrueData_.x.data(12:end), '-k', 'linewidth', 2);
|
||||
|
||||
figure(3)
|
||||
plot([TrueData_.x3.data(12:end)-TrueData.x3.data(12:end)], '-k', 'linewidth', 2);
|
||||
|
||||
end
|
||||
|
||||
fprintf('Max. abs. error is %s.\n', num2str(max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end))), 16));
|
||||
|
||||
if max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end)))>1e-5
|
||||
error('equation.evaluate() returned wrong values.')
|
||||
end
|
||||
|
|
@ -0,0 +1,147 @@
|
|||
// --+ options: json=compute, stochastic +--
|
||||
|
||||
PLOTS = true;
|
||||
|
||||
var x1 x2 x3 x1bar x2bar z y x u v s;
|
||||
|
||||
varexo ex1
|
||||
ex2
|
||||
ex1bar (used='estimationonly')
|
||||
ex2bar (used='estimationonly')
|
||||
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 cx cy beta
|
||||
lambda
|
||||
px3;
|
||||
|
||||
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;
|
||||
cx = 1.0;
|
||||
cy = 1.0;
|
||||
|
||||
|
||||
lambda = 0.5; // Share of optimizing agents.
|
||||
|
||||
px3 = -.1;
|
||||
|
||||
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:x3']
|
||||
diff(diff(x3)) = x ;
|
||||
|
||||
[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)*( cy*y + cx*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.01;
|
||||
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 for non logged variables, and one for logged variables
|
||||
init = ones(10, M_.endo_nbr+M_.exo_nbr);
|
||||
initialconditions = dseries(init, 2000Q1, vertcat(M_.endo_names,M_.exo_names));
|
||||
|
||||
// Simulate the model for 500 periods
|
||||
TrueData = simul_backward_model(initialconditions, 500);
|
||||
TrueData_ = copy(TrueData);
|
||||
|
||||
TrueData = equation.evaluate(TrueData, 'eq:x3', 2004Q1);
|
||||
|
||||
if PLOTS
|
||||
|
||||
figure(1)
|
||||
plot(TrueData_.x.data(12:end), TrueData_.x3.data(12:end)-TrueData.x3.data(12:end), '.k', 'linewidth', 2);
|
||||
|
||||
figure(2)
|
||||
plot(TrueData_.x.data(12:end), '-k', 'linewidth', 2);
|
||||
|
||||
figure(3)
|
||||
plot([TrueData_.x3.data(12:end)-TrueData.x3.data(12:end)], '-k', 'linewidth', 2);
|
||||
|
||||
end
|
||||
|
||||
fprintf('Max. abs. error is %s.\n', num2str(max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end))), 16));
|
||||
|
||||
if max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end)))>1e-12
|
||||
error('equation.evaluate() returned wrong values.')
|
||||
end
|
||||
|
|
@ -0,0 +1,147 @@
|
|||
// --+ options: json=compute, stochastic +--
|
||||
|
||||
PLOTS = true;
|
||||
|
||||
var x1 x2 x3 x1bar x2bar z y x u v s;
|
||||
|
||||
varexo ex1
|
||||
ex2
|
||||
ex1bar (used='estimationonly')
|
||||
ex2bar (used='estimationonly')
|
||||
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 cx cy beta
|
||||
lambda
|
||||
px3;
|
||||
|
||||
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;
|
||||
cx = 1.0;
|
||||
cy = 1.0;
|
||||
|
||||
|
||||
lambda = 0.5; // Share of optimizing agents.
|
||||
|
||||
px3 = -.1;
|
||||
|
||||
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:x3']
|
||||
diff(diff(log(x3))) = x ;
|
||||
|
||||
[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)*( cy*y + cx*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.01;
|
||||
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 for non logged variables, and one for logged variables
|
||||
init = ones(10, M_.endo_nbr+M_.exo_nbr);
|
||||
initialconditions = dseries(init, 2000Q1, vertcat(M_.endo_names,M_.exo_names));
|
||||
|
||||
// Simulate the model for 500 periods
|
||||
TrueData = simul_backward_model(initialconditions, 50);
|
||||
TrueData_ = copy(TrueData);
|
||||
|
||||
TrueData = equation.evaluate(TrueData, 'eq:x3', 2004Q1);
|
||||
|
||||
if PLOTS
|
||||
|
||||
figure(1)
|
||||
plot(TrueData_.x.data(12:end), TrueData_.x3.data(12:end)-TrueData.x3.data(12:end), '.k', 'linewidth', 2);
|
||||
|
||||
figure(2)
|
||||
plot(TrueData_.x.data(12:end), '-k', 'linewidth', 2);
|
||||
|
||||
figure(3)
|
||||
plot([TrueData_.x3.data(12:end)-TrueData.x3.data(12:end)], '-k', 'linewidth', 2);
|
||||
|
||||
end
|
||||
|
||||
fprintf('Max. abs. error is %s.\n', num2str(max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end))), 16));
|
||||
|
||||
if max(abs(TrueData.x3.data(12:end)-TrueData_.x3.data(12:end)))>1e-5
|
||||
error('equation.evaluate() returned wrong values.')
|
||||
end
|
||||
Loading…
Reference in New Issue