diff --git a/matlab/dsge_likelihood.m b/matlab/dsge_likelihood.m index 62bf68c1fa4a5ffab0e7650cf75642556bf0cf5c..8920413aa620e869cedcc84263a3b6790fb127fe 100644 --- a/matlab/dsge_likelihood.m +++ b/matlab/dsge_likelihood.m @@ -206,7 +206,7 @@ Q = Model.Sigma_e; H = Model.H; % Test if Q is positive definite. -if ~issquare(Q) && EstimatedParameters.ncx +if ~issquare(Q) || EstimatedParameters.ncx || isfield(EstimatedParameters,'calibrated_covariances') [Q_is_positive_definite, penalty] = ispd(Q); if ~Q_is_positive_definite fval = objective_function_penalty_base+penalty; @@ -214,10 +214,20 @@ if ~issquare(Q) && EstimatedParameters.ncx info = 43; return end + if isfield(EstimatedParameters,'calibrated_covariances') + correct_flag=check_consistency_covariances(Q); + if ~correct_flag + penalty = sum(Q(EstimatedParameters.calibrated_covariances.position).^2); + fval = objective_function_penalty_base+penalty; + exit_flag = 0; + info = 71; + return + end + end end % Test if H is positive definite. -if ~issquare(H) && EstimatedParameters.ncn +if ~issquare(H) || EstimatedParameters.ncn || isfield(EstimatedParameters,'calibrated_covariances_ME') [H_is_positive_definite, penalty] = ispd(H); if ~H_is_positive_definite fval = objective_function_penalty_base+penalty; @@ -225,6 +235,16 @@ if ~issquare(H) && EstimatedParameters.ncn info = 44; return end + if isfield(EstimatedParameters,'calibrated_covariances_ME') + correct_flag=check_consistency_covariances(H); + if ~correct_flag + penalty = sum(H(EstimatedParameters.calibrated_covariances_ME.position).^2); + fval = objective_function_penalty_base+penalty; + exit_flag = 0; + info = 72; + return + end + end end diff --git a/matlab/dynare_estimation_1.m b/matlab/dynare_estimation_1.m index 2c8783bfd2674bbceda38087da391909c004496e..da865af16000afd5be534462df74724a113d55ba 100644 --- a/matlab/dynare_estimation_1.m +++ b/matlab/dynare_estimation_1.m @@ -80,38 +80,16 @@ end [dataset_,xparam1, hh, M_, options_, oo_, estim_params_,bayestopt_] = dynare_estimation_init(var_list_, dname, [], M_, options_, oo_, estim_params_, bayestopt_); +%check for calibrated covariances before updating parameters +if ~isempty(estim_params_) + estim_params_=check_for_calibrated_covariances(xparam1,estim_params_,M_); +end % Set sigma_e_is_diagonal flag (needed if the shocks block is not declared in the mod file). M_.sigma_e_is_diagonal = 1; -if estim_params_.ncx || any(nnz(tril(M_.Sigma_e,-1))) +if estim_params_.ncx || any(nnz(tril(M_.Correlation_matrix,-1))) || isfield(estim_params_,'calibrated_covariances') M_.sigma_e_is_diagonal = 0; end -% Set the correlation matrix if necessary. -if ~isequal(estim_params_.ncx,nnz(tril(M_.Sigma_e,-1))) - M_.Correlation_matrix = diag(1./sqrt(diag(M_.Sigma_e)))*M_.Sigma_e*diag(1./sqrt(diag(M_.Sigma_e))); - % Remove NaNs appearing because of variances calibrated to zero. - if any(isnan(M_.Correlation_matrix)) - zero_variance_idx = find(~diag(M_.Sigma_e)); - for i=1:length(zero_variance_idx) - M_.Correlation_matrix(zero_variance_idx(i),:) = 0; - M_.Correlation_matrix(:,zero_variance_idx(i)) = 0; - end - end -end - -% Set the correlation matrix of measurement errors if necessary. -if ~isequal(estim_params_.ncn,nnz(tril(M_.H,-1))) - M_.Correlation_matrix_ME = diag(1./sqrt(diag(M_.H)))*M_.H*diag(1./sqrt(diag(M_.H))); - % Remove NaNs appearing because of variances calibrated to zero. - if any(isnan(M_.Correlation_matrix_ME)) - zero_variance_idx = find(~diag(M_.H)); - for i=1:length(zero_variance_idx) - M_.Correlation_matrix_ME(zero_variance_idx(i),:) = 0; - M_.Correlation_matrix_ME(:,zero_variance_idx(i)) = 0; - end - end -end - data = dataset_.data; rawdata = dataset_.rawdata; data_index = dataset_.missing.aindex; @@ -151,7 +129,7 @@ if exist(mh_scale_fname) end if ~isempty(estim_params_) - set_parameters(xparam1); + M_ = set_all_parameters(xparam1,estim_params_,M_); end % compute sample moments if needed (bvar-dsge) diff --git a/matlab/non_linear_dsge_likelihood.m b/matlab/non_linear_dsge_likelihood.m index 41ba13fdadd3030982c374d4a0ec59a9be8fb81d..474b8dbb1c7c9cfe570f16b7ea0e5914288eed8c 100644 --- a/matlab/non_linear_dsge_likelihood.m +++ b/matlab/non_linear_dsge_likelihood.m @@ -163,7 +163,7 @@ Model = set_all_parameters(xparam1,EstimatedParameters,Model); Q = Model.Sigma_e; H = Model.H; -if ~isscalar(Q) && EstimatedParameters.ncx +if ~issquare(Q) || EstimatedParameters.ncx || isfield(EstimatedParameters,'calibrated_covariances') [Q_is_positive_definite, penalty] = ispd(Q); if ~Q_is_positive_definite fval = objective_function_penalty_base+penalty; @@ -171,9 +171,20 @@ if ~isscalar(Q) && EstimatedParameters.ncx info = 43; return end + if isfield(EstimatedParameters,'calibrated_covariances') + correct_flag=check_consistency_covariances(Q); + if ~correct_flag + penalty = sum(Q(EstimatedParameters.calibrated_covariances.position).^2); + fval = objective_function_penalty_base+penalty; + exit_flag = 0; + info = 71; + return + end + end + end -if ~isscalar(H) && EstimatedParameters.ncn +if ~issquare(H) || EstimatedParameters.ncn || isfield(EstimatedParameters,'calibrated_covariances_ME') [H_is_positive_definite, penalty] = ispd(H); if ~H_is_positive_definite fval = objective_function_penalty_base+penalty; @@ -181,6 +192,17 @@ if ~isscalar(H) && EstimatedParameters.ncn info = 44; return end + if isfield(EstimatedParameters,'calibrated_covariances_ME') + correct_flag=check_consistency_covariances(H); + if ~correct_flag + penalty = sum(H(EstimatedParameters.calibrated_covariances_ME.position).^2); + fval = objective_function_penalty_base+penalty; + exit_flag = 0; + info = 72; + return + end + end + end %------------------------------------------------------------------------------ diff --git a/matlab/print_info.m b/matlab/print_info.m index 157275fb4b86f087d89a0ae5f9327f9e7e9863a2..691bda930031bba3e76039dee37d7ad85fa6c111 100644 --- a/matlab/print_info.m +++ b/matlab/print_info.m @@ -113,6 +113,10 @@ if ~noprint error(['Discretionary policy: some eigenvalues greater than options_.qz_criterium. Model potentially unstable.']); case 63 error(['Discretionary policy: NaN elements are present in the solution. Procedure failed.']); + case 71 + error(['Calibrated covariance of the structural errors implies correlation larger than +-1.']); + case 72 + error(['Calibrated covariance of the measurement errors implies correlation larger than +-1.']); % Aim Code Conversions by convertAimCodeToInfo.m case 102 error('Aim: roots not correctly computed by real_schur'); diff --git a/matlab/set_all_parameters.m b/matlab/set_all_parameters.m index 656ac34edfcc6524bfe00ea775ada69b556adc86..9fb3ef4bd603419ad240e0dede1b95166990d3da 100644 --- a/matlab/set_all_parameters.m +++ b/matlab/set_all_parameters.m @@ -56,9 +56,11 @@ nvn = estim_params.nvn; ncn = estim_params.ncn; np = estim_params.np; Sigma_e = M.Sigma_e; +Correlation_matrix = M.Correlation_matrix; H = M.H; - -% setting shocks variance +Correlation_matrix_ME = M.Correlation_matrix_ME; +% setting shocks variance on the diagonal of Covariance matrix; used later +% for updating covariances if nvx var_exo = estim_params.var_exo; for i=1:nvx @@ -69,7 +71,8 @@ end % update offset offset = nvx; -% setting measument error variance +% setting measument error variance; on the diagonal of Covariance matrix; used later +% for updating covariances if nvn for i=1:nvn k = estim_params.nvn_observable_correspondence(i,1); @@ -81,38 +84,42 @@ end offset = nvx+nvn; % setting shocks covariances -if ~isempty(M.Correlation_matrix) - Sigma_e = diag(sqrt(diag(Sigma_e)))*M.Correlation_matrix*diag(sqrt(diag(Sigma_e))); % use of old correlation matrix is correct due to the diagonal structure and later only using the hence correctly updated diagonal entries of Sigma_e -end if ncx corrx = estim_params.corrx; for i=1:ncx k1 = corrx(i,1); k2 = corrx(i,2); - M.Correlation_matrix(k1,k2) = xparam1(i+offset); - M.Correlation_matrix(k2,k1) = M.Correlation_matrix(k1,k2); - Sigma_e(k1,k2) = xparam1(i+offset)*sqrt(Sigma_e(k1,k1)*Sigma_e(k2,k2)); - Sigma_e(k2,k1) = Sigma_e(k1,k2); + Correlation_matrix(k1,k2) = xparam1(i+offset); + Correlation_matrix(k2,k1) = Correlation_matrix(k1,k2); end end - +%build covariance matrix from correlation matrix and variances already on +%diagonal +Sigma_e = diag(sqrt(diag(Sigma_e)))*Correlation_matrix*diag(sqrt(diag(Sigma_e))); +%if calibrated covariances, set them now to their stored value +if isfield(estim_params,'calibrated_covariances') + Sigma_e(estim_params.calibrated_covariances.position)=estim_params.calibrated_covariances.cov_value; +end % update offset offset = nvx+nvn+ncx; + % setting measurement error covariances -if ~isempty(M.Correlation_matrix_ME) - H = diag(sqrt(diag(H)))*M.Correlation_matrix_ME*diag(sqrt(diag(H))); -end if ncn corrn_observable_correspondence = estim_params.corrn_observable_correspondence; for i=1:ncn k1 = corrn_observable_correspondence(i,1); k2 = corrn_observable_correspondence(i,2); - M.Correlation_matrix_ME(k1,k2) = xparam1(i+offset); - M.Correlation_matrix_ME(k2,k1) = M.Correlation_matrix_ME(k1,k2); - H(k1,k2) = xparam1(i+offset)*sqrt(H(k1,k1)*H(k2,k2)); - H(k2,k1) = H(k1,k2); + Correlation_matrix_ME(k1,k2) = xparam1(i+offset); + Correlation_matrix_ME(k2,k1) = Correlation_matrix_ME(k1,k2); end end +%build covariance matrix from correlation matrix and variances already on +%diagonal +H = diag(sqrt(diag(H)))*Correlation_matrix_ME*diag(sqrt(diag(H))); +%if calibrated covariances, set them now to their stored value +if isfield(estim_params,'calibrated_covariances_ME') + H(estim_params.calibrated_covariances_ME.position)=estim_params.calibrated_covariances_ME.cov_value; +end % update offset offset = nvx+ncx+nvn+ncn; @@ -125,7 +132,9 @@ end % updating matrices in M if nvx || ncx M.Sigma_e = Sigma_e; + M.Correlation_matrix=Correlation_matrix; end if nvn || ncn M.H = H; + M.Correlation_matrix_ME=Correlation_matrix_ME; end \ No newline at end of file diff --git a/matlab/set_parameters.m b/matlab/set_parameters.m index d1e81afd81e9a841723ecbc72223c3f6d2666ce5..58e167d5544166bc44f75874aedf4fd8c0dec695 100644 --- a/matlab/set_parameters.m +++ b/matlab/set_parameters.m @@ -14,7 +14,7 @@ function set_parameters(xparam1) % SPECIAL REQUIREMENTS % none -% Copyright (C) 2003-2009 Dynare Team +% Copyright (C) 2003-2013 Dynare Team % % This file is part of Dynare. % @@ -39,9 +39,11 @@ nvn = estim_params_.nvn; ncn = estim_params_.ncn; np = estim_params_.np; Sigma_e = M_.Sigma_e; +Correlation_matrix = M_.Correlation_matrix; offset = 0; -% stderrs of the exogenous shocks +% setting shocks variance on the diagonal of Covariance matrix; used later +% for updating covariances if nvx var_exo = estim_params_.var_exo; for i=1:nvx @@ -58,10 +60,17 @@ if ncx for i=1:ncx k1 = corrx(i,1); k2 = corrx(i,2); - Sigma_e(k1,k2) = xparam1(i+offset)*sqrt(Sigma_e(k1,k1)*Sigma_e(k2,k2)); - Sigma_e(k2,k1) = Sigma_e(k1,k2); + Correlation_matrix(k1,k2) = xparam1(i+offset); + Correlation_matrix(k2,k1) = Correlation_matrix(k1,k2); end end +%build covariance matrix from correlation matrix and variances already on +%diagonal +Sigma_e = diag(sqrt(diag(Sigma_e)))*Correlation_matrix*diag(sqrt(diag(Sigma_e))); +if isfield(estim_params,'calibrated_covariances') + Sigma_e(estim_params.calibrated_covariances.position)=estim_params.calibrated_covariances.cov_value; +end + % and update offset offset = offset + ncx + ncn; @@ -70,4 +79,5 @@ if np M_.params(estim_params_.param_vals(:,1)) = xparam1(offset+1:end); end -M_.Sigma_e = Sigma_e; \ No newline at end of file +M_.Sigma_e = Sigma_e; +M_.Correlation_matrix=Correlation_matrix; diff --git a/matlab/set_prior.m b/matlab/set_prior.m index f82223ec853418e7aafe2e95d3616be6f50fe42d..9f2c40285be2580791e376a9079365f4333703ff 100644 --- a/matlab/set_prior.m +++ b/matlab/set_prior.m @@ -75,9 +75,10 @@ if nvx end if nvn estim_params_.nvn_observable_correspondence=NaN(nvn,1); % stores number of corresponding observable - if isequal(M_.H,0) + if isequal(M_.H,0) %if no previously set measurement error, initialize H nvarobs = size(options_.varobs,1); M_.H = zeros(nvarobs,nvarobs); + M_.Correlation_matrix_ME = eye(nvarobs); end for i=1:nvn obsi_ = strmatch(deblank(M_.endo_names(estim_params_.var_endo(i,1),:)),deblank(options_.varobs),'exact'); @@ -116,6 +117,7 @@ if ncn if isequal(M_.H,0) nvarobs = size(options_.varobs,1); M_.H = zeros(nvarobs,nvarobs); + M_.Correlation_matrix_ME = eye(nvarobs); end xparam1 = [xparam1; estim_params_.corrn(:,3)]; ub = [ub; max(min(estim_params_.corrn(:,5),1),-1)]; diff --git a/tests/estimation/fs2000_calibrated_covariance.mod b/tests/estimation/fs2000_calibrated_covariance.mod new file mode 100644 index 0000000000000000000000000000000000000000..20d34c84ad003789a503dacafccffb98b81006b3 --- /dev/null +++ b/tests/estimation/fs2000_calibrated_covariance.mod @@ -0,0 +1,86 @@ +//two covariances are calibrated and one covariance of the ME. One of the calibrated covariances is superseded by estimation + +var m P c e W R k d n l gy_obs gp_obs y dA; +varexo e_a e_m e_f; + +parameters alp bet gam mst rho psi del; + +alp = 0.33; +bet = 0.99; +gam = 0.003; +mst = 1.011; +rho = 0.7; +psi = 0.787; +del = 0.02; + +model; +dA = exp(gam+e_a); +log(m) = (1-rho)*log(mst) + rho*log(m(-1))+e_m; +-P/(c(+1)*P(+1)*m)+bet*P(+1)*(alp*exp(-alp*(gam+log(e(+1))))*k^(alp-1)*n(+1)^(1-alp)+(1-del)*exp(-(gam+log(e(+1)))))/(c(+2)*P(+2)*m(+1))=0; +W = l/n; +-(psi/(1-psi))*(c*P/(1-n))+l/n = 0; +R = P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(-alp)/W; +1/(c*P)-bet*P*(1-alp)*exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)/(m*l*c(+1)*P(+1)) = 0; +c+k = exp(-alp*(gam+e_a))*k(-1)^alp*n^(1-alp)+(1-del)*exp(-(gam+e_a))*k(-1); +P*c = m; +m-1+d = l; +e = exp(e_a); +y = k(-1)^alp*n^(1-alp)*exp(-alp*(gam+e_a)); +gy_obs = dA*y/y(-1); +gp_obs = (P/P(-1))*m(-1)/dA; +end; + +initval; +k = 6; +m = mst; +P = 2.25; +c = 0.45; +e = 1; +W = 4; +R = 1.02; +d = 0.85; +n = 0.19; +l = 0.86; +y = 0.6; +gy_obs = exp(gam); +gp_obs = exp(-gam); +dA = exp(gam); +end; + +varobs gp_obs gy_obs; + +shocks; +var e_a; stderr 0.01; +var e_f; stderr 0.01; +var gy_obs; stderr 0.01; +var gy_obs, gp_obs= 0.00005; +var e_a, e_m = 0.00005; +var e_m, e_f = 0.00001; +end; + +steady; +check; + +estimated_params; +stderr e_m, 0.008862; +corr e_a, e_m, 0.5; +stderr gp_obs, 0.035449; +end; + +estimated_params_init; +stderr e_m, 0.5; +corr e_a, e_m, 0.5; +stderr gp_obs, 0.5; +end; + +estimation(order=1,datafile=fsdat_simul,nobs=192, loglinear, mh_replic=2002, mh_nblocks=1, mh_jscale=0.8); + +if ~isequal(M_.Sigma_e(3,2),1e-5) || ~isequal(M_.Sigma_e(2,3),1e-5) + error('Problem in setting calibrated covariance of structural shocks') +end +if isequal(M_.Sigma_e(2,1),5e-5) || isequal(M_.Sigma_e(1,2),5e-5) + error('Problem in overriding calibrated covariance of structural shocks by estimated correlation') +end +if ~isequal(M_.H(2,1),5e-5) || ~isequal(M_.H(1,2),5e-5) + error('Problem in setting calibrated covariance of measurement errors') +end \ No newline at end of file