Added PEA-3 in rbcii-simple.

master
Stéphane Adjemian (Charybdis) 2016-04-29 22:51:00 +02:00
parent 5fbf944445
commit 6ed585a3a3
1 changed files with 205 additions and 0 deletions

View File

@ -0,0 +1,205 @@
locals;
dynare_config();
% Compile the original model.
dynare('rbc2.mod', 'console','noclearall');
truemodel.oo = oo_;
truemodel.M = M_;
truemodel.options = options_;
warning off all,
delete('rbc2.m', 'rbc2.log', 'rbc2_results.mat');
warning on all
clear('oo_', 'M_', 'options_')
% Define indices used for the approximation.
peaopt.idExpectedTerms = strmatch('ExpectedTerm', truemodel.M.endo_names);
peaopt.idExpectedEquations = 6;
peaopt.idEfficiency = strmatch('Efficiency', truemodel.M.endo_names);
peaopt.idCapital = strmatch('Capital', truemodel.M.endo_names);
peaopt.idOutput = strmatch('Output', truemodel.M.endo_names);
peaopt.idConsumption = strmatch('Consumption', truemodel.M.endo_names);
peaopt.idInvestment = strmatch('Investment', truemodel.M.endo_names);
peaopt.offsetparams = length(truemodel.M.params);
initialize_with_ep = false;
% Set the size of the simulated sample.
peaopt.samplesize = 100000;
% Remove the first simulations.
burnin = 1000;
peaopt.init = burnin+1;
% Simuation of the model with perturnation approach (2nd order).
truemodel.options.order = 2;
truemodel.options.periods = peaopt.samplesize+10;
check_model(truemodel.M);
truemodel.oo.dr = set_state_space(truemodel.oo.dr, truemodel.M, truemodel.options);
[truemodel.oo.dr, info, truemodel.M, truemodel.options, truemodel.oo] = resol(0, truemodel.M, truemodel.options, truemodel.oo);
set_dynare_seed(31415);
[simulatedvariables, truemodel.oo] = simult(truemodel.oo.dr.ys, truemodel.oo.dr, truemodel.M, truemodel.options, truemodel.oo);
if initialize_with_ep
ts = extended_path([], 5000, [], truemodel.options, truemodel.M, truemodel.oo);
ExpectedTerm = ts.ExpectedTerm.data(peaopt.init:end);
lCapital = log(ts.Capital.data(peaopt.init-1:end-1));
lEfficiency = log(ts.Efficiency.data(peaopt.init:end));
truemodel.oo.endo_simul = transpose(ts.data);
else
ExpectedTerm = simulatedvariables(peaopt.idExpectedTerms, peaopt.init:peaopt.samplesize)';
lCapital = log(simulatedvariables(peaopt.idCapital, peaopt.init-1:peaopt.samplesize-1))';
lEfficiency = log(simulatedvariables(peaopt.idEfficiency,peaopt.init:peaopt.samplesize))';
truemodel.oo.endo_simul = simulatedvariables;
end
% Set option for optimizer (NLS in PEA).
nlsopt = optimoptions('lsqcurvefit','Display','off');
% NLS: Estimate the PEA-2 parameters (used as a first guess)
X2 = [ones(rows(lEfficiency), 1), lCapital, lEfficiency, lCapital.*lCapital, lCapital.*lEfficiency, lEfficiency.*lEfficiency];
X3 = [X2, lCapital.*lCapital.*lCapital, lCapital.*lCapital.*lEfficiency, lCapital.*lEfficiency.*lEfficiency, lEfficiency.*lEfficiency.*lEfficiency];
old = zeros(10,1);
old(1) = log(truemodel.oo.steady_state(peaopt.idExpectedTerms));
[rparams_init, resnorm, residual] = lsqcurvefit(@pea_predictor, old, X3, ExpectedTerm, [], [], nlsopt); %%
rparams_init = [-2.288066265940801
6.133612871733019
-13.688559635587980
-5.446915204762234
24.874487636569398
-96.831528275604711
1.431335015127595
-11.058117331327608
69.010403456632815
-33.947595125598575];
rparams_init = rparams_init(:);
skipline(2)
disp('Initial guess for the PEA-2 parameters ')
skipline()
disp(sprintf(' Constant %2.6f',rparams_init(1)))
disp(sprintf(' Capital %2.6f',rparams_init(2)))
disp(sprintf(' Efficiency %2.6f',rparams_init(3)))
disp(sprintf(' Capital*Capital %2.6f',rparams_init(4)))
disp(sprintf(' Capital*Efficiency %2.6f',rparams_init(5)))
disp(sprintf(' Efficiency*Efficiency %2.6f',rparams_init(6)))
disp(sprintf(' Capital*Capital*Capital %2.6f',rparams_init(7)))
disp(sprintf(' Capital*Capital*Efficiency %2.6f',rparams_init(8)))
disp(sprintf(' Capital*Efficiency*Efficiency %2.6f',rparams_init(9)))
disp(sprintf(' Efficiency*Efficiency*Efficiency %2.6f',rparams_init(10)))
skipline()
s2old = (residual'*residual);
% Get values of the parameters.
BETA = truemodel.M.params(strmatch('beta', truemodel.M.param_names));
ALPHA = truemodel.M.params(strmatch('alpha', truemodel.M.param_names));
TAU = truemodel.M.params(strmatch('tau', truemodel.M.param_names));
DELTA = truemodel.M.params(strmatch('delta', truemodel.M.param_names));
RHO = truemodel.M.params(strmatch('rho', truemodel.M.param_names));
SIGMA = truemodel.M.params(strmatch('sigma', truemodel.M.param_names));
ASTAR = truemodel.M.params(strmatch('effstar', truemodel.M.param_names));
ISS = truemodel.oo.steady_state(strmatch('Investment', truemodel.M.endo_names));
lambda = 1;
tolerance = 1e-4; % Percentage change in the reduced form parameters.
maxiter = 50000;
% Generate productivity time series (once for all)
set_dynare_seed(31415);
efficiency = zeros(peaopt.samplesize, 1);
innovations = [0; SIGMA*randn(peaopt.samplesize-1, 1)];
for t=2:peaopt.samplesize
efficiency(t) = RHO*efficiency(t-1)+innovations(t);
end
Efficiency = ASTAR*exp(efficiency);
% Set options for NLS
nlsopt = optimoptions('lsqcurvefit','Display','off','Algorithm','levenberg-marquardt');
crits = Inf(2,1);
RPARAMS = NaN(10, 2);
RPARAMS(:,1) = rparams_init;
% Main loop ()
iteration = 1;
while 1
Capital = zeros(peaopt.samplesize, 1);
Capital(1) = truemodel.oo.steady_state(peaopt.idCapital);
ExpectedTerm = zeros(peaopt.samplesize, 1);
LagrangeMultiplier = zeros(peaopt.samplesize, 1);
X3 = zeros(peaopt.samplesize, 10);
for t=2:peaopt.samplesize
X3(t, 1:6) = [1, log(Capital(t-1)), log(Efficiency(t)), log(Capital(t-1))*log(Capital(t-1)), log(Capital(t-1))*log(Efficiency(t)), log(Efficiency(t))*log(Efficiency(t))];
X3(t, 7) = log(Capital(t-1))*log(Capital(t-1))*log(Capital(t-1));
X3(t, 8) = log(Capital(t-1))*log(Capital(t-1))*log(Efficiency(t));
X3(t, 9) = log(Capital(t-1))*log(Efficiency(t))*log(Efficiency(t));
X3(t, 10) = log(Efficiency(t))*log(Efficiency(t))*log(Efficiency(t));
ExpectedTerm(t) = exp(X3(t,:)*RPARAMS(:, 1));
Output = Efficiency(t)*Capital(t-1)^ALPHA;
Consumption = ExpectedTerm(t)^(-1/TAU);
Investment = Output - Consumption;
if Investment>.975*ISS
Capital(t) = Output+(1-DELTA)*Capital(t-1)-Consumption;
LagrangeMultiplier(t) = 0;
else
Capital(t) = (1-DELTA)*Capital(t-1)+.975*ISS;
LagrangeMultiplier(t) = Consumption^(-TAU)-ExpectedTerm(t);
end
end
ExpectedTerm_ = beta*ExpectedTerm(2:peaopt.samplesize).*(ALPHA*Efficiency(2:peaopt.samplesize).*Capital(2:peaopt.samplesize).^(ALPHA-1)+1-DELTA)-LagrangeMultiplier(2:peaopt.samplesize)*(1-DELTA)*BETA;
iparams = RPARAMS(:, 1); %X2(1:peaopt.samplesize-1,:)\ExpectedTerm_;
[rparams0, resnorm0,residual0,exitflag0] = lsqcurvefit(@pea_predictor, iparams, X3(peaopt.init:peaopt.samplesize-1,:), ExpectedTerm_(peaopt.init:end), [], [], nlsopt);%
s2 = residual0'*residual0;
RPARAMS(:,2) = RPARAMS(:,1);
RPARAMS(:,1) = lambda*rparams0 + (1-lambda)*RPARAMS(:,1);
crit = 100*max(abs(RPARAMS(:,1)-RPARAMS(:,2))./abs(RPARAMS(:,2)));
ds2 = 100*(s2-s2old)/s2old;
s2old = s2;
iteration = iteration + 1;
disp(sprintf('Iteration %s\t crit = %3.16f\t s2 = %1.8f\t D(s2) = %1.6f',int2str(iteration),crit,s2,ds2))
crits(1:end-1) = crits(2:end);
crits(end) = crit;
if all(crits<tolerance)
exitflag0
break
end
if iteration>maxiter
break
end
end
rparams = RPARAMS(:,1);
Capital = zeros(peaopt.samplesize+1, 1); Capital(1) = truemodel.oo.steady_state(peaopt.idCapital);
ExpectedTerm = zeros(peaopt.samplesize, 1);
Consumption = zeros(peaopt.samplesize, 1);
Output = zeros(peaopt.samplesize, 1);
Investment = zeros(peaopt.samplesize,1);
LagrangeMultiplier = zeros(peaopt.samplesize, 1);
X3 = zeros(peaopt.samplesize, 10);
for t=2:peaopt.samplesize
X3(t, 1:6) = [1, log(Capital(t-1)), log(Efficiency(t)), log(Capital(t-1))*log(Capital(t-1)), log(Capital(t-1))*log(Efficiency(t)), log(Efficiency(t))*log(Efficiency(t))];
X3(t, 7) = log(Capital(t-1))*log(Capital(t-1))*log(Capital(t-1));
X3(t, 8) = log(Capital(t-1))*log(Capital(t-1))*log(Efficiency(t));
X3(t, 9) = log(Capital(t-1))*log(Efficiency(t))*log(Efficiency(t));
X3(t, 10) = log(Efficiency(t))*log(Efficiency(t))*log(Efficiency(t));
ExpectedTerm(t) = exp(X3(t,:)*rparams);
Output(t) = Efficiency(t)*Capital(t-1)^ALPHA;
Consumption(t) = ExpectedTerm(t)^(-1/TAU);
Investissement = Output(t) - Consumption(t);
if Investissement>.975*ISS
Investment(t) = Investissement;
Capital(t) = (1-DELTA)*Capital(t-1) + Investment(t);
LagrangeMultiplier(t) = 0;
else
Investment(t) = .975*ISS;
Capital(t) = (1-DELTA)*Capital(t-1)+.975*ISS;
LagrangeMultiplier(t) = Consumption(t)^(-TAU)-ExpectedTerm(t);
end
end
save('rbc.pea3.mat','Output','Investment','Consumption','Capital','Efficiency','efficiency','LagrangeMultiplier','ExpectedTerm','innovations','rparams','rparams_init','peaopt','truemodel');
clearvars