mne-tools
ENH Enable common average reference projection for ECOG, DBS, and SEEG data
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.
Hello! 👋 Thanks for opening your first pull request here! ❤️ We will try to get back to you soon. 🚴🏽♂️
![welcome[bot] avatar](https://avatars.githubusercontent.com/in/4141?v=4 "Published by a pub.dev verified publisher"){kind=small}
@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.
@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.
================================== 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:
https://github.com/mne-tools/mne-python/blob/b74d00449e796f0af163c5cfb86386a7b8a98d6e/mne/io/reference.py#L315-L317
you'll need to figure out why that warning isn't triggered anymore when that test is run.
@moritz-gerster it was a bit more difficult than what I though. I pushed here. Let me know what you think.
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?
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
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.
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.
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?
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?
@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
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!
@agramfort this one will hopefully come back green or close to it. Ready for review/merge from my end
@moritz-gerster can you check if it works for you before we merge?
Yes, I'll check tomorrow!
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?
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: @.***>
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!
🎉 Congrats on merging your first pull request! 🥳 Looking forward to seeing more from you in the future! 💪
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()
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).