[Eeglablist] Rank deficiency and ASR+ICA

Mon Sep 23 11:58:30 PDT 2024

Hi Francisco,

   1. If HPF1.0Hz+ ICA changes results compared with HPF0.1Hz + ICA, you
   want to first identify what your data have between 0.1 and 1.0 Hz: Is it an
   artifact or brain signals?
   2. If you are convinced that it is an artifact, you do not want to use
   that freq band of the signals anyway.  Run HPF1.0Hz + ICA. You do not lose
   anything.
   3. If you are convinced that it is a brain signal, run HPF0.1Hz + ICA.
   You do not lose anything.

Personally I am not a big fan of using HPF1.0Hz + ICA and copying the data
to HPF0.1Hz. It is less straightforward and not a general solution. When it
is more desirable depends on your data after all.
Here is my analogical explanation of what happens when you do it.
https://sccn.ucsd.edu/wiki/Makoto's_preprocessing_pipeline#What_happens_to_the_.3C_1_Hz_data_if_ICA_is_calculated_on_.3E_1_Hz_data_and_applied_to_0.1_Hz_data.3F_.2805.2F18.2F2022_Updated.29

Makoto

On Mon, Sep 23, 2024 at 1:13 PM Francisco Zambrano <
panchozambrano84 at gmail.com> wrote:

> Hello again,
>
> Thanks a lot for your insights, the results from Cedric's reply gave me an
> answer that sounds more in line with what I've read.
>
> I've come with another question regarding ICA. If I'm trying to use ICA's
> trick of applying the >1 Hz data's unmixing matrix  to the <1 Hz data, but
> I don't know if It's better to do the whole processing pipeline(removing
> bad channels, interpolating, removing bad trials,...) first and then
> applying it back to the <1Hz data that is also been processed with the same
> steps; or if I'm supposed to do it in another way. Could you please help me
> with this?
>
> Thanks in advance.
> Francisco.
>
>
> El mié, 7 ago 2024 a la(s) 3:37 p.m., Makoto Miyakoshi via eeglablist (
> eeglablist at sccn.ucsd.edu) escribió:
>
>> Hi Pancho,
>>
>> > rank(EEG.data(:,:)) over the original data It says the rank is 11
>>
>> Try rank(double(EEG.data(:,:)) and you'll see more realistic values at
>> least. It is counterintuitive if you see it naively that single or double
>> makes such a difference here.
>>
>> For the complete answer, see Cedric's reply.
>>
>> Makoto
>>
>>
>> On Mon, Aug 5, 2024 at 9:40 PM Pancho Zambrano via eeglablist <
>> eeglablist at sccn.ucsd.edu> wrote:
>>
>> > Hi!,
>> >
>> > I'm currently working on a script for my thesis on an EEG channel
>> > reconstruction protocol for Universidad de Concepción, and I've found
>> > myself in some trouble regarding rank deficiency. The script it's about
>> > quantification of the positive yzing impact of ASR algorithm from
>> > cleanrawdata() on the data by analyzing the increase of  the amount of
>> > dipolar sources present before and after using this default toolbox
>> > function. This is accomplished by comparing the residual variance from
>> the
>> > dipfit model of the clean and original data. The problem is I use the
>> > matcorr function to find the best match between the original and clean
>> ICA
>> > mixing matrices, but they have different number of columns: the first
>> one
>> > is 64x64 and the clean one is 64x55. Since I've re-referenced the data
>> > after ASR, I thought this was the problem but I also noticed when I use
>> > rank(EEG.data(:,:)) over the original data It says the rank is 11, so
>> I'm
>> > lost here, please help. I got this method from the article:
>> >
>> > Evaluation of Artifact Subspace Reconstruction for Automatic Artifact
>> > Components Removal in Multi-Channel EEG Recordings
>> >
>> > from:
>> >
>> >
>> https://urldefense.com/v3/__https://ieeexplore.ieee.org/abstract/document/8768041?casa_token=IsR0oujDOEAAAAAA:T1NPgCRqT15sUYjG8KLtR8a9yW6PqQP7iaKzg5057R93Oef1yCmxDKt6ZtLM4IZk2t_xF4qWJ0I__;!!Mih3wA!ClbUhbvVFcY2DhNRP-wwMym8wBcZTXz2jeCl8yebysZiEmG01Z-cZm3OgV4T4R1GXwiDsWNEwyIFZcOglDZZkbBiLjHP0sFV$
>> > This is my code:
>> > % Initialize EEGLAB
>> > eeglab;
>> >
>> > % Load the dataset
>> > EEG = pop_loadset('filename', 'C110_filt_butt_zap2.set');
>> >
>> > % Perform ICA on the original data
>> > EEG = pop_runica(EEG, 'extended', 1);
>> >
>> > % Save the original ICA results
>> > ica_weights_orig = EEG.icaweights;
>> > ica_sphere_orig = EEG.icasphere;
>> > ica_act_orig = EEG.icaact;
>> > ica_mix_matrix_orig = pinv(ica_weights_orig * ica_sphere_orig);
>> >
>> > % ASR dataset
>> > EEG_clean = pop_loadset('filename', 'prueba_1hz_ica.set');
>> >
>> > % % Perform ICA on the ASR cleaned data
>> > % EEG_clean = pop_runica(EEG_clean, 'extended', 1);
>> >
>> > % Save the ASR ICA results
>> > ica_weights_clean = EEG_clean.icaweights;
>> > ica_sphere_clean = EEG_clean.icasphere;
>> > ica_act_clean = EEG_clean.icaact;
>> > ica_mix_matrix_clean = pinv(ica_weights_clean * ica_sphere_clean);
>> >
>> > % Calculate the correlations between the components
>> > % n_components_orig = size(ica_mix_matrix_orig, 2);
>> > % n_components_clean = size(ica_mix_matrix_clean, 2);
>> > %
>> > % corr_matrix = zeros(n_components_orig, n_components_clean);
>> > %
>> > % for i = 1:n_components_orig
>> > %     for j = 1:n_components_clean
>> > %         corr_matrix(i, j) = corr(ica_mix_matrix_orig(:, i),
>> > ica_mix_matrix_clean(:, j));
>> > %     end
>> > % end
>> >
>> > % Use the matcorr function to find the best matching
>> > rmmean = true; % Adjust as needed (remove mean or not)
>> > method = 0; % Correlation method
>> > weighting = []; % No specific weighting
>> >
>> > [corr, indx, indy, corrs] = matcorr(ica_mix_matrix_orig,
>> > ica_mix_matrix_clean, rmmean, method, weighting);
>> >
>> > % Display the matched correlations
>> > best_matches = corrs(sub2ind(size(corrs), indx, indy));
>> >
>> > disp('Correlations between the best matched components:');
>> > disp(best_matches);
>> >
>> > % Apply the spatial filter obtained from the ICA of the original data to
>> > the cleaned data
>> > W_orig = ica_weights_orig * ica_sphere_orig;
>> > Y_clean = W_orig * reshape(EEG_clean.data, size(EEG_clean.data, 1), []);
>> >
>> > % Calculate the activities of the ICs
>> > Y_clean = reshape(Y_clean, size(EEG_clean.icaact));
>> >
>> > % Calculate the mean power reduction for the IC activities
>> > Y_orig = W_orig * reshape(EEG.data, size(EEG.data, 1), []);
>> > Y_orig = reshape(Y_orig, size(EEG.icaact));
>> >
>> > power_reduction = mean(var(Y_orig, 0, 2)) - mean(var(Y_clean, 0, 2));
>> >
>> > disp('Mean power reduction:');
>> > disp(power_reduction);
>> >
>> > % Classify the ICs using iclabel
>> > EEG = iclabel(EEG, 'default');
>> > EEG_clean = iclabel(EEG_clean, 'default');
>> >
>> > % Get the IC classification labels
>> > ic_labels_orig = EEG.etc.ic_classification.ICLabel.classifications;
>> > ic_labels_clean =
>> EEG_clean.etc.ic_classification.ICLabel.classifications;
>> >
>> > % Set up DIPFIT for the original data
>> > EEG = pop_dipfit_settings(EEG, 'hdmfile',
>> 'standard_BEM/standard_vol.mat',
>> > ...
>> >     'coordformat', 'MNI', 'mrifile', 'standard_mri.mat', ...
>> >     'chanfile', 'standard_1005.elc', 'coord_transform', [0 0 0 0 0 0 1 1
>> > 1], ...
>> >     'chansel', 1:size(EEG.data, 1));
>> >
>> > % Fit dipoles to the original ICA components
>> > EEG = pop_multifit(EEG, 1:size(EEG.icaweights, 1), 'threshold', 100);
>> >
>> > % Get the residual variance of the fitted dipoles in the original data
>> > rv_orig = [EEG.dipfit.model.rv];
>> >
>> > % Set up DIPFIT for the ASR cleaned data
>> > EEG_clean = pop_dipfit_settings(EEG_clean, 'hdmfile',
>> > 'standard_BEM/standard_vol.mat', ...
>> >     'coordformat', 'MNI', 'mrifile', 'standard_mri.mat', ...
>> >     'chanfile', 'standard_1005.elc', 'coord_transform', [0 0 0 0 0 0 1 1
>> > 1], ...
>> >     'chansel', 1:size(EEG_clean.data, 1));
>> >
>> > % Fit dipoles to the ASR cleaned ICA components
>> > EEG_clean = pop_multifit(EEG_clean, 1:size(EEG_clean.icaweights, 1),
>> > 'threshold', 100);
>> >
>> > % Get the residual variance of the fitted dipoles in the ASR cleaned
>> data
>> > rv_clean = [EEG_clean.dipfit.model.rv];
>> >
>> > % Identify dipolar ICs (residual variance < 5%)
>> > dipolar_ICs_orig = find(rv_orig < 0.05);
>> > dipolar_ICs_clean = find(rv_clean < 0.05);
>> >
>> > % Compare the number of dipolar ICs before and after ASR
>> > num_dipolar_ICs_orig = length(dipolar_ICs_orig);
>> > num_dipolar_ICs_clean = length(dipolar_ICs_clean);
>> >
>> > % Display dipolar IC results
>> > fprintf('Number of dipolar ICs in original data: %d\n',
>> > num_dipolar_ICs_orig);
>> > fprintf('Number of dipolar ICs in ASR cleaned data: %d\n',
>> > num_dipolar_ICs_clean);
>> >
>> > % Save all results
>> > save('ica_evaluation_results.mat', 'ica_weights_orig',
>> 'ica_weights_clean',
>> > 'ica_sphere_orig', 'ica_sphere_clean', 'ica_act_orig', 'ica_act_clean',
>> > 'indx', 'indy', 'best_matches', 'Y_clean', 'power_reduction',
>> > 'ic_labels_orig', 'ic_labels_clean', 'rv_orig', 'rv_clean',
>> > 'dipolar_ICs_orig', 'dipolar_ICs_clean');
>> >
>> > Please, help.
>> >
>> > Francisco Alonso Zambrano Salamanca
>> > Biomedical Engineering
>> > Universidad de Concepción
>> > Chile
>> > _______________________________________________
>> > Eeglablist page: http://sccn.ucsd.edu/eeglab/eeglabmail.html
>> > To unsubscribe, send an empty email to
>> > eeglablist-unsubscribe at sccn.ucsd.edu
>> > For digest mode, send an email with the subject "set digest mime" to
>> > eeglablist-request at sccn.ucsd.edu
>> _______________________________________________
>> Eeglablist page: http://sccn.ucsd.edu/eeglab/eeglabmail.html
>> To unsubscribe, send an empty email to
>> eeglablist-unsubscribe at sccn.ucsd.edu
>> For digest mode, send an email with the subject "set digest mime" to
>> eeglablist-request at sccn.ucsd.edu
>
>