Possible issue with loglikelihood computation in kalman/kalman_filter.py update function

[sortega87]

The computation of the log-likelihood, on the kalman update function, is currently made after the state vector is updated with the last observation. That is, the log_likelihood is computed in line 1434 using the updated state of line 1420.

I believe the computation should be made with the prior state estimate; moving the log_likelihood computation just after line 1408 should solve the issue.

filterpy/filterpy/kalman/kalman_filter.py

Lines 1329 to 1436 in a437893

	def update(x, P, z, R, H=None, return_all=False):
	"""
	Add a new measurement (z) to the Kalman filter. If z is None, nothing
	is changed.
	This can handle either the multidimensional or unidimensional case. If
	all parameters are floats instead of arrays the filter will still work,
	and return floats for x, P as the result.
	update(1, 2, 1, 1, 1) # univariate
	update(x, P, 1
	Parameters
	----------
	x : numpy.array(dim_x, 1), or float
	State estimate vector
	P : numpy.array(dim_x, dim_x), or float
	Covariance matrix
	z : (dim_z, 1): array_like
	measurement for this update. z can be a scalar if dim_z is 1,
	otherwise it must be convertible to a column vector.
	R : numpy.array(dim_z, dim_z), or float
	Measurement noise matrix
	H : numpy.array(dim_x, dim_x), or float, optional
	Measurement function. If not provided, a value of 1 is assumed.
	return_all : bool, default False
	If true, y, K, S, and log_likelihood are returned, otherwise
	only x and P are returned.
	Returns
	-------
	x : numpy.array
	Posterior state estimate vector
	P : numpy.array
	Posterior covariance matrix
	y : numpy.array or scalar
	Residua. Difference between measurement and state in measurement space
	K : numpy.array
	Kalman gain
	S : numpy.array
	System uncertainty in measurement space
	log_likelihood : float
	log likelihood of the measurement
	"""
	#pylint: disable=bare-except
	if z is None:
	if return_all:
	return x, P, None, None, None, None
	return x, P
	if H is None:
	H = np.array([1])
	if np.isscalar(H):
	H = np.array([H])
	Hx = np.atleast_1d(dot(H, x))
	z = reshape_z(z, Hx.shape[0], x.ndim)
	# error (residual) between measurement and prediction
	y = z - Hx
	# project system uncertainty into measurement space
	S = dot(dot(H, P), H.T) + R
	# map system uncertainty into kalman gain
	try:
	K = dot(dot(P, H.T), linalg.inv(S))
	except:
	# can't invert a 1D array, annoyingly
	K = dot(dot(P, H.T), 1./S)
	# predict new x with residual scaled by the kalman gain
	x = x + dot(K, y)
	# P = (I-KH)P(I-KH)' + KRK'
	KH = dot(K, H)
	try:
	I_KH = np.eye(KH.shape[0]) - KH
	except:
	I_KH = np.array([1 - KH])
	P = dot(dot(I_KH, P), I_KH.T) + dot(dot(K, R), K.T)
	if return_all:
	# compute log likelihood
	log_likelihood = logpdf(z, dot(H, x), S)
	return x, P, y, K, S, log_likelihood
	return x, P

[Eheran1]

Have you compared how this affects the filter?

[sortega87]

@Eheran1 I had trouble with the current implementation while using the log likelihood to fit the parameters of a Dynamic Linear Model (DLM, as in here). I had no luck in fitting the model parameters with the current implementation. However, I was able to fit the model, and the resulting parameters were consistent with my data, with the suggested changes.

The log likelihood is computed using the posterior state estimate in the current implementation, I believe the log likelihood should be computed using the prior state estimate.