From 6584e2421012012ebc7d639d5d2d0dc21748deaa Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Adjemian=20=28Ry=C3=BBk=29?=
<stepan@adjemian.eu>
Date: Thu, 27 Mar 2025 18:47:18 +0100
Subject: [PATCH] Bug fix. The value of the likelihood returned by NLKF
(particle filter) was wrong.
---
matlab/nonlinear-filters/nonlinear_kalman_filter.m | 7 ++++++-
meson.build | 3 +++
2 files changed, 9 insertions(+), 1 deletion(-)
diff --git a/matlab/nonlinear-filters/nonlinear_kalman_filter.m b/matlab/nonlinear-filters/nonlinear_kalman_filter.m
index 354100c2f0..551345c1ac 100644
--- a/matlab/nonlinear-filters/nonlinear_kalman_filter.m
+++ b/matlab/nonlinear-filters/nonlinear_kalman_filter.m
@@ -123,6 +123,8 @@ StateVectorVarianceSquareRoot = chol(ReducedForm.StateVectorVariance)';
lik = NaN(sample_size,1);
LIK = NaN;
+n = rows(Y);
+
for t=1:sample_size
xbar = [StateVectorMean ; zeros(number_of_structural_innovations,1) ] ;
sqr_Px = [StateVectorVarianceSquareRoot zeros(number_of_state_variables,number_of_structural_innovations);
@@ -186,7 +188,10 @@ for t=1:sample_size
return
end
end
- lik(t) = log( sum(probability2(Y(:,t),H_lower_triangular_cholesky,tmp(mf1,:)).*weights,1) ) ;
+ PredictedObservedVarianceSquareRoot = chol(PredictedObservedVariance);
+ PredictedObservedInverseVarianceSquareRoot = PredictedObservedVarianceSquareRoot\eye(n); % Inverse of the Cholesky -> inv(PredictedObservedVariance) = A*A'
+ StandardPredictionError = PredictedObservedInverseVarianceSquareRoot'*PredictionError;
+ lik(t) = -.5*n*log(2*pi) + sum(log(diag(PredictedObservedInverseVarianceSquareRoot))) - .5*dot(StandardPredictionError, StandardPredictionError);
end
LIK = -sum(lik(start:end));
diff --git a/meson.build b/meson.build
index 89abf1f985..1f162f7627 100644
--- a/meson.build
+++ b/meson.build
@@ -1770,6 +1770,9 @@ mod_and_m_tests = [
{ 'test' : [ 'particle/noisy_ar1.mod',
'particle/noisy_ar1_nlf_sis.mod'],
'extra' : [ 'particle/noisy_ar1_data.m'] },
+ { 'test' : [ 'particle/noisy_ar1.mod',
+ 'particle/noisy_ar1_nlf_nlkf.mod'],
+ 'extra' : [ 'particle/noisy_ar1_data.m'] },
{ 'test' : [ 'particle/dsge_base2.mod' ],
'extra' : [ 'particle/risky.m',
'particle/extreme.m',
--
GitLab