show_ic_eeg

 Show ICs

 Copyright (C) 2011, ATR All Rights Reserved.
 License : New BSD License(see VBMEG_LICENSE.txt)

0001 function show_ic_eeg(p)
0002 % Show ICs
0003 %
0004 % Copyright (C) 2011, ATR All Rights Reserved.
0005 % License : New BSD License(see VBMEG_LICENSE.txt)
0006 
0007 disp(mfilename);
0008 
0009 % Classification program
0010 correct_ic_prog='correct_classify_ic_eeg';
0011 
0012 % Make directory to save figures
0013 save_dir = fullfile(p.fig_root, mfilename, p.sub);
0014 if exist(save_dir, 'dir') ~= 7
0015     vb_mkdir(save_dir);
0016 end
0017 
0018 % Load EEG data
0019 data_file = fullfile(p.proj_root, p.eeg_dirname, p.trial_dirname, 'br_s1.eeg.mat');
0020 [eeg, ~, time_info] = vb_load_meg_data(data_file);
0021 [Nch, Nt, Ntr] = size(eeg);
0022 samplingrate = time_info.sample_frequency;
0023 time = time_info.time;
0024 
0025 % Load positions of channels
0026 pos = vb_load_channel(data_file);
0027 
0028 % Load ICA result
0029 ica_dir = fullfile(p.proj_root, p.eeg_dirname, p.ica_dirname);
0030 ica_result = fullfile(ica_dir, [correct_ic_prog '.mat']);
0031 load(ica_result,'w','s','v','unmix','mix','brain_ic','noise_ic')
0032 
0033 % Calculate ICs
0034 ic = zeros(Nch, Nt, Ntr);
0035 for tr = 1:Ntr
0036     ic(:,:,tr) = unmix*eeg(:,:,tr);  % [component x sample x trial]
0037 end
0038 
0039 % Calculate stimulu-triggered average and power spectrum
0040 mic = mean(ic, 3);% Average each IC across trials
0041 [pow, f] = vb_fftpower(ic, samplingrate);% Calculate power spectrum
0042 mp = mean(log(pow), 3);% Average power spectrum across trials
0043 
0044 % Obtain classification
0045 classification = -ones(1, Nch);
0046 classification(brain_ic) = 1;% 1: brain_ic, -1: noise_ic
0047 
0048 % Show and save ICs
0049 h = figure;
0050 cc = 1;
0051 while cc <= Nch
0052     clf(h)
0053     set(0, 'CurrentFigure', h)
0054     k = 1;% Position of figure
0055     for c = cc:min([cc+4,Nch]);
0056         if classification(c) == 1
0057             color = 'b';
0058         else
0059             color = 'r';
0060         end
0061         
0062         % Plot spatial pattern of mixing vector
0063         subplot(5, 3, 1+3*(k-1))
0064         vb_plot_sensor_2d(pos, mix(:,c));
0065         vb_plot_sensor_2d_head_plot_add(gca);
0066         axis square off
0067         title(['IC ' num2str(c)])
0068         
0069         % Plot stimulus-triggered average of components
0070         subplot(5, 3, 2+3*(k-1))
0071         plot(time, mic(c,:), color);
0072         xlim([-0.5 1])
0073         if k == 1
0074             title('Stimulus-triggered average')
0075         elseif k == 5
0076             xlabel('Time [s]')
0077         end
0078         
0079         % Plot power spectrum
0080         subplot(5, 3, 3+3*(k-1))
0081         plot(f, mp(c,:), color);
0082         xlim([0 50])
0083         if k == 1
0084             title('Power spectrum (log)')
0085         elseif k == 5
0086             xlabel('Frequency [Hz]')
0087         end
0088         k = k+1;
0089     end
0090     ax1 = axes('Position', [0 0 1 1], 'Visible', 'off');
0091     text1 = {'Blue : Brain'; 'Red : Noise'};
0092     axes(ax1)
0093     text(0.93, 0.87, text1)
0094     
0095     % Save figure
0096     vb_savefig_as_shown(h, fullfile(save_dir, [num2str(cc) '-' num2str(cc+4)]))
0097     cc = cc+5;
0098 end
0099