Can't obtain graph when using detect_artifacts

Support & Discussions
,
1

MNE-Python version: 0.23
operating system: Mac OS 11.3.1, Python 3.9

I have a 10 minute recording with sampling frequency 200 hz and have dropped all but the first 9 channels for purposes of visualizing sources after performing ICA. I don’t have any information regarding electrode recording locations.

From the documentation I understand detect_artifacts is experimental, but I was under the impression that it plots all sources and indicates by highlighting those that can be dropped.

My code so far after dropping all but the first 9 channels from raw:

raw = raw.filter(l_freq=1, h_freq=None)  #high pass with 1 Hz cutoff
ica = ICA(method='fastica', random_state=180, max_iter="auto")
ica.fit(raw)

Using ica.plot_sources(raw, start=0,stop=180,show_scrollbars=True)
I wasn’t able to discern which of the 10 sources could be “excluded.”

I instead tried

ica.detect_artifacts(raw,start_find=0, stop_find=180)
print(ica.exclude)

The result indicated sources to be excluded were 0, 1, and 6,. This was for the first 180 seconds, and, not surprisingly, by changing start_find and stop_find , I found different values ofica.exclude. However, I was under the impression that ica.detect_artifacts should also plot all 10 sources, with the excluded ones not highlighted. Is that correct?

My goal is to look at the graphs produced by ica.plot_sources and ica.detect_sources over the same time intervals to gain a better sense of where the artifacts are located.

2

Hello @PaulF,

wow, how did you find out about ICA.detect_artifacts()? It hasn’t received much love over the past years, isn’t mentioned in any tutorial, and honestly I believe we should deprecate and eventually remove it (but this needs to be discussed with @cbrnr, @sappelhoff and @Denis)

Instead of using this method, I’d suggest you follow the procedures laid out in the ICA tutorial here:

https://mne.tools/stable/auto_tutorials/preprocessing/40_artifact_correction_ica.html#example-eog-and-ecg-artifact-repair

Specifically, you’ll want to use ICA.find_bads_eog() and potentially ICA.find_bads_ecg().

Please feel free to report back how it goes and if you start to hit any walls!

Best wishes,
Richard

1 Like
3

I found it here:

https://mne.tools/stable/generated/mne.preprocessing.ICA.html?highlight=ica%20plot_sources#mne.preprocessing.ICA.plot_sources

right below ICA.copy()

I’ll look at the tutorial. Thanks much.:slight_smile:

4

Ah wow, so you actually scrolled through the API docs :smiley: Thumbs up for your thorough research, in any case!

1 Like
5

Oh hell, this should have been deprecated years ago. The code comes from a very early phase of our exploration with support for the right API for ICA. git blame says I contributed this 9 years ago. +100 for deprecation. The core API for ICA as one can see it today consolidate ~5-6 years ago. We apparently forgot to remove detect_artifacts.

6

Well–it was a fun find!

With regard to ICA.find_bads_eog():

It returns eeg_idx and scores, which essentially tells us how to sort our EOG-related components by their corresponding correlation scores. (In decreasing order I presume?)

This reminds me of PCA, where I can order eigenvectors of a covariance matrix by their corresponding eigenvalues, listed in decreasing order. If N eigenvalues account for, say 90% of my total variance, then I might use the first N eigenvectors as my principal components.

Can a similar approach be used here?

7

Whoops, I’ll do that for the MNE-Python 0.24 release then. Will tag you on GitHub. Thanks for the quick response, @Denis!

It’s one score per ICA component, starting from the first component. It’s not ordered by the score, but by components, in ascending order.

I don’t think I fully understand your train of thoughts here, sorry :frowning: What you’ll typically want to do is set the threshold of find_bads_eog() such that EOG-related components are marked as such, while keeping a low false-alarm rate.

But I suppose you already know that and are actually asking about something else?

FWIW ICA‘s n_components parameter can also accept a floating-point number to determine which fraction of the total variance should be explained.

Best wishes,
Richard

8

I can possibly answer your question. Those scores refer to Pearson correlation coefficients between the signal in the EOG channel and the components. It is meant to identify the ICA components containing EOG artifacts. In other words, components with the highest correlation coeffients likely represent EOG artifacts (of different origin). After finding these, you can drop these components from your data and obtain processed data free of EOG artifacts.

So in that sense it is not related to PCA, which treats feature dimensionality.

1 Like
9

These comments have been extraordinarily helpful. Is this workflow the sort of thing I could be doing to repair my original signals, or am I overlooking some important concepts? (I’ll use just 12 of the original channels.)

channels=['Fp1', 'Fp2', 'F3', 'F4', 'C3', 'C4', 'P3', 'P4', 'O1', 'O2', 'F7', 'F8']
raw_filt = raw.filter(l_freq=1, h_freq=None)
ica = ICA(method='fastica',random_state=5,max_iter="auto")
ica.fit(raw_filt)   
eog_idx,eog_scores=ica.find_bads_eog(raw,ch_name=channels, measure='correlation',threshold=.90,reject_by_annotation=False)
raw_repaired=ica.apply(raw, exclude=eog_idx)

By lowering the threshold, I’m going to see more components I can exclude. For example, if the threshold is .9, I exclude ICA000. But if I lower it to .7, I get ICA000, ICA001, and ICA004.

1. Is there a reasonable threshold I could use in general, say .9?

Also, below is a plot of my components.

2. Is there a reason I would see ICA004 excluded before ICA005, which closely resembles ICA000?

It’s not exactly clear what I should be looking for.

3. If I also want to consider heartbeat artifacts, could I then use a process similar to that above, applying find_bads_eog to raw_repaired as defined above? Should one type of repair be performed before the other?

Screen Shot 2021-10-28 at 12.36.05 PM
Screen Shot 2021-10-28 at 12.36.05 PM1918×1472 151 KB

10

Well, I am only a user of MNE, but my advice would be to base your criterion on the visual plot of the components, as you show in the picture. EOG components have typcially the highest amplitude near the nose and are often symmetrical distributed. So based on the visuals, I would identify 00,04 an 05 as EOG components that you want to remove.You can also plot the power spectrum and confirm that the power peaks at frequencies between 1-13 Hz. frequency response

So 0.7 seems like a reasonable value.

ICA10 looks a bit weird btw, but no idea what is going on there, maybe a cable or something.

For the heartbeats, you could follow a similair procedure using ecg_events instead:https://mne.tools/stable/generated/mne.preprocessing.find_ecg_events.html

And you can use:

mne.preprocessing.ICA — MNE 0.23.4 documentation

for plotting the spectrum.

11

The “reasonable” threshold depends on the data you have and will vary between labs, EEG systems, and even within a study between participants. You need to find a threshold that seems to capture the EOG-related ICs while not producing false positives, as you don’t want to end up removing components that are associated with brain activity. If in doubt, choose a higher threshold and keep some artifacts to avoid accidentally removing physiological brain activity.

I assume that the the amplitudes of the component topographies are largely different. You can get colorbars by passing colorbar=True to plot_components().

At this stage I would recommend you start reading up on ICA in general and how components associated with common artifacts (EOG, ECG) are expected to look like. I don’t have a good literature reference at hand right now, but feel free to ask if you cannot find anything usable via Google.

FWIW, I found the best way to “learn” ICA is by sitting together with a coworker who’s experienced in ICA artifact removal and goes with you through the data together.

No, ideally you will want to run find_bads_eog and find_bads_ecg on EOG and ECG epochs, created via mne.preprocessing.create_eog_epochs() and mne.preprocessing.create_ecg_epochs(). Please work through the ICA tutorial to understand this approach:

https://mne.tools/stable/auto_tutorials/preprocessing/40_artifact_correction_ica.html#example-eog-and-ecg-artifact-repair

Once you have discovered a list of ICA components associated with EOG and another list of components associated with ECG activity, you will want to create the union of those components and add it to ica.exclude, e.g.:

excludes = sorted(set(excludes_eog + excludes_ecg))
ica.exclude = excludes
ica.apply(raw)

Best wishes,
Richard