Skip to content

Conversation

@mrcfon
Copy link

@mrcfon mrcfon commented May 28, 2024

This PR adds the following classes:

  • IPLSKalmanSmoother, iterated posterior linearisation smoother
  • AugmentedGaussianStatePrediction, the state prediction class that can host cross-covariance
  • LinearTransitionModel, a transition model that consumes explicitly specified transition matrix and bias
  • AugmentedKalmanPredictor, predictor that uses the bias vector in LinearTransitionModel
  • AugmentedUnscentedKalmanPredictor that generates the cross-covariance and is able to suppress prediction noise

The Iterated Posterior Linearization Smoother (IPLS) described in [1] performs statistical linear regression (SLR) of the nonlinear transition and measurement models w.r.t. to the current posterior approximation.

A running example can be found here.

[1] Á. F. García-Fernández, L. Svensson and S. Särkkä, "Iterated Posterior Linearization Smoother," in IEEE Transactions on Automatic Control, vol. 62, no. 4, pp. 2056-2063, April 2017, doi: 10.1109/TAC.2016.2592681.

@mrcfon mrcfon requested a review from a team as a code owner May 28, 2024 18:11
@mrcfon mrcfon requested review from hpritchett-dstl and nperree-dstl and removed request for a team May 28, 2024 18:11
@jswright-dstl jswright-dstl self-requested a review June 3, 2024 15:25
Copy link
Contributor

@csherman-dstl csherman-dstl left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There are a number of places where changes have been made to parts of the codebase that are not parts of this PR. Otherwise can I confirm that all changes that have been made in #1016 are included within this PR?

Comment on lines 126 to 128
def matrix(self, **kwargs): return self.meas_matrix

def bias(self, **kwargs): return self.bias_value
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

put returns on different lines

Comment on lines 726 to 727
from stonesoup.types.array import Matrix
from stonesoup.types.state import StateVector
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Imports should be at the top of the file


return self.noise_covar

class GeneralLinearGaussian(MeasurementModel, LinearModel, GaussianModel):
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This class can inherit most methods from LinearGaussian

self.noise_covar = CovarianceMatrix(self.noise_covar)

@property
def ndim_state(self):
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Identical to inherited method


from ..linear import KnownTurnRateSandwich
from ..linear import ConstantVelocity
from ..linear import LinearTransitionModel
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Import is not used

MeasurementPrediction.register(CompositeState) # noqa: E305


class AugmentedGaussianState(GaussianState):
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This class belongs in stonesoup.types.state

Comment on lines 503 to 519
transition_model: TransitionModel = Property(doc="The transition model to be used.")
control_model: ControlModel = Property(
default=None,
doc="The control model to be used. Default `None` where the predictor "
"will create a zero-effect linear :class:`~.ControlModel`.")
alpha: float = Property(
default=0.5,
doc="Primary sigma point spread scaling parameter. Default is 0.5.")
beta: float = Property(
default=2,
doc="Used to incorporate prior knowledge of the distribution. If the "
"true distribution is Gaussian, the value of 2 is optimal. "
"Default is 2")
kappa: float = Property(
default=0,
doc="Secondary spread scaling parameter. Default is calculated as "
"3-Ns")
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

All properties are inherited from UnscentedKalmanPredictor

Comment on lines 321 to 334
transition_model: TransitionModel = Property(doc="The transition model to be used.")
measurement_model: MeasurementModel = Property(default=None, doc="The measurement model to be used.")
alpha: float = Property(
default=0.5,
doc="Primary sigma point spread scaling parameter. Default is 0.5.")
beta: float = Property(
default=2,
doc="Used to incorporate prior knowledge of the distribution. If the "
"true distribution is Gaussian, the value of 2 is optimal. "
"Default is 2")
kappa: float = Property(
default=0,
doc="Secondary spread scaling parameter. Default is calculated as "
"3-Ns")
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Multiple properties are inherited from UnscentedKalamanSmoother

Comment on lines 382 to 386
while True:
# we have no test of convergence, but limited the number of iterations
if len(smoothed_tracks) >= self.n_iterations:
# warnings.warn("IPLS reached pre-specified number of iterations.")
break
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Does it make more sense for this to just be a for loop?

@narykov
Copy link

narykov commented May 23, 2025

@csherman-dstl, I’ve checked, and this PR does include the changes from #1016. I’ve applied the corrections you suggested and pushed the updates.

model = LinearTransitionModel(
transition_matrix=F,
bias_value=a_vector,
noise_covar=Q
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
noise_covar=Q
covariance_matrix=Q

The noise_covar property has been renamed to covariance_matrix but there are number of occurrences in the tests where this hasn't been changed.

Comment on lines 53 to 72
),
( # Augmented Kalman
AugmentedKalmanPredictor,
LinearTransitionModel(
transition_matrix=np.array([[0.8, 0.2], [0.3, 0.7]]),
bias_value=np.zeros([2, 1]),
noise_covar=np.diag([0.10961003, 0.88557178])
),
( # Augmented Unscented Kalman
AugmentedUnscentedKalmanPredictor,
LinearTransitionModel(
transition_matrix=np.array([[0.8, 0.2], [0.3, 0.7]]),
bias_value=np.zeros([2, 1]),
noise_covar=np.diag([0.10961003, 0.88557178])
),
],
ids=["standard", "extended", "unscented", "cubature", "stochasticIntegration"]
)
ids=["standard", "extended", "unscented", "cubature", "stochasticIntegration", "augmented_kalman", "augmented_unscented_kalman"]
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There's a mismatch in the number of brackets used here and also there are no given values for prior_mean and prior_covar

continue

""" Compute SLR parameters (transition). """
from stonesoup.types.prediction import Prediction
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
from stonesoup.types.prediction import Prediction

Prediction is already imported in this file

)
f_matrix, a_vector, lambda_cov_matrix = slr_definition(previous_state, trans_fun, force_symmetry=True)

"Perform linear time update"
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

# should be used for comments (this occurs multiple times in this method)

Suggested change
"Perform linear time update"
# Perform linear time update

Comment on lines 18 to 19
from functools import partial
from ..types.state import State
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

State is unused and from functools import partial should be near the top of the file

break

# SLR is wrt to tne approximated posterior in post_state, not the original prior in hypothesis.prediction
meas_fun = partial(super().predict_measurement, measurement_model=measurement_model,
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
meas_fun = partial(super().predict_measurement, measurement_model=measurement_model,
meas_fun = partial(self.predict_measurement, measurement_model=measurement_model,

As far as I can see predict_measurement has not been overwritten and so super() will call the same thing as just using self.

Comment on lines 464 to 477
smoothed_tracks = []

for iteration in range(self.n_iterations):

if iteration == 0:
# initialising by performing sigma-point smoothing via the UKF smoother
smoothed_track = UnscentedKalmanSmoother(transition_model=self.transition_model,
alpha=self.alpha,
beta=self.beta,
kappa=self.kappa).smooth(track)
smoothed_tracks.append(smoothed_track)
continue
else:
smoothed_track = smoothed_tracks[-1]
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
smoothed_tracks = []
for iteration in range(self.n_iterations):
if iteration == 0:
# initialising by performing sigma-point smoothing via the UKF smoother
smoothed_track = UnscentedKalmanSmoother(transition_model=self.transition_model,
alpha=self.alpha,
beta=self.beta,
kappa=self.kappa).smooth(track)
smoothed_tracks.append(smoothed_track)
continue
else:
smoothed_track = smoothed_tracks[-1]
smoothed_track = UnscentedKalmanSmoother(transition_model=self.transition_model,
alpha=self.alpha,
beta=self.beta,
kappa=self.kappa).smooth(track)
for _ in range(self.n_iterations):

I'm not sure of the need for the smoothed_tracks list, as far as I can tell only the latest version is ever used.
Also the smoothed_track can be initiated outside the loop instead of using the first iteration to only do this

@narykov
Copy link

narykov commented Jun 8, 2025

@csherman-dstl, thank you for your feedback. I've implemented your comments, and pushed the updates. I see that there are a number of test failures (311, 312, 313) following the check, but they appear to be related to test_road_network_shortest_path, which does not overlap with this PR.

Comment on lines 750 to 758
Parameters
----------
prior : :class:`~.State`
Prior state, :math:`\mathbf{x}_{k-1}`
timestamp : :class:`datetime.datetime`
Time to transit to (:math:`k`)
**kwargs : various, optional
These are passed to :meth:`~.TransitionModel.covar`
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

should control_input and transition_noise be included in the docstring?



class AugmentedGaussianState(GaussianState):
""" This is a new GaussianState class that can also store information on cross-covariance
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
""" This is a new GaussianState class that can also store information on cross-covariance
""" This is a GaussianState class that can also store information on cross-covariance

The class won't always be new

assert timestamp == prediction.timestamp


def test_augmentedgaussianstateprediction():
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should this test check the value of cross_covar

"KLDivergence between current and prior posterior state estimate.")
max_iterations: int = Property(
default=5,
doc="Number of iterations before while loop is exited and a non-convergence warning is "
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this docstring still accurate? As far as I can see the warning is commented out.

@narykov
Copy link

narykov commented Jun 30, 2025

@csherman-dstl, thanks again for the comments, all implemented now. The documentation build is currently failing due to the CI environment:

RuntimeError:
Kaleido requires Google Chrome to be installed.

As noted above, the other test failures are not related to this pull request.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants