diff --git a/matlab/+occbin/kalman_update_algo_1.m b/matlab/+occbin/kalman_update_algo_1.m index 11205e63e0c80159a2af1994e51d4ca4a5a5acd6..83cbd587ba8f3ae9e440bd46ce7e1743d54a5f1e 100644 --- a/matlab/+occbin/kalman_update_algo_1.m +++ b/matlab/+occbin/kalman_update_algo_1.m @@ -213,7 +213,12 @@ if any(myregime) || ~isequal(regimes_(1),regimes0(1)) if M_.occbin.constraint_nbr==1 regime_end(niter) = regimes_(1).regimestart(end); end - [a, a1, P, P1, v, alphahat, etahat, lik, V] = occbin_kalman_update0(a,a1,P,P1,data_index,Z,v,Y,H,QQQ,TT,RR,CC,iF,L,mm, options_.rescale_prediction_error_covariance, options_.occbin.likelihood.IF_likelihood, options_.occbin.filter.state_covariance); + [a, a1, P, P1, v, alphahat, etahat, lik, V, error_flag] = occbin_kalman_update0(a,a1,P,P1,data_index,Z,v,Y,H,QQQ,TT,RR,CC,iF,L,mm, options_.rescale_prediction_error_covariance, options_.occbin.likelihood.IF_likelihood, options_.occbin.filter.state_covariance); + if error_flag + etahat=NaN(size(QQQ,1),1); + warning(orig_warning_state); + return; + end etahat_hist(niter) = {etahat}; lik_hist(niter) = lik; opts_simul.SHOCKS(1,:) = etahat(:,2)'; @@ -278,7 +283,12 @@ if any(myregime) || ~isequal(regimes_(1),regimes0(1)) TT(:,:,2)=ss.T(my_order_var,my_order_var,1); RR(:,:,2)=ss.R(my_order_var,:,1); CC(:,2)=ss.C(my_order_var,1); - [a, a1, P, P1, v, alphahat, etahat, lik, V] = occbin_kalman_update0(a,a1,P,P1,data_index,Z,v,Y,H,QQQ,TT,RR,CC,iF,L,mm, options_.rescale_prediction_error_covariance, options_.occbin.likelihood.IF_likelihood, options_.occbin.filter.state_covariance); + [a, a1, P, P1, v, alphahat, etahat, lik, V, error_flag] = occbin_kalman_update0(a,a1,P,P1,data_index,Z,v,Y,H,QQQ,TT,RR,CC,iF,L,mm, options_.rescale_prediction_error_covariance, options_.occbin.likelihood.IF_likelihood, options_.occbin.filter.state_covariance); + if error_flag + etahat=NaN(size(QQQ,1),1); + warning(orig_warning_state); + return; + end end else error_flag = 330; @@ -294,7 +304,7 @@ end error_flag = out.error_flag; if ~error_flag && niter>options_.occbin.likelihood.max_number_of_iterations && ~isequal(regimes_(1),regimes0(1)) - error_flag = 1; + error_flag = 331; end if ~error_flag @@ -357,7 +367,7 @@ else F = ZZ*P(:,:,t)*ZZ' + H(di,di); sig=sqrt(diag(F)); if any(any(isnan(F))) - error_flag=1; + error_flag=325; warning(orig_warning_state); return; end diff --git a/matlab/+occbin/kalman_update_algo_3.m b/matlab/+occbin/kalman_update_algo_3.m index 457c47480666f57757ceeca79d35c83fe6883b5e..efffc9badbb50593472ba8681a4f0d2133733efc 100644 --- a/matlab/+occbin/kalman_update_algo_3.m +++ b/matlab/+occbin/kalman_update_algo_3.m @@ -282,7 +282,7 @@ end error_flag = out.error_flag; if ~error_flag && niter>options_.occbin.likelihood.max_number_of_iterations && ~isequal(regimes_(1),regimes0(1)) %fixed point algorithm did not converge - error_flag = 1; + error_flag = 331; end if ~error_flag diff --git a/matlab/+occbin/kalman_update_engine.m b/matlab/+occbin/kalman_update_engine.m index 647be1a95f627d57b58eff2bd70c0517df537cfc..3459a26712cdb504bbc9852ee227ee886914a9fd 100644 --- a/matlab/+occbin/kalman_update_engine.m +++ b/matlab/+occbin/kalman_update_engine.m @@ -413,4 +413,13 @@ if length(likvec)>1 % sum the likelihood of multiple solutions likx = -2*log(sum(exp(-likvec./2))); end + +if info(1)==0 + if isnan(likx) + info = 323; + elseif any(any(isnan(ax))) + info = 324; + end +end + end diff --git a/matlab/get_error_message.m b/matlab/get_error_message.m index 9680a675737846095a2609abd8d036494dc67189..7156e3c253c4e7d2dc1b60b35d92491e4d016227 100644 --- a/matlab/get_error_message.m +++ b/matlab/get_error_message.m @@ -197,6 +197,16 @@ switch info(1) message = 'Occbin: there was a problem in running the smoother. Simulation within smoother failed.'; case 322 message = 'Occbin: smoother did not converge.'; + case 323 + message = 'Piecewise linear Kalman filter: the likelihood is NaN.'; + case 324 + message = 'Piecewise linear Kalman filter: updated state vector is NaN.'; + case 325 + message = 'Piecewise linear Kalman filter: filter covariance NaN.'; + case 330 + message = 'Piecewise linear Kalman filter: update step did not reach a fixed point (periodic loop).'; + case 331 + message = 'Piecewise linear Kalman filter: update step did not reach a fixed point (max number of iterations reached).'; case 401 message = 'Cycle reduction reached the iteration limit. Try increasing maxit.'; case 402 diff --git a/matlab/kalman/likelihood/missing_observations_kalman_filter.m b/matlab/kalman/likelihood/missing_observations_kalman_filter.m index adfdda765267051bf38a262bddb952fa89409f4c..67b3adf6143d50501acb09497181ca29816850f8 100644 --- a/matlab/kalman/likelihood/missing_observations_kalman_filter.m +++ b/matlab/kalman/likelihood/missing_observations_kalman_filter.m @@ -276,12 +276,6 @@ while notsteady && t<=last [ax, a1x, Px, P1x, vx, Tx, Rx, Cx, regx, info, M_, likx] = occbin.kalman_update_engine(a0(:,t-1),a1(:,t-1:t),P0(:,:,t-1),P1(:,:,t-1:t),t,data_index(t-1:t),Z,vv(:,t-1:t),Y(:,t-1:t),H,Qt,T0,R0,TT(:,:,t-1:t),RR(:,:,t-1:t),CC(:,t-1:t),regimes_(t:t+1),base_regime,d_index,M_,dr, endo_steady_state,exo_steady_state,exo_det_steady_state,options_,occbin_options); % [ax, a1x, Px, P1x, vx, Tx, Rx, Cx, regimes_(t:t+2), info, M_, likx] = occbin.kalman_update_algo_1(a0(:,t-1),a1(:,t-1:t),P0(:,:,t-1),P1(:,:,t-1:t),data_index(t-1:t),Z,vv(:,t-1:t),Y(:,t-1:t),H,Qt,T0,R0,TT(:,:,t-1:t),RR(:,:,t-1:t),CC(:,t-1:t),regimes_(t:t+1),M_,dr,endo_steady_state,exo_steady_state,exo_det_steady_state,options_,occbin_options); end - if isnan(likx) || any(any(isnan(ax))) - if options_.debug - fprintf('\nmissing_observations_kalman_filter:PKF failed in period %u with: NaN in likelihood value\n', t); - end - return - end if info if options_.debug fprintf('\nmissing_observations_kalman_filter:PKF failed in period %u with: %s\n', t, get_error_message(info,options_));