Running ICA

What does this show?

Why do I need this step?

In this lesson you only run ICA. Classifying and removing components is the ICLabel lesson.

The reason is that slow drift confuses ICA, and the official tutorial says decompositions are notably higher quality with a high-pass above 1 Hz. But a 1 Hz filter distorts slow ERP components such as P3.

The solution is two copies of the same data:

1. The ICA copy, filtered at 1 Hz, used only to compute ICA.
2. The main dataset, filtered at 0.1 Hz, which receives the ICA result and carries on to the analysis.

The ICA result is a set of weights, a fixed way of combining channels. That is why it is valid to apply it to the main dataset, as the official tutorial explains.

How do I do it in EEGLAB?

Start with the main dataset loaded, as in the Starting point. We call it dataset 1.

First: create the ICA copy

1. Make sure dataset 1 is current in the Datasets menu.
2. Open Tools › then Filter the data › then Basic FIR filter (new, default).
3. In Lower edge of the frequency pass band (Hz), type 1. Leave the other fields as they are and press OK.
4. Name the new dataset D1 ica1hz and leave Overwrite it in memory (set=yes; unset=create a new dataset) unchecked. You now have dataset 2, and dataset 1 is unchanged.

Second: run ICA on the ICA copy

1. With dataset 2 current, open Tools › then Decompose data by ICA.
2. In ICA algorithm to use (click to select), keep the first option, Extended Infomax (runica.m; default). The options field below shows 'extended', 1.
3. Keep Reorder components by variance (if that’s not already the case) checked, leave Use only channel type(s) or indices empty, and press OK.
4. A window may say that the data rank is lower than the number of channels, with a Proposed rank: field. Compare it with the expected rank: the number of channels, minus the interpolated channels, minus 1 for the average reference. In D1 it is 31 out of 32. If the number differs, type the expected one, then press OK.
5. Wait. The command window prints the training steps, which can take from under a minute to a few minutes depending on your computer.

Third: display the component maps

1. Open Plot › then Component maps › then In 2-D and keep the defaults.
2. Look for the blink component: a strong map at the front of the head, symmetric between left and right. It is usually among the first components.

Fourth: copy the weights to the main dataset

1. Go back to dataset 1 from the Datasets menu.
2. Open Edit › then Dataset info.
3. Next to ICA weights array or text/binary file (if any):, press from other dataset, type 2, then OK. Three fields fill in automatically: ALLEEG(2).icaweights, ALLEEG(2).icasphere, and ALLEEG(2).icachansind.
4. Press OK to close the window.
5. Save the main dataset with File › then Save current dataset as. The ICLabel lesson starts from it.

How do I write it as code?

% Main dataset (high-pass 0.1 Hz, line noise removed, average reference) is in EEG
EEG_main = EEG;

% ICA copy: same data, high-pass at 1 Hz
EEG_ica = pop_eegfiltnew(EEG_main, 'locutoff', 1);
EEG_ica.setname = 'D1 ica1hz';

% Expected rank: channels - interpolated channels - 1 (average reference)
expectedRank = EEG_ica.nbchan - numel(badChans) - 1;
fprintf('Expected data rank: %d\n', expectedRank);

% Extended Infomax on the ICA copy
EEG_ica = pop_runica(EEG_ica, 'icatype', 'runica', 'extended', 1);
% If "Data rank (N)" is missing or N differs from expectedRank, run instead:
% EEG_ica = pop_runica(EEG_ica, 'icatype', 'runica', 'extended', 1, 'pca', expectedRank);

% Copy the decomposition back to the main dataset
EEG_main.icaweights  = EEG_ica.icaweights;
EEG_main.icasphere   = EEG_ica.icasphere;
EEG_main.icachansind = EEG_ica.icachansind;
EEG_main.icawinv     = [];   % recomputed by eeg_checkset
EEG_main.icaact      = [];
EEG = eeg_checkset(EEG_main);

% Scalp maps of the first 20 components
pop_topoplot(EEG, 0, 1:20);

EEG = pop_saveset(EEG, 'filename', 'D1_main_ica.set', 'filepath', pwd);

In a script, the rank window does not appear. pop_runica tries to detect the reduced rank, but its estimate can fail. When it detects it, it prints a line such as Data rank (31) is smaller than the number of channels (32). Compare that number with the Expected data rank the script printed. If the line is missing or the numbers differ, use the commented line that passes 'pca' with the expected value. If your data has bad channels, for example 32 channels with 2 interpolated, the expected rank is 29.

Which value should I choose?

How do I check it worked?

1. In the main EEGLAB window, with dataset 1 current, ICA weights should show Yes.
2. In the command window, type size(EEG.icawinv). In D1 without interpolated channels it shows 32 31: 32 channels and 31 ICA components. Each interpolated channel reduces the number of components by one.
3. Open the component maps with Plot › then Component maps › then In 2-D on dataset 1. You should see the same maps as on the ICA copy, including the blink component.
4. Reload the saved file with File › then Load existing dataset, and check that ICA weights still shows Yes.
5. The rank matches: the Data rank that pop_runica printed equals the Expected data rank. In D1 both are 31.

What mistakes do people make?

• Computing ICA on 0.1 Hz data, or continuing the ERP on the 1 Hz copy. The first gives a weaker decomposition, the second distorts P3 and N400. The two copies together are the solution.
• Trusting the automatic rank estimate. pop_runica tries to detect the rank, but the estimate can fail. Signs of a missed rank: several components that look almost identical or like pure noise, or ICA that does not converge. Always compare the printed Data rank with the expected rank (channels minus interpolated channels minus 1). If they differ, type 'pca', 31 in the options field, with your expected number instead of 31.
• Non-brain channels slipping in by accident. Exclude EOG and ECG channels according to your goal, and never let them in without a decision. Use Use only channel type(s) or indices to choose the channels. In D1 we decided to keep the two eye channels for a clear reason.
• Expecting the same result every time. ICA shuffles the data order at random during training, so the order and shape of components vary slightly between runs. Your component order may not match the picture in this lesson: the blink may be component 1 here and component 3 for you. Find a component by its shape, not its number.
• Changing the reference after ICA. The weights are tied to the reference used when ICA was computed. Re-reference before ICA and leave it unchanged afterward. EEGLAB itself warns you if you change it.

How do I write it in my paper?

A template: edit the data, settings and versions to match what you ran.

Independent component analysis (extended Infomax; Lee et al., 1999) was computed with runica ([31] components after PCA reduction to the data rank) on a copy of the data high-pass filtered at 1 Hz, and the resulting weights were applied to the 0.1 Hz high-pass filtered data.

Going further

• EEGLAB Tutorial: https://eeglab.org/tutorials/06_RejectArtifacts/RunICA.html (6. Reject artifacts > d. Indep. Comp. Analysis: Running ICA decompositions)
• EEGLAB Tutorial: https://eeglab.org/tutorials/06_RejectArtifacts/RunICA.html (How to deal with the aggressive high-pass filter applied before running ICA)
• EEGLAB Tutorial: https://eeglab.org/tutorials/06_RejectArtifacts/RunICA.html (Issues with data rank deficiencies)
• Luck, Applied Event-Related Potential Data Analysis (1st edition): https://socialsci.libretexts.org/Bookshelves/Psychology/Biological_Psychology/Applied_Event-Related_Potential_Data_Analysis_(Luck) (Chapter 9: Artifact Correction with ICA)
• ICA applied to EEG part 7: Running ICA in EEGLAB and visualizing components: https://www.youtube.com/watch?v=2hrYEYSycGI (Full video)