Commit d8fea3ff authored by Marco Ratto's avatar Marco Ratto
Browse files

Correlations and measurement errors in identification:

1) Filter out measurement errors with error message that suggests to explicitly write measurement errors in model definition
2) allow identification checks with correlations, by switching to numerical derivatives
parent 389f4301
......@@ -546,10 +546,10 @@ if analytic_derivation,
end
if full_Hess,
[dum, DT, DOm, DYss, dum2, D2T, D2Om, D2Yss] = getH(A, B, Model,DynareResults,DynareOptions,kron_flag,indparam,indexo,iv);
[dum, DT, DOm, DYss, dum2, D2T, D2Om, D2Yss] = getH(A, B, EstimatedParameters, Model,DynareResults,DynareOptions,kron_flag,indparam,indexo,iv);
clear dum dum2;
else
[dum, DT, DOm, DYss] = getH(A, B, Model,DynareResults,DynareOptions,kron_flag,indparam,indexo,iv);
[dum, DT, DOm, DYss] = getH(A, B, EstimatedParameters, Model,DynareResults,DynareOptions,kron_flag,indparam,indexo,iv);
end
else
DT = derivatives_info.DT(iv,iv,:);
......
......@@ -187,10 +187,18 @@ if prior_exist,
name = bayestopt_.name;
name_tex = char(M_.exo_names_tex(indexo,:),M_.param_names_tex(indx,:));
offset = estim_params_.nvx;
offset = offset + estim_params_.nvn;
offset = offset + estim_params_.ncx;
offset = offset + estim_params_.ncn;
if estim_params_.nvn || estim_params_.ncn,
error('Identification does not support measurement errors. Instead, define them explicitly in measurement equations in model definition.')
else
offset = estim_params_.nvx;
%offset = offset + estim_params_.nvn;
offset = offset + estim_params_.ncx;
if estim_params_.ncx
options_ident.analytic_derivation=0;
options_ident.analytic_derivation_mode=-1;
end
%offset = offset + estim_params_.ncn;
end
else
indx = [1:M_.param_nbr];
indexo = [1:M_.exo_nbr];
......
function [H, dA, dOm, Hss, gp, d2A, d2Om, H2ss] = getH(A, B, M_,oo_,options_,kronflag,indx,indexo,iv)
% function [H, dA, dOm, Hss, gp, d2A, d2Om, H2ss] = getH(A, B, M_,oo_,options_,kronflag,indx,indexo,iv)
function [H, dA, dOm, Hss, gp, d2A, d2Om, H2ss] = getH(A, B, estim_params_,M_,oo_,options_,kronflag,indx,indexo,iv)
% function [H, dA, dOm, Hss, gp, d2A, d2Om, H2ss] = getH(A, B, estim_params_,M_,oo_,options_,kronflag,indx,indexo,iv)
% computes derivative of reduced form linear model w.r.t. deep params
%
% Inputs:
......@@ -49,16 +49,16 @@ function [H, dA, dOm, Hss, gp, d2A, d2Om, H2ss] = getH(A, B, M_,oo_,options_,kro
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
if nargin<6 || isempty(kronflag)
if nargin<7 || isempty(kronflag)
kronflag = 0;
end
if nargin<7 || isempty(indx)
if nargin<8 || isempty(indx)
indx = [];
end
if nargin<8 || isempty(indexo)
if nargin<9 || isempty(indexo)
indexo = [];
end
if nargin<9 || isempty(iv)
if nargin<10 || isempty(iv)
iv = (1:length(A))';
end
......@@ -96,7 +96,7 @@ if kronflag==-1, % perturbation
end
if nargout>5,
H2 = hessian_sparse('thet2tau',[sqrt(diag(M_.Sigma_e(indexo,indexo))); M_.params(indx)], ...
options_.gstep,M_, oo_, indx,indexo,0,[],[],[],iv);
options_.gstep,estim_params_,M_, oo_, indx,indexo,0,[],[],[],iv);
H2ss = zeros(m1,tot_param_nbr,tot_param_nbr);
iax=find(triu(rand(tot_param_nbr,tot_param_nbr)));
H2 = H2(:,iax);
......@@ -126,7 +126,7 @@ if kronflag==-2,
if nargout>5,
[residual, g1, g2 ] = feval([M_.fname,'_dynamic'],yy0, oo_.exo_steady_state', ...
M_.params, oo_.dr.ys, 1);
g22 = hessian_sparse('thet2tau',[M_.params(indx)],options_.gstep,M_, oo_, indx,[],-1);
g22 = hessian_sparse('thet2tau',[M_.params(indx)],options_.gstep,estim_params_,M_, oo_, indx,[],-1);
H2ss=full(g22(1:M_.endo_nbr,:));
H2ss = reshape(H2ss,[M_.endo_nbr param_nbr param_nbr]);
for j=1:M_.endo_nbr,
......@@ -147,7 +147,7 @@ if kronflag==-2,
[residual, g1 ] = feval([M_.fname,'_dynamic'],yy0, oo_.exo_steady_state', ...
M_.params, oo_.dr.ys, 1);
end
gp = fjaco('thet2tau',[M_.params(indx)],M_, oo_, indx,[],-1);
gp = fjaco('thet2tau',[M_.params(indx)],estim_params_,M_, oo_, indx,[],-1);
Hss=gp(1:M_.endo_nbr,:);
gp=gp(M_.endo_nbr+1:end,:);
gp = reshape(gp,[size(g1) param_nbr]);
......
function [JJ, H, gam, gp, dA, dOm, dYss] = getJJ(A, B, M_,oo_,options_,kronflag,indx,indexo,mf,nlags,useautocorr)
% function [JJ, H, gam, gp, dA, dOm, dYss] = getJJ(A, B, M_,oo_,options_,kronflag,indx,indexo,mf,nlags,useautocorr)
function [JJ, H, gam, gp, dA, dOm, dYss] = getJJ(A, B, estim_params_, M_,oo_,options_,kronflag,indx,indexo,mf,nlags,useautocorr)
% function [JJ, H, gam, gp, dA, dOm, dYss] = getJJ(A, B, estim_params_, M_,oo_,options_,kronflag,indx,indexo,mf,nlags,useautocorr)
% computes derivatives of 1st and 2nd order moments of observables with
% respect to estimated parameters
%
......@@ -54,16 +54,16 @@ function [JJ, H, gam, gp, dA, dOm, dYss] = getJJ(A, B, M_,oo_,options_,kronflag,
% You should have received a copy of the GNU General Public License
% along with Dynare. If not, see <http://www.gnu.org/licenses/>.
if nargin<7 || isempty(indx)
if nargin<8 || isempty(indx)
% indx = [1:M_.param_nbr];
end,
if nargin<8 || isempty(indexo)
if nargin<9 || isempty(indexo)
indexo = [];
end,
if nargin<10 || isempty(nlags)
if nargin<11 || isempty(nlags)
nlags=3;
end
if nargin<11 || isempty(useautocorr)
if nargin<12 || isempty(useautocorr)
useautocorr=0;
end
......@@ -73,15 +73,16 @@ warning('off','MATLAB:divideByZero')
if kronflag == -1,
fun = 'thet2tau';
params0 = M_.params;
JJ = fjaco(fun,[sqrt(diag(M_.Sigma_e(indexo,indexo))); M_.params(indx)],M_, oo_, indx,indexo,1,mf,nlags,useautocorr);
para0 = get_all_parameters(estim_params_, M_);
JJ = fjaco(fun,para0,estim_params_,M_, oo_, indx,indexo,1,mf,nlags,useautocorr);
M_.params = params0;
params0 = M_.params;
H = fjaco(fun,[sqrt(diag(M_.Sigma_e(indexo,indexo))); M_.params(indx)],M_, oo_, indx,indexo,0,mf,nlags,useautocorr);
H = fjaco(fun,para0,estim_params_,M_, oo_, indx,indexo,0,mf,nlags,useautocorr);
M_.params = params0;
params0 = M_.params;
gp = fjaco(fun,[sqrt(diag(M_.Sigma_e(indexo,indexo))); M_.params(indx)],M_, oo_, indx,indexo,-1);
gp = fjaco(fun,para0,estim_params_,M_, oo_, indx,indexo,-1);
M_.params = params0;
offset = length(indexo);
offset = length(para0)-length(indx);
gp = gp(:,offset+1:end);
dYss = H(1:M_.endo_nbr,offset+1:end);
dA = reshape(H(M_.orig_endo_nbr+[1:numel(A)],:),[size(A),size(H,2)]);
......@@ -92,7 +93,7 @@ if kronflag == -1,
assignin('base','M_', M_);
assignin('base','oo_', oo_);
else
[H, dA, dOm, dYss, gp] = getH(A, B, M_,oo_,options_,kronflag,indx,indexo);
[H, dA, dOm, dYss, gp] = getH(A, B, estim_params_,M_,oo_,options_,kronflag,indx,indexo);
gp = reshape(gp,size(gp,1)*size(gp,2),size(gp,3));
gp = [dYss; gp];
% if isempty(H),
......@@ -194,4 +195,4 @@ gam = [oo_.dr.ys(oo_.dr.order_var(mf)); gam];
% if useautocorr,
warning('on','MATLAB:divideByZero')
% end
\ No newline at end of file
% end
......@@ -80,7 +80,7 @@ if info(1)==0,
oo_.dr.ys, 1);
vg1 = [oo_.dr.ys(oo_.dr.order_var); vec(g1)];
[JJ, H, gam, gp, dA, dOm, dYss] = getJJ(A, B, M_,oo0,options_,kron_flag,indx,indexo,bayestopt_.mf2,nlags,useautocorr);
[JJ, H, gam, gp, dA, dOm, dYss] = getJJ(A, B, estim_params_, M_,oo0,options_,kron_flag,indx,indexo,bayestopt_.mf2,nlags,useautocorr);
derivatives_info.DT=dA;
derivatives_info.DOm=dOm;
derivatives_info.DYss=dYss;
......@@ -95,7 +95,7 @@ if info(1)==0,
disp('The number of moments with non-zero derivative is smaller than the number of parameters')
disp(['Try increasing ar = ', int2str(nlags+1)])
nlags=nlags+1;
[JJ, H, gam, gp, dA, dOm, dYss] = getJJ(A, B, M_,oo0,options_,kron_flag,indx,indexo,bayestopt_.mf2,nlags,useautocorr);
[JJ, H, gam, gp, dA, dOm, dYss] = getJJ(A, B, estim_params_, M_,oo0,options_,kron_flag,indx,indexo,bayestopt_.mf2,nlags,useautocorr);
derivatives_info.DT=dA;
derivatives_info.DOm=dOm;
derivatives_info.DYss=dYss;
......@@ -127,6 +127,9 @@ if info(1)==0,
else
normaliz1=[];
end
if ~isempty(estim_params_.corrx),
normaliz1 = [normaliz1 estim_params_.corrx(:,8)']; % normalize with prior standard deviation
end
if ~isempty(estim_params_.param_vals),
normaliz1 = [normaliz1 estim_params_.param_vals(:,7)']; % normalize with prior standard deviation
end
......
function tau = thet2tau(params, M_, oo_, indx, indexo, flagmoments,mf,nlags,useautocorr,iv)
function tau = thet2tau(params, estim_params_, M_, oo_, indx, indexo, flagmoments,mf,nlags,useautocorr,iv)
%
% Copyright (C) 2011-2012 Dynare Team
......@@ -25,20 +25,24 @@ if nargin==1,
indexo = [];
end
if nargin<6,
if nargin<7,
flagmoments=0;
end
if nargin<9 || isempty(useautocorr),
if nargin<10 || isempty(useautocorr),
useautocorr=0;
end
if nargin<10 || isempty(iv),
if nargin<11 || isempty(iv),
iv=[1:M_.endo_nbr];
end
M_.params(indx) = params(length(indexo)+1:end);
if ~isempty(indexo)
M_.Sigma_e(indexo,indexo) = diag(params(1:length(indexo)).^2);
if length(params)>length(indx),
M_ = set_all_parameters(params,estim_params_,M_);
else
M_.params(indx) = params;
end
% if ~isempty(indexo)
% M_.Sigma_e(indexo,indexo) = diag(params(1:length(indexo)).^2);
% end
[A,B,tele,tubbies,M_,options_,oo_] = dynare_resolve(M_,options_,oo_);
if flagmoments==0,
ys=oo_.dr.ys(oo_.dr.order_var);
......@@ -71,4 +75,4 @@ else
tau = [tau;vec(dum(mf,mf))];
end
tau = [ oo_.dr.ys(oo_.dr.order_var(mf)); tau];
end
\ No newline at end of file
end
Supports Markdown
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