v4: update of sensitivity code [Marco]

git-svn-id: https://www.dynare.org/svn/dynare/dynare_v4@681 ac1d8469-bf42-47a9-8791-bf33cf982152
time-shift
michel 2006-03-14 13:21:10 +00:00
parent 1f3036649e
commit 0514df601f
2 changed files with 94 additions and 10 deletions

View File

@ -1,6 +1,6 @@
function x0 = stab_map_(Nsam, fload, alpha2, alpha, prepSA, pprior)
function x0 = stab_map_(Nsam, fload, alpha2, prepSA, pprior)
%
% function x0 = stab_map_(Nsam, fload, alpha2, alpha, prepSA, pprior)
% function x0 = stab_map_(Nsam, fload, alpha2, prepSA, pprior)
%
% Mapping of stability regions in the prior ranges applying
% Monte Carlo filtering techniques.
@ -13,8 +13,6 @@ function x0 = stab_map_(Nsam, fload, alpha2, alpha, prepSA, pprior)
% fload = 0 to run new MC; 1 to load prevoiusly generated analysis
% alpha2 = significance level for bivariate sensitivity analysis
% [abs(corrcoef) > alpha2]
% alpha = significance level for univariate sensitivity analysis
% (uses smirnov)
% prepSA = 1: save transition matrices for mapping reduced form
% = 0: no transition matrix saved (default)
% pprior = 1: sample from prior ranges (default): sample saved in
@ -31,7 +29,7 @@ function x0 = stab_map_(Nsam, fload, alpha2, alpha, prepSA, pprior)
% - stable subset (dotted lines)
% - unstable subset (solid lines)
% 3) Bivariate plots of significant correlation patterns
% ( abs(corrcoef) > alpha2) under the stable subset
% ( abs(corrcoef) > alpha2) under the stable and unstable subsets
%
% USES lptauSEQ,
% stab_map_1, stab_map_2
@ -77,13 +75,10 @@ end
if nargin<3,
alpha2=0.3;
end
if nargin<4 | isempty(alpha),
alpha=0.002;
end
if nargin<5,
if nargin<4,
prepSA=0;
end
if nargin<6,
if nargin<5,
pprior=1;
end

89
matlab/stab_map_1.m Normal file
View File

@ -0,0 +1,89 @@
function proba = stab_map_1(lpmat, ibehaviour, inonbehaviour, aname, ishock)
%function stab_map_1(lpmat, ibehaviour, inonbehaviour, aname, ishock)
%
% lpmat = Monte Carlo matrix
% ibehaviour = index of behavioural runs
% inonbehaviour = index of non-behavioural runs
% ishock = 1 estimated shocks included
% ishock = 0 estimated shocks excluded (default)
%
% Plots: dotted lines for BEHAVIOURAL
% solid lines for NON BEHAVIOURAL
% USES smirnov
global estim_params_ bayestopt_ M_
if nargin<5,
ishock=0;
end
fname_ = M_.fname;
nshock = estim_params_.nvx;
nshock = nshock + estim_params_.nvn;
nshock = nshock + estim_params_.ncx;
nshock = nshock + estim_params_.ncn;
number_of_grid_points = 2^9; % 2^9 = 512 !... Must be a power of two.
bandwidth = 0; % Rule of thumb optimal bandwidth parameter.
kernel_function = 'gaussian'; % Gaussian kernel for Fast Fourrier Transform approximaton.
%kernel_function = 'uniform'; % Gaussian kernel for Fast Fourrier Transform approximaton.
if ishock,
npar = nshock + estim_params_.np;
else
npar = estim_params_.np;
end
for i=1:ceil(npar/12),
figure,
for j=1+12*(i-1):min(npar,12*i),
subplot(3,4,j-12*(i-1))
optimal_bandwidth = mh_optimal_bandwidth(lpmat(ibehaviour,j),length(ibehaviour),bandwidth,kernel_function);
[x1,f1] = kernel_density_estimate(lpmat(ibehaviour,j),number_of_grid_points,...
optimal_bandwidth,kernel_function);
plot(x1, f1,':k','linewidth',2)
optimal_bandwidth = mh_optimal_bandwidth(lpmat(inonbehaviour,j),length(inonbehaviour),bandwidth,kernel_function);
[x1,f1] = kernel_density_estimate(lpmat(inonbehaviour,j),number_of_grid_points,...
optimal_bandwidth,kernel_function);
hold on, plot(x1, f1,'k','linewidth',2)
%hist(lpmat(ibehaviour,j),30)
if ishock,
title(bayestopt_.name{j},'interpreter','none')
else
title(bayestopt_.name{j+nshock},'interpreter','none')
end
end
saveas(gcf,[fname_,'_',aname,'_',int2str(i)])
end
% Smirnov test for Blanchard;
for i=1:ceil(npar/12),
figure,
for j=1+12*(i-1):min(npar,12*i),
subplot(3,4,j-12*(i-1))
if ~isempty(ibehaviour),
h=cumplot(lpmat(ibehaviour,j));
set(h,'color',[0 0 0], 'linestyle',':')
end
hold on,
if ~isempty(inonbehaviour),
h=cumplot(lpmat(inonbehaviour,j));
set(h,'color',[0 0 0])
end
% if exist('kstest2')==2 & length(inonbehaviour)>0 & length(inonbehaviour)<Nsam,
% [H,P,KSSTAT] = kstest2(lpmat(ibehaviour,j),lpmat(inonbehaviour,j));
% title([bayestopt_.name{j+nshock},'. K-S prob ', num2str(P)])
% else
[H,P,KSSTAT] = smirnov(lpmat(ibehaviour,j),lpmat(inonbehaviour,j));
if ishock,
tittxt = bayestopt_.name{j};
else
tittxt = bayestopt_.name{j+nshock};
end
title([tittxt,'. K-S prob ', num2str(P)],'interpreter','none')
proba(j)=P;
% end
end
saveas(gcf,[fname_,'_',aname,'_SA_',int2str(i)])
end