Better separation between array- and MNE object-accepting functions

[tsbinns]

Something that's been on my mind. Probably a 1.0 matter.

The problem

MNE-Python has a nice seperation of functions that operate on arrays, and those that operate on MNE objects, especially for spectral analyses. MNE-Connectivity doesn't. E.g., spectral_connectivity_epochs can accept array-like time series data, array-like source estimates, Epochs, and EpochsSpectrum/TFR objects.

It makes the code that bit messier and harder to maintain. Also for users, it leads to a very busy docstring (which I've been told from beginners to MNE is a bit tricky to understand):

Proposed solution

I would like there to be more consistency with MNE-Python, e.g., a spectral_connectivity_epochs function that accepts MNE objects, and a spectral_connectivity_epochs_array function that accepts arrays. I'd also propose this for non-spectral connectivivty analysis functions, but the spectral connectivity stuff is definitely the messiest right now.

Yes, there would be significant changes needed to existing code, but that's why I'd consider it a 1.0 matter.

Proposal to drop support for time series data from spectral connectivity functions

However, the question is what the accepted data should represent, especially for the spectral functions: do they expect spectral coefficients to be passed, or would time series data still be allowed?

Dropping support for time series data would allow us to clean up the code a lot. Would make maintenance much easier as it's really quite complicated to understand in it's current state. It also adds very little burden to users: if they have Epochs data, just pass in the result of Epochs.compute_psd/tfr(outputs="complex"); and likewise the results of the psd_array_... and tfr_array_... functions for array data.

However, there is the question of how source estimate data is handled if support for time series data is dropped. I've never worked with this before, so I don't know what the equivalent pipeline looks like using MNE-Python spectral functions.

Open to discussion on this, but I'd argue in the long-term it's a big improvement to code maintainability and understandability of the spectral connectivity functions (not accepting loads of different forms of data and the caveats of what other params need to be passed for each).

[larsoner]

That all sounds reasonable to me!

For source estimates, I think here is where we should figure out how to get from source space -> PSD -> connectivity rather than source space -> connectivity:

mne-connectivity/examples/mne_inverse_coherence_epochs.py

Lines 100 to 121 in 8820b52

	stcs = apply_inverse_epochs(
	epochs, inverse_operator, lambda2, method, pick_ori="normal", return_generator=True
	)
	###############################################################################
	# Compute the coherence between sources
	# -------------------------------------
	#
	# Now we are ready to compute the coherence in the alpha and beta band.
	# fmin and fmax specify the lower and upper freq. for each band, respectively.
	#
	# To speed things up, we use 2 parallel jobs and use mode='fourier', which
	# uses a FFT with a Hanning window to compute the spectra (instead of
	# a multitaper estimation, which has a lower variance but is slower).
	# By using faverage=True, we directly average the coherence in the alpha and
	# beta band, i.e., we will only get 2 frequency bins.
	fmin = (8.0, 13.0)
	fmax = (13.0, 30.0)
	sfreq = raw.info["sfreq"] # the sampling frequency
	coh = spectral_connectivity_epochs(

We would need a function that computed the spectral coefficients and returned a generator. Something like https://mne.tools/stable/generated/mne.time_frequency.psd_array_multitaper.html would work if it had a return_generator.

[ruuskas]

I agree that spectral_connectivity_epochs and spectral_connectivity_time should preferably take spectral coefficients instead of time series data. In my use cases I often compute the spectral coefficients using e.g. tfr_array_morlet anyway for other uses. If some help reorganizing the code and docs is needed, I could be able to find time. Also fully agree with Eric.

[tsbinns]

Thanks both for the input!

We would need a function that computed the spectral coefficients and returned a generator. Something like https://mne.tools/stable/generated/mne.time_frequency.psd_array_multitaper.html would work if it had a return_generator.

@larsoner I'll have a look at this. And if this MNE-Connectivity change to better separate the workflows came in a 1.0 release, it'd be fine to jump the min MNE-Python ver to whatever the much more recent one that supports returning generators, right?

If some help reorganizing the code and docs is needed, I could be able to find time.

Thanks @ruuskas! If this is a 1.0 change, I suppose it's still some time down the line. But any help on this would be very welcome.