[FieldTrip] Single trial time series of spectral data in source space

[Johan Liljefors]

Dear Fieldtrip users
I am working on calculating trial level time series of spectral power in source space, for a MEG dataset I have collected. I have pre-processed and cleaned my data, which for a single subject is stored in the variable data_meg. I have the following variables:

  *   data_meg: The cleaned and preprocessed dataset for a single subject
  *   Headmodel_mne_meg and sourcemodelT: aligned headmodels, and sourcemodel (which has been calculated using Freesurfer)

If anyone has experience they can share it would be very appreciated. I have made an attempt below, the code runs but I cant tell if the results makes sense or not.

Thank  you for your time,
Regards
Johan

My code is as follows:
% Load the data
[data_meg, headmodel_mne_meg, sourcemodelT] = load_data(wsl_path,par_i,false);

% Compute the forward solution
leadfield_mne = compute_leadfield(data_ meg, sourcemodelT, headmodel_mne_meg);

% Calculate evoked data and apply noise covariance matrix from empty room recording
cfg = [];
cfg.covariance = 'yes';
cfg.covariancewindow = [-1 7.5];
cfg.keeptrials = 'no';
data_ind_avg = ft_timelockanalysis(cfg, data_ meg);
data_ind_avg.cov = cov_noise{par_i};                        % cov_noise{par_i} is the subject covariance matrix from the empty room recording

% Calculate the spatial filter using the evoked dataset and the noise covariance
source_avg = source_recon(leadfield_mne, headmodel_mne_meg,  [],data_ind_avg);
filter = source_avg.avg.filter;

% Do trial level source reconstruction using the filter
source_trial = source_recon_trials(leadfield_mne, headmodel_mne_meg, filter, data_meg);

% Calculate trial level spectral power using superlets and save the power spectrum
% This is done by "faking" a fieldtrip structure by putting each single trial into the .avg field and then running ft_freqanalysis
foi = 8:26;
toi = -1:0.01:7.5;
powspctrm = zeros(length(data_meg.trialinfo),...
    size(sourcemodelT.pos,1),...
    length(foi),...
    length(toi));

parfor iTrial = 1:length(data_meg.trialinfo)
        source_single_trial = rmfield(source_trial,"trial");     % These two lines "fake" a single trial which is then processed by ft_freqanalysis
        source_single_trial.avg = source_trial.trial(iTrial);
        cfg = [];
        cfg.method              = 'superlet';
        cfg.output              = 'pow';
        cfg.pad                 = 'nextpow2';
        cfg.foi                 = foi;
        cfg.toi                 = toi;
        % cfg.order goes from 1 to 30 (https://www.biorxiv.org/content/10.1101/583732v4.full.pdf)
        cfg.order               = floor(linspace(1,30,numel(cfg.foi)));
        cfg.width               = 5;
        cfg.gwidth               = 3;
        TFR_induced = ft_freqanalysis(cfg,source_single_trial);
        powspctrm(iTrial,:,:,:) = TFR_induced.powspctrm;
end

% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%           LOCAL FUNCTIONS                        %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function leadfield_mne = compute_leadfield(data_meg, sourcemodelT, headmodel_mne_meg)
    cfg = [];
    cfg.grad             = data_meg.grad;              % sensor positions
    cfg.channel          = 'meggrad';                  % the used channels
    cfg.senstype         = 'meg';
    cfg.sourcemodel      = sourcemodelT;
    cfg.headmodel        = headmodel_mne_meg;          % volume conduction model
    leadfield_mne        = ft_prepare_leadfield(cfg,data_meg);
end

function source = source_recon(leadfield_mne, headmodel_mne_meg, filter, data)
    cfg                     = [];
    cfg.method              = 'mne';
    cfg.channel             = 'meggrad';
    cfg.senstype            = 'meg';
    cfg.latency             = [-1 7.5];
    cfg.sourcemodel         = leadfield_mne;
    cfg.headmodel           = headmodel_mne_meg;
    cfg.mne.prewhiten       = 'yes';
    cfg.mne.lambda          = 3;
    cfg.mne.scalesourcecov  = 'yes';
    if isempty(filter)
        cfg.mne.keepfilter      = 'yes';
    else
        cfg.sourcemodel.filter = filter;
    end
    source = ft_sourceanalysis(cfg,data);
end

function source_trial = source_recon_trials(leadfield_mne, headmodel_mne_meg, filter, data_meg)
    cfg                     = [];
    cfg.method              = 'mne';
    cfg.channel             = 'meggrad';
    cfg.senstype            = 'meg';
    cfg.latency             = [-1 7.5];
    cfg.sourcemodel         = leadfield_mne;
    cfg.headmodel           = headmodel_mne_meg;
    cfg.mne.prewhiten       = 'yes';
    cfg.mne.lambda          = 3;
    cfg.mne.scalesourcecov  = 'yes';
    if isempty(filter)
        cfg.mne.keepfilter      = 'yes';
    else
        cfg.sourcemodel.filter = filter;
    end
    cfg.keeptrials = 'yes';
    cfg.rawtrial = 'yes';
    source_trial = ft_sourceanalysis(cfg,data_meg);
end





N?r du skickar e-post till Karolinska Institutet (KI) inneb?r detta att KI kommer att behandla dina personuppgifter. H?r finns information om hur KI behandlar personuppgifter<https://ki.se/om-ki/integritetsskyddspolicy>.


Sending email to Karolinska Institutet (KI) will result in KI processing your personal data. You can read more about KI's processing of personal data here<https://staff.ki.se/data-protection-policy>.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.science.ru.nl/pipermail/fieldtrip/attachments/20240315/d8c22fac/attachment.htm>