estimate_source_current_eeg

 Estimate source currents from EEG data

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

0001 function estimate_source_current_eeg(p)
0002 % Estimate source currents from EEG data
0003 %
0004 % Copyright (C) 2011, ATR All Rights Reserved.
0005 % License : New BSD License(see VBMEG_LICENSE.txt)
0006 
0007 disp(mfilename);
0008 
0009 % Set parameters
0010 time_signal = [0 0.3]; % Period for stimulus-evoked activity [s]
0011 time_noise = [-0.3 0]; % Period for baseline [s]
0012 
0013 % Make directory to save current
0014 current_dir = fullfile(p.proj_root, p.eeg_dirname, p.current_dirname);
0015 if exist(current_dir, 'dir') ~= 7
0016     mkdir(current_dir);
0017 end
0018 
0019 for ta = 1:length(p.task_list)
0020     
0021     %% Estimate current variance
0022     clear bp
0023     bayes_parm = vb_set_bayes_default_parameters;
0024     
0025     % Set input files (relative path from proj_root)
0026     bayes_parm.brainfile = fullfile(p.brain_dirname, [p.struct_name, '.brain.mat']);
0027     bayes_parm.areafile = fullfile(p.brain_dirname, [p.struct_name, '.area.mat']);
0028     bayes_parm.actfile = fullfile(p.brain_dirname, [p.struct_name, '.act.mat']);
0029     bayes_parm.megfile{1} = fullfile(p.eeg_dirname, p.trial_dirname, ['crri_' p.task_list{ta} '.info.mat']);
0030     bayes_parm.basisfile = fullfile(p.eeg_dirname, p.leadfield_dirname, [p.struct_name '.basis.mat']);
0031     bayes_parm.megfile_baseline = bayes_parm.megfile;
0032     
0033     % Specify current area and fMRI activity
0034     bayes_parm.area_key  = 'Cortex';% Current area
0035     bayes_parm.act_key = p.prior{ta};% Prior for current variance
0036     bayes_parm.prior_weight = p.prior_weight(ta);% Relative influence of prior information (0-1)
0037     
0038     % Load time information
0039     [~, ~, time_info] = vb_load_meg_data(fullfile(p.proj_root, bayes_parm.megfile{1}));
0040     time = time_info.time;
0041     
0042     % Set time window for stimulus-evoked activity
0043     [~, from] = min(abs(time-time_signal(1)));
0044     [~, to] = min(abs(time-time_signal(2)));
0045     bayes_parm.twin_meg = [from to]; % Time window for analysis
0046     bayes_parm.Tperiod = to-from+1; % Time period for current estimation
0047     bayes_parm.Tnext = to-from+1; % Time step for next period
0048     
0049     % Set time window for baseline
0050     [~, from] = min(abs(time-time_noise(1)));
0051     [~, to] = min(abs(time-time_noise(2)));
0052     bayes_parm.twin_noise = [from to];
0053     bayes_parm.twin_baseline = [from to];
0054     
0055     bayes_parm.patch_norm = ON;% Patch area normalization is applied
0056     
0057     % Set output file
0058     bayes_parm.bayesfile = fullfile(p.eeg_dirname, p.current_dirname, [p.task_list{ta} '.bayes.mat']);
0059     
0060     % Set parameters for noise covariance estimation
0061     bayes_parm = vb_set_noise_estimation_model(bayes_parm, p.proj_root);
0062     
0063     % Estimate current variance
0064     vb_job_vb(p.proj_root, bayes_parm)
0065     
0066     
0067     %% Estimate source current
0068     
0069     % Set input files (relative path from proj_root)
0070     current_parm.megfile = bayes_parm.megfile;
0071     current_parm.bayesfile = bayes_parm.bayesfile;
0072     
0073     % Trial-averaged current will be estimated.
0074     current_parm.trial_average = ON;
0075     
0076     % Set time window for estimating current
0077     a = vb_load_meg_data(fullfile(p.proj_root, current_parm.megfile{1}));
0078     current_parm.twin_meg = [1 size(a, 2)];
0079     current_parm.Tperiod = size(a, 2);
0080     
0081     % Set output file (relative path from proj_root)
0082     current_parm.currfile = fullfile(p.eeg_dirname, p.current_dirname, [p.task_list{ta} '.curr.mat']);
0083     
0084     % Estimate source current
0085     vb_job_current(p.proj_root, current_parm);
0086 end
0087 
0088 
0089 
0090 
0091