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

MRG, ENH: Add DICS bias tests #8610

Merged
merged 13 commits into from Dec 15, 2020
Merged

MRG, ENH: Add DICS bias tests #8610

merged 13 commits into from Dec 15, 2020

Conversation

larsoner
Copy link
Member

@larsoner larsoner commented Dec 4, 2020

@wmvanvliet @britta-wstnr I am adding DICS bias tests along the lines of what we have for LCMV and minimum norm. After this PR is all good, we can add orientation tests to check the max-power-ori calculations.

However, it looks like noise_csd-based whitening is broken for DICS. I added a simple _cov_as_csd that just packs our data and noise cov into CSD objects, and the results should in principle match what we get for LCMV. I copied the code and limits from the LCMV tests. DICS performs within the same limits as LCMV as long as no whitening is used. Once things are whitened, things break. You can see that only the noise_cov=False versions of the bias tests are uncommented, except one -- the uncommented one went from the peak localization limit around 40 (meaning 40% of vertices were correctly localized exactly to the correct vertex) to being ~1%.

Two things I fixed so far:

  1. Bad channel handling -- if there are bads in info but CSDs have the bad channels, there is a picking problem.
  2. I added the same _restore_pos_semidef "trick" that we use for LCMV. Basically if EEG channels are used, the massive scale differences become a problem and the "trick" restores Cm to being positive semidefinite. (I'm not 100% sure it's correct but it seems to at least be okay for LCMV.) However, even without this hack workaround / fix, the localization limits have to be ~1%, so there is still something broken here.

@larsoner larsoner changed the title WIP: Add DICS bias tests MRG, ENH: Add DICS bias tests Dec 4, 2020
@larsoner larsoner added this to the 0.22 milestone Dec 4, 2020
@larsoner
Copy link
Member Author

larsoner commented Dec 4, 2020

Okay turns out the problem was that the whitener was not applied to the data. I still think there is something fishy with the ordering where real_filter happens before the whitener is applied, but the whitener can be complex (?). But we should fix that in another PR. @britta-wstnr @wmvanvliet can you review? We should get this in for 0.22

Copy link
Member

@britta-wstnr britta-wstnr left a comment

Choose a reason for hiding this comment

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

Mostly LGTM, @larsoner.
Just a question re the boundaries for the tests.

@@ -189,11 +193,14 @@ def make_dics(info, forward, csd, reg=0.05, noise_csd=None, label=None,
(freq, i + 1, n_freqs))

Cm = csd.get_data(index=i)

# XXX: Weird that real_filter happens *before* whitening, which could
# make things complex again...?
Copy link
Member

Choose a reason for hiding this comment

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

Agree, @wmvanvliet can you double check?

Copy link
Member

Choose a reason for hiding this comment

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

@larsoner what do we do with this? Keep as is and investigate what it would change in new PR?

mne/beamformer/tests/test_dics.py Show resolved Hide resolved
@larsoner
Copy link
Member Author

larsoner commented Dec 7, 2020

@britta-wstnr I added a comment, and while I was in there I added some orientation tests.

@larsoner
Copy link
Member Author

larsoner commented Dec 8, 2020

@britta-wstnr I also realized the problem with the no-reg solutions -- the bias test covariances being used were actually quite rank deficient because they were calculated using the product of the forward with its transpose for < 250 source locations but 366 channels. I added a tiny bit of reg directly to the covariance, then project with the projectors, which is probably much more realistic. Now those tests are much more stable (the ranges are only one or two percent instead of ~10)!

This also showed me that we could get rid of the positive-semidefinite "restoring" hack -- that was only needed for the bias tests because they were horribly rank deficient.

mne/beamformer/tests/test_lcmv.py Outdated Show resolved Hide resolved
ori = loc.copy() / np.linalg.norm(loc, axis=1, keepdims=True)
else:
# doesn't make sense for pooled (None) or max-power (can't be all 3)
ori = None
loc = np.linalg.norm(loc, axis=1) if pick_ori == 'vector' else np.abs(loc)
Copy link
Member

Choose a reason for hiding this comment

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

Wait, now I am confused. Didn't we set out to test whether the chosen orientation matches the simulated orientation? Does this just test the orientation the leadfield is giving?

Copy link
Member Author

Choose a reason for hiding this comment

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

The "simulated data" here is just the leadfield, and the (free ori) leadfield is oriented as X0, Y0, Z0, X1, Y1, Z1, ... , XN, YN, ZN for the N+1 sources, so the correct orientation of the simulated data is just np.eye(3) repeating over and over again

Copy link
Member

Choose a reason for hiding this comment

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

But we are not testing the actual selected orientation in the scalar beamformer cases, or are we?

Copy link
Member Author

Choose a reason for hiding this comment

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

The max_power_ori can't possibly work in this situation -- it selects a single orientation for each location, and we are simulating activity in each of the three cardinal directions. We probably should have a separate test for example that if we simulate data in the normal direction, the max_power_ori beamformer recovers the correct orientations. I'll see if it's easy to add, should be...

Copy link
Member Author

Choose a reason for hiding this comment

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

(can't possibly work to get the orientations right anyway)

Copy link
Member

Choose a reason for hiding this comment

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

Okay, I see. If we want to check for real in the orientation selection for the DICS beamformers, then we would indeed need a test with a single orientation.

@britta-wstnr
Copy link
Member

Ah, great that this is resolved and we don't need weirdly high regularization parameters anymore, @larsoner!

@larsoner
Copy link
Member Author

larsoner commented Dec 9, 2020

@britta-wstnr I say we fix the real_filter in another PR. It can be backported easily to 0.22 if we need to. Are you happy with this one otherwise?

@larsoner
Copy link
Member Author

@britta-wstnr okay with this one? Would be good to get this in since it has some useful fixes.

Copy link
Member

@britta-wstnr britta-wstnr left a comment

Choose a reason for hiding this comment

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

Just spelling errors, otherwise looks good to go!

@@ -189,11 +193,14 @@ def make_dics(info, forward, csd, reg=0.05, noise_csd=None, label=None,
(freq, i + 1, n_freqs))

Cm = csd.get_data(index=i)

# XXX: Weird that real_filter happens *before* whitening, which could
# make things complex again...?
Copy link
Member

Choose a reason for hiding this comment

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

@larsoner what do we do with this? Keep as is and investigate what it would change in new PR?

reg, weight_norm, use_cov, depth, lower, upper,
lower_ori, upper_ori):
"""Test orientation selection for bias for max-power LCMV."""
# we simulate data for the fixed orientation foward and beamform using
Copy link
Member

Choose a reason for hiding this comment

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

I have a neat GitHub Actions box saying foward ==> forward. Do you see those as well @larsoner or is that a new feature that is still being rolled out?

assert lower <= perc <= upper
# Comute the dot products of our forward normals and
Copy link
Member

Choose a reason for hiding this comment

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

Again, Github actions flagging Comute ==> Compute

@larsoner
Copy link
Member Author

Keep as is and investigate what it would change in new PR?

Yeah that's what I'm thinking. If it's wrong, it's wrong on master and this PR at least does not make it worse.

@larsoner larsoner merged commit 5ba77c8 into mne-tools:master Dec 15, 2020
2 checks passed
@larsoner larsoner deleted the dics-bias branch December 15, 2020 20:52
larsoner added a commit to agramfort/mne-python that referenced this pull request Dec 16, 2020
* upstream/master: (42 commits)
  MRG, ENH: Add DICS bias tests (mne-tools#8610)
  MRG, BUG, ENH: Add window option (mne-tools#8662)
  BUG: Fix alpha for volumes (mne-tools#8663)
  MRG, BUG: Fix bugs with envcorr (mne-tools#8658)
  MRG, ENH: Progressbar for csd_morlet (mne-tools#8608)
  Render is necessary now (mne-tools#8657)
  VIZ: Fix head size (mne-tools#8651)
  MRG, MAINT: bump sphinxcontrib-bitex version (mne-tools#8653)
  MRG, MAINT: Improve server env (mne-tools#8656)
  BUG: Mayavi center (mne-tools#8644)
  VIZ, ENH: allow show/hide annotations by label (mne-tools#8624)
  Add regression test for EEGLAB data with a chanlocs struct (mne-tools#8647)
  FIX: scalar_bar (mne-tools#8643)
  MRG: Small fix to tutorial; rename plot_events ordinate label to "Event id"; improve some SSP docstrings (mne-tools#8612)
  MRG, ENH: make plot alignment use defaults for colors (mne-tools#8553)
  BUG: Fix passing of channel type (mne-tools#8638)
  FIX: fixed loop over norm PSF/CTF options (mne-tools#8636)
  MRG, BUG: Pass kwargs (mne-tools#8630)
  DOC: Clearer error message (mne-tools#8631)
  BUG: Fix number of labels (mne-tools#8629)
  ...
larsoner added a commit to wmvanvliet/mne-python that referenced this pull request Dec 16, 2020
* upstream/master: (38 commits)
  MRG, ENH: Add DICS bias tests (mne-tools#8610)
  MRG, BUG, ENH: Add window option (mne-tools#8662)
  BUG: Fix alpha for volumes (mne-tools#8663)
  MRG, BUG: Fix bugs with envcorr (mne-tools#8658)
  MRG, ENH: Progressbar for csd_morlet (mne-tools#8608)
  Render is necessary now (mne-tools#8657)
  VIZ: Fix head size (mne-tools#8651)
  MRG, MAINT: bump sphinxcontrib-bitex version (mne-tools#8653)
  MRG, MAINT: Improve server env (mne-tools#8656)
  BUG: Mayavi center (mne-tools#8644)
  VIZ, ENH: allow show/hide annotations by label (mne-tools#8624)
  Add regression test for EEGLAB data with a chanlocs struct (mne-tools#8647)
  FIX: scalar_bar (mne-tools#8643)
  MRG: Small fix to tutorial; rename plot_events ordinate label to "Event id"; improve some SSP docstrings (mne-tools#8612)
  MRG, ENH: make plot alignment use defaults for colors (mne-tools#8553)
  BUG: Fix passing of channel type (mne-tools#8638)
  FIX: fixed loop over norm PSF/CTF options (mne-tools#8636)
  MRG, BUG: Pass kwargs (mne-tools#8630)
  DOC: Clearer error message (mne-tools#8631)
  BUG: Fix number of labels (mne-tools#8629)
  ...
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

Successfully merging this pull request may close these issues.

None yet

2 participants