Commit e07d97a1 authored by Stéphane Adjemian's avatar Stéphane Adjemian
Browse files

Merge pull request #1186 from rattoma/gsa

gsa and identification fixes
parents bebcfffe 1b769d1d
......@@ -238,6 +238,10 @@ if iload <=0,
disp('Testing ML Starting value')
else
switch parameters
case 'calibration'
parameters_TeX = 'Calibration';
disp('Testing calibration')
params(1,:) = get_all_parameters(estim_params_,M_);;
case 'posterior_mode'
parameters_TeX = 'Posterior mode';
disp('Testing posterior mode')
......
......@@ -48,6 +48,8 @@ if isfield(options_gsa,'graph_format'),
end
if isfield(options_gsa,'mode_file'),
options_.mode_file=options_gsa.mode_file;
elseif isfield(options_gsa,'neighborhood_width') && options_gsa_.neighborhood_width>0,
options_.mode_file='';
end
options_.order = 1;
......
......@@ -472,17 +472,24 @@ if iload==0,
mkdir(xdir)
end
nrun=length(y0);
nest=min(250,nrun);
nest=max(50,nrun/2);
nest=min(250,nest);
nfit=min(1000,nrun);
% dotheplots = (nfit<=nest);
% gsa_ = gsa_sdp(y0(1:nest), x0(1:nest,:), 2, [],[-1 -1 -1 -1 -1 0],[],0,[fname,'_est'], pnames);
[ys,is] = sort(y0);
istep = ceil(nrun/nest);
if istep>1,
iest = is(floor(istep/2):istep:end);
nest = length(iest);
irest = is(setdiff([1:nrun],[floor(istep/2):istep:nrun]));
istep = ceil(length(irest)/(nfit-nest));
ifit = union(iest, irest(1:istep:end));
else
warning('the number of samples is too small for ANOVA estimation')
si=nan(np,1);
return
end
if ~ismember(irest(end),ifit),
ifit = union(ifit, irest(end));
end
......
......@@ -241,7 +241,11 @@ if fload==0,
% end
% end
%load([fname_,'_mode'])
if neighborhood_width>0 && isempty(options_.mode_file),
xparam1 = get_all_parameters(estim_params_,M_);
else
eval(['load ' options_.mode_file '.mat;']);
end
if neighborhood_width>0,
for j=1:nshock,
lpmat0(:,j) = randperm(Nsam)'./(Nsam+1); %latin hypercube
......@@ -279,7 +283,7 @@ if fload==0,
iwrong=zeros(1,Nsam);
inorestriction=zeros(1,Nsam);
irestriction=zeros(1,Nsam);
infox=zeros(1,Nsam);
infox=zeros(Nsam,1);
for j=1:Nsam,
M_ = set_all_parameters([lpmat0(j,:) lpmat(j,:)]',estim_params_,M_);
%try stoch_simul([]);
......@@ -346,7 +350,6 @@ if fload==0,
info=endogenous_prior_restrictions(Tt,Rr,M_,options_,oo_);
infox(j,1)=info(1);
if info(1),
iwrong(j)=j;
inorestriction(j)=j;
else
iunstable(j)=0;
......@@ -389,7 +392,7 @@ if fload==0,
iunstable=iunstable(find(iunstable)); % violation of BK & restrictions & solution could not be found (whatever goes wrong)
iindeterm=iindeterm(find(iindeterm)); % indeterminacy
iwrong=iwrong(find(iwrong)); % dynare could not find solution
ixun=iunstable(find(~ismember(iunstable,[iindeterm,iwrong]))); % explosive roots
ixun=iunstable(find(~ismember(iunstable,[iindeterm,iwrong,inorestriction]))); % explosive roots
% % map stable samples
% istable=[1:Nsam];
......@@ -526,14 +529,21 @@ delete([OutputDirectoryName,filesep,fname_,'_',awrongname,'.*']);
if length(iunstable)>0 || length(iwrong)>0,
fprintf(['%4.1f%% of the prior support gives unique saddle-path solution.\n'],length(istable)/Nsam*100)
fprintf(['%4.1f%% of the prior support gives explosive dynamics.\n'],(length(iunstable)-length(iwrong)-length(iindeterm) )/Nsam*100)
fprintf(['%4.1f%% of the prior support gives explosive dynamics.\n'],(length(ixun) )/Nsam*100)
if ~isempty(iindeterm),
fprintf(['%4.1f%% of the prior support gives indeterminacy.'],length(iindeterm)/Nsam*100)
fprintf(['%4.1f%% of the prior support gives indeterminacy.\n'],length(iindeterm)/Nsam*100)
end
if ~isempty(iwrong),
skipline()
disp(['For ',num2str(length(iwrong)/Nsam*100,'%4.1f'),'% of the prior support dynare could not find a solution.'])
inorestriction = istable(find(~ismember(istable,irestriction))); % violation of prior restrictions
if ~isempty(iwrong) || ~isempty(inorestriction),
skipline()
if any(infox==49),
fprintf(['%4.1f%% of the prior support violates prior restrictions.\n'],(length(inorestriction) )/Nsam*100)
end
if ~isempty(iwrong),
skipline()
disp(['For ',num2str(length(iwrong)/Nsam*100,'%4.1f'),'% of the prior support dynare could not find a solution.'])
skipline()
end
if any(infox==1),
disp([' For ',num2str(length(find(infox==1))/Nsam*100,'%4.1f'),'% The model doesn''t determine the current variables uniquely.'])
end
......@@ -564,14 +574,12 @@ if length(iunstable)>0 || length(iwrong)>0,
if any(infox==30),
disp([' For ',num2str(length(find(infox==30))/Nsam*100,'%4.1f'),'% Ergodic variance can''t be computed.'])
end
if any(infox==49),
disp([' For ',num2str(length(find(infox==49))/Nsam*100,'%4.1f'),'% The model violates one (many) endogenous prior restriction(s).'])
end
end
skipline()
if length(iunstable)<Nsam || length(istable)>1
itot = [1:Nsam];
isolve = itot(find(~ismember(itot,iwrong))); % dynare could find a solution
% Blanchard Kahn
if neighborhood_width,
options_mcf.xparam1 = xparam1(nshock+1:end);
......@@ -585,7 +593,7 @@ if length(iunstable)>0 || length(iwrong)>0,
mcf_analysis(lpmat, istable, itmp, options_mcf, options_);
if ~isempty(iindeterm),
itmp = itot(find(~ismember(itot,iindeterm)));
itmp = isolve(find(~ismember(isolve,iindeterm)));
options_mcf.amcf_name = aindname;
options_mcf.amcf_title = aindtitle;
options_mcf.beha_title = 'NO indeterminacy';
......@@ -595,7 +603,7 @@ if length(iunstable)>0 || length(iwrong)>0,
end
if ~isempty(ixun),
itmp = itot(find(~ismember(itot,ixun)));
itmp = isolve(find(~ismember(isolve,ixun)));
options_mcf.amcf_name = aunstname;
options_mcf.amcf_title = aunsttitle;
options_mcf.beha_title = 'NO explosive solution';
......@@ -604,8 +612,8 @@ if length(iunstable)>0 || length(iwrong)>0,
mcf_analysis(lpmat, itmp, ixun, options_mcf, options_);
end
inorestriction = istable(find(~ismember(istable,irestriction))); % what went wrong beyong prior restrictions
iwrong = iwrong(find(~ismember(iwrong,inorestriction))); % what went wrong beyong prior restrictions
inorestriction = istable(find(~ismember(istable,irestriction))); % violation of prior restrictions
iwrong = iwrong(find(~ismember(iwrong,inorestriction))); % what went wrong beyond prior restrictions
if ~isempty(iwrong),
itmp = itot(find(~ismember(itot,iwrong)));
options_mcf.amcf_name = awrongname;
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment