MRG: Protect against bad projections #2838

larsoner · 2016-02-03T22:13:33Z

One possible solution to the projection issue. @agramfort can you take a look and see if you think this is a good approach? You can see that for e.g. our LCMV tests we had a lot of violations due to subselections. I'll need to make a lot of other cleanups to get rid of all the warnings in tests, but I wanted to make sure that's the right approach before doing it.

Related to #2760.

wmvanvliet · 2016-02-04T07:15:11Z

what's with the stacklevel parameter? Is that a temporary thing?

agramfort · 2016-02-04T08:25:17Z

for the beamformer test can't you just apply the projections before subselection? and raise a warning if proj have not been applied?

I don't like this new stacklevel parameter either ;)

larsoner · 2016-02-04T13:30:18Z

It's in private functions so don't worry about it. It makes it so toy actually see where you as a user made a call that created a problem. We should have probably been doing it with all our warnings. If that's not clear I'll give an example. I don't think applying first (e.g. to raw) will fully fix it. I think it had to do with the covariance that is loaded also having projectors that get partially applied. But I'll try again. Do you tentatively agree warnings are the way to go?

larsoner · 2016-02-04T14:56:21Z

FYI the stacklevel changes outputs from things like this (for the io module):

/home/larsoner/custombuilds/mne-python/mne/io/fiff/raw.py:186: UserWarning: This file contains raw Internal Active Shielding data. It may be distorted. Elekta recommends it be run through MaxFilter to produce reliable results. Consider closing the file and running MaxFilter on the data.
  warnings.warn(msg)  # , stacklevel=7)

i.e., not informative about user's line caused the warning, to:

/home/larsoner/Desktop/untitled3.py:3: UserWarning: This file contains raw Internal Active Shielding data. It may be distorted. Elekta recommends it be run through MaxFilter to produce reliable results. Consider closing the file and running MaxFilter on the data.
  allow_maxshield=True)

numpy does this with their warnings about slicing and the like. I think we should start doing it everywhere. It makes debugging warnings a lot easier.

larsoner · 2016-02-04T15:07:36Z

@agramfort if I do raw.apply_proj() and remove all the warning catches that I put in for the one test test_lcmv, I end up with 7 lines in the test that elicit warnings about too few channels during projection.

I think it's an interaction with the covariance code, where even if you apply proj, after doing a selection with epochs, you need to subselect the channels from the covariance that actually need to be used in beamforming. During that process, projectors get applied with a limited subset of channels. I suspect it's a bug -- it seems like we'd want to apply the projectors to the covariance matrix and then subselect -- but it looks like someone went through a lot of effort to let it work the other way (subselect, then apply re-normalized partial projection vectors). It might be related to needing to know the degrees of freedom. If you apply N vectors then subselect channels I don't think you can guarantee how many DOF you've lost, whereas if you subselect then apply N vectors, you know you've dropped your DOF by N.

agramfort · 2016-02-05T17:30:40Z

I suspect a bug too... I may have the time to look next time. leave me a broken PR I can take over if you want

larsoner · 2016-02-05T17:55:35Z

Once the stacklevel stuff is merged, this will be rebased and more or less ready to go. You can find the bugs by looking at where I had to add warning-catchers.

larsoner · 2016-02-09T16:33:37Z

@agramfort now you can see where all the proj options are problematic in our tests by the added warnings.catch_warnings calls. I'm fine to merge this as is, or have you take over, up to you.

agramfort · 2016-02-13T22:32:39Z

mne/beamformer/_lcmv.py


    # Handle whitening + data covariance
-    whitener, _ = compute_whitener(noise_cov, info, picks, rank=rank)
+    whitener = compute_whitener(noise_cov, info, picks, rank=rank)[0]


I prefer the older way. I find it less error prone as it makes the code break when the number of output args changes.

larsoner · 2016-02-26T15:33:03Z

@agramfort feel free to push directly to this branch if you want

agramfort · 2016-02-28T13:55:47Z

mne/io/proj.py

+                    orig_n = p['data']['data'].shape[1]
+                    if len(vecsel) < 0.9 * orig_n:
+                        warn('Projection vector "%s" magnitude %0.2f, '
+                             'applying projector with %s/%s of the original '


I don't get this sentence. Words are missing or something.

"vector %s has magnitude %f (should be unity)" make more sense?

yes it does

agramfort · 2016-02-28T14:04:38Z

about the "lot of efforts" are you talking about the _prepare_beamformer_input function?

larsoner · 2016-02-28T14:10:37Z

Not quite. The problem is that I'm not sure how you properly deal with using projections, but then later only subselecting a set of channels like is done in that example. Fortunately we don't do it in many places in the codebase (maybe only in tests to primarily save time?) but I suspect what we're doing isn't quite right. If you eliminate, say, half the channels than your previously-orthogonal vectors very will might not be orthogonal anymore. So things like rank estimates can get screwed up. But I suppose in that case we'd be conservative by saying that 2 projectors always decreases rank by 2, even with channel removal.

larsoner · 2016-02-28T14:11:38Z

Also see from above:

I think it's an interaction with the covariance code, where even if you apply proj, after doing a selection with epochs, you need to subselect the channels from the covariance that actually need to be used in beamforming. During that process, projectors get applied with a limited subset of channels. I suspect it's a bug -- it seems like we'd want to apply the projectors to the covariance matrix and then subselect -- but it looks like someone went through a lot of effort to let it work the other way (subselect, then apply re-normalized partial projection vectors). It might be related to needing to know the degrees of freedom. If you apply N vectors then subselect channels I don't think you can guarantee how many DOF you've lost, whereas if you subselect then apply N vectors, you know you've dropped your DOF by N.

larsoner · 2016-02-28T14:14:12Z

Changed warning text and rebased

agramfort · 2016-02-28T14:22:58Z

does this help:

agramfort@9bf4111

?

larsoner · 2016-02-28T14:27:43Z

The way to check would be to remove all the with warnings.catch_warnings(record=True) that I inserted in the tests. I'll try it.

larsoner · 2016-02-28T14:28:36Z

mne/beamformer/tests/test_lcmv.py

-                            picks=picks, baseline=(None, 0),
-                            preload=epochs_preload,
-                            reject=dict(grad=4000e-13, mag=4e-12, eog=150e-6))
+        with warnings.catch_warnings(record=True):


this one still needs to stay

larsoner · 2016-02-28T14:30:14Z

I stopped after the first three checks, but it looks like no, it didn't help

larsoner · 2016-02-28T14:31:30Z

I think it's because a lot of the warnings come from the guts of the covariance code, not the beamforming code

agramfort · 2016-02-28T14:32:38Z

ok I'll try to look later... I have to stop.

larsoner · 2016-02-28T14:33:48Z

Cool, thanks for looking. Feel free to push to my branch if you find a solution, I probably won't work on it in the meantime.

larsoner · 2016-03-14T16:54:56Z

@agramfort comments addressed, see what you think

larsoner · 2016-03-14T16:58:26Z

mne/cov.py

-    vals = np.abs(vals)  # avoid negative values (numerical errors)
+    # avoid negative (numerical errors) or zero (semi-definite matrix) values
+    tol = vals.max() * vals.size * np.finfo(np.float64).eps
+    vals = np.where(vals > tol, vals, tol)


We were getting some values that actually equalled zero, which threw a warning. This prevents that from happening hopefully by setting a hopefully reasonable minimum, but please double-check my logic @agramfort

larsoner · 2016-03-18T13:56:31Z

ping @agramfort

agramfort · 2016-03-24T21:15:00Z

mne/io/proj.py

+                             '%s/%s of the original channels available may '
+                             'be dangerous, consider recomputing and adding '
+                             'projection vectors for channels that are '
+                             'eventually used'


please update the warning to recommend to use normalize_proj method

larsoner · 2016-03-25T14:00:04Z

@agramfort warning updated. Anything else?

coveralls · 2016-03-25T14:52:01Z

Coverage decreased (-0.04%) to 90.731% when pulling 0628608 on Eric89GXL:proj-protect into 923bdc3 on mne-tools:master.

agramfort · 2016-03-26T13:42:06Z

@Eric89GXL I pushed a commit here. I think we're almost there...

coveralls · 2016-03-26T14:06:52Z

Coverage increased (+0.08%) to 90.844% when pulling 29021d8 on Eric89GXL:proj-protect into 923bdc3 on mne-tools:master.

larsoner · 2016-03-27T03:28:06Z

@agramfort anything left to do that you can see?

agramfort · 2016-03-27T20:58:30Z

I think we're good. Hopefully we did not break anything...

please rebase and add a clear message on what's new.

thanks heaps @Eric89GXL !

larsoner · 2016-03-28T11:57:41Z

Done and merged by rebase, thanks @agramfort

agramfort · 2016-03-28T12:02:05Z

thanks !

larsoner added this to the 0.12 milestone Feb 3, 2016

larsoner force-pushed the proj-protect branch from e50a95d to 50852af Compare February 9, 2016 16:31

larsoner changed the title ~~WIP: Protect against bad projections~~ ENH: Protect against bad projections Feb 9, 2016

agramfort reviewed Feb 13, 2016
View reviewed changes

larsoner force-pushed the proj-protect branch from 50852af to 63cd409 Compare February 15, 2016 14:43

larsoner mentioned this pull request Feb 26, 2016

ENH : better error message when unit is not known #2955

Merged

agramfort reviewed Feb 28, 2016
View reviewed changes

larsoner force-pushed the proj-protect branch from 63cd409 to 408e82e Compare February 28, 2016 14:13

larsoner reviewed Feb 28, 2016
View reviewed changes

larsoner and others added 9 commits March 14, 2016 12:18

WIP: Protect against bad projections

d7c0696

FIX: Fix docstring

cbbd754

FIX: Fix warning level

c5515ee

FIX: Fix warnings in tests

2c1aaf4

FIX: Fix stacklevel

4833cfe

FIX: Revert minor change

aadb92f

FIX: Fix warning text

59bd22d

pick after projection

814d0ea

FIX: Address comments

e124bb3

larsoner force-pushed the proj-protect branch from b04a50d to be96e6f Compare March 14, 2016 16:54

FIX: Better handling

2a4b805

larsoner force-pushed the proj-protect branch from be96e6f to 2a4b805 Compare March 14, 2016 16:56

larsoner reviewed Mar 14, 2016
View reviewed changes

agramfort reviewed Mar 24, 2016
View reviewed changes

FIX: Fix message

0628608

remove more warnings in examples

29021d8

larsoner closed this Mar 28, 2016

larsoner deleted the proj-protect branch March 28, 2016 11:57

larsoner mentioned this pull request Feb 1, 2022

confusing warning with applying projs #10286

Closed

Uh oh!

MRG: Protect against bad projections #2838

MRG: Protect against bad projections #2838

Uh oh!

Conversation

larsoner commented Feb 3, 2016

Uh oh!

wmvanvliet commented Feb 4, 2016

Uh oh!

agramfort commented Feb 4, 2016

Uh oh!

larsoner commented Feb 4, 2016 via email

Uh oh!

larsoner commented Feb 4, 2016

Uh oh!

larsoner commented Feb 4, 2016

Uh oh!

agramfort commented Feb 5, 2016 via email

Uh oh!

larsoner commented Feb 5, 2016

Uh oh!

larsoner commented Feb 9, 2016

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

larsoner commented Feb 26, 2016

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

agramfort commented Feb 28, 2016

Uh oh!

larsoner commented Feb 28, 2016

Uh oh!

larsoner commented Feb 28, 2016

Uh oh!

larsoner commented Feb 28, 2016

Uh oh!

agramfort commented Feb 28, 2016

Uh oh!

larsoner commented Feb 28, 2016

Uh oh!

Choose a reason for hiding this comment

Uh oh!

larsoner commented Feb 28, 2016

Uh oh!

larsoner commented Feb 28, 2016

Uh oh!

agramfort commented Feb 28, 2016 via email

Uh oh!

larsoner commented Feb 28, 2016

Uh oh!

larsoner commented Mar 14, 2016

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

larsoner commented Mar 18, 2016

Uh oh!

Choose a reason for hiding this comment

Uh oh!

larsoner commented Mar 25, 2016

Uh oh!

coveralls commented Mar 25, 2016

Uh oh!

agramfort commented Mar 26, 2016

Uh oh!

coveralls commented Mar 26, 2016

Uh oh!

larsoner commented Mar 27, 2016

Uh oh!

agramfort commented Mar 27, 2016

Uh oh!

larsoner commented Mar 28, 2016