Increasing L2 reg in STLSQ

[royalasdaircoding]

Hi,

I see in STLSQ there is now a line in the code which multiplies $\alpha$ by 2 if a LinAlgWarning is raised. Presumably this is to improve the condition number.

My question is, does the code inform you when this happens? Are there hypothetically cases where I fit a model with $\alpha=0.1$ but in reality my model is fit with $\alpha=0.4$, for example?

Thanks
Alasdair

[Jacob-Stevens-Haas]

Looks like you're right, it doesn't. Yes, that $\alpha=.4$ example is feasible. IDK if this was in the original STLSQ paper, so I'm hesitant to change behavior without more research, do you think the user ought to know, and if so, what's the right way of informing the user?

[royalasdaircoding]

As far as I know L2 regularisation and hence the unbiasing step weren't in the original paper at all (if you mean the 2016 SINDy paper). I'm just basing this off the code listed in the supplementary material of the original paper, I can't remember where the original MATLAB code can be found now. I don't remember this being in older releases of pysindy either.

For the behaviour, I assume that the number of cases this happens is small but I was thinking for parameter sensitivity studies it may be relevant? If you grid search $\lambda$ and $\alpha$, at higher $\lambda$ the thresholding of the initial least-squares may be sufficient that $\alpha$ isn't silently increased which may not be the case at lower $\lambda$.

I would say it would be better not filter the LinAlgWarning as you would want to know if the problem is ill-conditioned and maybe raise a warning the first time that $\alpha$ has been increased. The docstring could include a line about the behaviour as well. I also understand that feature libraries are often inherently correlated so I can understand why the filter is there.

[Jacob-Stevens-Haas]

I'm not even sure this behavior occurs, since LinAlgWarning is explicitly ignored. I think it was supposed to add an "error" filter instead of "ignore". Linking the relevant lines here. So I think we need to do some change there.

One challenge is that it's hard to detect a warning unless it's turned into an error, but after it's turned into an error, it can't be re-emitted as a warning:

>>> with warnings.catch_warnings():
...   warnings.filterwarnings("error", category=UserWarning)
...   try:
...     warnings.warn("hi")
...   except UserWarning:
...     warnings.filters = warnings.filters[1:]
...     print(warnings.filters)
...     raise
...     
[('default', None, <class 'DeprecationWarning'>, '__main__', 0), ('ignore', None, <class 'DeprecationWarning'>, None, 0), ('ignore', None, <class 'PendingDeprecationWarning'>, None, 0), ('ignore', None, <class 'ImportWarning'>, None, 0), ('ignore', None, <class 'ResourceWarning'>, None, 0)]
Traceback (most recent call last):
  File "<python-input-20>", line 5, in <module>
    warnings.warn("hi")
    ~~~~~~~~~~~~~^^^^^^
UserWarning: hi

I think we need to look at some options for how to elegantly capture this information. From the git log, it looks like it was added to support the paper referenced in these notebooks in commit f518f5c

[royalasdaircoding]

Ah good catch. Would something like this work?

with warnings.catch_warnings() as w: 
    warnings.filterwarnings("always", category=UserWarning)

    warnings.warn("hi")
    
    if w:
        print(w)

[Jacob-Stevens-Haas]

That is a good way to catch the warning, but rather than print we'd need to use warnings.warn_explicit() to re-create the caught warning. Something like:

# module constant
MOD_LOG = logging.getLogger(__name___)

...

with warnings.catch_warnings(record=True, category=LinAlgWarning) as ws:
    warnings.filterwarnings("always", category=LinAlgWarning)
    coef = ridge_regression(x, y, self.alpha, **kw)
for w in ws:
    warnings.warn_explicit(w.message, w.category, w.filename, w.lineno)
if ws:
    self.alpha = 2 * self.alpha
    MOD_LOG.info(f"Increasing regularization weight to {self.alpha} due to regression matrix ill-conditioning")

For tests, I think we'd need a small, rank-deficient regression problem, an STLSQ with a small alpha, and then verify that (a) warning was emitted, (b) log was recorded, and (c) optimizer.alpha changed.