rlabbe / filterpy

Python Kalman filtering and optimal estimation library. Implements Kalman filter, particle filter, Extended Kalman filter, Unscented Kalman filter, g-h (alpha-beta), least squares, H Infinity, smoothers, and more. Has companion book 'Kalman and Bayesian Filters in Python'.
MIT License
3.29k stars 617 forks source link

Estimation of covariance matricies #73

Closed hlbkin closed 4 years ago

hlbkin commented 6 years ago

Hello!

Does this package has any estimation techniques for estimating Q and R matricies? Also in case of time-varying Q_t and R_t?

Thanks!

rlabbe commented 6 years ago

filterpy.common includes a few functions. Q_discrete_white_noise if you have the noise variance, Q_continuous_white_noise if you know the spectral density.

van_loan_discretization will compute Q for a continuous model, and there is the unfortunately uncommented/undocumented linear_ode_discretation.

None of those seem to be getting at what you want, but that is what is there.

baogorek commented 6 years ago

It looks like StatsModels has some functions for estimating general state space models through maximum likelihood (see the "Custom state space models" section - http://www.statsmodels.org/dev/statespace.html). I'm just learning about all this, but it seems like in practice, people are just guessing these matrices, or maybe tuning a very simple form with cross validation. @rlabbe, do you agree? Are these models even identifiable in their most general case?

rlabbe commented 6 years ago

I think in general people guess and/or tune, though my only evidence of that is the utter lack of any rigorous math or well described heuristics in the literature for this sort of thing. It's not so bad with kinematic models, as we can more or less estimate the maximum acceleration possible, and then populate Q with that. To be honest I am not strong enough in statistics to fully grasp the area and significance. For example: with digital sensors, most sensors are discrete. Maybe I have a scale that weighs in 10 gram increments. How do you model that? It's not Gaussian, that's for sure. I just throw a UKF or particle filter at the problem and call it solved. But I'm not launching spacecraft into multiple year missions.

baogorek commented 6 years ago

Just finished a Kaggle competition using the Kalman filter (code here). I used statsmodels, which does allow estimation using maximum likelihood. For the model with the transition matrix below, I was able to estimate variance parameters for the moving mean (first column) but the other state parameters relating to day-of-week effects were mostly nan. Still, it didn't seem to hurt the algorithm. statsmodels also gives the ability to simulate data with various configurations, and I was able to recover the parameter values using the technique, so the good news is that these models can be estimated with data. The bad news is that there might be very little information about some of the variance parameters, especially if there are more than a couple of states. Also, the numerical optimizers failed constantly, with the exception of Nelder Mead.

        seasonal_design = np.array([[1,  0,  0,  0,  0,  0,  0],
                                    [0, -1, -1, -1, -1, -1, -1],
                                    [0,  1,  0,  0,  0,  0,  0],
                                    [0,  0,  1,  0,  0,  0,  0],
                                    [0,  0,  0,  1,  0,  0,  0],
                                    [0,  0,  0,  0,  1,  0,  0],
                                    [0,  0,  0,  0,  0,  1,  0]])

I don't think you need to be a stats genius to approach estimation, just an objective function and a numerical optimizer. Since the kalman filter is giving you the gaussian mean and covariance matrix at each time step, I suppose it could be as straightforward as sequentially computing the density of each time point and keeping a running sum of those logged values. You'd try to get that sum to be as large as possible by modulating the variance parameters.

For non-gaussian data, I'd stick with calling it an "approximation"!

BTW, (for what it's worth) this script where I used statsmodel, filterpy, etc. to get the same filtered result. I do hope you add estimation into filterpy. I like the way you designed the library better than the others.

Prokhozhijj commented 5 years ago

Perhaps these docs can be useful: