2010-01-11 00:27:43 +01:00
|
|
|
function [r,flag] = smm_objective(xparams,sample_moments,weighting_matrix,options,parallel)
|
2010-01-06 15:38:32 +01:00
|
|
|
% Evaluates the objective of the Simulated Moments Method.
|
|
|
|
%
|
|
|
|
% INPUTS:
|
2017-05-16 15:10:20 +02:00
|
|
|
% xparams [double] p*1 vector of estimated parameters.
|
2010-01-06 15:38:32 +01:00
|
|
|
% sample_moments [double] n*1 vector of sample moments (n>=p).
|
|
|
|
% weighting_matrix [double] n*n symetric, positive definite matrix.
|
|
|
|
% options [ ] Structure defining options for SMM.
|
2010-01-08 18:14:50 +01:00
|
|
|
% parallel [ ] Structure defining the parallel mode settings (optional).
|
2010-01-06 15:38:32 +01:00
|
|
|
%
|
2017-05-16 15:10:20 +02:00
|
|
|
% OUTPUTS:
|
2010-01-06 15:38:32 +01:00
|
|
|
% r [double] scalar, the value of the objective function.
|
|
|
|
% junk [ ] empty matrix.
|
|
|
|
%
|
|
|
|
% SPECIAL REQUIREMENTS
|
|
|
|
% The user has to provide a file where the moment conditions are defined.
|
|
|
|
|
2019-09-11 11:14:20 +02:00
|
|
|
% Copyright (C) 2010-2019 Dynare Team
|
2010-01-06 15:38:32 +01:00
|
|
|
%
|
|
|
|
% 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
|
2017-05-16 15:10:20 +02:00
|
|
|
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
|
2010-01-06 15:38:32 +01:00
|
|
|
|
2011-09-17 11:47:36 +02:00
|
|
|
global M_ options_ oo_
|
2010-01-06 15:38:32 +01:00
|
|
|
persistent mainStream mainState
|
|
|
|
persistent priorObjectiveValue
|
2010-01-11 00:27:43 +01:00
|
|
|
|
|
|
|
flag = 1;
|
|
|
|
|
2010-01-08 18:14:50 +01:00
|
|
|
if nargin<5
|
|
|
|
if isempty(mainStream)
|
2012-04-20 14:41:54 +02:00
|
|
|
if matlab_ver_less_than('7.12')
|
|
|
|
mainStream = RandStream.getDefaultStream;
|
|
|
|
else
|
|
|
|
mainStream = RandStream.getGlobalStream;
|
|
|
|
end
|
2010-01-08 18:14:50 +01:00
|
|
|
mainState = mainStream.State;
|
|
|
|
else
|
|
|
|
mainStream.State = mainState;
|
|
|
|
end
|
2010-01-06 15:38:32 +01:00
|
|
|
end
|
|
|
|
|
|
|
|
if isempty(priorObjectiveValue)
|
|
|
|
priorObjectiveValue = Inf;
|
|
|
|
end
|
|
|
|
|
2010-01-11 14:33:23 +01:00
|
|
|
penalty = 0;
|
|
|
|
for i=1:options.estimated_parameters.nb
|
|
|
|
if ~isnan(options.estimated_parameters.upper_bound(i)) && xparams(i)>options.estimated_parameters.upper_bound(i)
|
|
|
|
penalty = penalty + (xparams(i)-options.estimated_parameters.upper_bound(i))^2;
|
|
|
|
end
|
|
|
|
if ~isnan(options.estimated_parameters.lower_bound(i)) && xparams(i)<options.estimated_parameters.lower_bound(i)
|
|
|
|
penalty = penalty + (xparams(i)-options.estimated_parameters.lower_bound(i))^2;
|
|
|
|
end
|
|
|
|
end
|
|
|
|
|
|
|
|
if penalty>0
|
|
|
|
flag = 0;
|
2017-05-16 15:10:20 +02:00
|
|
|
r = priorObjectiveValue + penalty;
|
2017-05-16 12:42:01 +02:00
|
|
|
return
|
2010-01-11 14:33:23 +01:00
|
|
|
end
|
|
|
|
|
2010-01-10 14:49:18 +01:00
|
|
|
save('estimated_parameters.mat','xparams');
|
2010-01-06 15:38:32 +01:00
|
|
|
|
|
|
|
% Check for local determinacy of the deterministic steady state.
|
|
|
|
noprint = options_.noprint; options_.noprint = 1;
|
2019-09-11 11:14:20 +02:00
|
|
|
[~,local_determinacy_and_stability,info] = check(M_,options_,oo_); options_.noprint = noprint;
|
2010-01-06 15:38:32 +01:00
|
|
|
if ~local_determinacy_and_stability
|
|
|
|
r = priorObjectiveValue * (1+info(2));
|
2010-01-11 00:27:43 +01:00
|
|
|
flag = 0;
|
2010-01-06 15:38:32 +01:00
|
|
|
return
|
|
|
|
end
|
|
|
|
|
|
|
|
simulated_moments = zeros(size(sample_moments));
|
|
|
|
|
2010-01-11 00:27:43 +01:00
|
|
|
% Just to be sure that things don't mess up with persistent variables...
|
|
|
|
clear perfect_foresight_simulation;
|
|
|
|
|
2010-01-08 18:14:50 +01:00
|
|
|
if nargin<5
|
|
|
|
for s = 1:options.number_of_simulated_sample
|
|
|
|
time_series = extended_path([],options.simulated_sample_size,1);
|
|
|
|
data = time_series(options.observed_variables_idx,options.burn_in_periods+1:options.simulated_sample_size);
|
|
|
|
eval(['tmp = ' options.moments_file_name '(data);'])
|
|
|
|
simulated_moments = simulated_moments + tmp;
|
|
|
|
simulated_moments = simulated_moments / options.number_of_simulated_sample;
|
|
|
|
end
|
|
|
|
else% parallel mode.
|
2010-01-11 14:33:23 +01:00
|
|
|
if ~isunix
|
|
|
|
error('The parallel version of SMM estimation is not implemented for non unix platforms!')
|
|
|
|
end
|
2010-01-15 16:45:07 +01:00
|
|
|
job_number = 1;% Remark. First job is for the master.
|
2019-09-11 11:14:20 +02:00
|
|
|
[~,hostname] = unix('hostname --fqdn');
|
2010-01-09 01:40:10 +01:00
|
|
|
hostname = deblank(hostname);
|
2010-01-08 18:14:50 +01:00
|
|
|
for i=1:length(parallel)
|
2010-01-15 16:45:07 +01:00
|
|
|
machine = deblank(parallel(i).machine);
|
|
|
|
if ~strcmpi(hostname,machine)
|
|
|
|
% For the slaves on a remote computer.
|
|
|
|
unix(['scp estimated_parameters.mat ' , parallel(i).login , '@' , machine , ':' parallel(i).folder ' > /dev/null']);
|
|
|
|
else
|
|
|
|
if ~strcmpi(pwd,parallel(i).folder)
|
|
|
|
% For the slaves on this computer but not in the same directory as the master.
|
|
|
|
unix(['cp estimated_parameters.mat ' , parallel(i).folder]);
|
|
|
|
end
|
|
|
|
end
|
2010-01-08 18:33:31 +01:00
|
|
|
for j=1:parallel(i).number_of_jobs
|
2017-05-16 15:10:20 +02:00
|
|
|
if (strcmpi(hostname,machine) && j>1) || ~strcmpi(hostname,machine)
|
2010-01-15 16:45:07 +01:00
|
|
|
job_number = job_number + 1;
|
|
|
|
unix(['ssh -A ' parallel(i).login '@' machine ' ./call_matlab_session.sh job' int2str(job_number) '.m &']);
|
2010-01-08 18:14:50 +01:00
|
|
|
end
|
|
|
|
end
|
|
|
|
end
|
2010-01-15 16:45:07 +01:00
|
|
|
% Finally the Master do its job
|
|
|
|
tStartMasterJob = clock;
|
|
|
|
eval('job1;')
|
|
|
|
tElapsedMasterJob = etime(clock, tStartMasterJob);
|
|
|
|
TimeLimit = tElapsedMasterJob*1.2;
|
2017-05-16 15:10:20 +02:00
|
|
|
% Master waits for the slaves' output...
|
2010-01-08 18:14:50 +01:00
|
|
|
tStart = clock;
|
2010-01-15 16:45:07 +01:00
|
|
|
tElapsed = 0;
|
|
|
|
while tElapsed<TimeLimit
|
2010-01-10 14:49:18 +01:00
|
|
|
if ( length(dir('./intermediary_results_from_master_and_slaves/simulated_moments_slave_*.dat'))==job_number )
|
2010-01-15 16:45:07 +01:00
|
|
|
break
|
2010-01-08 18:14:50 +01:00
|
|
|
end
|
2010-01-15 16:45:07 +01:00
|
|
|
tElapsed = etime(clock, tStart);
|
|
|
|
end
|
|
|
|
try
|
|
|
|
tmp = zeros(length(sample_moments),1);
|
2010-01-08 18:14:50 +01:00
|
|
|
for i=1:job_number
|
2010-01-10 14:49:18 +01:00
|
|
|
simulated_moments = load(['./intermediary_results_from_master_and_slaves/simulated_moments_slave_' int2str(i) '.dat'],'-ascii');
|
2010-01-08 18:14:50 +01:00
|
|
|
tmp = tmp + simulated_moments;
|
|
|
|
end
|
2017-05-16 15:10:20 +02:00
|
|
|
simulated_moments = tmp / job_number;
|
2010-01-15 16:45:07 +01:00
|
|
|
catch
|
2010-01-08 18:14:50 +01:00
|
|
|
r = priorObjectiveValue*1.1;
|
2010-01-11 00:27:43 +01:00
|
|
|
flag = 0;
|
2010-01-08 18:14:50 +01:00
|
|
|
return
|
|
|
|
end
|
|
|
|
end
|
2010-01-06 15:38:32 +01:00
|
|
|
|
|
|
|
r = transpose(simulated_moments-sample_moments)*weighting_matrix*(simulated_moments-sample_moments);
|
2010-01-11 17:36:36 +01:00
|
|
|
priorObjectiveValue = r;
|
2010-01-06 15:38:32 +01:00
|
|
|
|
2010-02-17 16:20:37 +01:00
|
|
|
if (options.optimization_routine>0) && exist('optimization_path.mat')
|
2010-01-11 17:36:36 +01:00
|
|
|
load('optimization_path.mat');
|
|
|
|
new_state = [ r; xparams];
|
|
|
|
estimated_parameters_optimization_path = [ estimated_parameters_optimization_path , new_state ];
|
|
|
|
save('optimization_path.mat','estimated_parameters_optimization_path');
|
|
|
|
end
|