% REJSTATEPOCH - reject bad eeg trials based a statistical measure. Can
% be applied either to the raw eeg data or the ICA
% component activity. This is an interactive function.
%
% Usage:
% >> [ Irej, Irejdetails, n, threshold, thresholdg] = ...
% rejstatepoch( signal, rej, 'key1', value1...);
%
% Inputs:
% signal - 3 dimensional data array channel x points x trials
% (instead of channels, one might also use independent
% components).
% rej - rejection array, one value per trial and per channel or
% component (size channel x trials).
% By default, values are normalized by this
% function and trehshold is expressed in term of standard
% deviation of the mean.
%
% Optional inputs:
% 'plot' - ['on'|'off'] interactive mode or just rejection.
% In the interactive mode, it plots the normalized
% entropy of original signal (blue) and the limits
% (red) (default:'on')
% 'threshold' - percentage error threshold (default 1-0.25/nb_trials)
% for individual trials of individual channel/component.
% This threshold is expressed in term of standard
% deviation from the mean (default 5).
% 'global' - ['on'|'off'], also perform threshold on the global
% measure (by default, the mean over all channel or
% electrodes).
% 'rejglob' - rejection array, one value per trials. Use this
% argument when the global measure for all channels or
% does not correspond to their mean.
% 'thresholdg' - global threshold for the reunion of all channels
% or components. By default, it is equal to 'threshold'.
% 'normalize' - ['on'|'off'], normalize values before applying the
% threshold. Default is 'on'.
% 'plotcom' - string command to plot single trials. Default
% is none.
% 'title' - title of the graph. Default is none.
% 'labels' - labels for electrodes (array not cell).
%
% Outputs:
% Irej - indexes of trials to be rejected
% Irejdetails - array for rejected components or channel (nb_rejected x
% nb_channel or nb_rejected x nb_components)
% n - number of trials rejected
% threshold - percentage error threshold
% thresholdg - percentage error threshold for global rejection
%
% See also: EEGLAB
% Algorithm:
% normalise the measure given as input and reject trials based on
% an uniform distribution of the data
% Copyright (C) 2001 Arnaud Delorme, Salk Institute, arno@salk.edu
%
% This file is part of EEGLAB, see http://www.eeglab.org
% for the documentation and details.
%
% Redistribution and use in source and binary forms, with or without
% modification, are permitted provided that the following conditions are met:
%
% 1. Redistributions of source code must retain the above copyright notice,
% this list of conditions and the following disclaimer.
%
% 2. Redistributions in binary form must reproduce the above copyright notice,
% this list of conditions and the following disclaimer in the documentation
% and/or other materials provided with the distribution.
%
% THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
% AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
% IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
% ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
% LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
% CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
% SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
% INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
% CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
% ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF
% THE POSSIBILITY OF SUCH DAMAGE.
% userdata
% gcf : plotsig - signal to plot in the pop_out window
% pnts - number of points per epoch
% Irej - rejection of trials
% Irejelec - rejection of trials and electrodes
% labels - labels for curves
% plotwin : rej1 - rejection array (with all electrodes x trials)
% rej2 - global rejection array (size trials)
% thr1 - threshold (electrodes)
% thr2 - threshold global
function [ Irej, Irejdetails, n, threshold, thresholdg] = rejstatepoch( signal, ...
rej, varargin); % pnts, th_E, th_rejg, command, commandplot, typetitle, E, rejg);
if nargin < 1
help rejstatepoch;
return;
end
if ~ischar( signal )
if nargin < 2
help rejstatepoch;
return;
end
if ~isempty( varargin ), g=struct(varargin{:});
else g= []; end
try, g.plot; catch, g.plot='on'; end
try, g.threshold; catch, g.threshold=5; end
try, g.thresholdg; catch, g.thresholdg=5; end
try, g.global; catch, g.global='on'; end
try, g.rejglob; catch, g.rejglob=[]; end
try, g.normalize; catch, g.normalize='on'; end
try, g.plotcom; catch, g.plotcom=''; end
try, g.title; catch, g.title=''; end
try, g.labels; catch, g.labels=''; end
g.rej = rej;
clear rej
switch lower(g.plot)
case {'on', 'off'} ;
otherwise disp('Error: Plot must be either ''on'' or ''off'''); return;
end;
switch lower(g.global)
case {'on', 'off'} ;
otherwise disp('Error: Global must be either ''on'' or ''off'''); return;
end;
switch lower(g.normalize)
case {'on', 'off'} ;
otherwise disp('Error: Normalize must be either ''on'' or ''off'''); return;
end;
if ~ischar(g.plotcom)
disp('Error: Plotcom must be a string to evaluate'); return;
end;
if ~ischar(g.title)
disp('Error: Title must be a string'); return;
end;
try, g.threshold*2;
catch, disp('Error: Threhsold must be a number'); return;
end;
try, g.thresholdg*2;
catch, disp('Error: Threhsoldg must be a number'); return;
end;
if length(g.threshold(:)) > 1
disp('Error: Threhsold must be a single number'); return;
end;
if length(g.thresholdg(:)) > 1
disp('Error: Threhsoldg must be a single number'); return;
end;
if ~isempty(g.rejglob)
if length(g.rejglob) ~= size(g.rej,2)
disp('Error: Rejglob must be have the same length as rej columns'); return;
end
else
switch lower(g.global), case 'on', g.rejglob = sum(g.rej,1); end
end;
if size(signal,3) ~= size(g.rej,2)
disp('Error: Signal must be have the same number of element in 3rd dimension as rej have columns'); return;
end
if isempty(g.labels)
for index = 1:size(g.rej,1)
g.labels(index,:) = sprintf('%3d', index);
end
if ~isempty(g.rejglob)
g.labels(index+2,:) = 'g. ';
end
end;
switch lower(g.normalize),
case 'on',
g.rej = (g.rej-mean(g.rej,2)*ones(1, size(g.rej,2)))./ (std(g.rej, 0, 2)*ones(1, size(g.rej,2)));
switch lower(g.global), case 'on', g.rejglob = (g.rejglob(:)-mean(g.rejglob(:)))./ std(g.rejglob(:)); end
end;
switch lower(g.global), case 'off',g.rejglob = []; end
% plot the buttons
% ----------------
try, icadefs; catch, GUIBUTTONCOLOR = [0.8 0.8 0.8]; BACKCOLOR = [0.8 0.8 0.8]; end;
figure('color', BACKCOLOR);
set(gcf, 'name', 'Rejectrials');
pos = get(gca,'position'); % plot relative to current axes
set( gca, 'tag', 'mainaxis');
q = [pos(1) pos(2) 0 0];
s = [pos(3) pos(4) pos(3) pos(4)]./100; % allow to use normalized position [0 100] for x and y
axis('off');
plotsig = sum(abs(signal(:,:)),1);
set(gcf, 'userdata', { plotsig, size(signal, 2), [], [] }); % the two last arguments are the rejection
% Create axis
% -----------
h6 = axes('Units','Normalized', 'tag', 'Plotwin', 'Position',[-10 12 120 84].*s+q);
title(g.title);
set( h6, 'userdata', { g.rej g.rejglob g.threshold g.thresholdg g.labels }); % g.rej was put twice because it is used to compute the global entropy
% CANCEL button
% -------------
h = uicontrol(gcf, 'Style', 'pushbutton', 'string', 'Cancel', 'Units','Normalized','Position',...
[-15 -6 15 6].*s+q, 'callback', 'close(gcf);', 'backgroundcolor', GUIBUTTONCOLOR);
% Entropy component text and threshold
% ------------------------------------
makebutton( 'Single-channel', 'Plotwin' , [3 -6 27 6].*s+q, [30 -6 15 6].*s+q, 3, g.threshold, GUIBUTTONCOLOR);
makebutton( 'All-channel', 'Plotwin' , [50 -6 27 6].*s+q, [77 -6 15 6].*s+q, 4, g.thresholdg, GUIBUTTONCOLOR);
% EEGPLOT button
% --------------
if ~isempty(g.plotcom)
h = uicontrol(gcf, 'backgroundcolor', GUIBUTTONCOLOR, 'Style', 'pushbutton', 'string', 'EEGPLOT', 'Units','Normalized','Position',[95 -2 15 6].*s+q, ...
'callback',['TMPEEG = get(gcbf, ''userdata'');' ...
'TMPREJ = TMPEEG{3};' ...
'TMPREJELEC = TMPEEG{4};' ...
g.plotcom ] );
posaccept = [95 -10 15 6];
else posaccept = [95 -6 15 6];
end;
% ACCEPT button
% -------------
command = 'fprintf(''Rejection indices has been put in the matrix TMPREJ\n'')';
haccept = uicontrol(gcf, 'backgroundcolor', GUIBUTTONCOLOR, 'Style', 'pushbutton', 'string', 'UPDATE', 'Units','Normalized','Position', posaccept.*s+q, ...
'callback', [ 'set(gcbo, ''userdata'', 1);' ... %signal to signify termination
'TMPEEG = get(gcbf, ''userdata'');' ...
'TMPREJ = TMPEEG{3};' ...
'TMPREJELEC = TMPEEG{4};' ...
command ] );
command = [ 'entwin = get(gcbf, ''currentobject'');' ...
'tmptype = get(entwin,''type'');' ...
'if tmptype(1) == ''l'' entwin = get(entwin, ''parent''); end;' ...
'tagwin = get(entwin,''tag'');' ... % either entropy or kurtosis
'switch tagwin,' ...
' case ''Plotwin'',' ... % check if the user clicked on the right window
' alldata = get(gcf, ''userdata'');' ...
' plotsig = alldata{1};' ...
' pnts = alldata{2};' ...
' fig = figure(''position'', [100 300 600 400],''color'', [1 1 1]);' ...
' I = alldata{3};' ...
' sweeps = size(plotsig,2) /pnts;' ...
' h1 = axes(''parent'', fig, ''Units'',''Normalized'', ''Position'',[0.6 0.11 0.27 0.815]);' ...
' pos = get(entwin, ''currentpoint'');' ...
' component = round(pos(1) / sweeps + 0.5);' ... % determine the component number
' alldata = get(entwin, ''userdata'');' ...
' rej = alldata{1};' ...
' if component <= size(rej,1)' ... % component
' rej_threshold = alldata{3};' ...
' component = max(1,min(component, size(rej,1))); ' ...
' rej = rej(component, :);' ...
' titlegraph = sprintf(''' g.title ' #%d'', component);' ...
' colorgraph = ''b'';' ...
' else' ... % global
' rej = alldata{2};' ...
' rej_threshold = alldata{4};' ...
' titlegraph = sprintf(''' g.title ' global'');' ...
' colorgraph = ''g'';' ...
' end;' ...
' plot([1:length(rej)], rej, ''color'', colorgraph);' ...
' ss = get(h1, ''xlim'');' ...
' set(h1, ''view'', [90 90]);' ...
' set(h1, ''xdir'', ''reverse'');' ...
' set(h1, ''XLim'', ss);' ...
' hold on;' ... % plot component
' yl = get(h1, ''ylim'');' ...
' set(h1, ''xtickmode'', ''manual'', ''xtick'', sweeps/2, ''xticklabel'', component, ''xlim'', [ 0 sweeps ]);' ...
' title( titlegraph );' ...
' plot( get(h1, ''xlim''), [rej_threshold rej_threshold], ''r'');' ... % plot limit
' plot( get(h1, ''xlim''), [-rej_threshold -rej_threshold], ''r'');' ... % plot limit
' set(h1, ''xticklabel'', []);' ...
' hold on;' ...
' h2 = axes(''parent'', fig,''Units'',''Normalized'', ''Position'',[0.13 0.11 0.27 0.815]);' ...
' erpimage( plotsig, ones(1,size(plotsig,2)/pnts), [0:1/(pnts-1):1], '''', 3, 1, ''nosort'', ''noplot'');' ...
' title(''Currentset all chans''); xlabel(''''); ylabel(''''); ' ...
' set(gca, ''xticklabel'', []);' ...
' hold on;' ...
' h3 = axes(''parent'', fig,''Units'',''Normalized'', ''Position'',[0.45 0.11 0.1 0.815]);' ...
' if any(I ~= 0)' ...
' rejImage = (I'' * ones(1, 10))'';' ...
' imagesc( rejImage'' );' ...
' set(gca, ''ydir'', ''normal'');' ...
' end;' ...
' title(''Rejected (all)''); xlabel(''''); ylabel('''');' ...
' set(gca, ''xticklabel'', [], ''yticklabel'', []);' ...
'end;' ...
'clear fig tmptype tagwin alldata rej rejImage plotsig sweeps pnts rej_threshold ss q s h1 h2 h3 pos component yl;' ];
% ' erpimage( rejImage(:)'', ones(1,size(I,2)), [0:1/(10-1):1], '''', 1, 0, ''nosort'', ''noplot'');' ...
set(gcf, 'WindowButtonDownFcn', command);
rejstatepoch('draw');
switch g.plot,
case 'on', waitfor( haccept, 'userdata'); drawnow;
end
threshold = g.threshold;
thresholdg = g.thresholdg;
Irej = [];
Irejdetails = [];
n = 0;
try
TMPEEG = get(gcf, 'userdata');
Irej = TMPEEG{3};
Irejdetails = TMPEEG{4};
n = length(find(Irej == 1));
plothandler = findobj( 'parent', gcf, 'tag', 'Plotwin');
TMPEEG = get(plothandler, 'userdata');
threshold = TMPEEG{3};
thresholdg = TMPEEG{4};
close(gcf);
catch, end
else %if signal is a string draw everything
% retrieve data
% -------------
gcfdata = get(gcf, 'userdata');
plotsig = gcfdata {1};
pnts = gcfdata {2};
sweeps = size(plotsig,2)/pnts;
h6 = findobj('parent', gcf, 'tag', 'Plotwin');
alldata = get(h6, 'userdata');
g.rej = alldata {1};
g.rejg = alldata {2};
g.threshold = alldata {3};
g.thresholdg = alldata {4};
set(h6, 'userdata', alldata);
nbchans = size(g.rej,1);
% reject trials
% -------------
rejelec = abs(g.rej) > g.threshold;
rej = max(rejelec,[],1);
n1 = sum(rej(:));
if ~isempty( g.rejg )
rej2 = abs(g.rejg) > g.thresholdg;
n2 = sum(rej2(:));
rej = rej | rej2(:)';
end
fprintf('%d trials rejected (single:%d, all:%d)\n', sum(rej(:)), n1, n2);
gcfdata {3} = rej;
gcfdata {4} = rejelec;
set(gcf, 'userdata', gcfdata);
% plot the sorted entropy curve
% -----------------------------
plotstat( 'Plotwin');
end
return;
function plotstat( id );
% plot the sorted entropy curve
% -----------------------------
h6 = findobj('parent', gcf, 'tag', id);
axes(h6); cla;
ttmp = get(gca, 'title');
oldtitle = get(ttmp, 'string');
% get datas
% ---------
alldata = get(gca, 'userdata');
g.rej = alldata {1};
g.rejg = alldata {2};
g.threshold = alldata {3};
g.thresholdg = alldata {4};
g.labels = alldata {5};
nbchans = size(g.rej,1);
sweeps = size(g.rej,2);
% plot datas
% ----------
g.rej = g.rej'; plot(g.rej(:)); g.rej = g.rej';
hold on;
yl = get(gca, 'ylim');
% plot vertival bars to separate components and the trehsold
% ----------------------------------------------------------
set( gca, 'tag', id, 'ylimmode', 'manual');
set(gca, 'xtickmode', 'manual', 'xtick', [0:sweeps:(size(g.rej(:),1)-1+2*sweeps)] + sweeps/2, ...
'xticklabel', g.labels, 'xlim', [ 0 (size(g.rej(:),1)-1+2*sweeps)]);
plot( [1 size(g.rej(:),1)], [-g.threshold -g.threshold], 'r'); % plot threshold
plot( [1 size(g.rej(:),1)], [g.threshold g.threshold], 'r'); % plot threshold
if ~isempty(g.rejg) % plot global ?
plot([size(g.rej(:),1)+sweeps:size(g.rej(:),1)+2*sweeps-1], g.rejg(:), 'g');
pp = patch([size(g.rej(:),1) size(g.rej(:),1) size(g.rej(:),1)+sweeps size(g.rej(:),1)+sweeps], [yl(1)-1 yl(2)+1 yl(2)+1 yl(1)-1], get(gcf, 'color'), 'clipping', 'off');
set(pp, 'EdgeColor', get(gcf, 'color'));
plot( [size(g.rej(:),1)+sweeps length(g.rejg)+size(g.rej(:),1)+sweeps], [-g.thresholdg -g.thresholdg], 'r'); % plot threshold
plot( [size(g.rej(:),1)+sweeps length(g.rejg)+size(g.rej(:),1)+sweeps], [g.thresholdg g.thresholdg], 'r'); % plot threshold
plot([size(g.rej(:),1)+sweeps size(g.rej(:),1)+sweeps], yl, 'k');
else
pp = patch([size(g.rej(:),1) size(g.rej(:),1) size(g.rej(:),1)+2*sweeps size(g.rej(:),1)+2*sweeps], [yl(1)-1 yl(2)+1 yl(2)+1 yl(1)-1], get(gcf, 'color'), 'clipping', 'off');
set(pp, 'EdgeColor', get(gcf, 'color'));
end
for index = 0:sweeps:size(g.rej(:),1);
plot([index index], yl, 'k');
end
% restore properties
title(oldtitle);
set(h6, 'userdata', alldata);
return;
function makebutton( string, ax, pos1, pos2, userindex, init, GUIBUTTONCOLOR );
h = uicontrol(gcf , 'backgroundcolor', GUIBUTTONCOLOR, 'Style', 'radiobutton', 'string', string, 'value', fastif(init == 0, 0, 1), 'Units','Normalized', 'Position', pos1, ...
'callback', [ 'textresh = findobj(''parent'', gcbf, ''tag'', ''' string ''');' ...
'checkstatus = get(gcbo, ''value'');' ...
'ax = findobj(''parent'', gcbf, ''tag'', ''' ax ''');' ...
'tmpdat = get(ax, ''userdata'');' ...
'if checkstatus' ... % change status of the textbox
' set(textresh, ''enable'', ''on'');' ...
' tmpdat{' int2str(userindex) '} = str2num(get(textresh, ''string''));' ...
'else' ...
' set(textresh, ''enable'', ''off'');' ...
' tmpdat{' int2str(userindex) '} = 0;' ...
'end;' ...
'set(ax, ''userdata'' , tmpdat);' ...
'rejstatepoch(''draw'');' ...
'clear tmpdat ax checkstatus textresh;' ] );
h = uicontrol(gcf, 'Style', 'edit', 'backgroundcolor', [1 1 1], 'tag', string, 'string', num2str(init), 'enable', fastif(init == 0, 'off', 'on'), 'Units','Normalized', 'Position', pos2, ...
'callback', [ 'ax = findobj(''parent'', gcbf, ''tag'', ''' ax ''');' ...
'tmpdat = get(ax, ''userdata'');' ...
'tmpdat{' int2str(userindex) '} = str2num(get(gcbo, ''string''));' ...
'set(ax, ''userdata'' , tmpdat);' ...
'rejstatepoch(''draw'');' ...
'clear tmpdat ax;' ] );
return;
Datasets
Models
