-
Notifications
You must be signed in to change notification settings - Fork 5
/
recsFirstGuess.m
146 lines (130 loc) · 5.24 KB
/
recsFirstGuess.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
function [interp,x,z,exitflag,output] = recsFirstGuess(interp,model,s,sss,xss,options)
% RECSFIRSTGUESS finds a first guess using the perfect foresight solution
%
% RECSFIRSTGUESS tries to find a first-guess for the rational expectations model
% by solving the corresponding perfect foresight solution on all the grid points
% of the state variables. By default, it considers that the model goes back to
% its deterministic steady state in 50 periods.
%
% INTERP = RECSFIRSTGUESS(INTERP,MODEL,S,SSS,XSS) uses the interpolation structure
% defined in the structure INTERP to fit the perfect foresight solution of the
% model defined in the structure MODEL. The grid of state variables is provided in
% the n-by-d matrix S. A first guess for the steady state of the model is provided
% in SSS and XSS for state and response variables. RECSFIRSTGUESS returns an
% interpolation structure, INTERP, containing the first guess.
% INTERP is a structure, which has to include the following fields:
% fspace : a definition structure for the interpolation family (created
% by the function fundef)
% MODEL is an object created by recsmodel.
%
% INTERP = RECSFIRSTGUESS(INTERP,MODEL,S,SSS,XSS,OPTIONS) solves the problem
% with the parameters defined by the structure OPTIONS. The fields of the
% structure are
% eqsolver : 'fsolve', 'lmmcp' (default), 'ncpsolve' or 'path'
% eqsolveroptions : options structure to be passed to eqsolver (default:
% empty structure)
% fgmethod : 'auto' (default), 'perturbation', 'perfect-foresight' or
% 'steady-state'
% T : integer defining the time horizon at which the model is
% supposed to converge to its steady state (default: 50)
% UseParallel : 'always' (default) to use parallel calculation (require
% Parallel Computing Toolbox)' or never'
%
% [INTERP,X] = RECSFIRSTGUESS(INTERP,MODEL,S,SSS,XSS,...) returns X, n-by-m matrix,
% containing the value of the response variables in the first period.
%
% [INTERP,X,Z] = RECSFIRSTGUESS(INTERP,MODEL,S,SSS,XSS,...) returns Z, n-by-p matrix,
% containing the value of the expectations variables in the first period.
%
% [INTERP,X,Z,EXITFLAG] = RECSFIRSTGUESS(INTERP,MODEL,S,SSS,XSS,...)
%
% [INTERP,X,Z,EXITFLAG,OUTPUT] = RECSFIRSTGUESS(INTERP,MODEL,S,SSS,XSS,...)
%
% See also RECSSOLVEDETERMINISTICPB, RECSSOLVELOCAL, RECSSOLVEREE, RECSSS.
% Copyright (C) 2011-2016 Christophe Gouel
% Licensed under the Expat license, see LICENSE.txt
%% Initialization
% Get s from interp structure
if nargin<=2 || isempty(s), s = interp.s; end
if nargin <=3 || isempty(sss)
if isfield(model,'sss') || (isobject(model) && ~isempty(model.sss))
sss = model.sss;
else
sss = [];
end
end
if nargin <=4 || isempty(xss)
if isfield(model,'xss') || (isobject(model) && ~isempty(model.xss))
xss = model.xss;
else
xss = [];
end
end
defaultopt = struct('fgmethod' , 'auto',...
'T' , 50 ,...
'UseParallel' , 'always');
if nargin <=5
options = defaultopt;
else
options = catstruct(defaultopt,options);
end
fgmethod = options.fgmethod;
T = options.T;
switch lower(options.UseParallel)
case 'never'
UseParallel = 0;
case 'always'
UseParallel = size(s,1);
end
n = size(s,1);
p = model.dim{3};
params = model.params;
%% Solve for the deterministic steady state
[sss,xss,zss] = recsSS(model,sss,xss,catstruct(options,struct('display',0)));
%% Find first-guess
[LB,UB] = model.functions.b(s,params);
if strcmpi(fgmethod,'auto') && all(isinf([LB(:); UB(:)]))
fgmethod = 'perturbation';
elseif strcmpi(fgmethod,'auto')
fgmethod = 'perfect-foresight';
end
switch fgmethod
case 'perturbation'
%% First-order perturbation
model = recsSolveLocal(model);
x = min(max(model.LinearSolution.X(s),LB),UB);
z = model.LinearSolution.Z(s);
exitflag = 1;
output = [];
case 'perfect-foresight'
%% Solve the perfect foresight problem on each point of the grid
x = zeros(n,size(xss,2));
z = zeros(n,size(zss,2));
exitflag = zeros(n,1);
N = zeros(n,1);
parfor (i=1:n, UseParallel)
[X,~,Z,~,exitflag(i),N(i)] = recsSolveDeterministicPb(model,s(i,:),...
T,xss,zss,sss,options);
x(i,:) = X(1,:);
z(i,:) = Z(1,:);
end
output = struct('exitflag',exitflag,...
'N' ,N);
exitflag = ~any(exitflag~=1);
case 'steady-state'
%% Response variables equal to their steady-state values
x = min(max(repmat(xss,n,1),LB),UB);
z = NaN(n,p);
exitflag = 1;
output = [];
end % switch fgmethod
%% Prepare output
interp.cx = funfitxy(interp.fspace,interp.Phi,x);
interp.cz = funfitxy(interp.fspace,interp.Phi,z);
interp.x = x;
interp.z = z;
% Check if it is possible to approximate the expectations function
if strcmp(model.infos.model_type,'fgh1')
hv = model.functions.h(zeros(n,0),[],[],s,x,params);
interp.ch = funfitxy(interp.fspace,interp.Phi,hv);
end