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

[GSOC] Add EpochsTFR support to spectral connectivity functions #232

Open
wants to merge 22 commits into
base: main
Choose a base branch
from

Conversation

tsbinns
Copy link
Collaborator

@tsbinns tsbinns commented Sep 6, 2024

WIP follow up of #225 to finalise GSoC work.

So far adds support for coefficients from Epochs.compute_tfr(method="morlet", output="complex") to spectral_connectivity_epochs(), equivalent to the cwt_morlet mode. Still to be done is adding support for EpochsTFR objects in spectral_connectivity_time().

In addition to the Morlet approach, spectral_connectivity_time() supports the time-freq. multitaper mode. #126 could also be addressed by adding support for this to spectral_connectivity_epochs(), in the form of EpochsTFR objects.

However, when trying this I discovered a bug that prevents the time-freq. multitaper mode being computed from epoched data (mne-tools/mne-python#12831), but that seems like an easy fix.

Will continue to work on this next week.

@tsbinns
Copy link
Collaborator Author

tsbinns commented Sep 10, 2024

Have now added support for EpochsTFR to spectral_connectivity_time (using a custom MNE branch with the bug in #12831 fixed).

@tsbinns
Copy link
Collaborator Author

tsbinns commented Sep 17, 2024

Resetting tests to main mne branch now that mne-tools/mne-python#12842 is merged

mne_connectivity/spectral/epochs.py Outdated Show resolved Hide resolved
mne_connectivity/spectral/epochs.py Outdated Show resolved Hide resolved
mne_connectivity/spectral/epochs.py Outdated Show resolved Hide resolved
mne_connectivity/spectral/epochs.py Outdated Show resolved Hide resolved
mne_connectivity/spectral/epochs.py Outdated Show resolved Hide resolved
mne_connectivity/spectral/time.py Outdated Show resolved Hide resolved
mne_connectivity/spectral/time.py Outdated Show resolved Hide resolved
mne_connectivity/spectral/time.py Outdated Show resolved Hide resolved
mne_connectivity/spectral/time.py Outdated Show resolved Hide resolved
mne_connectivity/spectral/time.py Outdated Show resolved Hide resolved
@tsbinns
Copy link
Collaborator Author

tsbinns commented Sep 19, 2024

Thanks for the detailed review @drammock! Until those wider API questions are addressed I'll double check the unit tests and try to find where the deviation in Fourier/Welch pipelines.

@tsbinns
Copy link
Collaborator Author

tsbinns commented Sep 30, 2024

Conda test keeps freezing on the micromamba setup step, seems similar to this: mamba-org/setup-micromamba#225
Is it worth trying to fix now?

@larsoner
Copy link
Member

Yeah feel free! Sounds like there might be a workaround via pinning in some of those issue comments

@tsbinns
Copy link
Collaborator Author

tsbinns commented Dec 9, 2024

Haven't forgotten about this, just need to sort mne-tools/mne-python#12910 before I can finalise the changes and tests here.

@tsbinns
Copy link
Collaborator Author

tsbinns commented Jan 13, 2025

mne-tools/mne-python#12910 is now merged, so I can fix the outstanding issues here.

Comment on lines +448 to +476
def _tfr_csd_from_mt(x_mt, y_mt, weights_x, weights_y):
"""Compute time-frequency CSD from tapered spectra.

Parameters
----------
x_mt : array, shape (..., n_tapers, n_freqs, n_times)
The tapered time-frequency spectra for signals x.
y_mt : array, shape (..., n_tapers, n_freqs, n_times)
The tapered time-frequency spectra for signals y.
weights_x : array, shape (n_tapers, n_freqs)
Weights to use for combining the tapered spectra of x_mt.
weights_y : array, shape (n_tapers, n_freqs)
Weights to use for combining the tapered spectra of y_mt.

Returns
-------
csd : array, shape (..., n_freqs, n_times)
The CSD between x and y.
"""
# expand weights dims to match x_mt and y_mt
weights_x = np.expand_dims(weights_x, axis=(*np.arange(x_mt.ndim - 3), -1))
weights_y = np.expand_dims(weights_y, axis=(*np.arange(y_mt.ndim - 3), -1))
# compute CSD
csd = np.sum(weights_x * x_mt * (weights_y * y_mt).conj(), axis=-3)
denom = np.sqrt((weights_x * weights_x.conj()).real.sum(axis=-3)) * np.sqrt(
(weights_y * weights_y.conj()).real.sum(axis=-3)
)
csd *= 2 / denom
return csd
Copy link
Collaborator Author

@tsbinns tsbinns Jan 16, 2025

Choose a reason for hiding this comment

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

There isn't a _csd_from_mt equivalent for TFR data in MNE. Would people be open to me adding this to the main MNE package in a small PR?
Otherwise it can just remain here.

if not is_tfr_con: # normal spectra (multitaper or Fourier)
this_psd = _psd_from_mt(x_t, weights)
else: # TFR spectra (multitaper)
this_psd = np.array([_tfr_from_mt(epo_x, weights) for epo_x in x_t])
Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Realising now that the way I implemented _tfr_from_mt in mne-tools/mne-python#12910 works fine there, but isn't flexible to an epochs dimension (in contrast _psd_from_mt is).

Would people be open to a small PR in the main package to fix this?

Comment on lines +88 to +94
freqs : array_like | None
Array of frequencies of interest for time-frequency decomposition. Only the
frequencies within the range specified by ``fmin`` and ``fmax`` are used. If
``data`` is an :term:`array-like` or :class:`~mne.Epochs` object, the
frequencies must be specified. If ``data`` is an
:class:`~mne.time_frequency.EpochsTFR` object, ``data.freqs`` is used and this
parameter is ignored.
Copy link
Collaborator Author

Choose a reason for hiding this comment

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

For TFR data, freqs can be read from the object, but in the current state freqs must be passed.
If freqs could be None by default, is this an API change that would require a deprecation cycle?
Of course a check would be made that freqs != None when data isn't EpochsTFR.

Comment on lines +405 to +408
if not isinstance(data, EpochsTFR) and freqs is None:
raise TypeError(
"`freqs` must be specified when `data` is not an EpochsTFR object"
)
Copy link
Collaborator Author

Choose a reason for hiding this comment

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

This is the check mentioned above.

@tsbinns
Copy link
Collaborator Author

tsbinns commented Jan 16, 2025

I've addressed the previous issues.

There should now be full support for EpochsSpectrum and EpochsTFR data in spectral_connectivity_epochs (including time-resolved multitaper coeffs) and for EpochsTFR in spectral_connectivity_time, along with unit tests covering the changes.

@tsbinns
Copy link
Collaborator Author

tsbinns commented Jan 16, 2025

Given that the multitaper TFR requires the unreleased MNE 1.10, is it worth holding off on merging (once it's checked and ready to go) until after MNE-Conn 0.8 is released?

On that note, since the last major part of the GSoC project in #223 was merged, there's a decent amount for a new release: https://mne.tools/mne-connectivity/dev/whats_new.html
What do people think?

Opened an issue here -> #276

@tsbinns tsbinns mentioned this pull request Jan 16, 2025
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.

3 participants