ddaf2b546a20951c82bc4d6b5f904ba1542.
time-shift
Johannes Pfeifer 2013-07-01 22:19:30 +02:00
commit e77388af39
13 changed files with 434 additions and 140 deletions

@ -1 +1 @@
Subproject commit 35f6ad8d478cb7042d4523b9b636690757105a43
Subproject commit b642763676536c148ad827a0286862be87c5ae8e

View File

@ -2855,6 +2855,7 @@ steady;
@end deffn
@anchor{equation_tag_for_conditional_steady_state}
@node Replace some equations during steady state computations
@subsection Replace some equations during steady state computations
@ -4631,19 +4632,16 @@ Uses the diffuse Kalman filter (as described in
@cite{Durbin and Koopman (2001)} and @cite{Koopman and Durbin
(2003)}) to estimate models with non-stationary observed variables.
When @code{diffused_filter} is used the @code{lik_init} option of
When @code{diffuse_filter} is used the @code{lik_init} option of
@code{estimation} has no effect.
When there are nonstationary variables in a model, there is no unique
deterministic steady state. The user must supply a MATLAB/Octave
function that computes the steady state values of the stationary
variables in the model and returns dummy values for the nonstationary
ones. The function should be called with the name of the @file{.mod}
file followed by @file{_steadystate}. See @file{fs2000_steadystate.m}
in @file{examples} directory for an example.
When there are nonstationary exogenous variables in a model, there is no unique deterministic steady state. For instance, if productivity is a pure random walk:
Note that the nonstationary variables in the model must be integrated
processes (their first difference or k-difference must be stationary).
@math{a_t = a_{t-1} + e_t}
any value of @math{\bar a} of @math{a} is a deterministic steady state for productivity. Consequently, the model admits an infinity of steady states. In this situation, the user must help Dynare in selecting one steady state, except if zero is a trivial model's steady state, which happens when the @code{linear} option is used in the model declaration. The user can either provide the steady state to Dynare using a @code{steady_state_model} block (or writing a steady state file) if a closed form solution is available, @pxref{steady_state_model}, or specify some constraints on the steady state, @pxref{equation_tag_for_conditional_steady_state}, so that Dynare computes the steady state conditionally on some predefined levels for the non stationary variables. In both cases, the idea is to use dummy values for the steady state level of the exogenous non stationary variables.
Note that the nonstationary variables in the model must be integrated processes (their first difference or k-difference must be stationary).
@item selected_variables_only
Only run the smoother on the variables listed just after the

24
matlab/@dynDate/disp.m Normal file
View File

@ -0,0 +1,24 @@
function disp(d)
% Copyright (C) 2013 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/>.
if isempty(d)
fprintf('%s is an empty dynDate object.\n', inputname(1));
else
fprintf('%s = <dynDate: %s>\n', inputname(1), format(d));
end

View File

@ -17,4 +17,8 @@ function display(d)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
fprintf('%s = <dynDate: %s>\n', inputname(1), format(d));
if isempty(d)
fprintf('%s is an empty dynDate object.\n', inputname(1));
else
fprintf('%s = <dynDate: %s>\n', inputname(1), format(d));
end

49
matlab/@dynDate/horzcat.m Normal file
View File

@ -0,0 +1,49 @@
function a = horzcat(varargin)
%@info:
%! @deftypefn {Function file} {@var{a} =} horzcat (@var{b},@var{c}, ...)
%! @anchor{horzcat}
%! @sp 1
%! Method of the dynDate class.
%! @sp 1
%! Concatenate dynDate objects to form a dynDates object. This method overloads the horizontal concatenation operator, so that
%! two (or more) dynDate objects an be concatenated :
%!
%! a = [b, c, d];
%! @sp 2
%! @strong{Inputs}
%! @sp 1
%! @table @ @var
%! @item b
%! dynDate object, instantiated by @ref{dynDate}.
%! @item c
%! dynDate object, instantiated by @ref{dynDate}.
%! @end table
%! @sp 2
%! @strong{Outputs}
%! @sp 1
%! @table @var
%! @item a
%! dynDates object.
%! @end table
%! @end deftypefn
%@eod:
% Copyright (C) 2013 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/>.
a = dynDates(varargin{:});

49
matlab/@dynDates/disp.m Normal file
View File

@ -0,0 +1,49 @@
function disp(dd)
% Copyright (C) 2013 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/>.
if isempty(dd)
fprintf('%s is an empty dynDates object.\n', inputname(1));
return
end
max_displayed = 5;
first_displayed = 2;
fprintf('%s = <dynDates: ', inputname(1));
if dd.ndat<=max_displayed
for i=1:dd.ndat
fprintf(format(dynDate(dd.time(i,:),dd.freq)))
if i<dd.ndat
fprintf(', ')
else
fprintf('>\n')
end
end
else
for i=1:first_displayed
fprintf(format(dynDate(dd.time(i,:),dd.freq)))
fprintf(', ')
end
fprintf(' ..., ')
fprintf(format(dynDate(dd.time(dd.ndat-1,:),dd.freq)))
fprintf(', ')
fprintf(format(dynDate(dd.time(dd.ndat,:),dd.freq)))
fprintf('>\n')
end

View File

@ -17,6 +17,11 @@ function display(dd)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
if isempty(dd)
fprintf('%s is an empty dynDates object.\n', inputname(1));
return
end
max_displayed = 5;
first_displayed = 2;

View File

@ -29,10 +29,9 @@ end
switch format
case 'm'
if exist([basename, '.m'],'file')
fid = fopen([basename, '.new', '.m'],'w');
else
fid = fopen([basename, '.m'],'w');
copyfile([basename, '.m'],[basename, '.old.m'])
end
fid = fopen([basename, '.m'],'w');
fprintf(fid,'%% File created on %s.\n',datestr(now));
fprintf(fid,'\n');
fprintf(fid,'FREQ__ = %s;\n',num2str(A.freq));
@ -71,17 +70,14 @@ switch format
end
eval(str);
if exist([basename, '.mat'],'file')
save([basename '.new.mat'],'INIT__','FREQ__','NAMES__','TEX__',A.name{:});
else
save([basename '.mat'],'INIT__','FREQ__','NAMES__','TEX__',A.name{:});
copyfile([basename, '.mat'],[basename, '.old.mat'])
end
save([basename '.mat'],'INIT__','FREQ__','NAMES__','TEX__',A.name{:});
case 'csv'
if exist([basename, '.csv'],'file')
fid = fopen([basename, '.new', '.csv'],'w');
else
fid = fopen([basename, '.csv'],'w');
copyfile([basename, '.csv'],[basename, '.old.csv'])
end
fid = fopen([basename, '.csv'],'w');
fprintf(fid,',%s', A.name{:});
fprintf(fid,'\n');
for t=1:A.nobs

View File

@ -25,99 +25,168 @@ function A = subsasgn(A,S,B)
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
if length(S)>1
error('dynSeries::subsasgn: Wrong syntax!')
end
merge_dynSeries_objects = 1;
switch S.type
case '{}'
if ~isequal(numel(S.subs),numel(unique(S.subs)))
error('dynSeries::subsasgn: Wrong syntax!')
end
for i=1:numel(S.subs)
element = S.subs{i};
idArobase = strfind(element,'@');
if ~isempty(idArobase)
switch length(idArobase)
case 2
idComma = strfind(element(idArobase(1)+1:idArobase(2)-1),',');
if ~isempty(idComma)
elements = cell(1,numel(idComma)+1); j = 1;
expression = element(idArobase(1)+1:idArobase(2)-1);
while ~isempty(expression)
[token, expression] = strtok(expression,',');
elements(j) = {[element(1:idArobase(1)-1), token, element(idArobase(2)+1:end)]};
j = j + 1;
switch length(S)
case 1
switch S(1).type
case '{}' % Multiple variable selection.
if ~isequal(numel(S(1).subs),numel(unique(S(1).subs)))
error('dynSeries::subsasgn: Wrong syntax!')
end
for i=1:numel(S(1).subs)
element = S(1).subs{i};
idArobase = strfind(element,'@');
if ~isempty(idArobase)
switch length(idArobase)
case 2
idComma = strfind(element(idArobase(1)+1:idArobase(2)-1),',');
if ~isempty(idComma)
elements = cell(1,numel(idComma)+1); j = 1;
expression = element(idArobase(1)+1:idArobase(2)-1);
while ~isempty(expression)
[token, expression] = strtok(expression,',');
elements(j) = {[element(1:idArobase(1)-1), token, element(idArobase(2)+1:end)]};
j = j + 1;
end
S(1).subs = replace_object_in_a_one_dimensional_cell_array(S(1).subs, elements(:), i);
else
error('dynSeries::subsasgn: Wrong syntax, matlab''s regular expressions cannot be used here!')
end
case 4
idComma_1 = strfind(element(idArobase(1)+1:idArobase(2)-1),',');
idComma_2 = strfind(element(idArobase(3)+1:idArobase(4)-1),',');
if ~isempty(idComma_1)
elements = cell(1,(numel(idComma_1)+1)*(numel(idComma_2)+1)); j = 1;
expression_1 = element(idArobase(1)+1:idArobase(2)-1);
while ~isempty(expression_1)
[token_1, expression_1] = strtok(expression_1,',');
expression_2 = element(idArobase(3)+1:idArobase(4)-1);
while ~isempty(expression_2)
[token_2, expression_2] = strtok(expression_2,',');
elements(j) = {[element(1:idArobase(1)-1), token_1, element(idArobase(2)+1:idArobase(3)-1), token_2, element(idArobase(4)+1:end)]};
j = j+1;
end
end
S(1).subs = replace_object_in_a_one_dimensional_cell_array(S(1).subs, elements(:), i);
else
error('dynSeries::subsasgn: Wrong syntax, matlab''s regular expressions cannot be used here!')
end
otherwise
error('dynSeries::subsasgn: Wrong syntax!')
end
end
end
if ~isequal(length(S(1).subs),B.vobs)
error('dynSeries::subsasgn: Wrong syntax!')
end
if ~isequal(S(1).subs(:),B.name)
for i = 1:B.vobs
if ~isequal(S(1).subs{i},B.name{i})
% Rename a variable.
id = strmatch(S(1).subs{i},A.name);
if isempty(id)
% Add a new variable a change its name.
B.name(i) = {S(1).subs{i}};
B.tex(i) = {name2tex(S(1).subs{i})};
else
% Rename variable and change its content.
B.name(i) = A.name(id);
B.tex(i) = A.tex(id);
end
end
end
end
case '.' % Single variable selection.
if ~isequal(S(1).subs,B.name)
if ~isequal(S(1).subs,B.name{1})
% Rename a variable.
id = strmatch(S(1).subs,A.name);
if isempty(id)
% Add a new variable a change its name.
B.name(1) = {S(1).subs};
B.tex(1) = {name2tex(S(1).subs)};
else
% Rename variable and change its content.
B.name(1) = A.name(id);
B.tex(1) = A.tex(id);
end
end
end
case '()' % Date(s) selection
if isa(S(1).subs{1},'dynDates') || isa(S(1).subs{1},'dynDate')
[junk, tdx] = intersect(A.time.time,S(1).subs{1}.time,'rows');
if isa(B,'dynSeries')
[junk, tdy] = intersect(B.time.time,S(1).subs{1}.time,'rows');
if isempty(tdy)
error('dynSeries::subsasgn: Periods of the dynSeries objects on the left and right hand sides must intersect!')
end
if ~isequal(A.vobs,B.vobs)
error('dynSeries::subsasgn: Dimension error! The number of variables on the left and right hand side must match.')
end
A.data(tdx,:) = B.data(tdy,:);
elseif isnumeric(B)
merge_dynSeries_objects = 0;
if isequal(length(tdx),rows(B))
if isequal(columns(A.data),columns(B))
A.data(tdx,:) = B;
else
error('dynSeries::subsasgn: Dimension error! The number of variables on the left and right hand side must match.')
end
else
error('dynSeries::subsassgn: Dimension error! The number of periods on the left and right hand side must match.')
end
else
error('dynSeries::subsasgn: The object on the right hand side must be a dynSeries object or a numeric array!')
end
else
error('dynSeries::subsasgn: Wrong syntax!')
end
otherwise
error('dynSeries::subsasgn: Wrong syntax!')
end
case 2
merge_dynSeries_objects = 0;
if ((isequal(S(1).type,'{}') || isequal(S(1).type,'.')) && isequal(S(2).type,'()'))
sA = extract(A,S(1).subs{:});
if (isa(B,'dynSeries') && isequal(sA.vobs,B.vobs)) || (isnumeric(B) && isequal(sA.vobs,columns(B)))
if isa(S(2).subs{1},'dynDates') || isa(S(2).subs{1},'dynDate')
[junk, tdx] = intersect(sA.time.time,S(2).subs{1}.time,'rows');
if isa(B,'dynSeries')
[junk, tdy] = intersect(B.time.time,S(2).subs{1}.time,'rows');
if isempty(tdy)
error('dynSeries::subsasgn: Periods of the dynSeries objects on the left and right hand sides must intersect!')
end
S.subs = replace_object_in_a_one_dimensional_cell_array(S.subs, elements(:), i);
else
error('dynSeries::subsasgn: Wrong syntax, matlab''s regular expressions cannot be used here!')
end
case 4
idComma_1 = strfind(element(idArobase(1)+1:idArobase(2)-1),',');
idComma_2 = strfind(element(idArobase(3)+1:idArobase(4)-1),',');
if ~isempty(idComma_1)
elements = cell(1,(numel(idComma_1)+1)*(numel(idComma_2)+1)); j = 1;
expression_1 = element(idArobase(1)+1:idArobase(2)-1);
while ~isempty(expression_1)
[token_1, expression_1] = strtok(expression_1,',');
expression_2 = element(idArobase(3)+1:idArobase(4)-1);
while ~isempty(expression_2)
[token_2, expression_2] = strtok(expression_2,',');
elements(j) = {[element(1:idArobase(1)-1), token_1, element(idArobase(2)+1:idArobase(3)-1), token_2, element(idArobase(4)+1:end)]};
j = j+1;
sA.data(tdx,:) = B.data(tdy,:);
elseif isnumeric(B)
merge_dynSeries_objects = 0;
if isequal(length(tdx),rows(B))
if isequal(columns(sA.data),columns(B))
sA.data(tdx,:) = B;
else
error('dynSeries::subsasgn: Dimension error! The number of variables on the left and right hand side must match.')
end
else
error('dynSeries::subsassgn: Dimension error! The number of periods on the left and right hand side must match.')
end
S.subs = replace_object_in_a_one_dimensional_cell_array(S.subs, elements(:), i);
else
error('dynSeries::subsasgn: Wrong syntax, matlab''s regular expressions cannot be used here!')
error('dynSeries::subsasgn: The object on the right hand side must be a dynSeries object or a numeric array!')
end
otherwise
else
error('dynSeries::subsasgn: Wrong syntax!')
end
A = merge(A,sA);
else
error('dynSeries::subsasgn: Dimension error! The number of variables on the left and right hand side must match.')
end
end
if ~isequal(length(S.subs),B.vobs)
error('dynSeries::subsasgn: Wrong syntax!')
end
if ~isequal(S.subs(:),B.name)
for i = 1:B.vobs
if ~isequal(S.subs{i},B.name{i})
% Rename a variable.
id = strmatch(S.subs{i},A.name);
if isempty(id)
% Add a new variable a change its name.
B.name(i) = {S.subs{i}};
B.tex(i) = {name2tex(S.subs{i})};
else
% Rename variable and change its content.
B.name(i) = A.name(id);
B.tex(i) = A.tex(id);
end
end
end
end
case '.'
if ~isequal(S.subs,B.name)
if ~isequal(S.subs,B.name{1})
% Rename a variable.
id = strmatch(S.subs,A.name);
if isempty(id)
% Add a new variable a change its name.
B.name(1) = {S.subs};
B.tex(1) = {name2tex(S.subs)};
else
% Rename variable and change its content.
B.name(1) = A.name(id);
B.tex(1) = A.tex(id);
end
end
end
otherwise
error('dynSeries::subsasgn: Wrong syntax!')
end
A = merge(A,B);
if merge_dynSeries_objects
A = merge(A,B);
end
%@test:1
%$ % Define a datasets.
@ -369,7 +438,7 @@ A = merge(A,B);
%$ T = all(t);
%@eof:10
%@test:10
%@test:11
%$ % Define a datasets.
%$ A = rand(10,3); B = rand(10,5);
%$
@ -393,4 +462,66 @@ A = merge(A,B);
%$ %t(5) = dyn_assert(ts1.data,[B(:,1:2), A(:,3), B(:,3:4)],1e-15);
%$ end
%$ T = all(t);
%@eof:10
%@eof:11
%@test:12
%$ % Define a datasets.
%$ A = rand(40,3); B = rand(40,1);
%$
%$ % Instantiate two dynSeries object.
%$ ts1 = dynSeries(A,'1950Q1',{'A1';'A2';'A3'},[]);
%$ ts2 = dynSeries(B,'1950Q1',{'B1'},[]);
%$
%$ % modify first object.
%$ try
%$ d1 = dynDate('1950Q3');
%$ d2 = dynDate('1951Q3');
%$ rg = d1:d2;
%$ ts1{'A1'}(rg) = ts2{'B1'}(rg);
%$ t(1) = 1;
%$ catch
%$ t(1) = 0;
%$ end
%$
%$ % Instantiate a time series object.
%$ if t(1)
%$ t(2) = dyn_assert(ts1.vobs,3);
%$ t(3) = dyn_assert(ts1.nobs,40);
%$ t(4) = dyn_assert(ts1.name{2},'A2');
%$ t(5) = dyn_assert(ts1.name{1},'A1');
%$ t(6) = dyn_assert(ts1.name{3},'A3');
%$ t(7) = dyn_assert(ts1.data,[[A(1:2,1); B(3:7); A(8:end,1)], A(:,2:3)],1e-15);
%$ end
%$ T = all(t);
%@eof:12
%@test:13
%$ % Define a datasets.
%$ A = rand(40,3); B = rand(40,1);
%$
%$ % Instantiate two dynSeries object.
%$ ts1 = dynSeries(A,'1950Q1',{'A1';'A2';'A3'},[]);
%$ ts2 = dynSeries(B,'1950Q1',{'B1'},[]);
%$
%$ % modify first object.
%$ try
%$ d1 = dynDate('1950Q3');
%$ d2 = dynDate('1951Q3');
%$ rg = d1:d2;
%$ ts1{'A1'}(rg) = B(3:7);
%$ t(1) = 1;
%$ catch
%$ t(1) = 0;
%$ end
%$
%$ % Instantiate a time series object.
%$ if t(1)
%$ t(2) = dyn_assert(ts1.vobs,3);
%$ t(3) = dyn_assert(ts1.nobs,40);
%$ t(4) = dyn_assert(ts1.name{2},'A2');
%$ t(5) = dyn_assert(ts1.name{1},'A1');
%$ t(6) = dyn_assert(ts1.name{3},'A3');
%$ t(7) = dyn_assert(ts1.data,[[A(1:2,1); B(3:7); A(8:end,1)], A(:,2:3)],1e-15);
%$ end
%$ T = all(t);
%@eof:13

View File

@ -118,25 +118,50 @@ if options_.lik_init == 1 % Kalman filter
end
Pstar = lyapunov_symm(T,R*Q*transpose(R),options_.qz_criterium,options_.lyapunov_complex_threshold);
Pinf = [];
elseif options_.lik_init == 2 % Old Diffuse Kalman filter
elseif options_.lik_init == 2 % Old Diffuse Kalman filter
if kalman_algo ~= 2
kalman_algo = 1;
end
Pstar = options_.Harvey_scale_factor*eye(np);
Pinf = [];
elseif options_.lik_init == 3 % Diffuse Kalman filter
elseif options_.lik_init == 3 % Diffuse Kalman filter
if kalman_algo ~= 4
kalman_algo = 3;
end
[Z,ST,R1,QT,Pstar,Pinf] = schur_statespace_transformation(mf,T,R,Q,options_.qz_criterium);
elseif options_.lik_init == 4
% Start from the solution of the Riccati equation.
elseif options_.lik_init == 4 % Start from the solution of the Riccati equation.
[err, Pstar] = kalman_steady_state(transpose(T),R*Q*transpose(R),transpose(build_selection_matrix(mf,np,nobs)),H);
mexErrCheck('kalman_steady_state',err);
Pinf = [];
if kalman_algo~=2
kalman_algo = 1;
end
elseif options_.lik_init == 5 % Old diffuse Kalman filter only for the non stationary variables
[eigenvect, eigenv] = eig(T);
eigenv = diag(eigenv);
nstable = length(find(abs(abs(eigenv)-1) > 1e-7));
unstable = find(abs(abs(eigenv)-1) < 1e-7);
V = eigenvect(:,unstable);
indx_unstable = find(sum(abs(V),2)>1e-5);
stable = find(sum(abs(V),2)<1e-5);
nunit = length(eigenv) - nstable;
Pstar = options_.Harvey_scale_factor*eye(np);
if kalman_algo ~= 2
kalman_algo = 1;
end
R_tmp = R(stable, :);
T_tmp = T(stable,stable);
if options_.lyapunov_fp == 1
Pstar_tmp = lyapunov_symm(T_tmp,Q,options_.lyapunov_fixed_point_tol,options_.lyapunov_complex_threshold, 3, R_tmp);
elseif options_.lyapunov_db == 1
Pstar_tmp = disclyap_fast(T_tmp,R_tmp*Q*R_tmp',options_.lyapunov_doubling_tol);
elseif options_.lyapunov_srs == 1
Pstar_tmp = lyapunov_symm(T_tmp,Q,options_.lyapunov_fixed_point_tol,options_.lyapunov_complex_threshold, 4, R_tmp);
else
Pstar_tmp = lyapunov_symm(T_tmp,R_tmp*Q*R_tmp',options_.qz_criterium,options_.lyapunov_complex_threshold);
end
Pstar(stable, stable) = Pstar_tmp;
Pinf = [];
end
kalman_tol = options_.kalman_tol;
riccati_tol = options_.riccati_tol;

View File

@ -373,7 +373,6 @@ switch DynareOptions.lik_init
error(['diffuse filter: options_.kalman_algo can only be equal ' ...
'to 0 (default), 3 or 4'])
end
[Z,T,R,QT,Pstar,Pinf] = schur_statespace_transformation(Z,T,R,Q,DynareOptions.qz_criterium);
Zflag = 1;
% Run diffuse kalman filter on first periods.
@ -453,6 +452,32 @@ switch DynareOptions.lik_init
Pinf = [];
a = zeros(mm,1);
Zflag = 0;
case options_.lik_init == 5 % Old diffuse Kalman filter only for the non stationary variables
[eigenvect, eigenv] = eig(T);
eigenv = diag(eigenv);
nstable = length(find(abs(abs(eigenv)-1) > 1e-7));
unstable = find(abs(abs(eigenv)-1) < 1e-7);
V = eigenvect(:,unstable);
indx_unstable = find(sum(abs(V),2)>1e-5);
stable = find(sum(abs(V),2)<1e-5);
nunit = length(eigenv) - nstable;
Pstar = options_.Harvey_scale_factor*eye(np);
if kalman_algo ~= 2
kalman_algo = 1;
end
R_tmp = R(stable, :);
T_tmp = T(stable,stable);
if DynareOptions.lyapunov_fp == 1
Pstar_tmp = lyapunov_symm(T_tmp,Q,DynareOptions.lyapunov_fixed_point_tol,DynareOptions.lyapunov_complex_threshold, 3, R_tmp);
elseif DynareOptions.lyapunov_db == 1
Pstar_tmp = disclyap_fast(T_tmp,R_tmp*Q*R_tmp',DynareOptions.lyapunov_doubling_tol);
elseif DynareOptions.lyapunov_srs == 1
Pstar_tmp = lyapunov_symm(T_tmp,Q,DynareOptions.lyapunov_fixed_point_tol,DynareOptions.lyapunov_complex_threshold, 4, R_tmp);
else
Pstar_tmp = lyapunov_symm(T_tmp,R_tmp*Q*R_tmp',DynareOptions.qz_criterium,DynareOptions.lyapunov_complex_threshold);
end
Pstar(stable, stable) = Pstar_tmp;
Pinf = [];
otherwise
error('dsge_likelihood:: Unknown initialization approach for the Kalman filter!')
end

View File

@ -1,4 +1,4 @@
function oo = evaluate_smoother(parameters,var_list)
function evaluate_smoother(parameters,var_list)
% Evaluate the smoother at parameters.
%
% INPUTS
@ -46,11 +46,8 @@ global options_ M_ bayestopt_ oo_ estim_params_ % estim_params_ may be emty
persistent dataset_
if isempty(dataset_) || isempty(bayestopt_)
options = options_;
options.smoother = 1; % this is necessary because of a check in dynare_estimation_init()
[dataset_,xparam1, M_, options_, oo_, estim_params_,bayestopt_] = dynare_estimation_init(var_list, M_.fname, [], M_, options, oo_, estim_params_, bayestopt_);
[dataset_,xparam1, M_, options_, oo_, estim_params_,bayestopt_] = dynare_estimation_init(var_list, M_.fname, [], M_, options_, oo_, estim_params_, bayestopt_);
end
if nargin==0
@ -86,43 +83,35 @@ if ischar(parameters)
end
end
pshape_original = bayestopt_.pshape;
bayestopt_.pshape = Inf(size(bayestopt_.pshape));
clear('priordens')
[atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp] = ...
DsgeSmoother(parameters,dataset_.info.ntobs,dataset_.data,dataset_.missing.aindex,dataset_.missing.state);
oo.Smoother.SteadyState = ys;
oo.Smoother.TrendCoeffs = trend_coeff;
oo_.Smoother.SteadyState = ys;
oo_.Smoother.TrendCoeffs = trend_coeff;
if options_.filter_covariance
oo.Smoother.Variance = P;
oo_.Smoother.Variance = P;
end
i_endo = bayestopt_.smoother_saved_var_list;
if options_.nk ~= 0
oo.FilteredVariablesKStepAhead = ...
oo_.FilteredVariablesKStepAhead = ...
aK(options_.filter_step_ahead,i_endo,:);
if ~isempty(PK)
oo.FilteredVariablesKStepAheadVariances = ...
oo_.FilteredVariablesKStepAheadVariances = ...
PK(options_.filter_step_ahead,i_endo,i_endo,:);
end
if ~isempty(decomp)
oo.FilteredVariablesShockDecomposition = ...
oo_.FilteredVariablesShockDecomposition = ...
decomp(options_.filter_step_ahead,i_endo,:,:);
end
end
dr = oo_.dr;
order_var = oo_.dr.order_var;
for i=bayestopt_.smoother_saved_var_list'
i1 = order_var(bayestopt_.smoother_var_list(i));
eval(['oo.SmoothedVariables.' deblank(M_.endo_names(i1,:)) ' = atT(i,:)'';']);
eval(['oo.FilteredVariables.' deblank(M_.endo_names(i1,:)) ' = squeeze(aK(1,i,:));']);
eval(['oo.UpdatedVariables.' deblank(M_.endo_names(i1,:)) ' = updated_variables(i,:)'';']);
i1 = oo_.dr.order_var(bayestopt_.smoother_var_list(i));
eval(['oo_.SmoothedVariables.' deblank(M_.endo_names(i1,:)) ' = atT(i,:)'';']);
if options_.nk>0
eval(['oo_.FilteredVariables.' deblank(M_.endo_names(i1,:)) ' = squeeze(aK(1,i,:));']);
end
eval(['oo_.UpdatedVariables.' deblank(M_.endo_names(i1,:)) ' = updated_variables(i,:)'';']);
end
for i=1:M_.exo_nbr
eval(['oo.SmoothedShocks.' deblank(M_.exo_names(i,:)) ' = innov(i,:)'';']);
end
oo.dr = oo_.dr;
bayestopt_.pshape = pshape_original;
eval(['oo_.SmoothedShocks.' deblank(M_.exo_names(i,:)) ' = innov(i,:)'';']);
end

View File

@ -2412,10 +2412,9 @@ CalibSmootherStatement::writeOutput(ostream &output, const string &basename) con
{
options_list.writeOutput(output);
symbol_list.writeOutput("var_list_", output);
output << "options_.mode_compute = 0;" << endl
<< "options_.smoother = 1;" << endl
<< "options_.order = 1;" << endl
<< "dynare_estimation(var_list_);" << endl;
output << "options_.smoother = 1;" << endl;
output << "options_.order = 1;" << endl;
output << "evaluate_smoother('calibration',var_list_);" << endl;
}
ExtendedPathStatement::ExtendedPathStatement(const OptionsList &options_list_arg)