Skip to content
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.

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

Estimation of covariance matricies #73

Closed
hlbkin opened this issue Nov 15, 2017 · 5 comments
Closed

Estimation of covariance matricies #73

hlbkin opened this issue Nov 15, 2017 · 5 comments

Comments

@hlbkin
Copy link

hlbkin commented Nov 15, 2017

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
Copy link
Owner

rlabbe commented Dec 11, 2017

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
Copy link

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
Copy link
Owner

rlabbe commented Feb 7, 2018

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
Copy link

baogorek commented Feb 8, 2018

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
Copy link
Contributor

Prokhozhijj commented Oct 26, 2018

Perhaps these docs can be useful:

@rlabbe rlabbe closed this as completed May 4, 2020
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

No branches or pull requests

4 participants