Commit 188d2694 authored by Sébastien Villemot's avatar Sébastien Villemot
Browse files

Merge pull request #379 from JohannesPfeifer/master

Bugfix in initial_estimation_checks; making graphs more explicit
parents 5965c8e3 50f97ffc
/*
* This file implements the Baseline New Keynesian DSGE model described in
* much detail in Jess Fernndez-Villaverde (2006): "A Baseline DSGE
* Model", available at http://economics.sas.upenn.edu/~jesusfv/benchmark_DSGE.pdf
*
* The parametrization is based on the estimated version of this model in
* Jess Fernndez-Villaverde (2010): "The econometrics of DSGE models",
* SERIEs, Vol. 1, pp. 3-49, DOI 10.1007/s13209-009-0014-7
*
* This implementation was written by Benjamin Born and Johannes Pfeifer. In
* case you spot mistakes, email us at jpfeifer@gmx.de
*
* This mod-file implements a non-linearized version of the New Keynesian
* model based on a recursive formulation of the price and wage setting
* equations. Moreover, it makes use of a steady state file to i) set
* parameters that depend on other parameters that are potentially estimated
* and ii) solve a nonlinear equation using a numerical solver to find the steady
* state of labor.
* The model is written in the beginning of period stock notation. To make the model
* conform with Dynare's end of period stock notation, we use the
* predetermined_variables-command.
*
* Please note that the following copyright notice only applies to this Dynare
* implementation of the model.
*/
/*
* Copyright (C) 2013 Dynare Team
*
* 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
* along with Dynare. If not, see <http://www.gnu.org/licenses/>.
*/
var d //preference shock
c //consumption
mu_z //trend growth rate of the economy (from neutral and investment specific technology)
mu_I //growth rate of investment-specific technology growth
mu_A //growth rate of neutral technology
lambda //Lagrange multiplier
R //Nominal Interest rate
PI //Inflation
r //rental rate of capital
x //investment
u //capacity utilization
q //Tobin's marginal q
f //variable for recursive formulation of wage setting
ld //aggregate labor demand
w //real wage
wstar //optimal real wage
PIstarw //optimal wage inflation
PIstar //optimal price inflation
g1 //variable 1 for recursive formulation of price setting
g2 //variable 2 for recursive formulation of price setting
yd //aggregate output
mc //marginal costs
k //capital
vp //price dispersion term
vw //wage dispersion term
l //aggregate labor bundle
phi //labor disutility shock
F; //firm profits
varexo epsd epsphi epsmu_I epsA epsm;
predetermined_variables k;
parameters h //consumption habits
betta //discount factor
gammma1 //capital utilization, linear term
gammma2 //capital utilization, quadratic term
delta //depreciation rate
kappa //capital adjustment costs parameter
eta //elasticity of substitution between labor varieties
epsilon //elasticity of substitution between goods varieties
varpsi //labor disutility parameter
gammma //inverse Frisch elasticity
chiw //wage indexation parameter
chi //price indexation
thetap //Calvo parameter prices
thetaw //Calvo parameter wages
alppha //capital share
Rbar //steady state interest rate
PIbar //steady state inflation
gammmaR //interest smoothing coefficient Taylor rule
gammmaPI //feedback coefficient to inflation monetary policy rule
gammmay //feedback coefficient to output growth deviation in monetary policy rule
Phi //firms fixed costs
rhod //autocorrelation preference shock
rhophi //autocorrelation labor disutility shock
Lambdamu //steady state growth rate of investmentment-specific technology
LambdaA //steady state neutral technology growth
Lambdax //steady state growth rate of investment
LambdaYd //steady state growth rate of output
sigma_d //standard deviation preference shock
sigma_phi //standard deviation labor disutility shock
sigma_mu //standard deviation investment-specific technology
sigma_A //standard deviation neutral technology
sigma_m; //standard deviation preference shock
//Note that the parameter naming in FV(2010) differs from FV(2006)
//Fixed parameters, taken from FV(2010), Table 2, p. 37
delta=0.025;
epsilon=10;
eta= 10;
Phi=0;
gammma2=0.001;
//Estimated parameters, taken from FV(2010), Table 3, p. 38, median estimate parameters
betta =0.998;
h=0.97;
varpsi =8.92;
gammma = 1.17;
kappa =9.51;
alppha =0.21;
thetap =0.82;
chi = 0.63;
thetaw =0.68;
chiw =0.62;
gammmaR =0.77;
gammmay =0.19;
gammmaPI =1.29;
PIbar = 1.01;
rhod = 0.12;
rhophi = 0.93;
sigma_A = -3.97;
sigma_d = -1.51;
sigma_phi =-2.36;
sigma_mu =-5.43;
sigma_m =-5.85;
Lambdamu=3.4e-3;
LambdaA = 2.8e-3;
LambdaYd= (LambdaA+alppha*Lambdamu)/(1-alppha);
/*
The following parameters are set in the steady state file as they depend on other
deep parameters that were estimated in the original study. Setting them in the
steady state file means they are updated for every parameter draw in the MCMC
algorithm, while the parameters initialized here are only set once for the initial
values of the parameters they depend on:
gammma1 as it depends on LambdaA, alppha, Lambdamu, betta, and delta
Rbar =0 as it depends on PI, LambdaA, alppha, Lambdamu, and betta
Lambdax
*/
/*
The following model equations are the stationary model equations, taken from
FV(2006), p. 20, section 3.2.
*/
model;
//1. FOC consumption
d*(c-h*c(-1)*mu_z^(-1))^(-1)-h*betta*d(+1)*(c(+1)*mu_z(+1)-h*c)^(-1)=lambda;
//2. Euler equation
lambda=betta*lambda(+1)*mu_z(+1)^(-1)/PI(+1)*R;
//3. FOC capital utilization
r=gammma1+gammma2*(u-1);
//4. FOC capital
q=betta*lambda(+1)/lambda*mu_z(+1)^(-1)*mu_I(+1)^(-1)*((1-delta)*q(+1)+r(+1)*u(+1)-(gammma1*(u(+1)-1)+gammma2/2*(u(+1)-1)^2));
//5. FOC investment
1=q*(1-(kappa/2*(x/x(-1)*mu_z-Lambdax)^2)-(kappa*(x/x(-1)*mu_z-Lambdax)*x/x(-1)*mu_z))
+betta*q(+1)*lambda(+1)/lambda*mu_z(+1)^(-1)*kappa*(x(+1)/x*mu_z(+1)-Lambdax)*(x(+1)/x*mu_z(+1))^2;
//6-7. Wage setting
f=(eta-1)/eta*wstar^(1-eta)*lambda*w^eta*ld+betta*thetaw*(PI^chiw/PI(+1))^(1-eta)*(wstar(+1)/wstar*mu_z(+1))^(eta-1)*f(+1);
f=varpsi*d*phi*PIstarw^(-eta*(1+gammma))*ld^(1+gammma)+betta*thetaw*(PI^chiw/PI(+1))^(-eta*(1+gammma))*(wstar(+1)/wstar*mu_z(+1))^(eta*(1+gammma))*f(+1);
//8-10. firm's price setting
g1=lambda*mc*yd+betta*thetap*(PI^chi/PI(+1))^(-epsilon)*g1(+1);
g2=lambda*PIstar*yd+betta*thetap*(PI^chi/PI(+1))^(1-epsilon)*PIstar/PIstar(+1)*g2(+1);
epsilon*g1=(epsilon-1)*g2;
//11-12. optimal inputs
u*k/ld=alppha/(1-alppha)*w/r*mu_z*mu_I;
mc=(1/(1-alppha))^(1-alppha)*(1/alppha)^alppha*w^(1-alppha)*r^alppha;
//13. law of motion wages
1=thetaw*(PI(-1)^chiw/PI)^(1-eta)*(w(-1)/w*mu_z^(-1))^(1-eta)+(1-thetaw)*PIstarw^(1-eta);
//14. law of motion prices
1=thetap*(PI(-1)^chi/PI)^(1-epsilon)+(1-thetap)*PIstar^(1-epsilon);
//15. Taylor Rule
R/Rbar=(R(-1)/Rbar)^gammmaR*((PI/PIbar)^gammmaPI*((yd/yd(-1)*mu_z)/exp(LambdaYd))^gammmay)^(1-gammmaR)*exp(epsm);
//16-17. Market clearing
yd=c+x+mu_z^(-1)*mu_I^(-1)*(gammma1*(u-1)+gammma2/2*(u-1)^2)*k;
yd=(mu_A*mu_z^(-1)*(u*k)^alppha*ld^(1-alppha)-Phi)/vp;
//18-20. Price and wage dispersion terms
l=vw*ld;
vp=thetap*(PI(-1)^chi/PI)^(-epsilon)*vp(-1)+(1-thetap)*PIstar^(-epsilon);
vw=thetaw*(w(-1)/w*mu_z^(-1)*PI(-1)^chiw/PI)^(-eta)*vw(-1)+(1-thetaw)*(PIstarw)^(-eta);
//21. Law of motion for capital
k(+1)*mu_z*mu_I-(1-delta)*k-mu_z*mu_I*(1-kappa/2*(x/x(-1)*mu_z-Lambdax)^2)*x=0;
//22. Profits
F=yd-1/(1-alppha)*w*ld;
//23. definition optimal wage inflation
PIstarw=wstar/w;
//exogenous processes
//24. Preference Shock
log(d)=rhod*log(d(-1))+epsd;
//25. Labor disutility Shock
log(phi)=rhophi*log(phi(-1))+epsphi;
//26. Investment specific technology
log(mu_I)=Lambdamu+epsmu_I;
//27. Neutral technology
log(mu_A)=LambdaA+epsA;
//28. Defininition composite technology
mu_z=mu_A^(1/(1-alppha))*mu_I^(alppha/(1-alppha));
end;
shocks;
var epsd; stderr exp(sigma_d);
var epsphi; stderr exp(sigma_phi);
var epsmu_I; stderr exp(sigma_mu);
var epsA; stderr exp(sigma_A);
var epsm; stderr exp(sigma_m);
end;
steady;
check;
stoch_simul(order=1,irf=20) yd c R PI;
function [ys,check] = NK_baseline_steadystate(ys,exe)
global M_ lgy_
if isfield(M_,'param_nbr') == 1
NumberOfParameters = M_.param_nbr;
for i = 1:NumberOfParameters
paramname = deblank(M_.param_names(i,:));
eval([ paramname ' = M_.params(' int2str(i) ');']);
end
check = 0;
end
%% Enter model equations here
options=optimset(); % set options for numerical solver
% the steady state computation follows FV (2006), section 4.1
PI=PIbar;
u=1;
q=1;
d=1;
phi=1;
m=0;
zeta=1;
mu_z=exp(LambdaYd);
mu_I=exp(Lambdamu);
mu_A=exp(LambdaA);
%set the parameter Lambdax
Lambdax=mu_z;
%set the parameter gammma1
gammma1=mu_z*mu_I/betta-(1-delta);
r=1*gammma1;
R=1+(PI*mu_z/betta-1);
%set Rbar
Rbar=R;
PIstar=((1-thetap*PI^(-(1-epsilon)*(1-chi)))/(1-thetap))^(1/(1-epsilon));
PIstarw=((1-thetaw*PI^(-(1-chiw)*(1-eta))*mu_z^(-(1-eta)))/(1-thetaw))^(1/(1-eta));
mc=(epsilon-1)/epsilon*(1-betta*thetap*PI^((1-chi)*epsilon))/(1-betta*thetap*PI^(-(1-epsilon)*(1-chi)))*PIstar;
w=(1-alppha)*(mc*(alppha/r)^alppha)^(1/(1-alppha));
wstar=w*PIstarw;
vp=(1-thetap)/(1-thetap*PI^((1-chi)*epsilon))*PIstar^(-epsilon);
vw=(1-thetaw)/(1-thetaw*PI^((1-chiw)*eta)*mu_z^eta)*PIstarw^(-eta);
tempvaromega=alppha/(1-alppha)*w/r*mu_z*mu_I;
ld=fsolve(@(ld)(1-betta*thetaw*mu_z^(eta-1)*PI^(-(1-chiw)*(1-eta)))/(1-betta*thetaw*mu_z^(eta*(1+gammma))*PI^(eta*(1-chiw)*(1+gammma)))...
-(eta-1)/eta*wstar/(varpsi*PIstarw^(-eta*gammma)*ld^gammma)*((1-h*mu_z^(-1))^(-1)-betta*h*(mu_z-h)^(-1))*...
((mu_A*mu_z^(-1)*vp^(-1)*tempvaromega^alppha-tempvaromega*(1-(1-delta)*(mu_z*mu_I)^(-1)))*ld-vp^(-1)*Phi)^(-1),0.25,options);
l=vw*ld;
k=tempvaromega*ld;
x=(1-(1-delta)*(mu_z*mu_I)^(-1))*k;
yd=(mu_A/mu_z*k^alppha*ld^(1-alppha)-Phi)/vp;
c=(mu_A*mu_z^(-1)*vp^(-1)*tempvaromega^alppha-tempvaromega*(1-(1-delta)*(mu_z*mu_I)^(-1)))*ld-vp^(-1)*Phi;
lambda=(1-h*betta*mu_z^(-1))*(1-h/mu_z)^(-1)/c;
F=yd-1/(1-alppha)*w*ld;
f=(eta-1)/eta*wstar*PIstarw^(-eta)*lambda*ld/(1-betta*thetaw*mu_z^(eta-1)*PI^(-(1-chiw)*(1-eta)));
f2=varpsi*d*phi*PIstarw^(-eta*(1+gammma))*ld^(1+gammma)/(1-betta*thetaw*(PI^chiw/PI)^(-eta*(1+gammma))*(wstar/wstar*mu_z)^(eta*(1+gammma)));
g1=lambda*mc*yd/(1-betta*thetap*PI^((1-chi)*epsilon));
g2=epsilon/(epsilon-1)*g1;
%% end own model equations
for iter = 1:length(M_.params)
eval([ 'M_.params(' num2str(iter) ') = ' M_.param_names(iter,:) ';' ])
end
if isfield(M_,'param_nbr') == 1
if isfield(M_,'orig_endo_nbr') == 1
NumberOfEndogenousVariables = M_.orig_endo_nbr;
else
NumberOfEndogenousVariables = M_.endo_nbr;
end
ys = zeros(NumberOfEndogenousVariables,1);
for i = 1:NumberOfEndogenousVariables
varname = deblank(M_.endo_names(i,:));
eval(['ys(' int2str(i) ') = ' varname ';']);
end
else
ys=zeros(length(lgy_),1);
for i = 1:length(lgy_)
ys(i) = eval(lgy_(i,:));
end
check = 0;
end
function [omega,f] = UnivariateSpectralDensity(dr,var_list)
% This function computes the theoretical spectral density of each
% endogenous variable declared in var_list. Results are stored in
% oo_ and may be plotted.
% oo_ and may be plotted. Plots are saved into the graphs-folder.
%
% Adapted from th_autocovariances.m.
% Copyright (C) 2006-2012 Dynare Team
% Copyright (C) 2006-2013 Dynare Team
%
% This file is part of Dynare.
%
......@@ -167,12 +167,20 @@ else
end
if pltinfo
if ~exist(M_.fname, 'dir')
mkdir('.',M_.fname);
end
if ~exist([M_.fname '/graphs'])
mkdir(M_.fname,'graphs');
end
for i= 1:nvar
figure('Name',['Spectral Density of ' deblank(M_.endo_names(ivar(i),:)) '.'])
hh = dyn_figure(options_,'Name',['Spectral Density of ' deblank(M_.endo_names(ivar(i),:)) '.']);
plot(omega,f(i,:),'-k','linewidth',2)
xlabel('0 \leq \omega \leq \pi')
ylabel('f(\omega)')
box on
axis tight
axis tight
dyn_saveas(hh,[M_.fname ,filesep,'graphs', filesep, 'SpectralDensity_' deblank(M_.endo_names(ivar(i),:))],options_)
end
end
......@@ -79,8 +79,8 @@ if isempty(pdfinfo)
else
number_of_the_last_post_data_file = length(pdfinfo);
name_of_the_last_post_data_file = ...
[ './' M_.dname ...
'/metropolis/' ...
[ pwd filesep M_.dname ...
filesep 'metropolis' filesep ...
M_.fname '_' ...
generic_post_data_file_name ...
int2str(number_of_the_last_post_data_file) ...
......
......@@ -180,7 +180,7 @@ if ~isequal(DynareOptions.mode_compute,1) && any(xparam1<BayesInfo.lb)
k = find(xparam1<BayesInfo.lb);
fval = objective_function_penalty_base+sum((BayesInfo.lb(k)-xparam1(k)).^2);
exit_flag = 0;
info = 41;
info = [41 k'];
if analytic_derivation,
DLIK=ones(length(xparam1),1);
end
......@@ -192,7 +192,7 @@ if ~isequal(DynareOptions.mode_compute,1) && any(xparam1>BayesInfo.ub)
k = find(xparam1>BayesInfo.ub);
fval = objective_function_penalty_base+sum((xparam1(k)-BayesInfo.ub(k)).^2);
exit_flag = 0;
info = 42;
info = [42 k'];
if analytic_derivation,
DLIK=ones(length(xparam1),1);
end
......
......@@ -1127,6 +1127,7 @@ if (~((any(bayestopt_.pshape > 0) && options_.mh_replic) || (any(bayestopt_.psha
plot(1:gend,measurement_error(index(k),:),'-k','linewidth',1)
hold off
name = deblank(options_.varobs(index(k),:));
xlim([1 gend])
if isempty(NAMES)
NAMES = name;
else
......
......@@ -121,6 +121,15 @@ n_varobs = size(options_.varobs,1);
% Set priors over the estimated parameters.
if ~isempty(estim_params_)
[xparam1,estim_params_,bayestopt_,lb,ub,M_] = set_prior(estim_params_,M_,options_);
if ~isempty(options_.mode_file) && ~options_.mh_posterior_mode_estimation
junk=length(xparam1);
load(options_.mode_file,'xparam1');
if length(xparam1) ~= junk
error([ 'ESTIMATION: the posterior mode file ' options_.mode_file ' has been generated using another specification. Please delete it and recompute the posterior mode.'])
end
end
if any(bayestopt_.pshape > 0)
% Plot prior densities.
if ~options_.nograph && options_.plot_priors
......@@ -145,7 +154,7 @@ if ~isempty(estim_params_)
for ii=2:size(outside_bound_vars,1)
disp_string=[disp_string,', ',outside_bound_vars{ii,:}];
end
error(['Initial value(s) of ', disp_string ,' are outside parameter bounds. Potentially, you should set prior_trunc=0.'])
error(['Initial value(s) of ', disp_string ,' are outside parameter bounds. Potentially, you should set prior_trunc=0. If you used the mode_file-option, check whether your mode-file is consistent with the priors.'])
end
lb = bounds(:,1);
ub = bounds(:,2);
......
......@@ -267,13 +267,13 @@ if options_gsa.redform && ~isempty(options_gsa.namendo),% ...
x0 = stab_map_(OutputDirectoryName,options_gsa);
end
if strmatch(':',options_gsa.namendo,'exact'),
options_gsa.namendo=M_.endo_names;
options_gsa.namendo=M_.endo_names(1:M_.orig_endo_nbr,:);
end
if strmatch(':',options_gsa.namexo,'exact'),
options_gsa.namexo=M_.exo_names;
end
if strmatch(':',options_gsa.namlagendo,'exact'),
options_gsa.namlagendo=M_.endo_names;
options_gsa.namlagendo=M_.endo_names(1:M_.orig_endo_nbr,:);
end
% options_.opt_gsa = options_gsa;
if options_gsa.morris==1,
......
......@@ -67,7 +67,7 @@ for j= 1:nvar
dyn_saveas(hh,[ fname '/graphs/forcst' int2str(n_fig)],options_);
n_fig =n_fig+1;
eval(['hh=dyn_figure(options_,''Name'',''Forecast (' int2str(n_fig) ')'');']);
eval(['hh=dyn_figure(options_,''Name'',''Forecasts (' int2str(n_fig) ')'');']);
m = 1;
end
subplot(nr,nc,m);
......@@ -86,6 +86,7 @@ for j= 1:nvar
hold on
plot([NaN(obs,1); hpdsup.(vn)]);
title(vn,'Interpreter','none');
xlim([1 obs+length(hpdsup.(vn))])
hold off
m = m + 1;
end
......
......@@ -74,11 +74,11 @@ end
for j=1:length(a),
if strmatch([fname_,tmp],a(j).name),
disp(a(j).name)
delete([OutDir,'/',a(j).name])
delete([OutDir,filesep,a(j).name])
end,
if strmatch([fname_,tmp1],a(j).name),
disp(a(j).name)
delete([OutDir,'/',a(j).name])
delete([OutDir,filesep,a(j).name])
end,
end
disp('done !')
......@@ -128,7 +128,7 @@ if ~loadSA,
for jx=1:gend, data_indx(jx,data_index{jx})=true; end
%stock_gend=data_info.gend;
%stock_data = data_info.data;
load([DirectoryName '/' M_.fname '_data.mat']);
load([DirectoryName filesep M_.fname '_data.mat']);
filfilt = dir([DirectoryName filesep M_.fname '_filter_step_ahead*.mat']);
filsmooth = dir([DirectoryName filesep M_.fname '_smooth*.mat']);
filupdate = dir([DirectoryName filesep M_.fname '_update*.mat']);
......@@ -137,7 +137,7 @@ if ~loadSA,
logpo2=[];
sto_ys=[];
for j=1:length(filparam),
%load([DirectoryName '/' M_.fname '_param',int2str(j),'.mat']);
%load([DirectoryName filesep M_.fname '_param',int2str(j),'.mat']);
if isempty(strmatch([M_.fname '_param_irf'],filparam(j).name))
load([DirectoryName filesep filparam(j).name]);
x=[x; stock];
......@@ -294,9 +294,18 @@ else
if ~options_.nograph,
ifig=0;
for i=1:size(vvarvecm,1),
if options_.opt_gsa.ppost
temp_name='RMSE Posterior: Log Prior';
else
if options_.opt_gsa.pprior
temp_name='RMSE Prior: Log Prior';
else
temp_name='RMSE MC: Log Prior';
end
end
if mod(i,9)==1,
ifig=ifig+1;
hh=dyn_figure(options_,'name',['Prior ',int2str(ifig)]);
hh=dyn_figure(options_,'name',[temp_name,' ',int2str(ifig)]);
end
subplot(3,3,i-9*(ifig-1))
h=cumplot(lnprior(ixx(1:nfilt0(i),i)));
......@@ -307,21 +316,30 @@ else
title(vvarvecm(i,:),'interpreter','none')
if mod(i,9)==0 || i==size(vvarvecm,1)
if options_.opt_gsa.ppost
dyn_saveas(hh,[OutDir '/' fname_ '_rmse_post_lnprior',int2str(ifig)],options_);
dyn_saveas(hh,[OutDir filesep fname_ '_rmse_post_lnprior',int2str(ifig)],options_);
else
if options_.opt_gsa.pprior
dyn_saveas(hh,[OutDir '/' fname_ '_rmse_prior_lnprior',int2str(ifig) ],options_);
dyn_saveas(hh,[OutDir filesep fname_ '_rmse_prior_lnprior',int2str(ifig) ],options_);
else
dyn_saveas(hh,[OutDir '/' fname_ '_rmse_mc_lnprior',int2str(ifig) ],options_);
dyn_saveas(hh,[OutDir filesep fname_ '_rmse_mc_lnprior',int2str(ifig) ],options_);
end
end
end
end
ifig=0;
for i=1:size(vvarvecm,1),
if options_.opt_gsa.ppost
temp_name='RMSE Posterior: Log Likelihood';
else
if options_.opt_gsa.pprior
temp_name='RMSE Prior: Log Likelihood';
else
temp_name='RMSE MC: Log Likelihood';
end
end
if mod(i,9)==1,
ifig=ifig+1;
hh = dyn_figure(options_,'Name',['Likelihood ',int2str(ifig)]);
hh = dyn_figure(options_,'Name',[temp_name,' ',int2str(ifig)]);
end
subplot(3,3,i-9*(ifig-1))
h=cumplot(likelihood(ixx(1:nfilt0(i),i)));
......@@ -335,21 +353,30 @@ else
end
if mod(i,9)==0 || i==size(vvarvecm,1)
if options_.opt_gsa.ppost
dyn_saveas(hh,[OutDir '/' fname_ '_rmse_post_lnlik',int2str(ifig) ],options_);
dyn_saveas(hh,[OutDir filesep fname_ '_rmse_post_lnlik',int2str(ifig) ],options_);
else
if options_.opt_gsa.pprior
dyn_saveas(hh,[OutDir '/' fname_ '_rmse_prior_lnlik',int2str(ifig)],options_);
dyn_saveas(hh,[OutDir filesep fname_ '_rmse_prior_lnlik',int2str(ifig)],options_);
else
dyn_saveas(hh,[OutDir '/' fname_ '_rmse_mc_lnlik',int2str(ifig) ],options_);
dyn_saveas(hh,[OutDir filesep fname_ '_rmse_mc_lnlik',int2str(ifig) ],options_);
end
end
end
end
ifig=0;
for i=1:size(vvarvecm,1),
if options_.opt_gsa.ppost
temp_name='RMSE Posterior: Log Posterior';
else
if options_.opt_gsa.pprior
temp_name='RMSE Prior: Log Posterior';
else
temp_name='RMSE MC: Log Posterior';
end
end
if mod(i,9)==1,
ifig=ifig+1;
hh = dyn_figure(options_,'Name',['Posterior ',int2str(ifig)]);
hh = dyn_figure(options_,'Name',[temp_name,' ',int2str(ifig)]);
end
subplot(3,3,i-9*(ifig-1))
h=cumplot(logpo2(ixx(1:nfilt0(i),i)));
......@@ -363,12 +390,12 @@ else
end
if mod(i,9)==0 || i==size(vvarvecm,1)
if options_.opt_gsa.ppost
dyn_saveas(hh,[OutDir '/' fname_ '_rmse_post_lnpost',int2str(ifig) ],options_);
dyn_saveas(hh,[OutDir filesep fname_ '_rmse_post_lnpost',int2str(ifig) ],options_);
else
if options_.opt_gsa.pprior
dyn_saveas(hh,[OutDir '/' fname_ '_rmse_prior_lnpost',int2str(ifig)],options_);
dyn_saveas(hh,[OutDir filesep fname_ '_rmse_prior_lnpost',int2str(ifig)],options_);
else