Time-domain analysis (ERP)

Studying how a system reacts to perturbations is a very effective way to learn about that system. This is also true in brain research. The perturbations in the background EEG activity produced by a given stimulation are called evoked potentials (EPs). When EEG voltages are measured with respect to an eliciting stimulus, a good strategy to visualize even small EPs, which are not easily discernable in the continuous EEG recording, is to take the average over a number of repetitions of that stimulus. The result is referred to as an event-related potential (ERP).

The logic underlying the computation of an ERP is quite straightforward. Each trial (i.e., repetition) contains both signal (i.e., evoked potential) and noise (i.e., background EEG activity). If we assume that the signal is similar over trials, whereas the noise is random in each trial, then averaging over many trials will lead to a dramatic noise reduction and a consequent improvement in the ability to detect the ERP.

The signal-to-noise ratio (SNR) is a typical measure of the capacity to distinguish signal from noise. It is simply calculated as the ratio of the signal amplitude to the noise amplitude. It can be shown that, in the computation of an ERP, the SNR increases as a function of the square root of the number of trials, provided the noise is Gaussian distributed. This means, for instance, that we need to increase the number of trials by a factor of four to double the SNR.

To recap, studying an ERP is rather simple: we first align the time-domain EEG signals to the time = 0 (i.e., we align to the onset of an eliciting event), and then we average across trials at each time point.

In ERP terminology, the term component denotes a particular feature of an ERP waveform and, in most cases, the features of interest are peaks or troughs. Therefore, ERP components are generally identified by their polarity (i.e., peak or trough) and their temporal latency. For instance, the N100 component is a negative modulation of the waveform that reaches its minimum at about 100 ms.

It is worth emphasizing that ERP components (i.e., the voltage peaks and troughs in the ERP waveform) do not necessarily correspond with transient electrical potentials elicited by the activity of localized populations of neurons in the brain. In fact, because of volume conduction, the ERP voltage is the sum of the activity of many populations of neurons. In other words, ERP components do not represent independent components (ICs), but rather can be interpreted as a combination of ICs.

A common method to elicit ERPs is by using the so-called oddball paradigm. In an oddball paradigm, two different stimuli are presented following a pseudo-random sequence. One of them (standard stimulus) appears frequently, whereas the other one (deviant stimulus) appears less often and unexpectedly. We will analyze some datasets in which standard and deviant stimuli are tones of different frequencies. To this end, we will use FieldTrip, a Matlab toolbox for MEG and EEG analysis. The main advantage of FieldTrip is that it offers advanced non-parametric statistical testing, as well as advanced methods for frequency and time-frequency analysis. The main disadvantage (for some users) is the lack of any interactive GUI, which prompts the users to develop their own scripts.

In general, an analysis implemented in FieldTrip consists of a sequence of calls to specific functions that produce a data object as an output. The input to these functions is always a configuration (cfg) structure, which specifies the parameters used by the function. Usually, the data is passed as a second input argument.

Let us start analyzing the first dataset (click here to download all datasets). The datasets were preprocessed using EEGLAB and following the steps illustrated in the previous section.

To begin with, let us have a look at the different trigger codes present in the dataset S01OddballSession1. To do this, let us run the following snippet of code.

Return to the wiki
home page

% Looking at the different trigger codes present in the dataset
cfg = [];
cfg.dataset = 'C:\Users\CIMeC\Desktop\wikiDatasets\S01\S01Oddball.vhdr'; 
                % file path
cfg.trialfun = 'ft_trialfun_general'; % the default trial function is fine
                                      % when trials are defined according
                                      % to specified trigger values
cfg.trialdef.eventtype = '?';
ft_definetrial(cfg);

↑ Go up

As show in Figure 3.1, all event types and values are displayed in the Matlab command window. The most relevant events are those marked as ‘Stimulus’. In our dataset, the values ‘S 103’, ‘S 104’, ‘S 105’, ‘S 106’, ‘S 107’, ‘S 108’, ‘S 109’, ‘S 113’, ‘S 118’ represent standard tones, whereas the values ‘S 203’, ‘S 204’, ‘S 205’, ‘S 206’, ‘S 207’, ‘S 208’, ‘S 209’, ‘S 213’, ‘S 218’ represent deviant tones.

Figure 3.1: Events detection

Each tone lasted 80 ms and was separated from the following tone by a gap of 700 ms. We can segment our dataset into epochs from -250 ms to 750 ms, where the time-point 0 is defined by the trigger associated with a given tone. To this end, we will first use the ft_definetrial function, and then we will use the ft_preprocessing function.

% Defining trials using trigger values – Standard tones
cfg = [];
cfg.dataset = 'C:\Users\CIMeC\Desktop\wikiDatasets\S01\S01Oddball.set';
                          % file path. You may need to edit the path
cfg.trialfun = 'ft_trialfun_general'; % the default trial function is fine
                                      % when trials are defined according
                                      % to specified trigger values
cfg.trialdef.eventtype = 'Stimulus';
cfg.trialdef.eventvalue = {'S 103','S 104','S 105','S 106', ...
    'S 107','S 108','S 109','S 113','S 118'}; % standard tones
cfg.trialdef.prestim = 0.25; % in seconds
cfg.trialdef.poststim = 0.75; % in seconds
cfg = ft_definetrial(cfg);

↑ Go up

The output resulting from ft_definetrial contains the trial definition as a Nx4 matrix, as shown by the variable cfg.trl. The first column of the matrix is the begin sample of each trial, the second column is the end sample, the third column contains the offset that defines where the relative t = 0 point is for each trial, and the fourth column contains the corresponding event value.

It is desirable to discard trials containing discontinuities that have been created after continuous artifactual segments have been removed.

% Removing trials containing discontinuities
trialsToReject = zeros(length(cfg.trl),1);
for eventId = 1:length(cfg.event)
    if isempty(cfg.event(eventId).type)
        tmpSample = cfg.event(eventId).sample;
        for trialInfoId = 1:length(cfg.trl)
            if (tmpSample - cfg.trl(trialInfoId,1)) * ...
                    (tmpSample - cfg.trl(trialInfoId,2)) < 0
                trialsToReject(trialInfoId) = trialInfoId;
            end
        end
    end
endtrialsToReject = trialsToReject(trialsToReject ~= 0);
cfg.trl(trialsToReject,:) = [];

↑ Go up

We will use the cfg configuration to read the previously defined trials from the original data file on disk. We will also apply a baseline correction.

% Reading trials and applying baseline correction – Standard tones
cfg.demean = 'yes';
cfg.baselinewindow = [-0.1 0]; % in seconds
dataStd = ft_preprocessing(cfg);

 ↑ Go up

dataStd is a structure variable whose most important fields are trial, which contains the individual trials, and time, which contains the time vector for each trials.

To visualize a single trial (e.g., trial 1) from one channel (e.g., channel 5), run the following command:

% Plotting a single trial from one channel
figure;
plot(dataStd.time{1}, dataStd.trial{1}(5,:)) % trial 1 from channel 5 (Fz)

↑ Go up

The previous steps will be repeated for the second experimental condition (i.e., deviant tones).

% Defining trials using trigger values – Deviant tones
cfg = [];
cfg.dataset = 'C:\Users\CIMeC\Desktop\wikiDatasets\S01\S01Oddball.set'; 
                % file path
cfg.trialfun = 'ft_trialfun_general'; % the default trial function is fine
                                      % when trials are defined according
                                      % to specified trigger values
cfg.trialdef.eventtype = 'Stimulus';
cfg.trialdef.eventvalue = {'S 203','S 204','S 205','S 206','S 207',...
    'S 208','S 209','S 213','S 218'}; % deviant tones
cfg.trialdef.prestim = 0.25; % in seconds
cfg.trialdef.poststim = 0.75; % in seconds
cfg = ft_definetrial(cfg);
 
% Removing trials containing discontinuities
trialsToReject = zeros(length(cfg.trl),1);
for eventId = 1:length(cfg.event)
    if isempty(cfg.event(eventId).type)
        tmpSample = cfg.event(eventId).sample;
        for trialInfoId = 1:length(cfg.trl)
            if (tmpSample - cfg.trl(trialInfoId,1)) * ...
                    (tmpSample - cfg.trl(trialInfoId,2)) < 0
                trialsToReject(trialInfoId) = trialInfoId;
            end
        end
    end
end
    
trialsToReject = trialsToReject(trialsToReject ~= 0);
cfg.trl(trialsToReject,:) = [];

% Reading trials and applying baseline correction – Deviant tones
cfg.demean = 'yes'; 
cfg.baselinewindow = [-0.1 0]; % in seconds
dataDev = ft_preprocessing(cfg);

↑ Go up

A convenient way to apply the previous analyses to all subjects is to use a for-loop, as shown below:

%% Initializations

dataPath = 'C:\Users\CIMeC\Desktop\wikiDatasets\'; % It depends on where
                                                   % the wikiDatasets folder
                                                   % is stored
fileNameDescription = 'Oddball.set';
 
subjects = {'S01' 'S02' 'S03' ...
            'S04' 'S05'}; % the variable subjects contains the list
                          % of subject-specific directories
nSubjs = numel(subjects);

 
%% Defining trials (or epochs) using trigger values
 
for i = 1:nSubjs
    fileName = strcat(subjects{i},fileNameDescription);
    dataSet = fullfile(dataPath,subjects{i},fileName);
 
    % Defining trials using trigger values - Standard tones
    cfg = [];
    cfg.dataset = dataSet; % file path
    cfg.trialfun = 'ft_trialfun_general'; % the default trial function is fine
                                          % when trials are defined according
                                          % to specified trigger values
    cfg.trialdef.eventtype = 'Stimulus';
    cfg.trialdef.eventvalue = {'S 103','S 104','S 105','S 106','S 107',...
        'S 108','S 109','S 113','S 118'}; % standard tones
    cfg.trialdef.prestim = 0.25; % in seconds
    cfg.trialdef.poststim = 0.75; % in seconds
    cfg = ft_definetrial(cfg);
    
    % Removing trials containing discontinuities
    trialsToReject = zeros(length(cfg.trl),1);
    for eventId = 1:length(cfg.event)
        if isempty(cfg.event(eventId).type)
            tmpSample = cfg.event(eventId).sample;
            for trialInfoId = 1:length(cfg.trl)
                if (tmpSample - cfg.trl(trialInfoId,1)) * ...
                        (tmpSample - cfg.trl(trialInfoId,2)) < 0
                    trialsToReject(trialInfoId) = trialInfoId;
                end
            end
        end
    end
    
    trialsToReject = trialsToReject(trialsToReject ~= 0);
    cfg.trl(trialsToReject,:) = [];
    
    % Reading trials and applying baseline correction - Standard tones
    cfg.demean = 'yes';
    cfg.baselinewindow = [-0.1 0]; % in seconds
    cfg.outputfile = fullfile(dataPath,subjects{i},'dataStd');
    ft_preprocessing(cfg);
    
    % Defining trials using trigger values - Deviant tones
    cfg = [];
    cfg.dataset = dataSet; % file path
    cfg.trialfun = 'ft_trialfun_general'; % the default trial function is fine
                                          % when trials are defined according
                                          % to specified trigger values
    cfg.trialdef.eventtype = 'Stimulus';
    cfg.trialdef.eventvalue = {'S 203','S 204','S 205','S 206','S 207',...
        'S 208','S 209','S 213','S 218'}; % deviant tones
    cfg.trialdef.prestim = 0.25; % in seconds
    cfg.trialdef.poststim = 0.75; % in seconds
    cfg = ft_definetrial(cfg);
    
    % Removing trials containing discontinuities
    trialsToReject = zeros(length(cfg.trl),1);
    for eventId = 1:length(cfg.event)
        if isempty(cfg.event(eventId).type)
            tmpSample = cfg.event(eventId).sample;
            for trialInfoId = 1:length(cfg.trl)
                if (tmpSample - cfg.trl(trialInfoId,1)) * ...
                        (tmpSample - cfg.trl(trialInfoId,2)) < 0
                    trialsToReject(trialInfoId) = trialInfoId;
                end
            end
        end
    end
    
    trialsToReject = trialsToReject(trialsToReject ~= 0);
    cfg.trl(trialsToReject,:) = [];
    
    % Reading trials and applying baseline correction - Deviant tones
    cfg.demean = 'yes';
    cfg.baselinewindow = [-0.1 0]; % in seconds
    cfg.outputfile = fullfile(dataPath,subjects{i},'dataDev');
    ft_preprocessing(cfg);
end 

↑ Go up

To compute the ERPs for each subject, we will use the ft_timelockanalysis function, which averages the epoched time series.

%% Computing ERP for each subject

% Standard tones
for i = 1:nSubjs
    fileName = 'dataStd.mat';
    dataSet = fullfile(dataPath,subjects{i},fileName);
    load(dataSet); % loading the variable 'data'
    
    % Averaging across trials
    cfg = [];
    cfg.channel = {'all'};
    cfg.outputfile = fullfile(dataPath,subjects{i},'erpStd');
    ft_timelockanalysis(cfg,data);
    clear data
end
 
% Deviant tones
for i = 1:nSubjs
    fileName = 'dataDev.mat';
    dataSet = fullfile(dataPath,subjects{i},fileName);
    load(dataSet); % loading the variable 'data'
    
    % Averaging across trials
    cfg = [];
    cfg.channel = {'all'};
    cfg.outputfile = fullfile(dataPath,subjects{i},'erpDev');
    ft_timelockanalysis(cfg,data);
    clear data
end

↑ Go up

The ft_timelockgrandaverage allows computing the grand mean ERP, including all individual subjects’ ERP as an input argument.

%% Computing group average ERP

% Standard tones
% Loading single subject ERPs in a single cell array
allErpStd = cell(1,nSubjs);
for i = 1:nSubjs
    fileName = 'erpStd.mat';
    dataSet = fullfile(dataPath,subjects{i},fileName);
    load(dataSet); % loading the variable 'timelock'
    allErpStd{i} = timelock;
    clear timelock
end
 
% Computing grand average across subjects
cfg = [];
cfg.latency = [-0.1 0.6];
cfg.keepindividual = 'yes'; % it is sometimes useful to keep individual... 
                            % representations in the output
cfg.outputfile = fullfile(dataPath,'grandAvgErpStd');
ft_timelockgrandaverage(cfg,allErpStd{:});
 
% Deviant tones
% Loading single subject ERPs in a single cell array
allErpDev = cell(1,nSubjs);
for i = 1:nSubjs
    fileName = 'erpDev.mat';
    dataSet = fullfile(dataPath,subjects{i},fileName);
    load(dataSet); % loading the variable 'timelock'
    allErpDev{i} = timelock;
    clear timelock
end
 
% Computing grand average across subjects
cfg = [];
cfg.latency = [-0.1 0.6];
cfg.keepindividual = 'yes'; % it is sometimes useful to keep individual...
                            % representations in the output
cfg.outputfile = fullfile(dataPath,'grandAvgErpDev');
ft_timelockgrandaverage(cfg,allErpDev{:});

↑ Go up

FieldTrip offers several ways to visualize the spatiotemporal structure of EEG data. The ft_multiplotER allows displaying a set of ERP time courses on a channel layout (see below for further information about the channel layout); the ft_singleplotER allows displaying a single ERP time course or an average time course across a set of specified channels; the ft_topoplotER allows displaying a spatial topography in a specified latency window. All these functions require a layout representing a 2-dimensional projection of the electrode locations.

We will define a channel layout by using one of the template layouts included in FieldTrip. In the present tutorial, we will adopt the general template EEG1005.lay and then we will increase the layout width and height.

%% Defining channel layout

templateLayout = 'EEG1005.lay'; % one of the template layouts included
                                % in FieldTrip
cfg = [];
cfg.layout = which(templateLayout); 
layout = ft_prepare_layout(cfg);
 
figure
ft_plot_lay(layout);
title(templateLayout)

% Increasing layout width and height
layout.width = 0.0615 * ones(length(layout.width),1);
layout.height = 0.0356 * ones(length(layout.height),1);

↑ Go up

Now we are ready to plot the grand average ERP.

%% Visualization of the group average ERP

% Plotting a channel layout with the ERP time courses of standard and...
% deviant tones at each electrode
fileName = 'grandAvgErpStd';
dataSet = fullfile(dataPath,fileName);
standard = load(dataSet); % loading the variable 'grandavg'
 
fileName = 'grandAvgErpDev';
dataSet = fullfile(dataPath,fileName);
deviant = load(dataSet); % loading the variable 'grandavg'
 
cfg = [];
cfg.layout = layout;
cfg.graphcolor = 'br';
cfg.showcomment = 'no';
cfg.showscale = 'no';
%cfg.showlabels = 'yes';
figure;
ft_multiplotER(cfg,standard.grandavg,deviant.grandavg);

% Plotting mean ERPs over fronto-central electrodes
cfg = [];
%cfg.ylim = [-1 1];
cfg.channel = {'Fz','FCz','Cz'};
cfg.linewidth = 2;
cfg.graphcolor = 'br';
figure;
ft_singleplotER(cfg,standard.grandavg,deviant.grandavg);
legend({'Standard','Deviant'})
 
% Matlab commands to make the figure intelligible
set(gca,'Fontsize',14);
%title('Mean over central electrodes');
set(gca,'box','on');
xlabel('time (s)');
ylabel('amplitude (\muV)');
 
% Plotting mean ERP of deviant tones over fronto-central electrodes...
% toghether with standard deviation 
cfg = [];
cfg.channel = {'Fz','FCz','Cz'};
cfg.avgoverchan = 'yes'; % averaging over channels
grandavgSel = ft_selectdata(cfg,deviant.grandavg);
x = grandavgSel.time'; % defining x axis
y = mean(squeeze(grandavgSel.individual),1)'; % defining y axis
e = std(squeeze(grandavgSel.individual),1)'; % standard deviation 
low = y - e; % lower bound
high = y + e; % upper bound
figure;
hp = patch([x; x(end:-1:1); x(1)], [low; high(end:-1:1); low(1)], 'r');
hold on;
hl = line(x,y);
set(hp,'facecolor',[1 0.8 0.8],'edgecolor','none');
set(hl,'color','r','linewidth',2);
set(gca,'FontSize',12);
title ('Mean and Standard deviation over fronto-central electrodes');
set(gca,'box','on');
xlabel('time (s)');
ylabel('amplitude (\muV)');
legend({'Deviant (StD)','Deviant (Mean)'})
 
% Plotting topography of the N100 component
cfg = [];
cfg.layout = layout;
cfg.xlim = [0.1 0.12];
cfg.zlim = [-2 2];
cfg.colormap = parula;
cfg.marker = 'off';
cfg.style = 'fill';
cfg.comment = ' ';
cfg.colorbar = 'yes';
 
figure;
subplot(2,1,1);ft_topoplotER(cfg,standard.grandavg);
title('Standard');
c = colorbar;
c.LineWidth = 1;
c.FontSize = 14;
title(c,'\muV');
 
subplot(2,1,2);ft_topoplotER(cfg,deviant.grandavg);
title('Deviant');
c = colorbar;
c.LineWidth = 1;
c.FontSize = 14;
title(c,'\muV');
suptitle('N100 topography');

↑ Go up

As for statistical analysis, FieldTrip offers both parametric and non-parametric tests. To address the family-wise error rate, we suggest using non-parametric tests with spatiotemporal clustering. The general idea is that EEG signals have specific spatiotemporal structure that can be exploited to maximize statistical sensitivity. Because neighboring channels and neighboring time-points are likely to reflect the same neural phenomena, aggregating spatiotemporal neighbors into clusters is equivalent to aggregating the evidence for an effect to be significant. For further details on cluster-based statistical testing, we refer the reader to Maris and Oostenveld's paper (see References below).

We recommend reading the aforementioned reference before continuing with the tutorial.

To be able to define clusters across channels, we need a neighborhood structure.

%% Creating a neighbourhood structure to define how spatial clusters...
% are formed
fileName = 'dataStd.mat';
dataSet = fullfile(dataPath,subjects{1},fileName);
load(dataSet); % loading the variable 'data'
 
templateElec = 'standard_1005.elc';
elec = ft_read_sens(which(templateElec));
 
idx = ismember(elec.label,data.label);
elec.chanpos = elec.chanpos(idx,:);
elec.chantype = elec.chantype(idx);
elec.chanunit = elec.chanunit(idx);
elec.elecpos = elec.elecpos(idx,:);
elec.label = elec.label(idx);
 
cfg = [];
cfg.method = 'distance'; %'template' or 'distance' or 'triangulation'
cfg.neighbourdist = 80; % maximum distance between neighbouring sensors...
                        % (only for 'distance')
cfg.feedback = 'yes'; % show a neighbour plot
neighbours = ft_prepare_neighbours(cfg,elec);

↑ Go up

Now we can run the statistical analysis at the group level.

%% Statistical analysis of time-domain data at the group level

% Creating a FieldTrip design matrix
design = zeros(1,2*nSubjs); % for within-subjects analysis,...
                            % the design matrix contains two rows
for i = 1:nSubjs
  design(1,i) = i;
end
for i = 1:nSubjs
  design(1,nSubjs+i) = i;
end
design(2,1:nSubjs) = 1;
design(2,nSubjs+1:2*nSubjs) = 2;
 
% Computing the difference as a T-statistic and running an inferential test
cfg = [];
cfg.channel = {'all'}; 
%cfg.latency = [0.05 0.15];
cfg.method = 'montecarlo';
cfg.statistic = 'ft_statfun_depsamplesT';
cfg.clusterthreshold = 'nonparametric_common';
cfg.correctm = 'cluster';
cfg.alpha = 0.025; % 0.05/2 because two-sided
cfg.numrandomization = 1000;
cfg.neighbours = neighbours;
cfg.design = design;
cfg.uvar = 1;
cfg.ivar = 2;
statErp = ft_timelockstatistics(cfg,allErpStd{:},allErpDev{:});
 
% Plotting stat output
cfg = [];
cfg.layout = layout;
cfg.parameter = 'stat';
cfg.maskparameter = 'mask';
cfg.graphcolor = 'r';
cfg.showcomment = 'no';
cfg.showscale = 'no';
figure;
ft_multiplotER(cfg,statErp);
 
% Plotting stat difference ERP over fronto-central electrodes
cfg = [];
cfg.ylim = [-6 16];
cfg.channel = {'Fz','FCz','Cz'};
cfg.parameter = 'stat';
cfg.maskparameter = 'mask';
cfg.linewidth = 2;
cfg.graphcolor = 'r';
figure;
ft_singleplotER(cfg,statErp);
set(gca,'Fontsize',14);
title ('Mean over fronto-central electrodes');
set(gca,'box','on');
xlabel('time (s)');
ylabel('std vs. dev (t-val)');
 
% Plotting N100 topography in units of t-value
cfg = [];
cfg.layout = layout;
cfg.xlim = [0.1 0.12];
cfg.parameter = 'stat';
cfg.zlim = [-5 5];
cfg.colormap = parula;
cfg.marker = 'off';
cfg.style = 'fill';
cfg.comment = ' ';
cfg.colorbar = 'yes';
figure;
ft_topoplotER(cfg,statErp);
title('N100 topography')
c = colorbar;
c.LineWidth = 1;
c.FontSize = 14;
title(c,'t-val')

↑ Go up

It is also possible to analyze single-subject data. Here the aim is to test whether there is a systematic difference between trials in different experimental conditions.

%% Single-subject analysis of time-domain data

subjId = 1; % subject that we will analyze

% Loading the trials with standard tones
fileName = 'dataStd.mat';
dataSet = fullfile(dataPath,subjects{subjId},fileName);
load(dataSet); % loading the variable 'data'
 
% Averaging across trials
cfg = [];
cfg.channel = {'all'};
cfg.keeptrials = 'yes';
timeLockStd = ft_timelockanalysis(cfg,data);

% Loading the trials with deviant tones 
fileName = 'dataDev.mat';
dataSet = fullfile(dataPath,subjects{subjId},fileName);
load(dataSet); % loading the variable 'data'
 
% Averaging across trials
cfg = [];
cfg.channel = {'all'};
cfg.keeptrials = 'yes';
timeLockDev = ft_timelockanalysis(cfg,data);
 
% Creating a FieldTrip design matrix
design = zeros(1,size(timeLockStd.trial,1) + size(timeLockDev.trial,1)); % for between-trials analysis, the design matrix contains only one row
design(1,1:size(timeLockStd.trial,1)) = 1;
design(1,(size(timeLockStd.trial,1)+1):(size(timeLockStd.trial,1) + size(timeLockDev.trial,1)))= 2;
 
% Computing the difference as a T-statistic and running an inferential test
cfg = [];
cfg.channel = {'all'}; 
%cfg.latency = [0.05 0.15];
cfg.method = 'montecarlo';
cfg.statistic = 'ft_statfun_indepsamplesT'; % independent samples T-statistic for a between-trials analysis
cfg.clusterthreshold = 'nonparametric_common';
cfg.correctm = 'cluster';
cfg.alpha = 0.025; % 0.05/2 because two-sided
cfg.numrandomization = 1000;
cfg.neighbours = neighbours;
cfg.design = design;
cfg.ivar = 1;
statErp2 = ft_timelockstatistics(cfg,timeLockStd,timeLockDev);
 
% Plotting stat output
cfg = [];
cfg.layout = layout;
cfg.parameter = 'stat';
cfg.maskparameter = 'mask';
cfg.graphcolor = 'r';
cfg.showcomment = 'no';
cfg.showscale = 'no';
figure;
ft_multiplotER(cfg,statErp2);

↑ Go up

It is noteworthy to point out that the number of trials should be similar in each condition when analyzing the effects of different factors. If there are large differences in trial count (e.g., 93 trials in condition 1 and 28 trials in condition 2), potential biases may be introduced.

Finally, in FieldTrip we can analyze each experimental condition individually, by comparing the early ERP component with the pre-stimulus baseline window. In doing this, it is important that the length of pre-stimulus and post-stimulus windows be matched. By means of the ft_selectdata function, we split the single participant ERP data into two data structures: ‘act’ (i.e., activity) and ‘bl’ (i.e., baseline).

%% Group analysis of time-domain data (activity vs. baseline)

% Splitting data into activity and baseline
act = cell(1,nSubjs);
bl = cell(1,nSubjs);
for i = 1:nSubjs
    cfg = [];
    cfg.latency = [-0.1 0];
    bl{i} = ft_selectdata(cfg,allErpStd{i});
    cfg.latency = [0.05 0.15];
    act{i} = ft_selectdata(cfg,allErpStd{i});
    bl{i}.time = act{i}.time; % FieldTrip allows for a direct comparison only if the time-axes are equal
end
 
% Creating a FieldTrip design matrix, conditions: 1 = act, 2 = bl
design = zeros(2,2*nSubjs);
for i = 1:nSubjs
  design(1,i) = i;
end
for i = 1:nSubjs
  design(1,nSubjs+i) = i;
end
design(2,1:nSubjs) = 1;
design(2,nSubjs+1:2*nSubjs) = 2;
 
% Computing the difference as a T-statistic, and running an inferential test
cfg = [];
cfg.channel = {'all'};
cfg.method = 'montecarlo';
cfg.statistic = 'ft_statfun_depsamplesT'; % dependent samples T-statistics
cfg.clusterthreshold = 'nonparametric_common';
cfg.correctm = 'cluster';
cfg.alpha = 0.025;% .05/2 because two-sided
cfg.numrandomization = 1000;
cfg.neighbours = neighbours;
cfg.design = design;
cfg.uvar = 1;
cfg.ivar = 2;
statErp3 = ft_timelockstatistics(cfg,act{:},bl{:});
 
% Plotting stat output
cfg = [];
cfg.layout = layout;
cfg.parameter = 'stat';
cfg.maskparameter = 'mask';
cfg.graphcolor = 'r';
cfg.showcomment = 'no';
cfg.showscale = 'no';
figure;
ft_multiplotER(cfg,statErp3);

↑ Go up

References

• Maris, E, & Oostenveld, R. (2007). Nonparametric statistical testing of EEG- and MEG-data. J. Neurosci. Methods 164, 177-190.

↑ Go up