[MRG] Add optional different ways of averaging EpochsTFR

[adam2392]

Reference issue

Closes: #8878

What does this implement/fix?

Essentially same workflow as Epochs.average().

[adam2392]

Not sure why CI failed. Locally,

(mne) adam2392@Adams-MacBook-Pro mne-python % pytest ./mne/time_frequency/tests 
====================================================== test session starts =======================================================
platform darwin -- Python 3.8.5, pytest-6.2.1, py-1.9.0, pluggy-0.13.1
rootdir: /Users/adam2392/Documents/mne-python, configfile: setup.cfg
plugins: xdist-2.1.0, cov-2.10.1, profiling-1.7.0, timeout-1.4.2, forked-1.3.0
collected 87 items                                                                                                               

mne/time_frequency/tests/test_ar.py s.                                                                                     [  2%]
mne/time_frequency/tests/test_csd.py ..............                                                                        [ 18%]
mne/time_frequency/tests/test_multitaper.py ss                                                                             [ 20%]
mne/time_frequency/tests/test_psd.py ......                                                                                [ 27%]
mne/time_frequency/tests/test_stft.py ........................................                                             [ 73%]
mne/time_frequency/tests/test_stockwell.py .....                                                                           [ 79%]
mne/time_frequency/tests/test_tfr.py ..................                                                                    [100%]

Module-level timings require pytest-harvest


--------------------------- generated xml file: /Users/adam2392/Documents/mne-python/junit-results.xml ---------------------------
====================================================== slowest 20 durations ======================================================
4.02s call     mne/time_frequency/tests/test_tfr.py::test_getitem_epochsTFR
2.61s call     mne/time_frequency/tests/test_tfr.py::test_plot_joint
1.15s call     mne/time_frequency/tests/test_tfr.py::test_plot
0.81s call     mne/time_frequency/tests/test_stockwell.py::test_stockwell_ctf
0.61s call     mne/time_frequency/tests/test_psd.py::test_psd
0.61s call     mne/time_frequency/tests/test_tfr.py::test_io
0.61s call     mne/time_frequency/tests/test_tfr.py::test_compute_tfr
0.46s call     mne/time_frequency/tests/test_tfr.py::test_tfr_ctf
0.40s call     mne/time_frequency/tests/test_tfr.py::test_time_frequency
0.39s call     mne/time_frequency/tests/test_stockwell.py::test_stockwell_api
0.20s call     mne/time_frequency/tests/test_tfr.py::test_tfr_multitaper
0.20s call     mne/time_frequency/tests/test_psd.py::test_psd_welch_average_kwarg[raw]
0.18s call     mne/time_frequency/tests/test_psd.py::test_psd_welch_average_kwarg[epochs]
0.18s call     mne/time_frequency/tests/test_psd.py::test_psd_welch_average_kwarg[evoked]
0.18s call     mne/time_frequency/tests/test_csd.py::test_csd_multitaper
0.16s call     mne/time_frequency/tests/test_ar.py::test_ar_raw
0.13s call     mne/time_frequency/tests/test_csd.py::test_csd_fourier
0.04s call     mne/time_frequency/tests/test_psd.py::test_compares_psd
0.04s call     mne/time_frequency/tests/test_stockwell.py::test_stockwell_core
0.03s call     mne/time_frequency/tests/test_csd.py::test_csd_morlet
==================================================== short test summary info =====================================================
SKIPPED [1] mne/time_frequency/tests/test_ar.py:17: Requires statsmodels version >= 0.8
SKIPPED [1] mne/time_frequency/tests/test_multitaper.py:14: Test test_dpss_windows skipped, requires nitime. Got exception (No module named 'nitime')
SKIPPED [1] mne/time_frequency/tests/test_multitaper.py:39: Test test_multitaper_psd skipped, requires nitime. Got exception (No module named 'nitime')
================================================= 84 passed, 3 skipped in 15.60s =================================================

[larsoner]

It's pip pre which uses prerelease versions of libraries. It looks unrelated. @drammock can you look? It has to do with latest matplotlib breaking something about our drawing / clicking / interaction, not sure what. Probably some upstream change modified how things get drawn or something.

[agramfort]

@adam2392 can you add/update a test?

[larsoner]

We have to be a little bit careful here because of the output option -- median on output='complex' data will almost certainly not work properly as is, see numpy/numpy#12943 / scipy/scipy#12669 / scipy/scipy#12676. I think we can fix by:

1. In the docstring add a warning to Notes about not passing np.median as it's unsafe.
2. If 'median' string-option is used, internally we need to separate into real and imaginary parts, take medians, then recombine. We can do this by defining a _median_complex in mne.fixes (I think eventually a fix will land in NumPy along these lines) that triages based on the iscomplexobj like in the SciPy PR, and wrap to that call rather than directly to np.median.
3. Add to the notes that we use the marginal median in the complex case

[agramfort]

can you clarify this? should we even accept it then? What is marginal median, can you better explain
or provide a ref?

[larsoner]

Passing in the string 'median' should be fine, internally we should notice this and compute the marginal median. Passing in np.median is unsafe (and should be documented as unsafe here) because NumPy don't compute the marginal median, they sort the complex values by real part and return whatever value is computed.

The marginal median is the median of each component separately, which is easy to compute:

https://en.wikipedia.org/wiki/Median#Marginal_median

This is the way we are going for SciPy:

scipy/scipy#12676 (comment)

[adam2392]

Updated docstring in 9ad5741

[agramfort]

let's document this ! thx @larsoner

…

[agramfort]

are the remaining CI issues expected?

[agramfort]

have you tested trim_mean with complex numbers? it behaves as expected?

[adam2392]

No good point. I'll take that out for now since it was just a copy/paste of Epochs, which is only real data. The main plus here in this PR was to allow one to use np.median, for being "robust to outliers"

[agramfort]

@larsoner feel free to merge if CI failures are expected. I did not keep track of which ones to ignore. thx