From 531f962b14eebfc5f1cc4ecfb3dd88fdb9bee6a1 Mon Sep 17 00:00:00 2001 From: cmaumet Date: Thu, 12 Dec 2024 11:49:13 +0100 Subject: [PATCH] first-level stat analysis --- narps_open/pipelines/matlabbatch_R5K7.m | 122 ++++++++++++++++++++++++ 1 file changed, 122 insertions(+) diff --git a/narps_open/pipelines/matlabbatch_R5K7.m b/narps_open/pipelines/matlabbatch_R5K7.m index c22eaf97..fc55e23e 100644 --- a/narps_open/pipelines/matlabbatch_R5K7.m +++ b/narps_open/pipelines/matlabbatch_R5K7.m @@ -151,3 +151,125 @@ 'ABS_PATH/narps_open_pipelines/data/original/ds001734/sub-001/func/usub-001_task-MGT_run-04_bold.nii' }; matlabbatch{end}.spm.spatial.smooth.fwhm = [8 8 8]; + +% ##### 4) First-level statistical analysis + +% We'll reuse Python code from DC61 to generate the conditions with parametric +% modulation + +def get_subject_information(event_files: list): + from nipype.interfaces.base import Bunch + + subject_info = [] + + for run_id, event_file in enumerate(event_files): + + onsets = [] + durations = [] + gain_value = [] + loss_value = [] + reaction_time = [] + + # Parse the events file + with open(event_file, 'rt') as file: + next(file) # skip the header + + for line in file: + info = line.strip().split() + +% --> independent_vars_first_level : - event-related design with each trial +% modelled with a duration of 4 sec and 3 linear parametric modulators (PMs +% orthogonalized via de-meaning against task and preceding PMs, respectively) +% for gain, loss and reaction time (in that order) as given in the .tsv log +% files + + onsets.append(float(info[0])) + durations.append(float(info[1])) + gain_value.append(float(info[2])) + loss_value.append(float(info[3])) + reaction_time.append(float(info[4])) + + # Create a Bunch for the run + subject_info.append( + Bunch( + conditions = [f'gamble_run{run_id + 1}'], + onsets = [onsets], + durations = [durations], % duration is 4s + amplitudes = None, + tmod = None, + pmod = [ + Bunch( + name = ['gain_param', 'loss_param', 'rt_param'], + poly = [1, 1, 1], + param = [gain_value, loss_value, reaction_time] + ) + ], + regressor_names = None, + regressors = None + )) + + return subject_info + + +matlabbatch{1}.spm.stats.fmri_spec.dir = ''; +matlabbatch{1}.spm.stats.fmri_spec.timing.units = 'secs'; +matlabbatch{1}.spm.stats.fmri_spec.timing.RT = ''; % Same as usual = TR +matlabbatch{1}.spm.stats.fmri_spec.sess(1).scans = ''; % scans from session 1 +% Note that parametric modulation was done in the Python code above and is not repeated here +matlabbatch{1}.spm.stats.fmri_spec.sess(1).cond = struct('name', {}, 'onset', {}, 'duration', {}, 'tmod', {}, 'pmod', {}, 'orth', {}); +matlabbatch{1}.spm.stats.fmri_spec.sess(1).multi = {''}; +matlabbatch{1}.spm.stats.fmri_spec.sess(1).regress = struct('name', {}, 'val', {}); +% --> 6 motion regressors (1st-order only) reflecting the 6 realignment parameters for translation and rotation movements obtained during preprocessing +matlabbatch{1}.spm.stats.fmri_spec.sess(1).multi_reg = {''}; % link to parameter motion files created by realign +matlabbatch{1}.spm.stats.fmri_spec.sess(1).hpf = 128; +matlabbatch{1}.spm.stats.fmri_spec.sess(2).scans = ''; % scans from session 2 +matlabbatch{1}.spm.stats.fmri_spec.sess(2).cond = struct('name', {}, 'onset', {}, 'duration', {}, 'tmod', {}, 'pmod', {}, 'orth', {}); +matlabbatch{1}.spm.stats.fmri_spec.sess(2).multi = {''}; +matlabbatch{1}.spm.stats.fmri_spec.sess(2).regress = struct('name', {}, 'val', {}); +matlabbatch{1}.spm.stats.fmri_spec.sess(2).multi_reg = {''}; % link to parameter motion files created by realign +matlabbatch{1}.spm.stats.fmri_spec.sess(2).hpf = 128; +matlabbatch{1}.spm.stats.fmri_spec.sess(3).scans = ''; % scans from session 3 +matlabbatch{1}.spm.stats.fmri_spec.sess(3).cond = struct('name', {}, 'onset', {}, 'duration', {}, 'tmod', {}, 'pmod', {}, 'orth', {}); +matlabbatch{1}.spm.stats.fmri_spec.sess(3).multi = {''}; +matlabbatch{1}.spm.stats.fmri_spec.sess(3).regress = struct('name', {}, 'val', {}); +matlabbatch{1}.spm.stats.fmri_spec.sess(3).multi_reg = {''}; % link to parameter motion files created by realign +matlabbatch{1}.spm.stats.fmri_spec.sess(3).hpf = 128; +matlabbatch{1}.spm.stats.fmri_spec.sess(4).scans = ''; % scans from session 4 +matlabbatch{1}.spm.stats.fmri_spec.sess(4).cond = struct('name', {}, 'onset', {}, 'duration', {}, 'tmod', {}, 'pmod', {}, 'orth', {}); +matlabbatch{1}.spm.stats.fmri_spec.sess(4).multi = {''}; +matlabbatch{1}.spm.stats.fmri_spec.sess(4).regress = struct('name', {}, 'val', {}); +matlabbatch{1}.spm.stats.fmri_spec.sess(4).multi_reg = {''}; % link to parameter motion files created by realign +matlabbatch{1}.spm.stats.fmri_spec.sess(4).hpf = 128; +matlabbatch{1}.spm.stats.fmri_spec.fact = struct('name', {}, 'levels', {}); +% --> canonical HRF plus temporal derivative +matlabbatch{1}.spm.stats.fmri_spec.bases.hrf.derivs = [1 0]; + +% Those are just the default +% matlabbatch{1}.spm.stats.fmri_spec.volt = 1; +% matlabbatch{1}.spm.stats.fmri_spec.global = 'None'; +% matlabbatch{1}.spm.stats.fmri_spec.mthresh = 0.8; +% matlabbatch{1}.spm.stats.fmri_spec.mask = {''}; +% matlabbatch{1}.spm.stats.fmri_spec.cvi = 'AR(1)'; + + +% ##### 5) Contrast definition at the first-level +% --> After model estimation, sum contrast images for each regressor of +% interest [task, gain (PM1), loss (PM2) and RT (PM3)] were computed across +% the 4 sessions in each participant. +self.contrast_list = ['0001', '0002', '0003', '0004'] +self.subject_level_contrasts = [ + ('task', 'T', + [f'gamble_run{r}' for r in range(1, len(self.run_list) + 1)], + [1]*len(self.run_list)), + ('effect_of_gain', 'T', + [f'gamble_run{r}xgain_param^1' for r in range(1, len(self.run_list) + 1)], + [1]*len(self.run_list)), + ('effect_of_loss', 'T', + [f'gamble_run{r}xloss_param^1' for r in range(1, len(self.run_list) + 1)], + [1]*len(self.run_list)) + ('effect_of_RT', 'T', + [f'gamble_run{r}xRT_param^1' for r in range(1, len(self.run_list) + 1)], + [1]*len(self.run_list)) +] + +% ##### 6) Group-level statistical analysis \ No newline at end of file