Labels for Epoch objects

[sourestdeeds]

Describe the bug

Upon fitting with ICA on an epochs instance, i pass in the object (along with a fitted ica instance on the same epoch object) using

ic_labels = label_components(epochs, ica, method='iclabel')

The error generated is

/var/folders/ph/9w0w89cd5pq1d7z1_kn2mt8m0000gp/T/ipykernel_18124/129543192.py in process_eeg(cap_size, ppt_num, session)
     77     ica.fit(epochs)
     78     # ICA Labels
---> 79     ic_labels = label_components(epochs, ica, method='iclabel')
     80     labels = ic_labels["labels"]
     81     exclude_idx = [idx for idx, label in enumerate(labels) 

/usr/local/anaconda3/lib/python3.9/site-packages/mne_icalabel/label_components.py in label_components(inst, ica, method)
     54     _check_option("method", method, methods)
     55     _validate_inst_and_ica(inst, ica)
---> 56     labels_pred_proba = methods[method](inst, ica)
     57     labels_pred = np.argmax(labels_pred_proba, axis=1)
     58     labels = [ICLABEL_NUMERICAL_TO_STRING[label] for label in labels_pred]

/usr/local/anaconda3/lib/python3.9/site-packages/mne_icalabel/iclabel/label_components.py in iclabel_label_components(inst, ica)
     43     .. footbibliography::
     44     """
---> 45     features = get_iclabel_features(inst, ica)
     46     labels_pred_proba = run_iclabel(*features)
     47     return labels_pred_proba

/usr/local/anaconda3/lib/python3.9/site-packages/mne_icalabel/iclabel/features.py in get_iclabel_features(inst, ica)
     80 
     81     # compute psd feature (float32)
---> 82     psd = _eeg_rpsd(inst, ica, icaact)
     83 
     84     # compute autocorr feature (float32)

/usr/local/anaconda3/lib/python3.9/site-packages/mne_icalabel/iclabel/features.py in _eeg_rpsd(inst, ica, icaact)
    243     assert isinstance(inst, (BaseRaw, BaseEpochs))  # sanity-check
    244     constants = _eeg_rpsd_constants(inst, ica)
--> 245     psd = _eeg_rpsd_compute_psdmed(inst, icaact, *constants)
    246     psd = _eeg_rpsd_format(psd)
    247     return psd

/usr/local/anaconda3/lib/python3.9/site-packages/mne_icalabel/iclabel/features.py in _eeg_rpsd_compute_psdmed(inst, icaact, ncomp, nfreqs, n_points, nyquist, index, window, subset)
    314         elif isinstance(inst, BaseEpochs):
    315             temp = np.hstack([icaact[it, index[:, k], :] for k in range(index.shape[-1])])
--> 316             temp = temp.reshape(index.shape[0], len(inst), order="F")
    317         else:
    318             raise RuntimeError  # should never happen

ValueError: cannot reshape array of size 818048 into shape (128,913)

If I instead use the following raw instance to classify the epochs it works, is this intended? Essentially I'm performing ICA after epoching, so I'm not sure how using the fitted ica on epoch data and then using the raw to compare this with is correct.

ic_labels = label_components(raw, ica, method='iclabel')

Perhaps I have some funky rank deficiency going on somewhere from interpolation or something?

[mscheltienne]

Hello, thank you for reporting the issue! but I can not reproduce based on your description..

from mne import make_fixed_length_epochs
from mne.datasets import sample
from mne.io import read_raw
from mne.preprocessing import ICA

from mne_icalabel import label_components

# load raw instance
directory = sample.data_path() / "MEG" / "sample"
fname = "sample_audvis_filt-0-40_raw.fif"
raw = read_raw(directory / fname, preload=False)
raw.pick_types(eeg=True)
raw.crop(tmin=0, tmax=20)
raw.load_data()

# create epochs
epochs = make_fixed_length_epochs(raw, duration=1., preload=True)

# fit an ICA on epochs
ica = ICA(n_components=5, method="picard")
ica.fit(epochs)

# run ICLabel
label_components(epochs, ica, method="iclabel")

Outputs:

>>> {'y_pred_proba': array([0.96021754, 0.4062748 , 0.92534626, 0.95375717, 0.5524776 ],
       dtype=float32),
 'labels': ['eye blink', 'other', 'brain', 'brain', 'other']}

Could you try to reproduce the issue by using a dataset from the MNE sample dataset as I did above, or could you share a minimally reproducible example with the attached files that raises the error?

[mscheltienne]

Looks similar to the issue in #38

[sourestdeeds]

If it helps I'm having a look at the line this is caused by. The shape of my raw and epochs are

epochs.get_data().shape
(909, 122, 513)
raw.get_data().shape
(122, 584600)

The index variable calculated in _eeg_rpsd_constants for raw and epochs is

# raw
(128, 9133)
# epochs
(128, 7)

[mscheltienne]

It's difficult to debug if I can not reproduce the error. Could you share some script and files to reproduce it?

For the shapes, with the example from the MNE sample dataset, I have:

epochs.get_data().shape
>>> (20, 59, 150)
raw.get_data().shape
>>>  (59, 3004)

Which looks similar to you, assuming you do have 122 EEG channels and 909 epochs. And for the index variables in _eeg_rpsd_constants, breaking on line 274 (executing it): https://github.com/mne-tools/mne-icalabel/blob/ecc6a917e9540acdaf148845a145c882fcfc3b26/mne_icalabel/iclabel/features.py#L274

# with inst as Epochs
>>> (150, 1)
# with inst as Raw
>>> (150, 39)

It's difficult to say if it's normal or not as it depends on the number of samples in your data and on the sampling frequency.

---

What happens if you try to re-fit a fresh ICA instance on the epochs?

from mne.preprocessing import ICA
from mne_icalabel import label_components

ica = ICA(n_components=5, method="picard")  # just to speed-up
ica.fit(epochs)
label_components(raw, ica, method="iclabel")

[mscheltienne]

Am I correct in assuming the sampling frequency of your data is 128 Hz?

[sourestdeeds]

Yeah, I downsampled the raw to 128Hz and then went

raw.resample(128)
raw.notch_filter(freqs=None, method='spectrum_fit', n_jobs=12)
raw.filter(7, 35, n_jobs=12)
raw.set_eeg_reference()

events, event_ids = mne.events_from_annotations(
    raw, verbose = False)
epochs = mne.Epochs(
    raw=raw, 
    events=events, 
    event_id=event_ids,
    preload=True, 
    tmin=0,
    tmax=4,
    baseline=None, 
    event_repeated='merge',
)

ica = mne.preprocessing.ICA(
    method='picard', fit_params=dict(extended=True, ortho=False)
)
ica.fit(epochs)
# ICA Labels
ic_labels = label_components(raw, ica, method='iclabel') # epochs object causes error
labels = ic_labels["labels"]
exclude_idx = [idx for idx, label in enumerate(labels) 
                if label not in ["brain"]]
print(f"Excluding these ICA components: {exclude_idx}")
# Apply ICA
ica.plot_overlay(epochs.average(), exclude=exclude_idx)
ica.apply(epochs, exclude=exclude_idx)

[mscheltienne]

Thank you! I managed to reproduce it with this sample, changing the duration of epochs to 4 seconds:

from mne import make_fixed_length_epochs
from mne.datasets import sample
from mne.io import read_raw
from mne.preprocessing import ICA

from mne_icalabel import label_components

# load raw instance
directory = sample.data_path() / "MEG" / "sample"
fname = "sample_audvis_filt-0-40_raw.fif"
raw = read_raw(directory / fname, preload=False)
raw.pick_types(eeg=True)
raw.crop(tmin=0, tmax=20)
raw.load_data()

# create epochs
epochs = make_fixed_length_epochs(raw, duration=4., preload=True)

# fit an ICA on epochs
ica = ICA(n_components=5, method="picard")
ica.fit(epochs)

# run ICLabel
label_components(epochs, ica, method="iclabel")

I'll have a look.

[mscheltienne]

@sourestdeeds It's fixed on main and will be fixed in release >= 0.2.0.

[sourestdeeds]

Thats amazing! Thanks so much, my automated workflow is working super well now.

Hopefully i can use this to get some classification results on my super noisy data