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

[BUG] Cross-spectral density missing frequencies #12633

Merged
merged 2 commits into from
May 29, 2024

Conversation

tsbinns
Copy link
Contributor

@tsbinns tsbinns commented May 29, 2024

Reference issue

Discussed in call with @larsoner & @wmvanvliet.

What does this implement/fix?

Currently, computing the CSD with the csd_multitaper, csd_fourier, csd_array_multitaper, or csd_array_fourier functions of the time_frequency module will return a CrossSpectralDensity object whose frequencies do not match the fmin and fmax parameters supplied to the functions.

E.g. calling the functions with fmin=5, fmax=10 could be expected to return a CSD object with frequencies [5, 6, 7, 8, 9, 10] (assuming a 1 Hz resolution). Instead, the following frequencies are returned [6, 7, 8, 9].

This occurs because of the logic used to generate the internal frequency mask:

orig_frequencies = rfftfreq(n_fft, 1.0 / sfreq)
freq_mask = (orig_frequencies > fmin) & (orig_frequencies < fmax)

As such, the fmin and fmax entries of 5 and 10 Hz may be included in the list of possible frequencies, but they are ignored when generating the frequency mask (e.g. as 5 is not greater than itself, and 10 is not less than itself).

A simple solution if to allow entries of orig_frequencies to also match fmin and fmax, e.g.:

freq_mask = (
(orig_frequencies > 0) & (orig_frequencies >= fmin) & (orig_frequencies <= fmax)
)

This change ensures that fmin=5, fmax=10 will return a CSD object with frequencies [5, 6, 7, 8, 9, 10]. Note that (orig_frequencies > 0) is important to prevent the CSD being returned for the zero frequency, which would not match the behaviour of other spectral methods (e.g. PSD computation).

@larsoner This does not include logic for cases that can be considered as floating point errors, but if I recall you said this was optional.

This requires a small change to the unit tests for catching failing cases where no frequency bins are captured due to fmin and fmax being too close, e.g.:
https://github.com/tsbinns/mne-python/blob/74e7348be7536b42bca2fa3ce48ecd4f0536a89b/mne/time_frequency/tests/test_csd.py#L397

Copy link
Member

@mscheltienne mscheltienne left a comment

Choose a reason for hiding this comment

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

Great explanation in this PR, thank you!

@larsoner larsoner enabled auto-merge (squash) May 29, 2024 19:24
@larsoner larsoner merged commit 9a222ba into mne-tools:main May 29, 2024
32 checks passed
@larsoner
Copy link
Member

Just tacking on a little test here...

@meeseeksdev Hello

Copy link

lumberbot-app bot commented May 29, 2024

Look at me, @larsoner, I'm Mr. Meeseeks!

@larsoner
Copy link
Member

@meeseeksdev backport to maint/1.7

meeseeksmachine pushed a commit to meeseeksmachine/mne-python that referenced this pull request May 29, 2024
larsoner pushed a commit that referenced this pull request May 29, 2024
…missing frequencies) (#12634)

Co-authored-by: Thomas S. Binns <[email protected]>
@larsoner larsoner mentioned this pull request May 29, 2024
@tsbinns tsbinns deleted the fix_csd_freq_mask branch May 30, 2024 17:24
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