ENH Enable common average reference projection for ECOG, DBS, and SEEG data #10656

moritz-gerster · 2022-05-20T14:47:20Z

Summary

Using a common average projection currently only works for EEG data. It should also be available for intracranial data. This PR tries to solve it.

What is the problem?

The method set_eeg_reference works, when using a projection, only for EEG data. When working with intracranial data such as ECoG it breaks because it does not find any EEG channels:

Exception has occurred: ValueError
Cannot create EEG average reference projector (no EEG data found)

What is the fix/enhancement?

With this PR, it is possible to specify the type of channels used for the reference projection. It is straightforward to implement since set_eeg_reference has the argument "ch_type" implemented already.

Using a common average projection currently only works for EEG data. It should also be available for intracranial data. This commit tries to solve it.

welcome · 2022-05-20T14:47:22Z

Hello! 👋 Thanks for opening your first pull request here! ❤️ We will try to get back to you soon. 🚴🏽‍♂️

agramfort · 2022-05-21T13:56:41Z

mne/io/proj.py


 @verbose
-def make_eeg_average_ref_proj(info, activate=True, verbose=None):
+def make_eeg_average_ref_proj(info, activate=True, ch_dict=None,


I would have gone with ch_type=('eeg',) as API and
default parameter

would that do the job for you?

Yes, thank you for the quick review @agramfort! I changed it in this commit 😊

@agramfort

Suggestion by @agramfort.

agramfort · 2022-05-21T15:00:12Z

mne/io/proj.py

+
+    ch_dict = {**{type_: True for type_ in ch_type},
+               'meg': False, 'ref_meg': False}
+    eeg_sel = pick_types(info, **ch_dict, exclude='bads')


I think you can reuse the private function _picks_to_idx

Thanks @agramfort. I don't think I understand this private function though. How can I include it?

eeg_sel = _picks_to_idx(info, ch_type, none='data', exclude='bads', allow_empty=False, with_ref_meg=False)???

I would go for none='eeg'

drammock · 2022-05-21T15:00:29Z

@moritz-gerster there are some failing tests:

style test:

=================================== FAILURES ===================================
__________________________ test_docstring_parameters ___________________________
mne/tests/test_docstring_parameters.py:164: in test_docstring_parameters
    raise AssertionError(msg)
E   AssertionError: 
E   mne.io.proj.make_eeg_average_ref_proj : PR01 : Parameters {'ch_type'} not documented
E   mne.io.proj.make_eeg_average_ref_proj : PR02 : Unknown parameters {'ch_dict'}
E   2 errors

functionality tests:

================================== FAILURES ===================================
___________________________ test_set_eeg_reference ____________________________
mne\io\tests\test_reference.py:172: in test_set_eeg_reference
    reref, ref_data = set_eeg_reference(raw, copy=False, projection=True)
E   Failed: DID NOT WARN. No warnings of type (<class 'RuntimeWarning'>,) were emitted. The list of emitted warnings is: [].

Let us know if you need help understanding the problems or finding the right fix.

moritz-gerster · 2022-05-26T08:14:27Z

@drammock thank you for your help! 😊

@moritz-gerster there are some failing tests:

style test:

=================================== FAILURES ===================================
__________________________ test_docstring_parameters ___________________________
mne/tests/test_docstring_parameters.py:164: in test_docstring_parameters
    raise AssertionError(msg)
E   AssertionError: 
E   mne.io.proj.make_eeg_average_ref_proj : PR01 : Parameters {'ch_type'} not documented
E   mne.io.proj.make_eeg_average_ref_proj : PR02 : Unknown parameters {'ch_dict'}
E   2 errors

I did update the docstring now but it seems there are still some doc tests failing. I don't understand why.

functionality tests:

================================== FAILURES ===================================
___________________________ test_set_eeg_reference ____________________________
mne\io\tests\test_reference.py:172: in test_set_eeg_reference
    reref, ref_data = set_eeg_reference(raw, copy=False, projection=True)
E   Failed: DID NOT WARN. No warnings of type (<class 'RuntimeWarning'>,) were emitted. The list of emitted warnings is: [].

I don't know how to fix this...

Any help would be appreciated! In the meantime, I will keep trying to find out myself.

…eg_ave_ref_proj

agramfort · 2022-05-26T16:56:24Z

mne/io/proj.py

    activate : bool
        If True projections are activated.
+    ch_type : tuple
+        List of channel types to use for reference projection.


@moritz-gerster can you list the valid channel types. I would say this makes sense for EEG, iEEG and ECoG?

Yes, @agramfort but also for DBS! We use the DBS leads to measure local field potentials. So it should be possible for the mne types eeg, seeg, dbs, and ecog.

drammock · 2022-05-27T21:26:11Z

================================== FAILURES ===================================
___________________________ test_set_eeg_reference ____________________________
mne\io\tests\test_reference.py:172: in test_set_eeg_reference
    reref, ref_data = set_eeg_reference(raw, copy=False, projection=True)
E   Failed: DID NOT WARN. No warnings of type (<class 'RuntimeWarning'>,) were emitted. The list of emitted warnings is: [].

This means that on line 172 of the file mne/io/tests/test_reference.py there is code that explicitly expects a warning to be raised, but the warning wasn't raised. It's probably this:

mne-python/mne/io/reference.py

Lines 315 to 317 in b74d004

    
           if _has_eeg_average_ref_proj(inst.info['projs']): 
        
               warn('An average reference projection was already added. The data ' 
        
                    'has been left untouched.')

you'll need to figure out why that warning isn't triggered anymore when that test is run.

agramfort · 2022-05-28T13:41:39Z

@moritz-gerster it was a bit more difficult than what I though. I pushed here. Let me know what you think.

mne/io/reference.py

…mne-python into eeg_ave_ref_proj

agramfort

ok so we still need to sort out the referencing with different channel types. One way I see to test for now is to get a proj for each channel type and use add_proj method to add both. Would this do the job for you @moritz-gerster ?

agramfort · 2022-07-09T10:17:00Z

mne/__init__.py

@@ -1,3 +1,4 @@
+import bullshit


I don't know this package ;)

Sorry for not coming back earlier, I am currently busy but I think I can come back to this next week!

In the meantime, I took a programming workshop and the instructor helped me to get my local version of mne running in my own virtual environment. To test that it works we added an error by importing bullshit :D I'll remove that of course.

mne/io/proj.py

moritz-gerster · 2022-07-09T11:31:38Z

ok so we still need to sort out the referencing with different channel types. One way I see to test for now is to get a proj for each channel type and use add_proj method to add both. Would this do the job for you @moritz-gerster ?

Dear @agramfort, Yes that sounds good. I am a little lost now though: What exactly should I change in the code?

larsoner · 2022-07-13T21:00:31Z

I think the idea would be to allow ch_type to be a list | tuple and iterate over it, creating and adding a separate projector for each type in the list

moritz-gerster · 2022-07-14T07:24:10Z

I think the idea would be to allow ch_type to be a list | tuple and iterate over it, creating and adding a separate projector for each type in the list

@larsoner thanks for the hint. I am not sure though if this would solve my problem. I would like to apply an average reference projector to my combined LFP+ECOG electrodes. I don't want a separate reference projector for the ECOG and LFP electrodes respectively. I need a common average reference across both electrode types.

larsoner · 2022-07-14T11:32:31Z

Ahh okay.

The separate average reference case can actually be solved easily by the user looping over N types and passing them as arguments one by one to get N separate projectors. So how about if they pass a list of N types, it makes one projector that jointly average references all of those types.

In other words, you still add support for list of types here, but in that case it produces a single projector.

moritz-gerster · 2022-07-15T15:29:55Z

So how about if they pass a list of N types, it makes one projector that jointly average references all of those types.

@larsoner how do I average all these ch_type specific projections that I added to yield a single projection?

larsoner · 2022-07-15T19:03:58Z

how do I average all these ch_type specific projections that I added to yield a single projection?

You would not average the projectors together, you would create a single projector that includes all the channels (of all of the given types together). Make sense?

larsoner · 2022-08-23T08:08:04Z

@moritz-gerster do you plan to work on this soon-ish? If not, then someone else might be able to push this over the finish line

moritz-gerster · 2022-08-23T08:13:35Z

Dear @larsoner, I am not sure. It is on my list but I cannot start working on it in the next four weeks. Also, even if I would come back I am not sure if I am proficient enough to finish it. If someone else would start working on that I'd be very happy!

larsoner · 2022-09-15T16:50:10Z

@agramfort this one will hopefully come back green or close to it. Ready for review/merge from my end

agramfort

just a nitpick

@moritz-gerster can you check if it works for you before we merge?

mne/io/reference.py

Co-authored-by: Alexandre Gramfort <alexandre.gramfort@m4x.org>

moritz-gerster · 2022-09-15T20:22:15Z

@moritz-gerster can you check if it works for you before we merge?

Yes, I'll check tomorrow!

moritz-gerster · 2022-09-16T08:19:20Z

Dear @larsoner and @agramfort,

I looked into the code and the projection now works! Thank you very much!!

There was some behavior of the code I found surprising though and it took me quite some time to understand it. I am not sure if this is a bug or on purpose. I guess this is not limited to intracranial data but to reference projections in general:

To test whether I can easily switch between the CAR projection and the monopolar reference, I plotted ECoG spectra side by side. I used

raw_CAR, _ = mne.set_eeg_reference(raw, ref_channels='average', copy=True, projection=True, ch_type='ecog')
psd_CAR, freqs = mne.time_frequency.psd_welch(raw_CAR, proj=True)
psd, freqs = mne.time_frequency.psd_welch(raw_CAR, proj=False)

However, the spectra looked identical. Only if I used

raw_CAR, _ = mne.set_eeg_reference(raw, ref_channels='average', copy=True, projection=True, ch_type='ecog')
raw_CAR.apply_proj()
psd_CAR, freqs = mne.time_frequency.psd_welch(raw_CAR)
psd, freqs = mne.time_frequency.psd_welch(raw)

the expected behavior occurred. The documentation of psd_welch (I know its deprecated) states

proj [bool]
Apply SSP projection vectors. If inst is ndarray this is not used.

The same is written in the new method compute_psd.

Does that mean only SSP projections should be used for plotting spectra but not reference projections?

agramfort · 2022-09-16T09:27:03Z

since you used projection=True you need raw_CAR.apply_proj() as otherwise the projection is not immediately applied. If you use projection=False then it's immediate. It's independent from the plotting code that takes the data as they are. ok?

…

Message ID: ***@***.***>

moritz-gerster · 2022-09-16T12:36:25Z

ok @agramfort.

For me this was confusing because I thought the great benefit of the projection is that I can easily turn it on- and off during plotting etc. while just having a single raw instance in my code. But if this is intended the code is good to go!

welcome · 2022-09-16T13:22:40Z

🎉 Congrats on merging your first pull request! 🥳 Looking forward to seeing more from you in the future! 💪

agramfort · 2022-09-16T13:22:54Z

Thx @moritz-gerster @larsoner 🙏🙏

drammock · 2022-09-16T13:23:18Z

For me this was confusing because I thought the great benefit of the projection is that I can easily turn it on- and off during plotting etc.

I agree that what you suggest should have worked, although since psd_welch is deprecated and we're about to release 1.2 we may not fix it. Can you try this instead? If this is broken we will definitely fix.

raw_CAR, _ = mne.set_eeg_reference(raw, ref_channels='average', copy=True,
                                   projection=True, ch_type='ecog')
for proj in (False, True):
    raw_CAR.compute_psd(proj=proj).plot()

moritz-gerster · 2022-09-16T19:14:27Z

Dear @drammock,

I tested it for compute_psd(). It shows the exact same behavior as psd_welch(): The projection is not applied during plotting, even if the projection kwarg is set to True raw_CAR.compute_psd(proj=True).

moritz-gerster added 3 commits May 20, 2022 15:52

Enable common average reference projection for ECOG, DBS, and SEEG

9d0ca3c

Using a common average projection currently only works for EEG data. It should also be available for intracranial data. This commit tries to solve it.

Update documentation of make_eeg_average_ref_proj

7232c10

Update changelog and contributors

fa9c864

Correct PR number

bf62bb2

agramfort reviewed May 21, 2022

View reviewed changes

Change default argument to ch_type=('eeg',)

7e369d8

Suggestion by @agramfort.

agramfort reviewed May 21, 2022

View reviewed changes

Update docstring

0aa47e6

agramfort added 2 commits May 26, 2022 18:45

Merge branch 'main' of https://github.com/mne-tools/mne-python into e…

ae834fd

…eg_ave_ref_proj

comist in docstring

8720de1

agramfort reviewed May 26, 2022

View reviewed changes

agramfort added 2 commits May 28, 2022 14:50

fix tests?

c6c560d

adding test + update API

df38254

moritz-gerster commented Jun 1, 2022

View reviewed changes

mne/io/reference.py Outdated Show resolved Hide resolved

moritz-gerster added 2 commits July 6, 2022 15:25

commit bullshit

21cb770

Merge branch 'eeg_ave_ref_proj' of https://github.com/moritz-gerster/…

0c3d885

…mne-python into eeg_ave_ref_proj

agramfort reviewed Jul 9, 2022

View reviewed changes

remove implemented error

f91559c

Loop over ch_types and add reference projections

edb85c4

larsoner added this to the 1.2 milestone Sep 13, 2022

larsoner added 5 commits September 14, 2022 06:46

Merge remote-tracking branch 'upstream/main' into eeg_ave_ref_proj

f497851

WIP: All at once

e127872

FIX: Add tests

2345ba4

FIX: ZMQ

ba9b977

FIX: Doc

9652e2b

larsoner added 4 commits September 15, 2022 10:14

FIX: Joint

6000b1e

FIX: Doc

82da781

FIX: More

b995d01

FIX: Flake

f81c5cd

larsoner self-assigned this Sep 15, 2022

FIX: Test

33adefa

agramfort reviewed Sep 15, 2022

View reviewed changes

mne/io/reference.py Outdated Show resolved Hide resolved

Update mne/io/reference.py

8d0a71e

Co-authored-by: Alexandre Gramfort <alexandre.gramfort@m4x.org>

agramfort merged commit fefa7db into mne-tools:main Sep 16, 2022

drammock mentioned this pull request Sep 16, 2022

compute_psd does not apply projs correctly #11177

Closed

Uh oh!

ENH Enable common average reference projection for ECOG, DBS, and SEEG data #10656

ENH Enable common average reference projection for ECOG, DBS, and SEEG data #10656

Uh oh!

Conversation

moritz-gerster commented May 20, 2022 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Summary

What is the problem?

What is the fix/enhancement?

Uh oh!

welcome bot commented May 20, 2022

Uh oh!

agramfort May 21, 2022

Choose a reason for hiding this comment

Uh oh!

moritz-gerster May 21, 2022

Choose a reason for hiding this comment

Uh oh!

agramfort May 21, 2022

Choose a reason for hiding this comment

Uh oh!

moritz-gerster May 21, 2022

Choose a reason for hiding this comment

Uh oh!

wmvanvliet Aug 23, 2022

Choose a reason for hiding this comment

Uh oh!

drammock commented May 21, 2022

Uh oh!

moritz-gerster commented May 26, 2022

Uh oh!

agramfort May 26, 2022

Choose a reason for hiding this comment

Uh oh!

moritz-gerster May 26, 2022 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

drammock commented May 27, 2022

Uh oh!

agramfort commented May 28, 2022

Uh oh!

Uh oh!

agramfort left a comment

Choose a reason for hiding this comment

Uh oh!

agramfort Jul 9, 2022

Choose a reason for hiding this comment

Uh oh!

moritz-gerster Jul 9, 2022 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Uh oh!

moritz-gerster commented Jul 9, 2022

Uh oh!

larsoner commented Jul 13, 2022

Uh oh!

moritz-gerster commented Jul 14, 2022

Uh oh!

larsoner commented Jul 14, 2022

Uh oh!

moritz-gerster commented Jul 15, 2022

Uh oh!

larsoner commented Jul 15, 2022

Uh oh!

larsoner commented Aug 23, 2022

Uh oh!

moritz-gerster commented Aug 23, 2022

Uh oh!

larsoner commented Sep 15, 2022

Uh oh!

agramfort left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

moritz-gerster commented Sep 15, 2022

Uh oh!

moritz-gerster commented Sep 16, 2022 • edited by drammock Loading Uh oh! There was an error while loading. Please reload this page.

moritz-gerster commented May 20, 2022 •

edited

Loading

moritz-gerster May 26, 2022 •

edited

Loading

moritz-gerster Jul 9, 2022 •

edited

Loading

moritz-gerster commented Sep 16, 2022 •

edited by drammock

Loading

moritz-gerster commented Sep 16, 2022 •

edited

Loading