sig_wave() - Illustrates results of a mass univariate t-test on a
butterfly plot of grand average ERPs, difference waves, or
spectra in voltage or t-scores. Time windows used in the test
are represented with vertical dashed lines. If plotting the
waves in units of t-scores, the critical t-score for
significance is represented with a horizontal dashed line.
Note that you can click on a wave to see the name of the
electrode where the wave was measured and the precise time
corresponding to that data point. This information will
appear in a small box. Click on the box again to make it
disappear.
Usage:
>> sig_wave(GND_GRP_specGND_or_fname,test_id,varargin);
Required Inputs:
GND_GRP_specGND_or_fname - A GND/GRP/specGND structure variable or the
filename of such a variable that has been saved to
disk. To create a GND variable from Kutaslab ERP
files (e.g., *.nrm files) use avgs2GND.m. To do the
same from EEGLAB *.set files use sets2GND.m.
To create a GRP structure use GNDs2GRP.m. To create
a specGND variable use sets2specGND.m. See Mass
Univariate ERP Toolbox documentation for detailed
information about the format of GND and GRP
variables. If you specifiy a filename be sure to
include the file's path, unless the file is in the
current working directory.
test_id - [integer] The index # of the t-test results
stored in the GND/GRP/specGND variable that you wish
to visualize. To see what test results are
available, look at the t_tests field of the variable
(e.g., GND.t_tests) or use headinfo.m or
headinfo_spec.m.
Optional Inputs:
x_ticks - [vector] The times/frequencies you want marked with ticks
on the x-axis. Note, because the EEG has been sampled
at discrete time points, not all times between the
minimum and maximum analyzed times/frequencies can be
used. When a tick is requested that is not available,
the closest possible value is used. If not
specified, x_ticks are automatically chosen based
on the time/frequency range covered by the diagram.
This option has no effect if the t-test was
executed based on mean voltages/power within specified
time windows/frequency bands.
wave_ticks - [vector] The waveform values you want marked with ticks
on the y-axis. {default: MATLAB default ticks}
fig_id - [integer] The index # of the MATLAB figure in which
the diagram will be produced. Useful for overwriting
old figures. {default: lowest unused index}
units - ['uV','dB', or 't'] If 'uV' topographies will be
visualized in units of microvolts. If 'dB' spectral
power will be shown in decibels (10*log10((uV^2)/Hz)).
If 't', they will be shown in t-scores. {default:
't'}
use_color - [integer] If non-zero, non-gray scale colors will be
used in the diagram. {default: 0}
x_limits - [min max] Limits for the x-axis in units of ms for
ERPs/difference waves (e.g., [-100 500]) and Hz for
power spectra (e.g., [9 13]). {default: first and
last time point or frequency}
wave_limits - [min max] Limits for the y-axis in units of t-scores
or uV (e.g., [-10 15]. {default: slightly bigger than
min and max of waveform}
ydir - [1 or -1] If 1, positive is plotted up for ERPs and
difference waves. Otherwise, positive is represented
as down. This argument has no effect if plotting
spectra. {default: -1}
test_line_width - [integer] Width of the lines indicating test window
boundaries and critical t-scores. {default: 2}
title - [string] The axis title. {default: descriptor for the
bin being visualized}
crit_val_box - ['yes' or 'no'] If 'yes', a box indicating the alpha
or q level for significance will be created next to
the horizontal line indicating the critical t-score
for significance.
verblevel - An integer specifiying the amount of information you want
this function to provide about what it is doing during runtime.
Options are:
0 - quiet, only show errors, warnings, and EEGLAB reports
1 - stuff anyone should probably know
2 - stuff you should know the first time you start working
with a data set {default value}
3 - stuff that might help you debug (show all
reports)
Global Variable:
VERBLEVEL - Mass Univariate ERP Toolbox level of verbosity (i.e., tells functions
how much to report about what they're doing during
runtime) set by the optional function argument 'verblevel'
Notes:
-The printed/exported figure will have the same dimensions as the figure
on the display. Thus you can undo figure clutter by re-sizing it.
Author:
David Groppe
Kutaslab, 12/20100001 % sig_wave() - Illustrates results of a mass univariate t-test on a 0002 % butterfly plot of grand average ERPs, difference waves, or 0003 % spectra in voltage or t-scores. Time windows used in the test 0004 % are represented with vertical dashed lines. If plotting the 0005 % waves in units of t-scores, the critical t-score for 0006 % significance is represented with a horizontal dashed line. 0007 % Note that you can click on a wave to see the name of the 0008 % electrode where the wave was measured and the precise time 0009 % corresponding to that data point. This information will 0010 % appear in a small box. Click on the box again to make it 0011 % disappear. 0012 % 0013 % Usage: 0014 % >> sig_wave(GND_GRP_specGND_or_fname,test_id,varargin); 0015 % 0016 % Required Inputs: 0017 % GND_GRP_specGND_or_fname - A GND/GRP/specGND structure variable or the 0018 % filename of such a variable that has been saved to 0019 % disk. To create a GND variable from Kutaslab ERP 0020 % files (e.g., *.nrm files) use avgs2GND.m. To do the 0021 % same from EEGLAB *.set files use sets2GND.m. 0022 % To create a GRP structure use GNDs2GRP.m. To create 0023 % a specGND variable use sets2specGND.m. See Mass 0024 % Univariate ERP Toolbox documentation for detailed 0025 % information about the format of GND and GRP 0026 % variables. If you specifiy a filename be sure to 0027 % include the file's path, unless the file is in the 0028 % current working directory. 0029 % test_id - [integer] The index # of the t-test results 0030 % stored in the GND/GRP/specGND variable that you wish 0031 % to visualize. To see what test results are 0032 % available, look at the t_tests field of the variable 0033 % (e.g., GND.t_tests) or use headinfo.m or 0034 % headinfo_spec.m. 0035 % 0036 % Optional Inputs: 0037 % x_ticks - [vector] The times/frequencies you want marked with ticks 0038 % on the x-axis. Note, because the EEG has been sampled 0039 % at discrete time points, not all times between the 0040 % minimum and maximum analyzed times/frequencies can be 0041 % used. When a tick is requested that is not available, 0042 % the closest possible value is used. If not 0043 % specified, x_ticks are automatically chosen based 0044 % on the time/frequency range covered by the diagram. 0045 % This option has no effect if the t-test was 0046 % executed based on mean voltages/power within specified 0047 % time windows/frequency bands. 0048 % wave_ticks - [vector] The waveform values you want marked with ticks 0049 % on the y-axis. {default: MATLAB default ticks} 0050 % fig_id - [integer] The index # of the MATLAB figure in which 0051 % the diagram will be produced. Useful for overwriting 0052 % old figures. {default: lowest unused index} 0053 % units - ['uV','dB', or 't'] If 'uV' topographies will be 0054 % visualized in units of microvolts. If 'dB' spectral 0055 % power will be shown in decibels (10*log10((uV^2)/Hz)). 0056 % If 't', they will be shown in t-scores. {default: 0057 % 't'} 0058 % use_color - [integer] If non-zero, non-gray scale colors will be 0059 % used in the diagram. {default: 0} 0060 % x_limits - [min max] Limits for the x-axis in units of ms for 0061 % ERPs/difference waves (e.g., [-100 500]) and Hz for 0062 % power spectra (e.g., [9 13]). {default: first and 0063 % last time point or frequency} 0064 % wave_limits - [min max] Limits for the y-axis in units of t-scores 0065 % or uV (e.g., [-10 15]. {default: slightly bigger than 0066 % min and max of waveform} 0067 % ydir - [1 or -1] If 1, positive is plotted up for ERPs and 0068 % difference waves. Otherwise, positive is represented 0069 % as down. This argument has no effect if plotting 0070 % spectra. {default: -1} 0071 % test_line_width - [integer] Width of the lines indicating test window 0072 % boundaries and critical t-scores. {default: 2} 0073 % title - [string] The axis title. {default: descriptor for the 0074 % bin being visualized} 0075 % crit_val_box - ['yes' or 'no'] If 'yes', a box indicating the alpha 0076 % or q level for significance will be created next to 0077 % the horizontal line indicating the critical t-score 0078 % for significance. 0079 % verblevel - An integer specifiying the amount of information you want 0080 % this function to provide about what it is doing during runtime. 0081 % Options are: 0082 % 0 - quiet, only show errors, warnings, and EEGLAB reports 0083 % 1 - stuff anyone should probably know 0084 % 2 - stuff you should know the first time you start working 0085 % with a data set {default value} 0086 % 3 - stuff that might help you debug (show all 0087 % reports) 0088 % 0089 % Global Variable: 0090 % VERBLEVEL - Mass Univariate ERP Toolbox level of verbosity (i.e., tells functions 0091 % how much to report about what they're doing during 0092 % runtime) set by the optional function argument 'verblevel' 0093 % 0094 % Notes: 0095 % -The printed/exported figure will have the same dimensions as the figure 0096 % on the display. Thus you can undo figure clutter by re-sizing it. 0097 % 0098 % 0099 % Author: 0100 % David Groppe 0101 % Kutaslab, 12/2010 0102 0103 %%%%%%%%%%%%%%%% REVISION LOG %%%%%%%%%%%%%%%%% 0104 % 12/15/2010-Compatible with specGND variables now too 0105 % 0106 % 6/4/2011-Throws error when user attemps to visualize cluster-based tests 0107 % 0108 0109 %%%%%%%%%%%%%%%% FUTURE WORK %%%%%%%%%%%%%%%%% 0110 % 0111 0112 0113 function sig_wave(GND_GRP_specGND_or_fname,test_id,varargin) 0114 0115 p=inputParser; 0116 p.addRequired('GND_GRP_specGND_or_fname',@(x) isstruct(x) || ischar(x)); 0117 p.addRequired('test_id',@(x) isnumeric(x) && (length(x)==1)); 0118 p.addParamValue('x_ticks',[],@isnumeric); 0119 p.addParamValue('wave_ticks',[],@isnumeric); 0120 p.addParamValue('fig_id',[],@(x) isnumeric(x) && (length(x)==1)); 0121 p.addParamValue('use_color',1,@(x) isnumeric(x) || ischar(x)); 0122 p.addParamValue('verblevel',[],@(x) isnumeric(x) && (length(x)==1)); 0123 p.addParamValue('x_limits',[],@(x) isnumeric(x) && (length(x)==2)); 0124 p.addParamValue('wave_limits',[],@(x) isnumeric(x) && (length(x)==2)); 0125 p.addParamValue('ydir',-1,@(x) isnumeric(x) && (length(x)==1)); 0126 p.addParamValue('test_line_width',2,@(x) isnumeric(x) && (length(x)==1) && x>0); 0127 p.addParamValue('title',[],@(x) ischar(x)); 0128 p.addParamValue('units','t',@ischar); 0129 p.addParamValue('crit_val_box','yes',@(x) ischar(x) && (strcmpi('yes',x) || strcmpi('no',x))); 0130 0131 p.parse(GND_GRP_specGND_or_fname,test_id,varargin{:}); 0132 0133 % Manage VERBLEVEL 0134 if isempty(p.Results.verblevel), 0135 VERBLEVEL=2; %not global, just local 0136 else 0137 VERBLEVEL=p.Results.verblevel; 0138 end 0139 0140 if ~(strcmpi(p.Results.units,'t') || strcmpi(p.Results.units,'uV') || strcmpi(p.Results.units,'dB')) 0141 error('Value of optional input ''units'' must be ''t'', ''uV'', or ''dB''.'); 0142 end 0143 0144 use_color=str2bool(p.Results.use_color); 0145 0146 0147 %Load GND, GRP, or specGND struct 0148 if ischar(GND_GRP_specGND_or_fname), 0149 fprintf('Loading GND, GRP, or specGND struct variable from file %s.\n',GND_GRP_specGND_or_fname); 0150 load(GND_GRP_specGND_or_fname,'-MAT'); 0151 if ~exist('GND','var') && ~exist('GRP','var') && ~exist('specGND','var') 0152 error('File %s does not contain a GND, GRP, or specGND variable.',GND_GRP_specGND_or_fname); 0153 end 0154 if exist('GRP','var'), 0155 GND=GRP; %for simplicity GRP variables are re-named GND 0156 clear GRP; 0157 elseif exist('specGND','var'), 0158 GND=specGND; %for simplicity specGND variables are re-named GND 0159 clear specGND; 0160 end 0161 VerbReport(sprintf('Experiment: %s',GND.exp_desc),2,VERBLEVEL); 0162 else 0163 GND=GND_GRP_specGND_or_fname; %for simplicity GRP and specGND variables are re-named GND 0164 clear GND_GRP_specGND_or_fname; 0165 end 0166 0167 n_test=length(GND.t_tests); 0168 if test_id>n_test, 0169 error('There are only %d t-tests stored with these data, but you requested Test #%d',n_test,test_id); 0170 end 0171 0172 %FDR or permutation test correction for multiple comparisons 0173 if isnan(GND.t_tests(test_id).n_perm) 0174 fdr_crct=1; 0175 if VERBLEVEL, 0176 fprintf('Correcting for multiple comparisons via FDR procedure: %s\n', ... 0177 GND.t_tests(test_id).mult_comp_method); 0178 end 0179 else 0180 fdr_crct=0; 0181 if VERBLEVEL, 0182 fprintf('Correcting for multiple comparisons via permutation test.\n'); 0183 end 0184 end 0185 0186 if isfield(GND.t_tests(1),'freq_band') 0187 %create temporary fields to make frequency domain plotting 0188 % easily compatible with time domain plot 0189 GND.t_tests(test_id).time_wind=GND.t_tests(test_id).freq_band; 0190 GND.t_tests(test_id).mean_wind=GND.t_tests(test_id).mean_band; 0191 GND.t_tests(test_id).used_tpt_ids=GND.t_tests(test_id).used_freq_ids; 0192 GND.grands_t=GND.grands_pow_dB_t; 0193 GND.grands=GND.grands_pow_dB; 0194 freq_step=GND.freqs(2)-GND.freqs(1); 0195 freq_ord=orderofmag(freq_step); 0196 ord_pow=log10(freq_ord); %useful for formatting x-tick labels 0197 if ord_pow>=0, 0198 n_dig_past_dot=0; 0199 else 0200 n_dig_past_dot=-ord_pow; 0201 end 0202 GND.time_pts=rnd_orderofmag(GND.freqs); 0203 freq_domain=1; 0204 else 0205 freq_domain=0; 0206 end 0207 0208 %Clean up code by copying some GND.t_tests variables to their own 0209 %variable 0210 use_bin=GND.t_tests(test_id).bin; 0211 if strcmpi(GND.t_tests(test_id).mean_wind,'yes'), 0212 error('You can''t use sig_wave.m to visualize the results of t-tests performed on potentials averaged across a time window or frequency band. Use sig_raster.m or sig_topo.m instead.'); 0213 elseif strcmpi(GND.t_tests(test_id).mult_comp_method,'cluster mass perm test'), 0214 error('You can''t use sig_wave.m to visualize the results of cluster based permutation tests. Use sig_raster.m, gui_erp.m, or sig_topo.m instead.'); 0215 else 0216 if strcmpi(p.Results.units,'uV') || strcmpi(p.Results.units,'dB') , 0217 grands_wav=squeeze(GND.grands(GND.t_tests(test_id).used_chan_ids, ... 0218 :,use_bin)); 0219 else 0220 if GND.t_tests(test_id).null_mean, 0221 %recompute t-scores because mean of null hypothesis is not zero 0222 grands_wav=squeeze( (GND.grands(GND.t_tests(test_id).used_chan_ids, ... 0223 :,use_bin)-GND.t_tests(test_id).null_mean)./ ... 0224 GND.grands_stder(GND.t_tests(test_id).used_chan_ids, ... 0225 :,use_bin)); 0226 else 0227 grands_wav=squeeze(GND.grands_t(GND.t_tests(test_id).used_chan_ids, ... 0228 :,use_bin)); 0229 end 0230 end 0231 end 0232 0233 0234 if VERBLEVEL, 0235 if fdr_crct, 0236 fprintf('q level of critical t-scores: %f.\n',GND.t_tests(test_id).desired_alphaORq); 0237 else 0238 fprintf('Estimated alpha level of critical t-scores: %f.\n',GND.t_tests(test_id).estimated_alpha); 0239 end 0240 end 0241 0242 if isempty(p.Results.fig_id) 0243 fig_h=figure; 0244 else 0245 fig_h=figure(p.Results.fig_id); clf; 0246 end 0247 if freq_domain, 0248 set(fig_h,'name',['Bin ' int2str(use_bin) ' [' GND.bindesc{use_bin} ']'],'paperpositionmode','auto'); 0249 else 0250 set(fig_h,'name',['Bin ' int2str(use_bin) ' [' GND.bin_info(use_bin).bindesc ']'],'paperpositionmode','auto'); 0251 end 0252 %setting paperpositionmode to 'auto' means that if the figure is manually 0253 %resized, the printed version of the figure will reflect the whatever the 0254 %shown size was (at the time of printing) 0255 0256 0257 %% Plot waveforms 0258 %First plot time=0 (if time domain), wave=0 lines 0259 if isempty(p.Results.x_limits), 0260 time_lim=[GND.time_pts(1) GND.time_pts(end)]; 0261 else 0262 time_lim=p.Results.x_limits; 0263 end 0264 hd=plot(time_lim,[0 0],'k'); hold on; %time=0 line 0265 set(hd,'linewidth',p.Results.test_line_width); 0266 mx_abs_t=max(max(grands_wav)); 0267 mn_abs_t=min(min(grands_wav)); 0268 t_mx=mx_abs_t*1.1; 0269 t_mn=mn_abs_t*1.1; 0270 if isempty(p.Results.wave_limits), 0271 t_lim=[t_mn t_mx]; 0272 else 0273 t_lim=p.Results.wave_limits; 0274 end 0275 if freq_domain 0276 hd=plot([0 0],t_lim,'k'); % time=0 line 0277 set(hd,'linewidth',p.Results.test_line_width); 0278 end 0279 0280 0281 %plot waveforms 0282 h_wav=plot(GND.time_pts,grands_wav'); 0283 axis([time_lim t_lim]); 0284 chan_ct=0; 0285 for w=h_wav', 0286 chan_ct=chan_ct+1; 0287 dat.lab=GND.chanlocs(GND.t_tests(test_id).used_chan_ids(chan_ct)).labels; 0288 dat.times=GND.time_pts; 0289 if freq_domain, 0290 bdf_code=['Cp = get(gca,''CurrentPoint''); ' ... 0291 'Xp=Cp(2,1);', ... 0292 'Yp=Cp(2,2);', ... 0293 'dat=get(gcbo,''userdata'');', ... 0294 'ht=text(Xp,Yp,[int2str(rnd_orderofmag(Xp)) '' Hz, '' dat.lab]);' ... 0295 'set(ht,''backgroundcolor'',''w'',''horizontalalignment'',''center'',''verticalalignment'',''middle'',''buttondownfcn'',''delete(gcbo);'',''edgecolor'',''k'');']; 0296 else 0297 bdf_code=['Cp = get(gca,''CurrentPoint''); ' ... 0298 'Xp=Cp(2,1);', ... 0299 'Yp=Cp(2,2);', ... 0300 'dat=get(gcbo,''userdata'');', ... 0301 'ht=text(Xp,Yp,[int2str(round(Xp)) '' ms, '' dat.lab]);' ... 0302 'set(ht,''backgroundcolor'',''w'',''horizontalalignment'',''center'',''verticalalignment'',''middle'',''buttondownfcn'',''delete(gcbo);'',''edgecolor'',''k'');']; 0303 end 0304 set(w,'ButtonDownFcn',bdf_code,'userdata',dat); 0305 if ~use_color, 0306 set(w,'color','k'); 0307 end 0308 end 0309 0310 0311 %% plot test windows with critical vals 0312 if use_color, 0313 crit_color='r'; 0314 else 0315 crit_color='k'; 0316 end 0317 0318 n_wind=size(GND.t_tests(test_id).time_wind,1); 0319 for w=1:n_wind, 0320 hd=plot([1 1]*GND.t_tests(test_id).time_wind(w,1),t_lim,'k--'); 0321 set(hd,'linewidth',p.Results.test_line_width); 0322 hd=plot([1 1]*GND.t_tests(test_id).time_wind(w,2),t_lim,'k--'); 0323 set(hd,'linewidth',p.Results.test_line_width); 0324 if ~strcmpi(p.Results.units,'uV'), 0325 %don't plot critical values if waveform is in uV 0326 hd=plot(GND.t_tests(test_id).time_wind(w,:),[1 1]*GND.t_tests(test_id).crit_t(1),'r--'); 0327 set(hd,'linewidth',p.Results.test_line_width,'color',crit_color); 0328 crit_t_for_text=GND.t_tests(test_id).crit_t(1); %t coordinate for writing alpha=.05 (or whatever) 0329 if length(GND.t_tests(test_id).crit_t)>1, 0330 hd=plot(GND.t_tests(test_id).time_wind(w,:),[1 1]*GND.t_tests(test_id).crit_t(2),'r--'); 0331 set(hd,'linewidth',p.Results.test_line_width,'color',crit_color); 0332 if abs(t_lim(2))>abs(t_lim(1)), 0333 crit_t_for_text=GND.t_tests(test_id).crit_t(2); 0334 end 0335 end 0336 end 0337 end 0338 0339 %% Axis labels 0340 if strcmpi(p.Results.units,'t'), 0341 ht=ylabel('t-score'); 0342 else 0343 if freq_domain, 0344 ht=ylabel('10*log10(\muV^2/Hz)'); 0345 else 0346 ht=ylabel('\muV'); 0347 set(ht,'rotation',0); 0348 end 0349 end 0350 set(ht,'fontsize',14); 0351 if freq_domain, 0352 ht=xlabel('Hz'); 0353 else 0354 ht=xlabel('Time (ms)'); 0355 end 0356 set(ht,'fontsize',14); 0357 if p.Results.ydir==-1 && ~freq_domain, 0358 set(gca,'ydir','reverse'); 0359 end 0360 0361 %% Title 0362 if isempty(p.Results.title), 0363 if freq_domain, 0364 ht=title(GND.bindesc{use_bin}); 0365 else 0366 ht=title(GND.bin_info(use_bin).bindesc); 0367 end 0368 set(ht,'fontsize',16); 0369 else 0370 ht=title(p.Results.title); 0371 set(ht,'fontsize',16); 0372 end 0373 0374 if strcmpi(p.Results.units,'t') && strcmpi(p.Results.crit_val_box,'yes'), 0375 %Add box denoting alpha level/FDR level 0376 0377 %Figure out where to plot box 0378 tm_rng=time_lim(2)-time_lim(1); 0379 min_wind=min(min(GND.t_tests(test_id).time_wind)); 0380 max_wind=max(max(GND.t_tests(test_id).time_wind)); 0381 min_space=min_wind-time_lim(1); 0382 max_space=time_lim(2)-max_wind; 0383 if min_space>max_space, 0384 box_xloc=min_wind-tm_rng*.02; 0385 halign='right'; 0386 else 0387 box_xloc=max_wind+tm_rng*.02; 0388 halign='left'; 0389 end 0390 0391 if strcmpi(GND.t_tests(test_id).mult_comp_method(1),'b') 0392 %FDR control used 0393 h_alph1=text(box_xloc,crit_t_for_text, ... 0394 ['q=' num2str(rnd_orderofmag(GND.t_tests(test_id).desired_alphaORq))]); 0395 else 0396 %FWER control used 0397 h_alph1=text(box_xloc,crit_t_for_text, ... 0398 ['\alpha=' num2str(rnd_orderofmag(GND.t_tests(test_id).estimated_alpha))]); 0399 end 0400 set(h_alph1,'color',crit_color,'fontweight','normal','fontsize',12, ... 0401 'horizontalalignment',halign,'backgroundcolor',[1 1 1], ... 0402 'edgecolor',[1 1 1]*.3,'clipping','on','fontname','fixedwidth'); 0403 end 0404 0405 %% Tick marks specified? 0406 if ~isempty(p.Results.x_ticks), 0407 set(gca,'xtick',p.Results.x_ticks); 0408 end 0409 if ~isempty(p.Results.wave_ticks), 0410 set(gca,'ytick',p.Results.wave_ticks); 0411 end 0412 0413 0414 % 0415 %% %%%%%%%%%%%%%%%%%%%%% function str2bool() %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0416 % 0417 function bool=str2bool(str) 0418 %function bool=str2bool(str) 0419 0420 if ischar(str), 0421 if strcmpi(str,'yes') || strcmpi(str,'y') 0422 bool=1; 0423 else 0424 bool=0; 0425 end 0426 else 0427 bool=str; 0428 end