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

New optimizer for stable linear models of arbitrary dimension #269

Merged
merged 11 commits into from
Dec 28, 2022

Conversation

akaptano
Copy link
Collaborator

@akaptano akaptano commented Dec 6, 2022

Hi all, this branch isn't fully polished yet, but a few users have been asking for linear models of relatively high dimension, so I implemented an optimizer that does this and promotes linear matrices that are negative definite (and therefore globally stable). The optimization trick is similar to the trapping optimizer. I also added functionality for the ConstrainedSR3 and new StableLinearSR3 optimizers to use both equality and inequality constraints. This and the new optimizer are now illustrated with simple examples in Notebook 1. There are three more linear model examples in examples/17_linear_stable_models/. There is a bit more work to do adding new unit tests and comments. Let me know what you think!

Alan Kaptanoglu added 8 commits November 11, 2022 15:18
…unctionality. Doesnt quite work yet I think for Zachs problem because only a subset of the variables are observable.
…blem. Having weird error from ConstrainedSR3, which is affecting StableLinearSR3. Going to merge with master first.
@akaptano akaptano added the enhancement New feature or request label Dec 6, 2022
Copy link
Collaborator

@znicolaou znicolaou left a comment

Choose a reason for hiding this comment

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

This looks mostly good to me--it is indeed useful to enforce a stable linear model. Just a couple comments:

  1. In the future, I think we should also include some functionality for a constrained STLSQ optimizer. The scipy lsq_linear routine (https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.lsq_linear.html) should make it pretty trivial for some inequality constraints, and the equality constraints can be enforce pretty easily with projections, e.g. https://en.wikipedia.org/wiki/Constrained_generalized_inverse. This would be useful for enforcing symmetries in PDEs, for example, and it would be nice to be able to compare the SR3 and the STLSQ cases.
  2. I think some data files data/fileWi3_5.npy, data/fileWi4.npy are missing for viscoelastic_pod_models.ipynb, so I couldn't run that notebook yet.
  3. You should probably change the example directory to 18_linear_stable_models to merge well with PR Parametric library #273

@akaptano
Copy link
Collaborator Author

@znicolaou Thanks for the notes -- hopefully will get this merge-ready within the week. Cassio also reports that the StableLinearSR3 optimizer gets pretty slow for large models (100 or 200-dimensional + ), and since this is exactly the kind of scenario we would like such stable models, I am planning to address this.

Fixing this is straightforward -- don't use CVXPY for the solve over the coefficients at each iteration. Not using CVXPY means that users cannot input inequality constraints, but actually the L1-constrained problem with equality constraints is exactly the same SR3 problem as in Champion et al "Unified..." so it has an analytic solution (which is also used in the TrappingSR3 solve for speed). I'll try to put this in soon and check the speedup, and this description is mostly so that I have it written somewhere.

@akaptano
Copy link
Collaborator Author

@znicolaou Regarding your first point, I'm actually not sure I understand how to make the Constrained STLSQ optimizer work (for the L0 norm, which is what it is good for). Generically, the greedy truncation of STLSQ will be done in such a way that breaks the constraints at each iteration, and I don't think I've seen anyone publish on this. Maybe you could do something at each iteration like truncate and then project that solution onto nearest solution that satisfies the constraints? Or make it so that your truncation choices always respect the constraints? My guess is the algorithm performance would take a hit in this scenario. Happy to discuss more in person next quarter.

…the viscoelastic examples, which turned out to be much too large. will add some smaller files later.
@codecov-commenter
Copy link

codecov-commenter commented Dec 21, 2022

Codecov Report

Base: 92.69% // Head: 92.39% // Decreases project coverage by -0.29% ⚠️

Coverage data is based on head (1400a55) compared to base (1e1b9bd).
Patch coverage: 84.90% of modified lines in pull request are covered.

Additional details and impacted files
@@            Coverage Diff             @@
##           master     #269      +/-   ##
==========================================
- Coverage   92.69%   92.39%   -0.30%     
==========================================
  Files          34       36       +2     
  Lines        3475     3682     +207     
==========================================
+ Hits         3221     3402     +181     
- Misses        254      280      +26     
Impacted Files Coverage Δ
pysindy/__init__.py 82.75% <50.00%> (ø)
pysindy/optimizers/__init__.py 68.00% <50.00%> (-3.43%) ⬇️
pysindy/optimizers/stable_linear_sr3.py 86.42% <86.42%> (ø)
pysindy/optimizers/constrained_sr3.py 91.30% <90.90%> (+1.56%) ⬆️
pysindy/optimizers/trapping_sr3.py 90.88% <0.00%> (+0.62%) ⬆️
pysindy/optimizers/sr3.py 94.92% <0.00%> (+0.72%) ⬆️

Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here.

☔ View full report at Codecov.
📢 Do you have feedback about the report comment? Let us know in this issue.

@znicolaou
Copy link
Collaborator

@akaptano agreed that generically there are challenges enforcing constraints during thresholding. For specific equality constraints that are linear in the coordinates, I think that adding a threshold would just amount to finding a lower-dimensional linear subspace to project onto. Those constraints would be useful for, e.g., enforcing rotationally-invariant PDEs with x_11-x_22=0, but I think the stability constraints involving eigenvalues would not be linear, and would be more complicated indeed...

@akaptano akaptano merged commit a299cd8 into master Dec 28, 2022
@akaptano akaptano deleted the cassio_branch branch December 28, 2022 02:23
jpcurbelo pushed a commit to jpcurbelo/pysindy_fork that referenced this pull request Apr 30, 2024
New optimizer for stable linear models of arbitrary dimension
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement New feature or request
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants