Skip to content

Kalman variable dt #151

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
pavelkomarov merged 4 commits into from
Aug 23, 2025
Merged

Kalman variable dt #151

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
5f360a9
work in progress
pavelkomarov Aug 23, 2025
9143fa7
think I sorted it out
pavelkomarov Aug 23, 2025
63b4f2a
added variable dt test
pavelkomarov Aug 23, 2025
d06a564
fine touches to docstrings
pavelkomarov Aug 23, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 2 additions & 1 deletion docs/source/kalman_smooth.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,4 +8,5 @@ kalman_smooth
.. autofunction:: pynumdiff.kalman_smooth.constant_velocity
.. autofunction:: pynumdiff.kalman_smooth.constant_acceleration
.. autofunction:: pynumdiff.kalman_smooth.constant_jerk
.. autofunction:: pynumdiff.kalman_smooth.known_dynamics
Copy link
Collaborator Author

@pavelkomarov pavelkomarov Aug 23, 2025 •
edited
Loading

Choose a reason for hiding this comment

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

I've decided to make my kalman filter function public, which then becomes kind of redundant with known_dynamics. I'm pretty sure no one was using known_dynamics, so I've removed it.

.. autofunction:: pynumdiff.kalman_smooth.kalman_filter
.. autofunction:: pynumdiff.kalman_smooth.rts_smooth
200 changes: 28 additions & 172 deletions examples/1_basic_tutorial.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

78 changes: 26 additions & 52 deletions examples/2a_optimizing_parameters_with_dxdt_known.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

78 changes: 26 additions & 52 deletions examples/2b_optimizing_parameters_with_dxdt_unknown.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

23 changes: 13 additions & 10 deletions examples/3_automatic_method_suggestion.ipynb
View file
Open in desktop

Large diffs are not rendered by default.

2 changes: 1 addition & 1 deletion pynumdiff/__init__.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,5 +6,5 @@
from .smooth_finite_difference import meandiff, mediandiff, gaussiandiff, friedrichsdiff, butterdiff
from .polynomial_fit import splinediff, polydiff, savgoldiff
from .total_variation_regularization import tvrdiff, velocity, acceleration, jerk, iterative_velocity, smooth_acceleration, jerk_sliding
from .kalman_smooth import rtsdiff, constant_velocity, constant_acceleration, constant_jerk, known_dynamics
from .kalman_smooth import kalman_filter, rts_smooth, rtsdiff, constant_velocity, constant_acceleration, constant_jerk
from .linear_model import spectraldiff, lineardiff, rbfdiff
2 changes: 1 addition & 1 deletion pynumdiff/kalman_smooth/__init__.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,3 +1,3 @@
"""This module implements Kalman filters
"""
from ._kalman_smooth import rtsdiff, constant_velocity, constant_acceleration, constant_jerk, known_dynamics
from ._kalman_smooth import kalman_filter, rts_smooth, rtsdiff, constant_velocity, constant_acceleration, constant_jerk
231 changes: 104 additions & 127 deletions pynumdiff/kalman_smooth/_kalman_smooth.py
View file
Open in desktop

Large diffs are not rendered by default.

8 changes: 4 additions & 4 deletions pynumdiff/linear_model/_linear_model.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -238,19 +238,19 @@ def rbfdiff(x, _t, sigma=1, lmbd=0.01):

:param np.array[float] x: data to differentiate
:param float or array[float] _t: This function supports variable step size. This parameter is either the constant
step size if given as a single float, or data locations if given as an array of same length as :code:`x`.
:math:`\\Delta t` if given as a single float, or data locations if given as an array of same length as :code:`x`.
:param float sigma: controls width of radial basis function
:param float lmbd: controls strength of bias toward data

:return: tuple[np.array, np.array] of\n
- **x_hat** -- estimated (smoothed) x
- **dxdt_hat** -- estimated derivative of x
"""
if isinstance(_t, (np.ndarray, list)): # support variable step size for this function
if np.isscalar(_t):
t = np.arange(len(x))*_t
else: # support variable step size for this function
if len(x) != len(_t): raise ValueError("If `_t` is given as array-like, must have same length as `x`.")
t = _t
else:
t = np.arange(len(x))*_t

# The below does the approximate equivalent of this code, but sparsely in O(N sigma), since the rbf falls off rapidly
# t_i, t_j = np.meshgrid(t,t)
Expand Down
8 changes: 4 additions & 4 deletions pynumdiff/polynomial_fit/_polynomial_fit.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ def splinediff(x, _t, params=None, options={}, degree=3, s=None, num_iterations=

:param np.array[float] x: data to differentiate
:param float or array[float] _t: This function supports variable step size. This parameter is either the constant
step size if given as a single float, or data locations if given as an array of same length as :code:`x`.
:math:`\\Delta t` if given as a single float, or data locations if given as an array of same length as :code:`x`.
:param list params: (**deprecated**, prefer :code:`degree`, :code:`cutoff_freq`, and :code:`num_iterations`)
:param dict options: (**deprecated**, prefer :code:`num_iterations`) an empty dictionary or {'iterate': (bool)}
:param int degree: polynomial degree of the spline. A kth degree spline can be differentiated k times.
Expand All @@ -30,11 +30,11 @@ def splinediff(x, _t, params=None, options={}, degree=3, s=None, num_iterations=
if 'iterate' in options and options['iterate']:
num_iterations = params[2]

if isinstance(_t, (np.ndarray, list)): # support variable step size for this function
if np.isscalar(_t):
t = np.arange(len(x))*_t
else: # support variable step size for this function
if len(x) != len(_t): raise ValueError("If `_t` is given as array-like, must have same length as `x`.")
t = _t
else:
t = np.arange(len(x))*_t

x_hat = x
for _ in range(num_iterations):
Expand Down
55 changes: 31 additions & 24 deletions pynumdiff/tests/test_diff_methods.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
from ..linear_model import lineardiff, spectraldiff, rbfdiff
from ..polynomial_fit import polydiff, savgoldiff, splinediff
from ..total_variation_regularization import velocity, acceleration, jerk, iterative_velocity, smooth_acceleration, jerk_sliding
from ..kalman_smooth import constant_velocity, constant_acceleration, constant_jerk
from ..kalman_smooth import rtsdiff, constant_velocity, constant_acceleration, constant_jerk
from ..smooth_finite_difference import mediandiff, meandiff, gaussiandiff, friedrichsdiff, butterdiff
# Function aliases for testing cases where parameters change the behavior in a big way, so error limits can be indexed in dict
def iterated_second_order(*args, **kwargs): return second_order(*args, **kwargs)
Expand Down Expand Up @@ -44,9 +44,10 @@ def spline_irreg_step(*args, **kwargs): return splinediff(*args, **kwargs)
(spline_irreg_step, {'degree':5, 's':2}),
(polydiff, {'degree':2, 'window_size':3}), (polydiff, [2, 3]),
(savgoldiff, {'degree':2, 'window_size':5, 'smoothing_win':5}), (savgoldiff, [2, 5, 5]),
(constant_velocity, {'r':1e-4, 'q':1e-2}), (constant_velocity, [1e-4, 1e-2]),
(constant_acceleration, {'r':1e-4, 'q':1e-1}), (constant_acceleration, [1e-4, 1e-1]),
(constant_jerk, {'r':1e-4, 'q':10}), (constant_jerk, [1e-4, 10]),
(constant_velocity, {'r':1e-2, 'q':1e3}), (constant_velocity, [1e-2, 1e3]),
(constant_acceleration, {'r':1e-3, 'q':1e4}), (constant_acceleration, [1e-3, 1e4]),
(constant_jerk, {'r':1e-4, 'q':1e5}), (constant_jerk, [1e-4, 1e5]),
(rtsdiff, {'order':2, 'qr_ratio':1e7, 'forwardbackward':True}),
(velocity, {'gamma':0.5}), (velocity, [0.5]),
(acceleration, {'gamma':1}), (acceleration, [1]),
(jerk, {'gamma':10}), (jerk, [10]),
Expand Down Expand Up @@ -184,24 +185,30 @@ def spline_irreg_step(*args, **kwargs): return splinediff(*args, **kwargs)
[(-1, -1), (0, 0), (0, -1), (1, 0)],
[(0, 0), (2, 2), (0, 0), (2, 2)],
[(1, 1), (3, 3), (1, 1), (3, 3)]],
constant_velocity: [[(-25, -25), (-25, -25), (0, -1), (0, 0)],
[(-3, -3), (-2, -2), (0, -1), (0, 0)],
[(-1, -2), (0, 0), (0, -1), (0, 0)],
[(-1, -1), (1, 0), (0, -1), (1, 0)],
[(1, 1), (2, 2), (1, 1), (2, 2)],
[(1, 1), (3, 3), (1, 1), (3, 3)]],
constant_acceleration: [[(-25, -25), (-25, -25), (0, -1), (0, 0)],
[(-3, -4), (-3, -3), (0, -1), (0, 0)],
[(-3, -3), (-2, -2), (0, -1), (0, 0)],
[(-1, -1), (0, 0), (0, -1), (0, 0)],
[(1, 1), (2, 2), (1, 1), (2, 2)],
[(1, 1), (3, 3), (1, 1), (3, 3)]],
constant_jerk: [[(-25, -25), (-25, -25), (0, -1), (0, 0)],
[(-4, -5), (-3, -4), (0, -1), (0, 0)],
[(-3, -4), (-2, -3), (0, -1), (0, 0)],
[(-1, -2), (0, 0), (0, -1), (0, 0)],
[(1, 0), (2, 1), (1, 0), (2, 1)],
[(1, 1), (3, 3), (1, 1), (3, 3)]],
constant_velocity: [[(-25, -25), (-25, -25), (0, -1), (1, 1)],
Copy link
Collaborator Author

@pavelkomarov pavelkomarov Aug 23, 2025

Choose a reason for hiding this comment

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

Mostly these improved.

[(-4, -5), (-3, -3), (0, -1), (1, 1)],
[(-3, -3), (0, 0), (0, -1), (1, 1)],
[(-3, -3), (1, 0), (0, -1), (1, 1)],
[(-1, -1), (2, 2), (0, 0), (2, 2)],
[(-1, -1), (3, 3), (0, 0), (3, 3)]],
constant_acceleration: [[(-25, -25), (-25, -25), (0, -1), (1, 1)],
[(-5, -5), (-4, -4), (0, -1), (1, 1)],
[(-4, -5), (-3, -3), (0, -1), (1, 1)],
[(-3, -3), (0, 0), (0, -1), (1, 1)],
[(-1, -1), (1, 1), (0, -1), (1, 1)],
[(0, 0), (3, 3), (0, 0), (3, 3)]],
constant_jerk: [[(-25, -25), (-25, -25), (0, -1), (1, 1)],
[(-6, -6), (-5, -5), (0, -1), (1, 1)],
[(-5, -5), (-4, -4), (0, -1), (1, 1)],
[(-3, -3), (-1, -1), (0, -1), (1, 1)],
[(-1, -1), (1, 1), (0, -1), (1, 1)],
[(0, 0), (3, 3), (0, 0), (3, 3)]],
rtsdiff: [[(-25, -25), (-25, -25), (0, -1), (1, 1)],
Copy link
Collaborator Author

@pavelkomarov pavelkomarov Aug 23, 2025

Choose a reason for hiding this comment

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

This one is tested in the irregular dt list, and it works great.

[(-5, -5), (-4, -4), (0, -1), (1, 1)],
[(-4, -4), (-3, -3), (0, -1), (1, 1)],
[(-2, -3), (0, 0), (0, -1), (1, 1)],
[(-1, -2), (1, 1), (0, -1), (1, 1)],
[(0, 0), (3, 3), (0, 0), (3, 3)]],
spectraldiff: [[(-9, -10), (-14, -15), (-1, -1), (0, 0)],
[(0, 0), (1, 1), (0, 0), (1, 1)],
[(1, 1), (1, 1), (1, 1), (1, 1)],
Expand Down Expand Up @@ -238,7 +245,7 @@ def test_diff_method(diff_method_and_params, test_func_and_deriv, request): # re
except: warn(f"Cannot import cvxpy, skipping {diff_method} test."); return

# sample the true function and true derivative, and make noisy samples
if diff_method in [spline_irreg_step, rbfdiff]: # list that can handle variable dt
if diff_method in [spline_irreg_step, rtsdiff, rbfdiff]: # list that can handle variable dt
x = f(t_irreg)
dxdt = df(t_irreg)
_t = t_irreg
Expand All @@ -259,7 +266,7 @@ def test_diff_method(diff_method_and_params, test_func_and_deriv, request): # re
# plotting code
if request.config.getoption("--plot") and not isinstance(params, list): # Get the plot flag from pytest configuration
fig, axes = request.config.plots[diff_method] # get the appropriate plot, set up by the store_plots fixture in conftest.py
t_ = t_irreg if diff_method in [spline_irreg_step, rbfdiff] else t
t_ = t_irreg if diff_method in [spline_irreg_step, rtsdiff, rbfdiff] else t
axes[i, 0].plot(t_, f(t_))
axes[i, 0].plot(t_, x, 'C0+')
axes[i, 0].plot(t_, x_hat, 'C2.', ms=4)
Expand Down