[Eeglablist] fft analysis

[goharbak at uni-bremen.de]

Hi dear community,

I have one question regarding my fft analysis. I preprocessed my data in eeglab lab and epoched them according to my two blocks including inclusion and exclusion). Now I want to run fft with hanning window, but I am not sure about my coding. Do you have any feedback on my coding?

1.I extracted the data from some channels and got the average of them to have one channel.

2.I extracted my data for inlcusion block and exlcusion block according to my markers separately.

3.I created a hanning window and run fft for each epoch of inlcusion and exclsuion blcoks separately

4.I made a loop through alpha and theta frequencies and extract the right numbers
You can see my coding below:


Subject_Index = {'mother1.set'};

%% Defining parameters
Condition_Marker = {'S 200', 'S 201'}; % enter markers
Condition_Names={'inclusion', 'exclusion'}; % enter name
% Frequency settings
cFreqNames = {'Alphamother', 'Thetamother'}; % Names for Alpha, Theta, and their peaks

% Frequency ranges for Alpha and Theta for mothers
cFreqs{1} = [8 13];  % Alpha
cFreqs{2} = [4 8];  % Theta


cNFreqs = 2;  % Number of frequency ranges (2 bands and their peaks)


cWinSec = 1.0;
srate = 500;      % sampling rate
cWinSize = srate * cWinSec; % FFT resolution will be .5 Hz
cNSubjects = length(Subject_Index);
mChans    = [1, 2, 9, 10, 11]; % average of midfrontal channels
oChans    = [3, 4, 5, 6, 7, 8];  % average of occipital channels
cNChans= 1

%%

for subj=1:length(Subject_Index)
    %% START LOOP

    EEG = pop_loadset('filename', Subject_Index{subj}, 'filepath', rawdata_location)
    EEG = eeg_checkset(EEG);



    % Set ROI channels and create data matrix with the data in the corresponding channels
    fields = {EEG.chanlocs.labels};
    channel1 = 'F7';
    channel2 = 'Fz';
    channel3 = 'F4';
    channel4 = 'F8';
    channel5 = 'Cz';

    % Check if channels exist and collect their indices
    channels = [];
    if ismember(channel1, fields), channels = [channels, find(strcmp(channel1, fields))]; end
    if ismember(channel2, fields), channels = [channels, find(strcmp(channel2, fields))]; end
    if ismember(channel3, fields), channels = [channels, find(strcmp(channel3, fields))]; end
    if ismember(channel4, fields), channels = [channels, find(strcmp(channel4, fields))]; end
    if ismember(channel5, fields), channels = [channels, find(strcmp(channel5, fields))]; end

    Frequency.ROI = channels;

    % Extract data for the selected channels while preserving the epochs
    data_raw = EEG.data(Frequency.ROI, :, :);  % data_raw is now 5 x samples x epochs

    % Average the data across channels while keeping the epoch dimension
    data_power = mean(data_raw, 1);  % Average across channels (1st dimension)


%% extracting epoched data

    selected_epochs.inclusion = [];
    selected_epochs.exclusion = [];
    for cond = 1:length(Condition_Marker)
        current_marker = Condition_Marker{cond};  % Current event marker (S 200 or S 201)
        current_name = Condition_Names{cond};     % Corresponding condition name (inclusion or exclusion)

        % Find the event indices corresponding to the current marker
        event_indices = find(strcmp({EEG.event.type}, current_marker)); % Find epochs for the current marker

        % Initialize a matrix to store the selected epochs
        % Dimensions: [1 x timepoints x trials]
        selected_epochs_power = zeros(1, 600, length(event_indices));  % Adjust 600 to your actual timepoints

        % Loop through each event index and fill in the selected_epochs_power
        for trial = 1:length(event_indices)
            % Select the corresponding trial's power data
            selected_epochs_power(:, :, trial) = data_power(:, :, event_indices(trial));
        end

        % Store the data for each condition separately
        selected_epochs.(current_name) = selected_epochs_power; % Save for inclusion or exclusion
    end


    inclusion_data=selected_epochs.inclusion
    exclusion_data=selected_epochs.exclusion





    numTrials_in = size(inclusion_data, 3);  % Number of trials for inclusion
    numTrials_ex = size(exclusion_data, 3);  % Number of trials for exclusion

    chanpowr_in = zeros(cWinSize / 2 + 1, numTrials_in);  % Power for inclusion
    chanpowr_ex = zeros(cWinSize / 2 + 1, numTrials_ex);  % Power for exclusion

    % Create Hanning window
    hanningWindow = hanning(cWinSize);  % Hanning window of size cWinSize

    % Step 2: Compute the FFT for each epoch
    for trial = 1:numTrials_in
        % Apply Hanning window to the inclusion epoch
        windowedData_in = inclusion_data(:, :, trial) .* hanningWindow;

        % Compute FFT and power spectrum
        F_in = fft(windowedData_in);
        chanpowr_in(:, trial) = (2 * abs(F_in(1:cWinSize / 2 + 1)) / cWinSize).^2;  % Keep only positive frequencies
    end

    for trial = 1:numTrials_ex
        % Apply Hanning window to the exclusion epoch
        windowedData_ex = exclusion_data(:, :, trial) .* hanningWindow;

        % Compute FFT and power spectrum
        F_ex = fft(windowedData_ex);
        chanpowr_ex(:, trial) = (2 * abs(F_ex(1:cWinSize / 2 + 1)) / cWinSize).^2;  % Keep only positive frequencies
    end

    % Step 3: Average power over trials for each condition
    chanpowr_avg_in = mean(chanpowr_in, 2);  % Average over trials for inclusion
    chanpowr_avg_ex = mean(chanpowr_ex, 2);  % Average over trials for exclusion
    fs = linspace(0,srate./2,cWinSize./2+1);




    AllParticipants_power_in = cell(length(Subject_Index), length(cFreqs)); % Initialize before the loop

    % Now loop through freqs and extract the right numbers
    for freq = 1:length(cFreqs)
        % which to select and average?
        freqndx = fs >= cFreqs{freq}(1) & fs <= cFreqs{freq}(2); % find indices of frequencies that we are interested in
        AllParticipants_power_in{subj, freq} = mean(chanpowr_avg_in(freqndx, :)); % Average power for the frequency range
    end




     AllParticipants_power_ex = cell(length(Subject_Index), length(cFreqs)); % Initialize before the loop

    % Now loop through freqs and extract the right numbers
    for freq = 1:length(cFreqs)
        % which to select and average?
        freqndx = fs >= cFreqs{freq}(1) & fs <= cFreqs{freq}(2); % find indices of frequencies that we are interested in
        AllParticipants_power_ex{subj, freq} = mean(chanpowr_avg_ex(freqndx, :)); % Average power for the frequency range
    end



end




Thank you very much for your help in advance.


Kind regards,
Niloofar




----

Niloofar Goharbakhsh


PhD student in Developmental and Educational Psychology,
University of Bremen, Germany