From 15f09acb394851a19a96bc0d18a364eea469b701 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 29 Mar 2023 11:34:09 +0200 Subject: [PATCH 01/25] [REFAC] pipeline J7F9 to conform with the Pipeline class --- narps_open/pipelines/team_J7F9.py | 1115 ++++++++++++++++------------- tests/pipelines/test_team_J7F9.py | 55 ++ 2 files changed, 667 insertions(+), 503 deletions(-) create mode 100644 tests/pipelines/test_team_J7F9.py diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index b96f2db6..6023a367 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -1,508 +1,617 @@ -from nipype.interfaces.spm import (Coregister, Smooth, OneSampleTTestDesign, EstimateModel, EstimateContrast, - Level1Design, TwoSampleTTestDesign, RealignUnwarp, - Normalize12, NewSegment, FieldMap) -from nipype.interfaces.fsl import ExtractROI -from nipype.interfaces.spm import Threshold -from nipype.algorithms.modelgen import SpecifySPMModel +#!/usr/bin/python +# coding: utf-8 + +""" Write the work of NARPS' team J7F9 using Nipype """ + +from os.path import join + +from nipype import Workflow, Node, MapNode from nipype.interfaces.utility import IdentityInterface, Function from nipype.interfaces.io import SelectFiles, DataSink +from nipype.interfaces.spm import ( + Smooth, Level1Design, OneSampleTTestDesign, TwoSampleTTestDesign, + EstimateModel, EstimateContrast, Threshold + ) +from nipype.algorithms.modelgen import SpecifySPMModel from nipype.algorithms.misc import Gunzip -from nipype import Workflow, Node, MapNode, JoinNode -from nipype.interfaces.base import Bunch - -from os.path import join as opj -import os -import json - -def get_subject_infos(event_files, runs): - ''' - Create Bunchs for specifySPMModel. - - Parameters : - - event_files: list of str, list of events files (one per run) for the subject - - runs: list of str, list of runs to use - - Returns : - - subject_info : list of Bunch for 1st level analysis. - ''' - from nipype.interfaces.base import Bunch - import numpy as np - - cond_names = ['trial', 'missed'] - onset = {} - duration = {} - weights_gain = {} - weights_loss = {} - - for r in range(len(runs)): # Loop over number of runs. - onset.update({s + '_run' + str(r+1) : [] for s in cond_names}) # creates dictionary items with empty lists - duration.update({s + '_run' + str(r+1) : [] for s in cond_names}) - weights_gain.update({'gain_run' + str(r+1) : []}) - weights_loss.update({'loss_run' + str(r+1) : []}) - - for r, run in enumerate(runs): - - f_events = event_files[r] - - with open(f_events, 'rt') as f: - next(f) # skip the header - - for line in f: - info = line.strip().split() - - for cond in cond_names: - val = cond + '_run' + str(r+1) # trial_run1 - val_gain = 'gain_run' + str(r+1) # gain_run1 - val_loss = 'loss_run' + str(r+1) # loss_run1 - if cond =='trial': - onset[val].append(float(info[0])) # onsets for trial_run1 - duration[val].append(float(0)) - weights_gain[val_gain].append(float(info[2])) # weights gain for trial_run1 - weights_loss[val_loss].append(float(info[3])) # weights loss for trial_run1 - elif cond=='missed': - if float(info[4]) < 0.1 or str(info[5]) == 'NoResp': - onset[val].append(float(info[0])) + +from narps_open.pipelines import Pipeline + +class PipelineTeamJ7F9(Pipeline): + """ A class that defines the pipeline of team J7F9. """ + + def __init__(self): + super().__init__() + self.fwhm = 8.0 + self.team_id = 'J7F9' + self.contrast_list = ['0001', '0002', '0003', '0004'] + + def get_preprocessing(self): + """ No preprocessing has been done by team J7F9 """ + return None + + def get_run_level_analysis(self): + """ No run level analysis has been done by team J7F9 """ + return None + + def get_subject_infos(event_files, runs): + """ + Create Bunchs for specifySPMModel. + + Parameters : + - event_files: list of str, list of events files (one per run) for the subject + - runs: list of str, list of runs to use + + Returns : + - subject_info : list of Bunch for 1st level analysis. + """ + from numpy import mean + from nipype.interfaces.base import Bunch + + condition_names = ['trial', 'missed'] + onset = {} + duration = {} + weights_gain = {} + weights_loss = {} + + for run_id in range(len(runs)): # Loop over number of runs. + # creates dictionary items with empty lists + onset.update({s + '_run' + str(run_id + 1) : [] for s in condition_names}) + duration.update({s + '_run' + str(run_id + 1) : [] for s in condition_names}) + weights_gain.update({'gain_run' + str(run_id + 1) : []}) + weights_loss.update({'loss_run' + str(run_id + 1) : []}) + + for run_id, _ in enumerate(runs): + f_events = event_files[run_id] + + with open(f_events, 'rt') as file: + next(file) # skip the header + + for line in file: + info = line.strip().split() + + for condition in condition_names: + val = condition + '_run' + str(run_id + 1) # trial_run1 + val_gain = 'gain_run' + str(run_id + 1) # gain_run1 + val_loss = 'loss_run' + str(run_id + 1) # loss_run1 + if condition == 'trial': + onset[val].append(float(info[0])) # onsets for trial_run1 duration[val].append(float(0)) - - - for gain_val in weights_gain.keys(): - weights_gain[gain_val] = weights_gain[gain_val] - np.mean(weights_gain[gain_val]) - weights_gain[gain_val] = weights_gain[gain_val].tolist() - - for loss_val in weights_loss.keys(): - weights_loss[loss_val] = weights_loss[loss_val] - np.mean(weights_loss[loss_val]) - weights_loss[loss_val] = weights_loss[loss_val].tolist() - - # Bunching is done per run, i.e. trial_run1, trial_run2, etc. - # But names must not have '_run1' etc because we concatenate runs - subject_info = [] - for r in range(len(runs)): - - if len(onset['missed_run' + str(r+1)]) ==0: - cond_names = ['trial'] - - cond = [c + '_run' + str(r+1) for c in cond_names] - gain = 'gain_run' + str(r+1) - loss = 'loss_run' + str(r+1) - - subject_info.insert(r, - Bunch(conditions=cond_names, - onsets=[onset[k] for k in cond], - durations=[duration[k] for k in cond], - amplitudes=None, - tmod=None, - pmod=[Bunch(name=['gain', 'loss'], - poly=[1, 1], - param=[weights_gain[gain], - weights_loss[loss]])], - regressor_names=None, - regressors=None)) - - return subject_info - -def get_contrasts(subject_id): - ''' - Create the list of tuples that represents contrasts. - Each contrast is in the form : - (Name,Stat,[list of condition names],[weights on those conditions]) - ''' - # list of condition names - conditions = ['trial', 'trialxgain^1', 'trialxloss^1'] - - # create contrasts - trial = ('trial', 'T', conditions, [1, 0, 0]) - - effect_gain = ('effect_of_gain', 'T', conditions, [0, 1, 0]) - - effect_loss = ('effect_of_loss', 'T', conditions, [0, 0, 1]) - - # contrast list - contrasts = [trial, effect_gain, effect_loss] - - return contrasts - - -def get_parameters_file(filepaths, subject_id, result_dir, working_dir): - ''' - Create new tsv files with only desired parameters per subject per run. - - Parameters : - - filepaths : paths to subject parameters file (i.e. one per run) - - subject_id : subject for whom the 1st level analysis is made - - result_dir: str, directory where results will be stored - - working_dir: str, name of the sub-directory for intermediate results - - Return : - - parameters_file : paths to new files containing only desired parameters. - ''' - import pandas as pd - import numpy as np - from os.path import join as opj - import os - - if not isinstance(filepaths, list): - filepaths = [filepaths] - parameters_file = [] - - for i, file in enumerate(filepaths): - df = pd.read_csv(file, sep = '\t', header=0) - temp_list = np.array([df['X'], df['Y'], df['Z'], - df['RotX'], df['RotY'], df['RotZ'], - df['CSF'], df['WhiteMatter'], df['GlobalSignal']]) # Parameters we want to use for the model - retained_parameters = pd.DataFrame(np.transpose(temp_list)) - new_path =opj(result_dir, working_dir, 'parameters_file', f"parameters_file_sub-{subject_id}_run0{str(i+1)}.tsv") - if not os.path.isdir(opj(result_dir, working_dir, 'parameters_file')): - os.mkdir(opj(result_dir, working_dir, 'parameters_file')) - writer = open(new_path, "w") - writer.write(retained_parameters.to_csv(sep = '\t', index = False, header = False, na_rep = '0.0')) - writer.close() - - parameters_file.append(new_path) - - return parameters_file - -def rm_gunzip_files(files, subject_id, result_dir, working_dir): - import shutil - from os.path import join as opj - - gunzip_dir = opj(result_dir, working_dir, 'l1_analysis', f"_subject_id_{subject_id}", 'gunzip_func') - - try: - shutil.rmtree(gunzip_dir) - except OSError as e: - print(e) - else: - print("The directory is deleted successfully") - - return files - -def rm_smoothed_files(files, subject_id, result_dir, working_dir): - import shutil - from os.path import join as opj - - smooth_dir = opj(result_dir, working_dir, 'l1_analysis', f"_subject_id_{subject_id}", 'smooth') - - try: - shutil.rmtree(smooth_dir) - except OSError as e: - print(e) - else: - print("The directory is deleted successfully") - - return files - - -def get_l1_analysis(subject_list, TR, fwhm, run_list, exp_dir, result_dir, working_dir, output_dir): - """ - Returns the first level analysis workflow. - - Parameters: - - exp_dir: str, directory where raw data are stored - - result_dir: str, directory where results will be stored - - working_dir: str, name of the sub-directory for intermediate results - - output_dir: str, name of the sub-directory for final results - - subject_list: list of str, list of subject for which you want to do the analysis - - run_list: list of str, list of runs for which you want to do the analysis - - fwhm: float, fwhm for smoothing step - - TR: float, time repetition used during acquisition - - Returns: - - l1_analysis : Nipype WorkFlow - """ - # Infosource Node - To iterate on subjects - infosource = Node(IdentityInterface(fields = ['subject_id', 'exp_dir', 'result_dir', - 'working_dir', 'run_list'], - exp_dir = exp_dir, result_dir = result_dir, working_dir = working_dir, - run_list = run_list), - name = 'infosource') - - infosource.iterables = [('subject_id', subject_list)] - - # Templates to select files node - param_file = opj('derivatives', 'fmriprep', 'sub-{subject_id}', 'func', - 'sub-{subject_id}_task-MGT_run-*_bold_confounds.tsv') - - func_file = opj('derivatives', 'fmriprep', 'sub-{subject_id}', 'func', - 'sub-{subject_id}_task-MGT_run-*_bold_space-MNI152NLin2009cAsym_preproc.nii.gz') - - event_files = opj('sub-{subject_id}', 'func', 'sub-{subject_id}_task-MGT_run-*_events.tsv') - - template = {'param' : param_file, 'func' : func_file, 'event' : event_files} - - # SelectFiles node - to select necessary files - selectfiles = Node(SelectFiles(template, base_directory=exp_dir), name = 'selectfiles') - - # DataSink Node - store the wanted results in the wanted repository - datasink = Node(DataSink(base_directory=result_dir, container=output_dir), name='datasink') - - # GUNZIP NODE : SPM do not use .nii.gz files - gunzip_func = MapNode(Gunzip(), name = 'gunzip_func', iterfield = ['in_file']) - - ## Smoothing node - smooth = Node(Smooth(fwhm = fwhm), name = 'smooth') - - # Get Subject Info - get subject specific condition information - subject_infos = Node(Function(input_names=['event_files', 'runs'], - output_names=['subject_info'], - function=get_subject_infos), - name='subject_infos') - - # SpecifyModel - Generates SPM-specific Model - specify_model = Node(SpecifySPMModel(concatenate_runs = True, input_units = 'secs', output_units = 'secs', - time_repetition = TR, high_pass_filter_cutoff = 128), name='specify_model') - - # Level1Design - Generates an SPM design matrix - l1_design = Node(Level1Design(bases = {'hrf': {'derivs': [0, 0]}}, timing_units = 'secs', - interscan_interval = TR), name='l1_design') - - # EstimateModel - estimate the parameters of the model - l1_estimate = Node(EstimateModel(estimation_method={'Classical': 1}), - name="l1_estimate") - - # Node contrasts to get contrasts - contrasts = Node(Function(function=get_contrasts, - input_names=['subject_id'], - output_names=['contrasts']), - name='contrasts') - - # Node parameters to get parameters files - parameters = Node(Function(function=get_parameters_file, - input_names=['filepaths', 'subject_id', 'result_dir', 'working_dir'], - output_names=['parameters_file']), - name='parameters') - - # EstimateContrast - estimates contrasts - contrast_estimate = Node(EstimateContrast(), name="contrast_estimate") - - remove_gunzip_files = Node(Function(input_names = ['files', 'subject_id', 'result_dir', 'working_dir'], - output_names = ['files'], - function = rm_gunzip_files), name = 'remove_gunzip_files') - - remove_smoothed_files = Node(Function(input_names = ['files', 'subject_id', 'result_dir', 'working_dir'], - output_names = ['files'], - function = rm_smoothed_files), name = 'remove_smoothed_files') - - # Create l1 analysis workflow and connect its nodes - l1_analysis = Workflow(base_dir = opj(result_dir, working_dir), name = "l1_analysis") - - l1_analysis.connect([(infosource, selectfiles, [('subject_id', 'subject_id')]), - (infosource, subject_infos, [('exp_dir', 'exp_dir'), ('run_list', 'runs')]), - (infosource, contrasts, [('subject_id', 'subject_id')]), - (infosource, remove_gunzip_files, [('subject_id', 'subject_id'), ('result_dir', 'result_dir'), - ('working_dir', 'working_dir')]), - (infosource, remove_smoothed_files, [('subject_id', 'subject_id'), - ('result_dir', 'result_dir'), - ('working_dir', 'working_dir')]), - (subject_infos, specify_model, [('subject_info', 'subject_info')]), - (contrasts, contrast_estimate, [('contrasts', 'contrasts')]), - (selectfiles, parameters, [('param', 'filepaths')]), - (selectfiles, subject_infos, [('event', 'event_files')]), - (infosource, parameters, [('subject_id', 'subject_id'), ('result_dir', 'result_dir'), - ('working_dir', 'working_dir')]), - (selectfiles, gunzip_func, [('func', 'in_file')]), - (gunzip_func, smooth, [('out_file', 'in_files')]), - (smooth, remove_gunzip_files, [('smoothed_files', 'files')]), - (remove_gunzip_files, specify_model, [('files', 'functional_runs')]), - (parameters, specify_model, [('parameters_file', 'realignment_parameters')]), - (specify_model, l1_design, [('session_info', 'session_info')]), - (l1_design, l1_estimate, [('spm_mat_file', 'spm_mat_file')]), - (l1_estimate, contrast_estimate, [('spm_mat_file', 'spm_mat_file'), - ('beta_images', 'beta_images'), - ('residual_image', 'residual_image')]), - (contrast_estimate, datasink, [('con_images', 'l1_analysis.@con_images'), - ('spmT_images', 'l1_analysis.@spmT_images'), - ('spm_mat_file', 'l1_analysis.@spm_mat_file')]), - (contrast_estimate, remove_smoothed_files, [('spmT_images', 'files')]) - ]) - - return l1_analysis - - -def get_subset_contrasts(file_list, method, subject_list, participants_file): - ''' - Parameters : - - file_list : original file list selected by selectfiles node - - subject_list : list of subject IDs that are in the wanted group for the analysis - - participants_file: str, file containing participants caracteristics - - method: str, one of "equalRange", "equalIndifference" or "groupComp" - - This function return the file list containing only the files belonging to subject in the wanted group. - ''' - equalIndifference_id = [] - equalRange_id = [] - equalIndifference_files = [] - equalRange_files = [] - - with open(participants_file, 'rt') as f: - next(f) # skip the header - - for line in f: + # weights gain for trial_run1 + weights_gain[val_gain].append(float(info[2])) + # weights loss for trial_run1 + weights_loss[val_loss].append(float(info[3])) + elif condition == 'missed': + if float(info[4]) < 0.1 or str(info[5]) == 'NoResp': + onset[val].append(float(info[0])) + duration[val].append(float(0)) + + for gain_key, gain_value in weights_gain: + gain_value = gain_value - mean(gain_value) + weights_gain[gain_key] = gain_value.tolist() + + for loss_key, loss_value in weights_loss: + loss_value = loss_value - mean(loss_value) + weights_loss[loss_key] = loss_value.tolist() + + # Bunching is done per run, i.e. trial_run1, trial_run2, etc. + # But names must not have '_run1' etc because we concatenate runs + subject_info = [] + for run_id in range(len(runs)): + + if len(onset['missed_run' + str(run_id + 1)]) ==0: + condition_names = ['trial'] + + conditions = [c + '_run' + str(run_id + 1) for c in condition_names] + gain = 'gain_run' + str(run_id + 1) + loss = 'loss_run' + str(run_id + 1) + + subject_info.insert( + run_id, + Bunch( + conditions = condition_names, + onsets = [onset[k] for k in conditions], + durations = [duration[k] for k in conditions], + amplitudes = None, + tmod = None, + pmod = [ + Bunch( + name = ['gain', 'loss'], + poly = [1, 1], + param = [weights_gain[gain], weights_loss[loss]] + ) + ], + regressor_names = None, + regressors=None + ) + ) + + return subject_info + + def get_contrasts(subject_id): + """ + Create the list of tuples that represents contrasts. + Each contrast is in the form : + (Name,Stat,[list of condition names],[weights on those conditions]) + """ + # List of condition names + conditions = ['trial', 'trialxgain^1', 'trialxloss^1'] + + # Create contrasts + trial = ('trial', 'T', conditions, [1, 0, 0]) + effect_gain = ('effect_of_gain', 'T', conditions, [0, 1, 0]) + effect_loss = ('effect_of_loss', 'T', conditions, [0, 0, 1]) + + # Return contrast list + return [trial, effect_gain, effect_loss] + + + def get_parameters_file(filepaths, subject_id, working_dir): + """ + Create new tsv files with only desired parameters per subject per run. + + Parameters : + - filepaths : paths to subject parameters file (i.e. one per run) + - subject_id : subject for whom the 1st level analysis is made + - working_dir: str, name of the directory for intermediate results + + Return : + - parameters_file : paths to new files containing only desired parameters. + """ + from os import mkdir + from os.path import join, isdir + from pandas import DataFrame, read_csv + from numpy import array, transpose + + if not isinstance(filepaths, list): + filepaths = [filepaths] + + parameters_file = [] + for file_id, file in enumerate(filepaths): + data_frame = read_csv(file, sep = '\t', header=0) + + # Parameters we want to use for the model + temp_list = array([ + data_frame['X'], data_frame['Y'], data_frame['Z'], + data_frame['RotX'], data_frame['RotY'], data_frame['RotZ'], + data_frame['CSF'], data_frame['WhiteMatter'], data_frame['GlobalSignal']]) + retained_parameters = DataFrame(transpose(temp_list)) + + # Write parameters to a parameters file + # TODO : warning !!! filepaths must be ordered (1,2,3,4) for the following code to work + new_path = join(working_dir, 'parameters_file', + f'parameters_file_sub-{subject_id}_run-{str(file_id + 1).zfill(2)}.tsv') + + if not isdir(join(working_dir, 'parameters_file')): + mkdir(join(working_dir, 'parameters_file')) + + with open(new_path, 'w') as writer: + writer.write(retained_parameters.to_csv( + sep = '\t', index = False, header = False, na_rep = '0.0')) + + parameters_file.append(new_path) + + return parameters_file + + def remove_gunzip_files(_, subject_id, working_dir): + """ + This method is used in a Function node to fully remove + the files generated by the gunzip node, once they aren't needed anymore. + + Parameters: + - _: Node input only used for triggering the Node + - subject_id: str, TODO + - working_id: str, TODO + """ + from shutil import rmtree + from os.path import join + + try: + rmtree(join(working_dir, 'l1_analysis', f'_subject_id_{subject_id}', 'gunzip_func')) + except OSError as error: + print(error) + else: + print('The directory is deleted successfully') + + def remove_smoothed_files(_, subject_id, working_dir): + """ + This method is used in a Function node to fully remove + the files generated by the smoothing node, once they aren't needed anymore. + + Parameters: + - _: Node input only used for triggering the Node + - subject_id: str, TODO + - working_id: str, TODO + """ + from shutil import rmtree + from os.path import join + + try: + rmtree(join(working_dir, 'l1_analysis', f'_subject_id_{subject_id}', 'smooth')) + except OSError as error: + print(error) + else: + print('The directory is deleted successfully') + + def get_subject_level_analysis(self): + """ + Create the subject level analysis workflow. + + Returns: + - l1_analysis : nipype.WorkFlow + """ + # Infosource Node - To iterate on subjects + infosource = Node(IdentityInterface( + fields = ['subject_id', 'exp_dir', 'result_dir', 'working_dir', 'run_list'], + dataset_dir = self.directories.dataset_dir, + results_dir = self.directories.results_dir, + working_dir = self.directories.working_dir, + run_list = self.run_list), + name = 'infosource') + infosource.iterables = [('subject_id', self.subject_list)] + + # Templates to select files node + template = { + # Parameter file + 'param' : join('derivatives', 'fmriprep', 'sub-{subject_id}', 'func', + 'sub-{subject_id}_task-MGT_run-*_bold_confounds.tsv'), + 'func' : join('derivatives', 'fmriprep', 'sub-{subject_id}', 'func', + 'sub-{subject_id}_task-MGT_run-*_bold_space-MNI152NLin2009cAsym_preproc.nii.gz'), + 'event' : join('sub-{subject_id}', 'func', + 'sub-{subject_id}_task-MGT_run-*_events.tsv') + } + + # SelectFiles - to select necessary files + selectfiles = Node(SelectFiles(template, base_directory = self.directories.dataset_dir), + name = 'selectfiles') + + # DataSink - store the wanted results in the wanted repository + datasink = Node(DataSink(base_directory = self.directories.output_dir), + name='datasink') + + # Gunzip - gunzip files because SPM do not use .nii.gz files + gunzip_func = MapNode(Gunzip(), + name = 'gunzip_func', + iterfield = ['in_file']) + + # Smooth - smoothing node + smooth = Node(Smooth(fwhm = self.fwhm), + name = 'smooth') + + # Funcion node get_subject_infos - get subject specific condition information + subject_infos = Node(Function( + function = self.get_subject_infos, + input_names = ['event_files', 'runs'], + output_names = ['subject_info']), + name = 'subject_infos') + + # SpecifyModel - generates SPM-specific Model + specify_model = Node(SpecifySPMModel( + concatenate_runs = True, input_units = 'secs', output_units = 'secs', + time_repetition = self.tr, high_pass_filter_cutoff = 128), + name='specify_model') + + # Level1Design - Generates an SPM design matrix + l1_design = Node(Level1Design( + bases = {'hrf': {'derivs': [0, 0]}}, timing_units = 'secs', + interscan_interval = self.tr), name='l1_design') + + # EstimateModel - estimate the parameters of the model + l1_estimate = Node(EstimateModel( + estimation_method={'Classical': 1}), + name='l1_estimate') + + # Function node get_contrasts - get the contrasts + contrasts = Node(Function( + function = self.get_contrasts, + input_names = [], + output_names = ['contrasts']), + name = 'contrasts') + + # Function node get_parameters_file - get parameters files + parameters = Node(Function( + function = self.get_parameters_file, + input_names = ['filepaths', 'subject_id', 'working_dir'], + output_names = ['parameters_file']), + name = 'parameters') + parameters.inputs.working_dir = self.directories.working_dir + + # EstimateContrast - estimates contrasts + contrast_estimate = Node(EstimateContrast(), + name = 'contrast_estimate') + + # Function node remove_gunzip_files - remove output of the gunzip node + remove_gunzip_files = Node(Function( + function = self.remove_gunzip_files, + input_names = ['_', 'subject_id', 'working_dir'], + output_names = []), + name = 'remove_gunzip_files') + remove_gunzip_files.inputs.working_dir = self.directories.working_dir + + # Function node remove_smoothed_files - remove output of the smooth node + remove_smoothed_files = Node(Function( + function = self.remove_smoothed_files, + input_names = ['_', 'subject_id', 'working_dir'], + output_names = []), + name = 'remove_smoothed_files') + remove_smoothed_files.inputs.working_dir = self.directories.working_dir + + # Create l1 analysis workflow and connect its nodes + l1_analysis = Workflow(base_dir = self.directories.working_dir, name = 'l1_analysis') + l1_analysis.connect([ + (infosource, selectfiles, [('subject_id', 'subject_id')]), + (infosource, subject_infos, [('exp_dir', 'exp_dir'), ('run_list', 'runs')]), + (infosource, contrasts, [('subject_id', 'subject_id')]), + (infosource, remove_gunzip_files, [('subject_id', 'subject_id')]), + (infosource, remove_smoothed_files, [('subject_id', 'subject_id')]), + (subject_infos, specify_model, [('subject_info', 'subject_info')]), + (contrasts, contrast_estimate, [('contrasts', 'contrasts')]), + (selectfiles, parameters, [('param', 'filepaths')]), + (selectfiles, subject_infos, [('event', 'event_files')]), + (infosource, parameters, [('subject_id', 'subject_id')]), + (selectfiles, gunzip_func, [('func', 'in_file')]), + (gunzip_func, smooth, [('out_file', 'in_files')]), + (smooth, remove_gunzip_files, [('smoothed_files', '_')]), + (smooth, specify_model, [('smoothed_files', 'functional_runs')]), + (parameters, specify_model, [('parameters_file', 'realignment_parameters')]), + (specify_model, l1_design, [('session_info', 'session_info')]), + (l1_design, l1_estimate, [('spm_mat_file', 'spm_mat_file')]), + (l1_estimate, contrast_estimate, [ + ('spm_mat_file', 'spm_mat_file'), + ('beta_images', 'beta_images'), + ('residual_image', 'residual_image')]), + (contrast_estimate, datasink, [ + ('con_images', 'l1_analysis.@con_images'), + ('spmT_images', 'l1_analysis.@spmT_images'), + ('spm_mat_file', 'l1_analysis.@spm_mat_file')]), + (contrast_estimate, remove_smoothed_files, [('spmT_images', '_')]) + ]) + + return l1_analysis + + def get_subset_contrasts(file_list, method, subject_list, participants_file): + """ + Parameters : + - file_list : original file list selected by selectfiles node + - subject_list : list of subject IDs that are in the wanted group for the analysis + - participants_file: str, file containing participants caracteristics + - method: str, one of 'equalRange', 'equalIndifference' or 'groupComp' + + This function return the file list containing only the files belonging + to subject in the wanted group. + """ + equal_indifference_id = [] + equal_range_id = [] + equal_indifference_files = [] + equal_range_files = [] + + with open(participants_file, 'rt') as file: + next(file) # skip the header + for line in file: info = line.strip().split() - - if info[0][-3:] in subject_list and info[1] == "equalIndifference": - equalIndifference_id.append(info[0][-3:]) - elif info[0][-3:] in subject_list and info[1] == "equalRange": - equalRange_id.append(info[0][-3:]) - - for file in file_list: - sub_id = file.split('/') - if sub_id[-1][-7:-4] in equalIndifference_id: - equalIndifference_files.append(file) - elif sub_id[-1][-7:-4] in equalRange_id: - equalRange_files.append(file) - - return equalIndifference_id, equalRange_id, equalIndifference_files, equalRange_files - - -def get_l2_analysis(subject_list, n_sub, contrast_list, method, exp_dir, output_dir, working_dir, result_dir, data_dir): - """ - Returns the 2nd level of analysis workflow. - - Parameters: - - exp_dir: str, directory where raw data are stored - - result_dir: str, directory where results will be stored - - working_dir: str, name of the sub-directory for intermediate results - - output_dir: str, name of the sub-directory for final results - - subject_list: list of str, list of subject for which you want to do the preprocessing - - contrast_list: list of str, list of contrasts to analyze - - n_sub: float, number of subjects used to do the analysis - - method: one of "equalRange", "equalIndifference" or "groupComp" - - Returns: - - l2_analysis: Nipype WorkFlow - """ - # Infosource - a function free node to iterate over the list of subject names - infosource_groupanalysis = Node(IdentityInterface(fields=['contrast_id', 'subjects'], - subjects = subject_list), - name="infosource_groupanalysis") - - infosource_groupanalysis.iterables = [('contrast_id', contrast_list)] - - # SelectFiles - contrast_file = opj(data_dir, 'NARPS-J7F9', "hypo_*_{contrast_id}_con_sub-*.nii") - - participants_file = opj(exp_dir, 'participants.tsv') - - templates = {'contrast' : contrast_file, 'participants' : participants_file} - - selectfiles_groupanalysis = Node(SelectFiles(templates, base_directory=result_dir, force_list= True), - name="selectfiles_groupanalysis") - - # Datasink node : to save important files - datasink_groupanalysis = Node(DataSink(base_directory = result_dir, container = output_dir), - name = 'datasink_groupanalysis') - - # Node to select subset of contrasts - sub_contrasts = Node(Function(input_names = ['file_list', 'method', 'subject_list', 'participants_file'], - output_names = ['equalIndifference_id', 'equalRange_id', 'equalIndifference_files', 'equalRange_files'], - function = get_subset_contrasts), - name = 'sub_contrasts') - - sub_contrasts.inputs.method = method - - ## Estimate model - estimate_model = Node(EstimateModel(estimation_method={'Classical':1}), name = "estimate_model") - - ## Estimate contrasts - estimate_contrast = Node(EstimateContrast(group_contrast=True), - name = "estimate_contrast") - - ## Create thresholded maps - threshold = MapNode(Threshold(use_fwe_correction=False, height_threshold = 0.001, - extent_fdr_p_threshold = 0.05, use_topo_fdr = False, force_activation = True), name = "threshold", iterfield = ["stat_image", "contrast_index"]) - - l2_analysis = Workflow(base_dir = opj(result_dir, working_dir), name = f"l2_analysis_{method}_nsub_{n_sub}") - - l2_analysis.connect([(infosource_groupanalysis, selectfiles_groupanalysis, [('contrast_id', 'contrast_id')]), - (infosource_groupanalysis, sub_contrasts, [('subjects', 'subject_list')]), - (selectfiles_groupanalysis, sub_contrasts, [('contrast', 'file_list'), ('participants', 'participants_file')]), - (estimate_model, estimate_contrast, [('spm_mat_file', 'spm_mat_file'), - ('residual_image', 'residual_image'), - ('beta_images', 'beta_images')]), - (estimate_contrast, threshold, [('spm_mat_file', 'spm_mat_file'), - ('spmT_images', 'stat_image')]), - (estimate_model, datasink_groupanalysis, [('mask_image', f"l2_analysis_{method}_nsub_{n_sub}.@mask")]), - (estimate_contrast, datasink_groupanalysis, [('spm_mat_file', f"l2_analysis_{method}_nsub_{n_sub}.@spm_mat"), - ('spmT_images', f"l2_analysis_{method}_nsub_{n_sub}.@T"), - ('con_images', f"l2_analysis_{method}_nsub_{n_sub}.@con")]), - (threshold, datasink_groupanalysis, [('thresholded_map', f"l2_analysis_{method}_nsub_{n_sub}.@thresh")])]) - - if method=='equalRange' or method=='equalIndifference': - contrasts = [('Group', 'T', ['mean'], [1]), ('Group', 'T', ['mean'], [-1])] - ## Specify design matrix - one_sample_t_test_design = Node(OneSampleTTestDesign(), name = "one_sample_t_test_design") - - l2_analysis.connect([(sub_contrasts, one_sample_t_test_design, [(f"{method}_files", 'in_files')]), - (one_sample_t_test_design, estimate_model, [('spm_mat_file', 'spm_mat_file')])]) - - threshold.inputs.contrast_index = [1, 2] - threshold.synchronize = True - - elif method == 'groupComp': - contrasts = [('Eq range vs Eq indiff in loss', 'T', ['Group_{1}', 'Group_{2}'], [1, -1])] - # Node for the design matrix - two_sample_t_test_design = Node(TwoSampleTTestDesign(unequal_variance=True), name = 'two_sample_t_test_design') - - l2_analysis.connect([(sub_contrasts, two_sample_t_test_design, [('equalRange_files', "group1_files"), - ('equalIndifference_files', 'group2_files')]), - (two_sample_t_test_design, estimate_model, [("spm_mat_file", "spm_mat_file")])]) - - threshold.inputs.contrast_index = [1] - threshold.synchronize = True - - estimate_contrast.inputs.contrasts = contrasts - - return l2_analysis - - -def reorganize_results(result_dir, output_dir, n_sub, team_ID): - """ - Reorganize the results to analyze them. - - Parameters: - - result_dir: str, directory where results will be stored - - output_dir: str, name of the sub-directory for final results - - n_sub: float, number of subject used for the analysis - - team_ID: str, ID of the team to reorganize results - - """ - from os.path import join as opj - import os - import shutil - import gzip - - h1 = opj(result_dir, output_dir, f"l2_analysis_equalIndifference_nsub_{n_sub}", '_contrast_id_gain') - h2 = opj(result_dir, output_dir, f"l2_analysis_equalRange_nsub_{n_sub}", '_contrast_id_gain') - h3 = opj(result_dir, output_dir, f"l2_analysis_equalIndifference_nsub_{n_sub}", '_contrast_id_gain') - h4 = opj(result_dir, output_dir, f"l2_analysis_equalRange_nsub_{n_sub}", '_contrast_id_gain') - h5 = opj(result_dir, output_dir, f"l2_analysis_equalIndifference_nsub_{n_sub}", '_contrast_id_loss') - h6 = opj(result_dir, output_dir, f"l2_analysis_equalRange_nsub_{n_sub}", '_contrast_id_loss') - h7 = opj(result_dir, output_dir, f"l2_analysis_equalIndifference_nsub_{n_sub}", '_contrast_id_loss') - h8 = opj(result_dir, output_dir, f"l2_analysis_equalRange_nsub_{n_sub}", '_contrast_id_loss') - h9 = opj(result_dir, output_dir, f"l2_analysis_groupComp_nsub_{n_sub}", '_contrast_id_loss') - - h = [h1, h2, h3, h4, h5, h6, h7, h8, h9] - - repro_unthresh = [opj(filename, "spmT_0002.nii") if i in [4, 5] else opj(filename, - "spmT_0001.nii") for i, filename in enumerate(h)] - - repro_thresh = [opj(filename, "_threshold1", - "spmT_0002_thr.nii") if i in [4, 5] else opj(filename, - "_threshold0", "spmT_0001_thr.nii") for i, filename in enumerate(h)] - - if not os.path.isdir(opj(result_dir, "NARPS-reproduction")): - os.mkdir(opj(result_dir, "NARPS-reproduction")) - - for i, filename in enumerate(repro_unthresh): - f_in = filename - f_out = opj(result_dir, "NARPS-reproduction", f"team_{team_ID}_nsub_{n_sub}_hypo{i+1}_unthresholded.nii") - shutil.copyfile(f_in, f_out) - - for i, filename in enumerate(repro_thresh): - f_in = filename - f_out = opj(result_dir, "NARPS-reproduction", f"team_{team_ID}_nsub_{n_sub}_hypo{i+1}_thresholded.nii") - shutil.copyfile(f_in, f_out) - - print(f"Results files of team {team_ID} reorganized.") - - + if info[0][-3:] in subject_list and info[1] == 'equalIndifference': + equal_indifference_id.append(info[0][-3:]) + elif info[0][-3:] in subject_list and info[1] == 'equalRange': + equal_range_id.append(info[0][-3:]) + + for file in file_list: + sub_id = file.split('/') + if sub_id[-1][-7:-4] in equal_indifference_id: + equal_indifference_files.append(file) + elif sub_id[-1][-7:-4] in equal_range_id: + equal_range_files.append(file) + + return equal_indifference_id, equal_range_id, equal_indifference_files, equal_range_files + + def reorganize_results(team_id, nb_sub, output_dir, results_dir): + """ + Reorganize the results to analyze them. + + Parameters: + - result_dir: str, directory where results will be stored + - output_dir: str, name of the sub-directory for final results + - nb_sub: float, number of subject used for the analysis + - team_id: str, ID of the team to reorganize results + + """ + from os import mkdir + from os.path import join, isdir + from shutil import copyfile + + hypotheses = [ + join(output_dir, f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_gain'), + join(output_dir, f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_gain'), + join(output_dir, f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_gain'), + join(output_dir, f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_gain'), + join(output_dir, f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_loss'), + join(output_dir, f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_loss'), + join(output_dir, f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_loss'), + join(output_dir, f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_loss'), + join(output_dir, f'l2_analysis_groupComp_nsub_{nb_sub}', '_contrast_id_loss') + ] + + # Build lists of files for unthresholded and thresholded maps + repro_unthresh = [] + repro_thresh = [] + for file_id, filename in enumerate(hypotheses): + if file_id in [4,5]: + repro_unthresh.append(join(filename, 'spmT_0002.nii')) + repro_thresh.append(join(filename, '_threshold1', 'spmT_0002_thr.nii')) + else: + repro_unthresh.append(join(filename, 'spmT_0001.nii')) + repro_thresh.append(join(filename, '_threshold0', 'spmT_0001_thr.nii')) + + if not isdir(join(results_dir, 'NARPS-reproduction')): + mkdir(join(results_dir, 'NARPS-reproduction')) + + for file_id, filename in enumerate(repro_unthresh): + f_in = filename + f_out = join(results_dir, + 'NARPS-reproduction', + f'team_{team_id}_nsub_{nb_sub}_hypo{file_id + 1}_unthresholded.nii') + copyfile(f_in, f_out) + + for file_id, filename in enumerate(repro_thresh): + f_in = filename + f_out = join(results_dir, + 'NARPS-reproduction', + f'team_{team_id}_nsub_{nb_sub}_hypo{file_id + 1}_thresholded.nii') + copyfile(f_in, f_out) + + print(f'Results files of team {team_id} reorganized.') + + def get_group_level_analysis(self): + """ + Return all workflows for the group level analysis. + + Returns; + - a list of nipype.WorkFlow + """ + + methods = ['equalRange', 'equalIndifference', 'groupComp'] + return [self.get_group_level_analysis_sub_workflow(method) for method in methods] + + def get_group_level_analysis_sub_workflow(self, method): + """ + Return a workflow for the group level analysis. + + Parameters: + - method: one of 'equalRange', 'equalIndifference' or 'groupComp' + + Returns: + - l2_analysis: nipype.WorkFlow + """ + # Compute the number of participants used to do the analysis + nb_subjects = len(self.subject_list) + + # Infosource - a function free node to iterate over the list of subject names + infosource_groupanalysis = Node( + IdentityInterface( + fields=['contrast_id', 'subjects'], + subjects = self.subject_list), + name='infosource_groupanalysis') + infosource_groupanalysis.iterables = [('contrast_id', self.contrast_list)] + + # SelectFiles + templates = { + # Contrasts for all participants + 'contrast' : join(self.directories.output_dir, + 'l1_analysis', '_subject_id_*', 'con_{contrast_id}.nii'), + # Participants file + 'participants' : join(self.directories.dataset_dir, 'participants.tsv') + } + + selectfiles_groupanalysis = Node(SelectFiles( + templates, + base_directory = self.directories.results_dir, + force_list = True), + name="selectfiles_groupanalysis") + + # Datasink - save important files + datasink_groupanalysis = Node(DataSink( + base_directory = str(self.directories.output_dir) + ), + name = 'datasink_groupanalysis') + + # Function node reorganize_results - organize results once computed + reorganize_res = Node(Function( + function = self.reorganize_results, + input_names = ['team_id', 'nb_subjects', 'results_dir', 'output_dir']), + name = 'reorganize_res') + reorganize_res.inputs.team_id = self.team_id + reorganize_res.inputs.nb_subjects = nb_subjects + reorganize_res.inputs.results_dir = self.directories.results_dir + reorganize_res.inputs.output_dir = self.directories.output_dir + + # Node to select subset of contrasts + sub_contrasts = Node(Function( + function = self.get_subset_contrasts, + input_names = ['file_list', 'method', 'subject_list', 'participants_file'], + output_names = [ + 'equalIndifference_id', + 'equalRange_id', + 'equalIndifference_files', + 'equalRange_files']), + name = 'sub_contrasts') + sub_contrasts.inputs.method = method + + # Estimate model + estimate_model = Node(EstimateModel( + estimation_method = {'Classical':1}), + name = 'estimate_model') + + # Estimate contrasts + estimate_contrast = Node(EstimateContrast( + group_contrast = True), + name = 'estimate_contrast') + + ## Create thresholded maps + threshold = MapNode(Threshold( + use_fwe_correction=False, height_threshold = 0.001, extent_fdr_p_threshold = 0.05, + use_topo_fdr = False, force_activation = True), + name = 'threshold', iterfield = ['stat_image', 'contrast_index']) + + l2_analysis = Workflow( + base_dir = self.directories.working_dir, + name = f'l2_analysis_{method}_nsub_{nb_subjects}') + l2_analysis.connect([ + (infosource_groupanalysis, selectfiles_groupanalysis, [ + ('contrast_id', 'contrast_id')]), + (infosource_groupanalysis, sub_contrasts, [ + ('subjects', 'subject_list')]), + (selectfiles_groupanalysis, sub_contrasts, [ + ('contrast', 'file_list'), + ('participants', 'participants_file')]), + (estimate_model, estimate_contrast, [ + ('spm_mat_file', 'spm_mat_file'), + ('residual_image', 'residual_image'), + ('beta_images', 'beta_images')]), + (estimate_contrast, threshold, [ + ('spm_mat_file', 'spm_mat_file'), + ('spmT_images', 'stat_image')]), + (estimate_model, datasink_groupanalysis, [ + ('mask_image', f'l2_analysis_{method}_nsub_{nb_subjects}.@mask')]), + (estimate_contrast, datasink_groupanalysis, [ + ('spm_mat_file', f'l2_analysis_{method}_nsub_{nb_subjects}.@spm_mat'), + ('spmT_images', f'l2_analysis_{method}_nsub_{nb_subjects}.@T'), + ('con_images', f'l2_analysis_{method}_nsub_{nb_subjects}.@con')]), + (threshold, datasink_groupanalysis, [ + ('thresholded_map', f'l2_analysis_{method}_nsub_{nb_subjects}.@thresh')])]) + + if method in ('equalRange', 'equalIndifference'): + contrasts = [('Group', 'T', ['mean'], [1]), ('Group', 'T', ['mean'], [-1])] + + # Specify design matrix + one_sample_t_test_design = Node(OneSampleTTestDesign(), + name = 'one_sample_t_test_design') + + threshold.inputs.contrast_index = [1, 2] + threshold.synchronize = True + + l2_analysis.connect([ + (sub_contrasts, one_sample_t_test_design, [(f'{method}_files', 'in_files')]), + (one_sample_t_test_design, estimate_model, [('spm_mat_file', 'spm_mat_file')])]) + + elif method == 'groupComp': + contrasts = [ + ('Eq range vs Eq indiff in loss', 'T', ['Group_{1}', 'Group_{2}'], [1, -1])] + + threshold.inputs.contrast_index = [1] + threshold.synchronize = True + + # Node for the design matrix + two_sample_t_test_design = Node(TwoSampleTTestDesign( + unequal_variance=True), + name = 'two_sample_t_test_design') + + l2_analysis.connect([ + (sub_contrasts, two_sample_t_test_design, [ + ('equalRange_files', 'group1_files'), + ('equalIndifference_files', 'group2_files')]), + (two_sample_t_test_design, estimate_model, [ + ('spm_mat_file', 'spm_mat_file')]) + ]) + + estimate_contrast.inputs.contrasts = contrasts + + return l2_analysis diff --git a/tests/pipelines/test_team_J7F9.py b/tests/pipelines/test_team_J7F9.py new file mode 100644 index 00000000..78819c6a --- /dev/null +++ b/tests/pipelines/test_team_J7F9.py @@ -0,0 +1,55 @@ +#!/usr/bin/python +# coding: utf-8 + +""" Tests of the 'narps_open.pipelines.team_J7F9' module. + +Launch this test with PyTest + +Usage: +====== + pytest -q test_team_J7F9.py + pytest -q test_team_J7F9.py -k +""" + +from statistics import mean + +from pytest import raises, helpers, mark +from nipype import Workflow + +from narps_open.pipelines.team_J7F9 import PipelineTeamJ7F9 + +class TestPipelinesTeamJ7F9: + """ A class that contains all the unit tests for the PipelineTeamJ7F9 class.""" + + @staticmethod + @mark.unit_test + def test_create(): + """ Test the creation of a PipelineTeamJ7F9 object """ + + pipeline = PipelineTeamJ7F9() + + # 1 - check the parameters + assert pipeline.fwhm == 8.0 + assert pipeline.team_id == 'J7F9' + + # 2 - check workflows + assert pipeline.get_preprocessing() is None + assert pipeline.get_run_level_analysis() is None + assert isinstance(pipeline.get_subject_level_analysis(), Workflow) + group_level = pipeline.get_group_level_analysis() + + assert len(group_level) == 3 + for sub_workflow in group_level: + assert isinstance(sub_workflow, Workflow) + + @staticmethod + @mark.pipeline_test + def test_execution(): + """ Test the execution of a PipelineTeamJ7F9 and compare results """ + results_4_subjects = helpers.test_pipeline( + 'J7F9', + '/references/', + '/data/', + '/output/', + 4) + assert mean(results_4_subjects) > .003 From 360dcd09d625995b79fc22499e652a22dc26a399 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 30 Aug 2023 11:27:30 +0200 Subject: [PATCH 02/25] [TEST] init a test for conftest.py --- tests/test_conftest.py | 139 +++++++++++++++++++++++++++++++++++++++++ 1 file changed, 139 insertions(+) create mode 100644 tests/test_conftest.py diff --git a/tests/test_conftest.py b/tests/test_conftest.py new file mode 100644 index 00000000..5b1a1f4f --- /dev/null +++ b/tests/test_conftest.py @@ -0,0 +1,139 @@ +#!/usr/bin/python +# coding: utf-8 + +""" Tests of the 'conftest.py' module. + +Launch this test with PyTest + +Usage: +====== + pytest -q test_conftest.py + pytest -q test_conftest.py -k +""" + +from os import remove +from os.path import join, isfile, abspath +from pathlib import Path + +from datetime import datetime + +from pytest import raises, mark + +from nipype import Node, Workflow +from nipype.interfaces.utility import Function + +from narps_open.utils.configuration import Configuration +from narps_open.runner import PipelineRunner +from narps_open.pipelines import Pipeline +from narps_open.pipelines.team_2T6S import PipelineTeam2T6S + +class MockupPipeline(Pipeline): + """ A simple Pipeline class for test purposes """ + + def __init__(self): + super().__init__() + self.test_file = abspath( + join(Configuration()['directories']['test_runs'], 'test_conftest.txt')) + if isfile(self.test_file): + remove(self.test_file) + + def __del__(self): + if isfile(self.test_file): + remove(self.test_file) + + # @staticmethod + def write_to_file(_, text_to_write: str, file_path: str): + """ Method used inside a nipype Node, to write a line in a test file """ + with open(file_path, 'a', encoding = 'utf-8') as file: + file.write(text_to_write) + + def create_workflow(self, workflow_name: str): + """ Return a nipype workflow with two nodes writing in a file """ + node_1 = Node(Function( + input_names = ['_', 'text_to_write', 'file_path'], + output_names = ['_'], + function = self.write_to_file), + name = 'node_1' + ) + # this input is set to now(), so that it changes at every run, thus preventing + # nipype's cache to work + node_1.inputs._ = datetime.now() + node_1.inputs.text_to_write = 'MockupPipeline : '+workflow_name+' node_1\n' + node_1.inputs.file_path = self.test_file + + node_2 = Node(Function( + input_names = ['_', 'text_to_write', 'file_path'], + output_names = [], + function = self.write_to_file), + name = 'node_2' + ) + node_2.inputs.text_to_write = 'MockupPipeline : '+workflow_name+' node_2\n' + node_2.inputs.file_path = self.test_file + + workflow = Workflow( + base_dir = Configuration()['directories']['test_runs'], + name = workflow_name + ) + workflow.add_nodes([node_1, node_2]) + workflow.connect(node_1, '_', node_2, '_') + + return workflow + + def get_preprocessing(self): + """ Return a fake preprocessing workflow """ + return self.create_workflow('TestPipelineRunner_preprocessing_workflow') + + def get_run_level_analysis(self): + """ Return a fake run level workflow """ + return self.create_workflow('TestPipelineRunner_run_level_workflow') + + def get_subject_level_analysis(self): + """ Return a fake subject level workflow """ + return self.create_workflow('TestPipelineRunner_subject_level_workflow') + + def get_group_level_analysis(self): + """ Return a fake subject level workflow """ + return self.create_workflow('TestPipelineRunner_group_level_workflow') + + def get_preprocessing_outputs(self): + """ Return a list of templates of the output files generated by the preprocessing """ + return [join(Configuration()['directories']['test_runs'], 'preprocessing_output.md')] + + def get_run_level_outputs(self): + """ Return a list of templates of the output files generated by the run level analysis. + Templates are expressed relatively to the self.directories.output_dir. + """ + return [join(Configuration()['directories']['test_runs'], 'run_output.md')] + + def get_subject_level_outputs(self): + """ Return a list of templates of the output files generated by the subject level analysis. + Templates are expressed relatively to the self.directories.output_dir. + """ + templates = [ + join(Configuration()['directories']['test_runs'], 'subject_{subject_id}_output_1.md'), + join(Configuration()['directories']['test_runs'], 'subject_{subject_id}_output_2.md') + ] + return_list = [] + for subject_id in self.subject_list: + return_list += [t.format(subject_id = subject_id) for t in templates] + + return return_list + + def get_group_level_outputs(self): + """ Return a list of templates of the output files generated by the group level analysis. + Templates are expressed relatively to the self.directories.output_dir. + """ + templates = [ + join(Configuration()['directories']['test_runs'], 'group_{nb_subjects}_output_a.md'), + join(Configuration()['directories']['test_runs'], 'group_{nb_subjects}_output_b.md') + ] + return [t.format(nb_subjects = len(self.subject_list)) for t in templates] + + def get_hypotheses_outputs(self): + """ Return the names of the files used by the team to answer the hypotheses of NARPS. + """ + template = join(Configuration()['directories']['test_runs'], 'hypothesis_{id}.md') + return [template.format(id = i) for i in range(1,18)] + +class TestPipelineRunner: + """ A class that contains all the unit tests for the PipelineRunner class.""" From 2b45a4808d2b70e15daf5be53a2ddb92cf35b80b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 30 Aug 2023 14:28:43 +0200 Subject: [PATCH 03/25] [REFAC] adding analyses outputs --- narps_open/pipelines/team_J7F9.py | 180 ++++++++++++++++++------------ 1 file changed, 108 insertions(+), 72 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 6023a367..314ab29e 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -24,7 +24,7 @@ def __init__(self): super().__init__() self.fwhm = 8.0 self.team_id = 'J7F9' - self.contrast_list = ['0001', '0002', '0003', '0004'] + self.contrast_list = ['0001', '0002', '0003'] def get_preprocessing(self): """ No preprocessing has been done by team J7F9 """ @@ -200,8 +200,8 @@ def remove_gunzip_files(_, subject_id, working_dir): Parameters: - _: Node input only used for triggering the Node - - subject_id: str, TODO - - working_id: str, TODO + - subject_id: str, id of th subject for which to remove the unzipped file + - working_dir: str, path to the working directory """ from shutil import rmtree from os.path import join @@ -220,8 +220,8 @@ def remove_smoothed_files(_, subject_id, working_dir): Parameters: - _: Node input only used for triggering the Node - - subject_id: str, TODO - - working_id: str, TODO + - subject_id: str, id of th subject for which to remove the smoothed file + - working_dir: str, path to the working directory """ from shutil import rmtree from os.path import join @@ -369,6 +369,33 @@ def get_subject_level_analysis(self): return l1_analysis + def get_subject_level_outputs(self): + """ Return the names of the files the subject level analysis is supposed to generate. """ + + # Contrat maps + templates = [join( + self.directories.output_dir, + 'l1_analysis', '_subject_id_{subject_id}', f'con_{contrast_id}.nii')\ + for contrast_id in self.contrast_list] + + # SPM.mat file + templates += [join( + self.directories.output_dir, + 'l1_analysis', '_subject_id_{subject_id}', 'SPM.mat')] + + # spmT maps + templates += [join( + self.directories.output_dir, + 'l1_analysis', '_subject_id_{subject_id}', f'spmT_{contrast_id}.nii')\ + for contrast_id in self.contrast_list] + + # Format with subject_ids + return_list = [] + for template in templates: + return_list += [template.format(subject_id = s) for s in self.subject_list] + + return return_list + def get_subset_contrasts(file_list, method, subject_list, participants_file): """ Parameters : @@ -403,63 +430,6 @@ def get_subset_contrasts(file_list, method, subject_list, participants_file): return equal_indifference_id, equal_range_id, equal_indifference_files, equal_range_files - def reorganize_results(team_id, nb_sub, output_dir, results_dir): - """ - Reorganize the results to analyze them. - - Parameters: - - result_dir: str, directory where results will be stored - - output_dir: str, name of the sub-directory for final results - - nb_sub: float, number of subject used for the analysis - - team_id: str, ID of the team to reorganize results - - """ - from os import mkdir - from os.path import join, isdir - from shutil import copyfile - - hypotheses = [ - join(output_dir, f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_gain'), - join(output_dir, f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_gain'), - join(output_dir, f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_gain'), - join(output_dir, f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_gain'), - join(output_dir, f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_loss'), - join(output_dir, f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_loss'), - join(output_dir, f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_loss'), - join(output_dir, f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_loss'), - join(output_dir, f'l2_analysis_groupComp_nsub_{nb_sub}', '_contrast_id_loss') - ] - - # Build lists of files for unthresholded and thresholded maps - repro_unthresh = [] - repro_thresh = [] - for file_id, filename in enumerate(hypotheses): - if file_id in [4,5]: - repro_unthresh.append(join(filename, 'spmT_0002.nii')) - repro_thresh.append(join(filename, '_threshold1', 'spmT_0002_thr.nii')) - else: - repro_unthresh.append(join(filename, 'spmT_0001.nii')) - repro_thresh.append(join(filename, '_threshold0', 'spmT_0001_thr.nii')) - - if not isdir(join(results_dir, 'NARPS-reproduction')): - mkdir(join(results_dir, 'NARPS-reproduction')) - - for file_id, filename in enumerate(repro_unthresh): - f_in = filename - f_out = join(results_dir, - 'NARPS-reproduction', - f'team_{team_id}_nsub_{nb_sub}_hypo{file_id + 1}_unthresholded.nii') - copyfile(f_in, f_out) - - for file_id, filename in enumerate(repro_thresh): - f_in = filename - f_out = join(results_dir, - 'NARPS-reproduction', - f'team_{team_id}_nsub_{nb_sub}_hypo{file_id + 1}_thresholded.nii') - copyfile(f_in, f_out) - - print(f'Results files of team {team_id} reorganized.') - def get_group_level_analysis(self): """ Return all workflows for the group level analysis. @@ -513,16 +483,6 @@ def get_group_level_analysis_sub_workflow(self, method): ), name = 'datasink_groupanalysis') - # Function node reorganize_results - organize results once computed - reorganize_res = Node(Function( - function = self.reorganize_results, - input_names = ['team_id', 'nb_subjects', 'results_dir', 'output_dir']), - name = 'reorganize_res') - reorganize_res.inputs.team_id = self.team_id - reorganize_res.inputs.nb_subjects = nb_subjects - reorganize_res.inputs.results_dir = self.directories.results_dir - reorganize_res.inputs.output_dir = self.directories.output_dir - # Node to select subset of contrasts sub_contrasts = Node(Function( function = self.get_subset_contrasts, @@ -615,3 +575,79 @@ def get_group_level_analysis_sub_workflow(self, method): estimate_contrast.inputs.contrasts = contrasts return l2_analysis + + def get_group_level_outputs(self): + """ Return all names for the files the group level analysis is supposed to generate. """ + + # Handle equalRange and equalIndifference + parameters = { + 'contrast_id': self.contrast_list, + 'method': ['equalRange', 'equalIndifference'], + 'file': [ + 'con_0001.nii', 'con_0002.nii', 'mask.nii', 'SPM.mat', + 'spmT_0001.nii', 'spmT_0002.nii', + join('_threshold0', 'spmT_0001_thr.nii'), join('_threshold1', 'spmT_0002_thr.nii') + ], + 'nb_subjects' : [str(len(self.subject_list))] + } + parameter_sets = product(*parameters.values()) + template = join( + self.directories.output_dir, + 'l2_analysis_{method}_nsub_{nb_subjects}', + '_contrast_id_{contrast_id}', + '{file}' + ) + + return_list = [template.format(**dict(zip(parameters.keys(), parameter_values)))\ + for parameter_values in parameter_sets] + + # Handle groupComp + parameters = { + 'contrast_id': self.contrast_list, + 'method': ['groupComp'], + 'file': [ + 'con_0001.nii', 'mask.nii', 'SPM.mat', 'spmT_0001.nii', + join('_threshold0', 'spmT_0001_thr.nii') + ], + 'nb_subjects' : [str(len(self.subject_list))] + } + parameter_sets = product(*parameters.values()) + template = join( + self.directories.output_dir, + 'l2_analysis_{method}_nsub_{nb_subjects}', + '_contrast_id_{contrast_id}', + '{file}' + ) + + return_list += [template.format(**dict(zip(parameters.keys(), parameter_values)))\ + for parameter_values in parameter_sets] + + return return_list + + def get_hypotheses_outputs(self): + """ Return all hypotheses output file names. + Note that hypotheses 5 to 8 correspond to the maps given by the team in their results ; + but they are not fully consistent with the hypotheses definitions as expected by NARPS. + """ + nb_sub = len(self.subject_list) + files = [ + join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), + join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0002', 'spmT_0001.nii'), + join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), + join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0002', 'spmT_0001.nii'), + join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), + join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0002', 'spmT_0001.nii'), + join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), + join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0002', 'spmT_0001.nii'), + join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0003', '_threshold1', 'spmT_0002_thr.nii'), + join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0002.nii'), + join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0003', '_threshold1', 'spmT_0001_thr.nii'), + join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0001.nii'), + join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0003', '_threshold0', 'spmT_0001_thr.nii'), + join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0001.nii'), + join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0003', '_threshold0', 'spmT_0002_thr.nii'), + join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0002.nii'), + join(f'l2_analysis_groupComp_nsub_{nb_sub}', '_contrast_id_0003', '_threshold0', 'spmT_0001_thr.nii'), + join(f'l2_analysis_groupComp_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0001.nii') + ] + return [join(self.directories.output_dir, f) for f in files] From ae59600a9154b60b08ff736854b9b8287734477c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 30 Aug 2023 14:39:51 +0200 Subject: [PATCH 04/25] [TEST] updating tests for J7F9 --- narps_open/pipelines/__init__.py | 2 +- tests/pipelines/test_pipelines.py | 9 ++++----- tests/pipelines/test_team_J7F9.py | 33 +++++++++++++++++++++---------- 3 files changed, 28 insertions(+), 16 deletions(-) diff --git a/narps_open/pipelines/__init__.py b/narps_open/pipelines/__init__.py index c3834fb6..e4b27f00 100644 --- a/narps_open/pipelines/__init__.py +++ b/narps_open/pipelines/__init__.py @@ -51,7 +51,7 @@ 'I52Y': None, 'I9D6': None, 'IZ20': None, - 'J7F9': None, + 'J7F9': 'PipelineTeamJ7F9', 'K9P0': None, 'L1A8': None, 'L3V8': None, diff --git a/tests/pipelines/test_pipelines.py b/tests/pipelines/test_pipelines.py index 9016aeb7..089f7e4c 100644 --- a/tests/pipelines/test_pipelines.py +++ b/tests/pipelines/test_pipelines.py @@ -131,12 +131,11 @@ def test_create(): class TestUtils: """ A class that contains all the unit tests for the utils in module pipelines.""" - @staticmethod @mark.unit_test def test_utils(): """ Test the utils methods of PipelineRunner """ - # 1 - Get number of not implemented pipelines - assert len(get_not_implemented_pipelines()) == 69 + # 1 - Get not implemented pipelines + assert '1K0E' in get_not_implemented_pipelines() - # 2 - Get number of implemented pipelines - assert len(get_implemented_pipelines()) == 1 + # 2 - Get implemented pipelines + assert '2T6S' in get_implemented_pipelines() diff --git a/tests/pipelines/test_team_J7F9.py b/tests/pipelines/test_team_J7F9.py index 78819c6a..da9023c4 100644 --- a/tests/pipelines/test_team_J7F9.py +++ b/tests/pipelines/test_team_J7F9.py @@ -11,9 +11,7 @@ pytest -q test_team_J7F9.py -k """ -from statistics import mean - -from pytest import raises, helpers, mark +from pytest import helpers, mark from nipype import Workflow from narps_open.pipelines.team_J7F9 import PipelineTeamJ7F9 @@ -42,14 +40,29 @@ def test_create(): for sub_workflow in group_level: assert isinstance(sub_workflow, Workflow) + @staticmethod + @mark.unit_test + def test_outputs(): + """ Test the expected outputs of a PipelineTeamJ7F9 object """ + pipeline = PipelineTeamJ7F9() + # 1 - 1 subject outputs + pipeline.subject_list = ['001'] + assert len(pipeline.get_preprocessing_outputs()) == 0 + assert len(pipeline.get_run_level_outputs()) == 0 + assert len(pipeline.get_subject_level_outputs()) == 7 + assert len(pipeline.get_group_level_outputs()) == 63 + assert len(pipeline.get_hypotheses_outputs()) == 18 + + # 2 - 4 subjects outputs + pipeline.subject_list = ['001', '002', '003', '004'] + assert len(pipeline.get_preprocessing_outputs()) == 0 + assert len(pipeline.get_run_level_outputs()) == 0 + assert len(pipeline.get_subject_level_outputs()) == 28 + assert len(pipeline.get_group_level_outputs()) == 63 + assert len(pipeline.get_hypotheses_outputs()) == 18 + @staticmethod @mark.pipeline_test def test_execution(): """ Test the execution of a PipelineTeamJ7F9 and compare results """ - results_4_subjects = helpers.test_pipeline( - 'J7F9', - '/references/', - '/data/', - '/output/', - 4) - assert mean(results_4_subjects) > .003 + helpers.test_pipeline_evaluation('J7F9') From 4b7e8314151ee43036cf4202904ca333c926f157 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 30 Aug 2023 14:55:17 +0200 Subject: [PATCH 05/25] Cleaning connections --- narps_open/pipelines/team_J7F9.py | 17 ++++++----------- 1 file changed, 6 insertions(+), 11 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 314ab29e..4c5bda6b 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -4,6 +4,7 @@ """ Write the work of NARPS' team J7F9 using Nipype """ from os.path import join +from itertools import product from nipype import Workflow, Node, MapNode from nipype.interfaces.utility import IdentityInterface, Function @@ -128,7 +129,7 @@ def get_subject_infos(event_files, runs): return subject_info - def get_contrasts(subject_id): + def get_contrasts(): """ Create the list of tuples that represents contrasts. Each contrast is in the form : @@ -145,7 +146,6 @@ def get_contrasts(subject_id): # Return contrast list return [trial, effect_gain, effect_loss] - def get_parameters_file(filepaths, subject_id, working_dir): """ Create new tsv files with only desired parameters per subject per run. @@ -242,11 +242,7 @@ def get_subject_level_analysis(self): """ # Infosource Node - To iterate on subjects infosource = Node(IdentityInterface( - fields = ['subject_id', 'exp_dir', 'result_dir', 'working_dir', 'run_list'], - dataset_dir = self.directories.dataset_dir, - results_dir = self.directories.results_dir, - working_dir = self.directories.working_dir, - run_list = self.run_list), + fields = ['subject_id']), name = 'infosource') infosource.iterables = [('subject_id', self.subject_list)] @@ -278,12 +274,13 @@ def get_subject_level_analysis(self): smooth = Node(Smooth(fwhm = self.fwhm), name = 'smooth') - # Funcion node get_subject_infos - get subject specific condition information + # Function node get_subject_infos - get subject specific condition information subject_infos = Node(Function( function = self.get_subject_infos, input_names = ['event_files', 'runs'], output_names = ['subject_info']), name = 'subject_infos') + subject_infos.inputs.runs = self.run_list # SpecifyModel - generates SPM-specific Model specify_model = Node(SpecifySPMModel( @@ -340,8 +337,6 @@ def get_subject_level_analysis(self): l1_analysis = Workflow(base_dir = self.directories.working_dir, name = 'l1_analysis') l1_analysis.connect([ (infosource, selectfiles, [('subject_id', 'subject_id')]), - (infosource, subject_infos, [('exp_dir', 'exp_dir'), ('run_list', 'runs')]), - (infosource, contrasts, [('subject_id', 'subject_id')]), (infosource, remove_gunzip_files, [('subject_id', 'subject_id')]), (infosource, remove_smoothed_files, [('subject_id', 'subject_id')]), (subject_infos, specify_model, [('subject_info', 'subject_info')]), @@ -401,7 +396,7 @@ def get_subset_contrasts(file_list, method, subject_list, participants_file): Parameters : - file_list : original file list selected by selectfiles node - subject_list : list of subject IDs that are in the wanted group for the analysis - - participants_file: str, file containing participants caracteristics + - participants_file: str, file containing participants characteristics - method: str, one of 'equalRange', 'equalIndifference' or 'groupComp' This function return the file list containing only the files belonging From 2d1739ad21c3067e4e456487b1b20dd61287c127 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 30 Aug 2023 15:01:43 +0200 Subject: [PATCH 06/25] [TEST] bug with pipelines/test_utils --- tests/pipelines/test_pipelines.py | 1 + 1 file changed, 1 insertion(+) diff --git a/tests/pipelines/test_pipelines.py b/tests/pipelines/test_pipelines.py index 089f7e4c..c38cf36a 100644 --- a/tests/pipelines/test_pipelines.py +++ b/tests/pipelines/test_pipelines.py @@ -131,6 +131,7 @@ def test_create(): class TestUtils: """ A class that contains all the unit tests for the utils in module pipelines.""" + @staticmethod @mark.unit_test def test_utils(): """ Test the utils methods of PipelineRunner """ From 678adfe3bbc355ada9ef047b9233723d5b78f455 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Thu, 31 Aug 2023 09:22:26 +0200 Subject: [PATCH 07/25] [BUG] issue with iterating over dicts --- narps_open/pipelines/team_J7F9.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 4c5bda6b..071d235e 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -87,11 +87,11 @@ def get_subject_infos(event_files, runs): onset[val].append(float(info[0])) duration[val].append(float(0)) - for gain_key, gain_value in weights_gain: + for gain_key, gain_value in weights_gain.items(): gain_value = gain_value - mean(gain_value) weights_gain[gain_key] = gain_value.tolist() - for loss_key, loss_value in weights_loss: + for loss_key, loss_value in weights_loss.items(): loss_value = loss_value - mean(loss_value) weights_loss[loss_key] = loss_value.tolist() From 4bd3ab53a6bf0a03ec30f76dcd1c109852a7479d Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Thu, 31 Aug 2023 12:11:56 +0200 Subject: [PATCH 08/25] [BUG] replace mkdir by makedirs in get_parameters --- narps_open/pipelines/team_J7F9.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 071d235e..e4e9fcc4 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -158,7 +158,7 @@ def get_parameters_file(filepaths, subject_id, working_dir): Return : - parameters_file : paths to new files containing only desired parameters. """ - from os import mkdir + from os import makedirs from os.path import join, isdir from pandas import DataFrame, read_csv from numpy import array, transpose @@ -182,8 +182,7 @@ def get_parameters_file(filepaths, subject_id, working_dir): new_path = join(working_dir, 'parameters_file', f'parameters_file_sub-{subject_id}_run-{str(file_id + 1).zfill(2)}.tsv') - if not isdir(join(working_dir, 'parameters_file')): - mkdir(join(working_dir, 'parameters_file')) + makedirs(join(working_dir, 'parameters_file'), exist_ok = True) with open(new_path, 'w') as writer: writer.write(retained_parameters.to_csv( From 6986e90351bca4071ea61387bee1f52f1125aceb Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Thu, 31 Aug 2023 14:35:10 +0200 Subject: [PATCH 09/25] [BUG] inside unit_tests workflow --- .github/workflows/unit_tests.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/unit_tests.yml b/.github/workflows/unit_tests.yml index 20f20ea3..d0097882 100644 --- a/.github/workflows/unit_tests.yml +++ b/.github/workflows/unit_tests.yml @@ -34,7 +34,7 @@ jobs: - name: Checkout repository uses: actions/checkout@v3 - - name: Load configuration for self-hosted runner + - name: Load configuration for self-hosted runner run: cp /home/neuro/local_testing_config.toml narps_open/utils/configuration/testing_config.toml - name: Install dependencies From 0f8365a1e846b3b3a671fc866b50a7c6066f5e71 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Fri, 1 Sep 2023 18:03:41 +0200 Subject: [PATCH 10/25] [TEST] testing the conftest module --- narps_open/data/results/__init__.py | 6 +- tests/conftest.py | 13 +- tests/data/test_results.py | 4 +- tests/test_conftest.py | 254 +++++++++++++++++++++------- 4 files changed, 205 insertions(+), 72 deletions(-) diff --git a/narps_open/data/results/__init__.py b/narps_open/data/results/__init__.py index a1108034..c36e6106 100644 --- a/narps_open/data/results/__init__.py +++ b/narps_open/data/results/__init__.py @@ -56,7 +56,7 @@ def get_uid(self): def get_file_urls(self): """ Return a dict containing the download url for each file of the collection. - * dict key is the file base name (with no extension) + * dict key is the file base name (with extension) * dict value is the download url for the file on Neurovault """ @@ -69,7 +69,7 @@ def get_file_urls(self): file_urls = {} for result in json['results']: # Get data for a file in the collection - file_urls[result['name']] = result['file'] + file_urls[result['name']+'.nii.gz'] = result['file'] return file_urls @@ -84,7 +84,7 @@ def download(self): for file_name, file_url in self.files.items(): urlretrieve( file_url, - join(self.directory, file_name+'.nii.gz'), + join(self.directory, file_name), show_download_progress ) diff --git a/tests/conftest.py b/tests/conftest.py index d740eeac..5ca91a07 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -58,7 +58,9 @@ def test_pipeline_execution( # Get missing subjects missing_subjects = set() for file in runner.get_missing_first_level_outputs(): - missing_subjects.add(get_subject_id(file)) + subject_id = get_subject_id(file) + if subject_id is not None: + missing_subjects.add(subject_id) # Leave if no missing subjects if not missing_subjects: @@ -66,7 +68,10 @@ def test_pipeline_execution( # Start pipeline runner.subjects = missing_subjects - runner.start(True, False) + try: # This avoids errors in the workflow to make the test fail + runner.start(True, False) + except(RuntimeError) as err: + print('RuntimeError: ', err) # Check missing files for the last time missing_files = runner.get_missing_first_level_outputs() @@ -80,7 +85,6 @@ def test_pipeline_execution( # Indices and keys to the unthresholded maps indices = list(range(1, 18, 2)) - keys = [f'hypo{i}_unthresh.nii.gz' for i in range(1, 10)] # Retrieve the paths to the reproduced files reproduced_files = runner.pipeline.get_hypotheses_outputs() @@ -88,7 +92,8 @@ def test_pipeline_execution( # Retrieve the paths to the results files collection = ResultsCollection(team_id) - results_files = [join(collection.directory, collection.files[k]) for k in keys] + results_files = [join(collection.directory, f) for f in collection.files.keys()] + results_files = [results_files[i] for i in indices] # Compute the correlation coefficients return [ diff --git a/tests/data/test_results.py b/tests/data/test_results.py index b32abeed..46d2d1f5 100644 --- a/tests/data/test_results.py +++ b/tests/data/test_results.py @@ -48,7 +48,7 @@ def test_create(): assert collection.uid == '4881' assert 'results/orig/4881_2T6S' in collection.directory test_str = 'http://neurovault.org/media/images/4881/hypo1_thresholded_revised.nii.gz' - assert collection.files['hypo1_thresh'] == test_str + assert collection.files['hypo1_thresh.nii.gz'] == test_str collection = ResultsCollection('43FJ') assert collection.team_id == '43FJ' @@ -56,7 +56,7 @@ def test_create(): assert 'results/orig/4824_43FJ' in collection.directory test_str = 'http://neurovault.org/media/images/4824/' test_str += 'Zstat_Thresholded_Negative_Effect_of_Loss_Equal_Indifference.nii.gz' - assert collection.files['hypo5_thresh'] == test_str + assert collection.files['hypo5_thresh.nii.gz'] == test_str @staticmethod @mark.unit_test diff --git a/tests/test_conftest.py b/tests/test_conftest.py index 5b1a1f4f..4dfce4d2 100644 --- a/tests/test_conftest.py +++ b/tests/test_conftest.py @@ -11,13 +11,13 @@ pytest -q test_conftest.py -k """ -from os import remove -from os.path import join, isfile, abspath -from pathlib import Path +from os import makedirs +from os.path import join, abspath, isdir, isfile +from shutil import rmtree from datetime import datetime -from pytest import raises, mark +from pytest import raises, mark, helpers, fixture from nipype import Node, Workflow from nipype.interfaces.utility import Function @@ -25,115 +25,243 @@ from narps_open.utils.configuration import Configuration from narps_open.runner import PipelineRunner from narps_open.pipelines import Pipeline -from narps_open.pipelines.team_2T6S import PipelineTeam2T6S +from narps_open.data.results import ResultsCollection + +TEST_DIR = abspath(join(Configuration()['directories']['test_runs'], 'test_conftest')) + +@fixture +def set_test_directory(scope = 'function'): + rmtree(TEST_DIR, ignore_errors = True) + makedirs(TEST_DIR, exist_ok = True) + + yield + + # Comment this line for debugging + #rmtree(TEST_DIR, ignore_errors = True) class MockupPipeline(Pipeline): """ A simple Pipeline class for test purposes """ def __init__(self): super().__init__() - self.test_file = abspath( - join(Configuration()['directories']['test_runs'], 'test_conftest.txt')) - if isfile(self.test_file): - remove(self.test_file) - - def __del__(self): - if isfile(self.test_file): - remove(self.test_file) - - # @staticmethod - def write_to_file(_, text_to_write: str, file_path: str): - """ Method used inside a nipype Node, to write a line in a test file """ - with open(file_path, 'a', encoding = 'utf-8') as file: - file.write(text_to_write) - - def create_workflow(self, workflow_name: str): - """ Return a nipype workflow with two nodes writing in a file """ - node_1 = Node(Function( - input_names = ['_', 'text_to_write', 'file_path'], - output_names = ['_'], - function = self.write_to_file), - name = 'node_1' + self.test_file = join(TEST_DIR, 'test_conftest.txt') + + # Init the test_file : write a number of execution set to zero + with open(self.test_file, 'w', encoding = 'utf-8') as file: + file.write(str(0)) + + def update_execution_count(_, file_path: str): + """ Method used inside a nipype Node, to update the execution count inside the file. + Arguments: + - file_path:str, path to the execution count file + + Return: the updated number of executions + """ + + # Read file counter + execution_counter = 0 + with open(file_path, 'r', encoding = 'utf-8') as file: + execution_counter = int(file.readline()) + + execution_counter += 1 + + # Write execution count back + with open(file_path, 'w', encoding = 'utf-8') as file: + file.write(str(execution_counter)) + + return execution_counter + + def decide_exception(execution_counter: int): + """ Method used inside a nipype Node, to simulate an exception during one run """ + if execution_counter == 1: + raise AttributeError + + def write_files(_, file_list: list, execution_counter: int): + """ Method used inside a nipype Node, to create a set of files """ + from pathlib import Path + + if execution_counter != 0: + for file_path in file_list: + Path(file_path).touch() + + def create_workflow(self, workflow_name: str, file_list: list): + """ Return a nipype workflow with two nodes. + First node updates the number of executions of the workflow? + Second node produces an exception, every two executions of the workflow. + Third node writes the output files, except once every three executions + of the workflow. + + Arguments: + - workflow_name: str, the name of the workflow to create + - file_list: list, list of the files that the workflow is supposed to generate + """ + node_count = Node(Function( + input_names = ['_', 'file_path'], + output_names = ['execution_counter'], + function = self.update_execution_count), + name = 'node_count' ) # this input is set to now(), so that it changes at every run, thus preventing # nipype's cache to work - node_1.inputs._ = datetime.now() - node_1.inputs.text_to_write = 'MockupPipeline : '+workflow_name+' node_1\n' - node_1.inputs.file_path = self.test_file + node_count.inputs._ = datetime.now() + node_count.inputs.file_path = self.test_file + + node_decide = Node(Function( + input_names = ['execution_counter'], + output_names = ['_'], + function = self.decide_exception), + name = 'node_decide' + ) - node_2 = Node(Function( - input_names = ['_', 'text_to_write', 'file_path'], + node_files = Node(Function( + input_names = ['_', 'file_list', 'execution_counter'], output_names = [], - function = self.write_to_file), - name = 'node_2' + function = self.write_files), + name = 'node_files' ) - node_2.inputs.text_to_write = 'MockupPipeline : '+workflow_name+' node_2\n' - node_2.inputs.file_path = self.test_file + node_files.inputs.file_list = file_list workflow = Workflow( - base_dir = Configuration()['directories']['test_runs'], + base_dir = TEST_DIR, name = workflow_name ) - workflow.add_nodes([node_1, node_2]) - workflow.connect(node_1, '_', node_2, '_') + workflow.add_nodes([node_count, node_decide, node_files]) + workflow.connect(node_count, 'execution_counter', node_decide, 'execution_counter') + workflow.connect(node_count, 'execution_counter', node_files, 'execution_counter') + workflow.connect(node_decide, '_', node_files, '_') return workflow def get_preprocessing(self): """ Return a fake preprocessing workflow """ - return self.create_workflow('TestPipelineRunner_preprocessing_workflow') + return self.create_workflow( + 'TestConftest_preprocessing_workflow', + self.get_preprocessing_outputs() + ) def get_run_level_analysis(self): """ Return a fake run level workflow """ - return self.create_workflow('TestPipelineRunner_run_level_workflow') + return self.create_workflow( + 'TestConftest_run_level_workflow', + self.get_run_level_outputs() + ) def get_subject_level_analysis(self): """ Return a fake subject level workflow """ - return self.create_workflow('TestPipelineRunner_subject_level_workflow') + return self.create_workflow( + 'TestConftest_subject_level_workflow', + self.get_subject_level_outputs() + ) def get_group_level_analysis(self): """ Return a fake subject level workflow """ - return self.create_workflow('TestPipelineRunner_group_level_workflow') + return self.create_workflow( + 'TestConftest_group_level_workflow', + self.get_group_level_outputs() + ) def get_preprocessing_outputs(self): """ Return a list of templates of the output files generated by the preprocessing """ - return [join(Configuration()['directories']['test_runs'], 'preprocessing_output.md')] + template = join(TEST_DIR, 'subject_id_{subject_id}_output_preprocessing_1.md') + return [template.format(subject_id = s) for s in self.subject_list] def get_run_level_outputs(self): """ Return a list of templates of the output files generated by the run level analysis. Templates are expressed relatively to the self.directories.output_dir. """ - return [join(Configuration()['directories']['test_runs'], 'run_output.md')] + template = join(TEST_DIR, 'subject_id_{subject_id}_output_run_1.md') + return [template.format(subject_id = s) for s in self.subject_list] def get_subject_level_outputs(self): """ Return a list of templates of the output files generated by the subject level analysis. Templates are expressed relatively to the self.directories.output_dir. """ - templates = [ - join(Configuration()['directories']['test_runs'], 'subject_{subject_id}_output_1.md'), - join(Configuration()['directories']['test_runs'], 'subject_{subject_id}_output_2.md') - ] - return_list = [] - for subject_id in self.subject_list: - return_list += [t.format(subject_id = subject_id) for t in templates] - - return return_list + template = join(TEST_DIR, 'subject_id_{subject_id}_output_analysis_1.md') + return [template.format(subject_id = s) for s in self.subject_list] def get_group_level_outputs(self): """ Return a list of templates of the output files generated by the group level analysis. Templates are expressed relatively to the self.directories.output_dir. """ templates = [ - join(Configuration()['directories']['test_runs'], 'group_{nb_subjects}_output_a.md'), - join(Configuration()['directories']['test_runs'], 'group_{nb_subjects}_output_b.md') + join(TEST_DIR, 'group_{nb_subjects}_output_a.md'), + join(TEST_DIR, 'group_{nb_subjects}_output_b.md') ] - return [t.format(nb_subjects = len(self.subject_list)) for t in templates] + return_list = [t.format(nb_subjects = len(self.subject_list)) for t in templates] + + template = join(TEST_DIR, 'hypothesis_{id}.md') + return_list += [template.format(id = i) for i in range(1,19)] + + return return_list def get_hypotheses_outputs(self): - """ Return the names of the files used by the team to answer the hypotheses of NARPS. - """ - template = join(Configuration()['directories']['test_runs'], 'hypothesis_{id}.md') - return [template.format(id = i) for i in range(1,18)] + """ Return the names of the files used by the team to answer the hypotheses of NARPS. """ + template = join(TEST_DIR, 'hypothesis_{id}.md') + return [template.format(id = i) for i in range(1,19)] + +class TestConftest: + """ A class that contains all the unit tests for the conftest module.""" + + @staticmethod + @mark.unit_test + def test_test_correlation_results(mocker): + """ Test the test_correlation_result helper """ + + mocker.patch( + 'conftest.Configuration', + return_value = { + 'testing': { + 'pipelines': { + 'correlation_thresholds' : [0.30, 0.70, 0.80, 0.85, 0.93] + } + } + } + ) + + assert helpers.test_correlation_results( + [0.35, 0.35, 0.36, 0.37, 0.38, 0.39, 0.99, 0.82, 0.40], 20) + assert helpers.test_correlation_results( + [0.75, 0.75, 0.76, 0.77, 0.88, 0.79, 0.99, 0.82, 0.80], 40) + assert helpers.test_correlation_results( + [0.85, 0.85, 0.86, 0.87, 0.98, 0.99, 0.99, 0.82, 0.81], 60) + assert helpers.test_correlation_results( + [0.86, 0.86, 0.86, 0.87, 0.88, 0.89, 0.99, 0.92, 0.95], 80) + assert helpers.test_correlation_results( + [0.95, 0.95, 0.96, 0.97, 0.98, 0.99, 0.99, 0.95, 1.0], 108) + assert not helpers.test_correlation_results( + [0.3, 0.29, 0.3, 0.3, 0.3, 0.39, 0.99, 0.82, 0.40], 20) + assert not helpers.test_correlation_results( + [0.60, 0.7, 0.7, 0.7, 0.8, 0.79, 0.99, 0.82, 0.80], 40) + assert not helpers.test_correlation_results( + [0.8, 0.79, 0.8, 0.8, 0.9, 0.99, 0.99, 0.82, 0.81], 60) + assert not helpers.test_correlation_results( + [0.8, 0.8, 0.8, 0.8, 0.88, 0.89, 0.99, 0.92, 0.95], 80) + assert not helpers.test_correlation_results( + [0.99, 0.99, 0.99, 1., 1., 1., -1., 0.95, 1.0], 108) + + @staticmethod + @mark.unit_test + def test_test_pipeline_execution(mocker, set_test_directory): + """ Test the test_pipeline_execution helper """ + + # Create mocks + mocker.patch('conftest.get_correlation_coefficient', return_value = 1.0) + fake_runner = PipelineRunner('2T6S') + fake_runner._pipeline = MockupPipeline() + mocker.patch('conftest.PipelineRunner', return_value = fake_runner) + fake_collection = ResultsCollection('2T6S') + mocker.patch('conftest.ResultsCollection', return_value = fake_collection) + + # Run pipeline + helpers.test_pipeline_execution('test_conftest', 20) + + # Check outputs + assert isdir(join(TEST_DIR, 'TestConftest_preprocessing_workflow')) + assert isdir(join(TEST_DIR, 'TestConftest_run_level_workflow')) + assert isdir(join(TEST_DIR, 'TestConftest_subject_level_workflow')) + assert isdir(join(TEST_DIR, 'TestConftest_group_level_workflow')) -class TestPipelineRunner: - """ A class that contains all the unit tests for the PipelineRunner class.""" + @staticmethod + @mark.unit_test + def test_test_pipeline_evaluation(): + """ Test the test_pipeline_evaluation helper """ From afae43ed416ab0be97e3878e7b5db7714fbc3deb Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 3 Jan 2024 14:20:26 +0100 Subject: [PATCH 11/25] [REFAC] [skip ci] --- narps_open/pipelines/team_J7F9.py | 599 ++++++++---------- .../utils/configuration/testing_config.toml | 4 +- tests/pipelines/test_team_J7F9.py | 48 ++ tests/test_data/pipelines/events.tsv | 6 + 4 files changed, 334 insertions(+), 323 deletions(-) create mode 100644 tests/test_data/pipelines/events.tsv diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index e4e9fcc4..b3e85996 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -17,6 +17,11 @@ from nipype.algorithms.misc import Gunzip from narps_open.pipelines import Pipeline +from narps_open.data.task import TaskInformation +from narps_open.data.participants import get_group +from narps_open.core.common import ( + remove_file, list_intersection, elements_in_string, clean_list + ) class PipelineTeamJ7F9(Pipeline): """ A class that defines the pipeline of team J7F9. """ @@ -26,6 +31,11 @@ def __init__(self): self.fwhm = 8.0 self.team_id = 'J7F9' self.contrast_list = ['0001', '0002', '0003'] + self.subject_level_contrasts = [ + ['trial', 'T', ['trial', 'trialxgain^1', 'trialxloss^1'], [1, 0, 0]], + ['effect_of_gain', 'T', ['trial', 'trialxgain^1', 'trialxloss^1'], [0, 1, 0]], + ['effect_of_loss', 'T', ['trial', 'trialxgain^1', 'trialxloss^1'], [0, 0, 1]] + ] def get_preprocessing(self): """ No preprocessing has been done by team J7F9 """ @@ -35,91 +45,71 @@ def get_run_level_analysis(self): """ No run level analysis has been done by team J7F9 """ return None - def get_subject_infos(event_files, runs): + def get_subject_information(event_files): """ Create Bunchs for specifySPMModel. Parameters : - event_files: list of str, list of events files (one per run) for the subject - - runs: list of str, list of runs to use Returns : - - subject_info : list of Bunch for 1st level analysis. + - subject_information : list of Bunch for subject level analysis. """ from numpy import mean from nipype.interfaces.base import Bunch - condition_names = ['trial', 'missed'] - onset = {} - duration = {} - weights_gain = {} - weights_loss = {} + subject_information = [] - for run_id in range(len(runs)): # Loop over number of runs. - # creates dictionary items with empty lists - onset.update({s + '_run' + str(run_id + 1) : [] for s in condition_names}) - duration.update({s + '_run' + str(run_id + 1) : [] for s in condition_names}) - weights_gain.update({'gain_run' + str(run_id + 1) : []}) - weights_loss.update({'loss_run' + str(run_id + 1) : []}) + # Create one Bunch per run + for event_file in event_files: - for run_id, _ in enumerate(runs): - f_events = event_files[run_id] + # Create empty lists + onsets = [] + durations = [] + weights_gain = [] + weights_loss = [] + onsets_missed = [] + durations_missed = [] - with open(f_events, 'rt') as file: + # Parse event file + with open(event_file, 'rt') as file: next(file) # skip the header for line in file: info = line.strip().split() - for condition in condition_names: - val = condition + '_run' + str(run_id + 1) # trial_run1 - val_gain = 'gain_run' + str(run_id + 1) # gain_run1 - val_loss = 'loss_run' + str(run_id + 1) # loss_run1 - if condition == 'trial': - onset[val].append(float(info[0])) # onsets for trial_run1 - duration[val].append(float(0)) - # weights gain for trial_run1 - weights_gain[val_gain].append(float(info[2])) - # weights loss for trial_run1 - weights_loss[val_loss].append(float(info[3])) - elif condition == 'missed': - if float(info[4]) < 0.1 or str(info[5]) == 'NoResp': - onset[val].append(float(info[0])) - duration[val].append(float(0)) - - for gain_key, gain_value in weights_gain.items(): - gain_value = gain_value - mean(gain_value) - weights_gain[gain_key] = gain_value.tolist() - - for loss_key, loss_value in weights_loss.items(): - loss_value = loss_value - mean(loss_value) - weights_loss[loss_key] = loss_value.tolist() - - # Bunching is done per run, i.e. trial_run1, trial_run2, etc. - # But names must not have '_run1' etc because we concatenate runs - subject_info = [] - for run_id in range(len(runs)): - - if len(onset['missed_run' + str(run_id + 1)]) ==0: - condition_names = ['trial'] - - conditions = [c + '_run' + str(run_id + 1) for c in condition_names] - gain = 'gain_run' + str(run_id + 1) - loss = 'loss_run' + str(run_id + 1) - - subject_info.insert( - run_id, + # Trials + onsets.append(float(info[0])) + durations.append(0.0) + weights_gain.append(float(info[2])) + weights_loss.append(float(info[3])) + + # Missed trials + if float(info[4]) < 0.1 or str(info[5]) == 'NoResp': + onsets_missed.append(float(info[0])) + durations_missed.append(0.0) + + # Mean center gain and loss weights + mean_weight = mean(weights_gain) + for element_id, element in enumerate(weights_gain): + weights_gain[element_id] = element - mean_weight + mean_weight = mean(weights_loss) + for element_id, element in enumerate(weights_loss): + weights_loss[element_id] = element - mean_weight + + # Fill Bunch + subject_information.append( Bunch( - conditions = condition_names, - onsets = [onset[k] for k in conditions], - durations = [duration[k] for k in conditions], + conditions = ['trial'] if not onsets_missed else ['trial', 'missed'], + onsets = [onsets] if not onsets_missed else [onsets, onsets_missed], + durations = [durations] if not onsets_missed else [durations, durations_missed], amplitudes = None, tmod = None, pmod = [ Bunch( name = ['gain', 'loss'], poly = [1, 1], - param = [weights_gain[gain], weights_loss[loss]] + param = [weights_gain, weights_loss] ) ], regressor_names = None, @@ -127,24 +117,7 @@ def get_subject_infos(event_files, runs): ) ) - return subject_info - - def get_contrasts(): - """ - Create the list of tuples that represents contrasts. - Each contrast is in the form : - (Name,Stat,[list of condition names],[weights on those conditions]) - """ - # List of condition names - conditions = ['trial', 'trialxgain^1', 'trialxloss^1'] - - # Create contrasts - trial = ('trial', 'T', conditions, [1, 0, 0]) - effect_gain = ('effect_of_gain', 'T', conditions, [0, 1, 0]) - effect_loss = ('effect_of_loss', 'T', conditions, [0, 0, 1]) - - # Return contrast list - return [trial, effect_gain, effect_loss] + return subject_information def get_parameters_file(filepaths, subject_id, working_dir): """ @@ -192,58 +165,18 @@ def get_parameters_file(filepaths, subject_id, working_dir): return parameters_file - def remove_gunzip_files(_, subject_id, working_dir): - """ - This method is used in a Function node to fully remove - the files generated by the gunzip node, once they aren't needed anymore. - - Parameters: - - _: Node input only used for triggering the Node - - subject_id: str, id of th subject for which to remove the unzipped file - - working_dir: str, path to the working directory - """ - from shutil import rmtree - from os.path import join - - try: - rmtree(join(working_dir, 'l1_analysis', f'_subject_id_{subject_id}', 'gunzip_func')) - except OSError as error: - print(error) - else: - print('The directory is deleted successfully') - - def remove_smoothed_files(_, subject_id, working_dir): - """ - This method is used in a Function node to fully remove - the files generated by the smoothing node, once they aren't needed anymore. - - Parameters: - - _: Node input only used for triggering the Node - - subject_id: str, id of th subject for which to remove the smoothed file - - working_dir: str, path to the working directory - """ - from shutil import rmtree - from os.path import join - - try: - rmtree(join(working_dir, 'l1_analysis', f'_subject_id_{subject_id}', 'smooth')) - except OSError as error: - print(error) - else: - print('The directory is deleted successfully') - def get_subject_level_analysis(self): """ Create the subject level analysis workflow. Returns: - - l1_analysis : nipype.WorkFlow + - subject_level_analysis : nipype.WorkFlow """ # Infosource Node - To iterate on subjects - infosource = Node(IdentityInterface( + information_source = Node(IdentityInterface( fields = ['subject_id']), - name = 'infosource') - infosource.iterables = [('subject_id', self.subject_list)] + name = 'information_source') + information_source.iterables = [('subject_id', self.subject_list)] # Templates to select files node template = { @@ -257,12 +190,12 @@ def get_subject_level_analysis(self): } # SelectFiles - to select necessary files - selectfiles = Node(SelectFiles(template, base_directory = self.directories.dataset_dir), - name = 'selectfiles') + select_files = Node(SelectFiles(template), name = 'select_files') + select_files.inputs.base_directory = self.directories.dataset_dir # DataSink - store the wanted results in the wanted repository - datasink = Node(DataSink(base_directory = self.directories.output_dir), - name='datasink') + data_sink = Node(DataSink(), name = 'data_sink') + data_sink.inputs.base_directory = self.directories.output_dir # Gunzip - gunzip files because SPM do not use .nii.gz files gunzip_func = MapNode(Gunzip(), @@ -270,39 +203,34 @@ def get_subject_level_analysis(self): iterfield = ['in_file']) # Smooth - smoothing node - smooth = Node(Smooth(fwhm = self.fwhm), - name = 'smooth') + smoothing = Node(Smooth(), name = 'smoothing') + smoothing.inputs.fwhm = self.fwhm - # Function node get_subject_infos - get subject specific condition information - subject_infos = Node(Function( - function = self.get_subject_infos, + # Function node get_subject_information - get subject specific condition information + subject_information = Node(Function( + function = self.get_subject_information, input_names = ['event_files', 'runs'], output_names = ['subject_info']), - name = 'subject_infos') - subject_infos.inputs.runs = self.run_list + name = 'subject_information') + subject_information.inputs.runs = self.run_list # SpecifyModel - generates SPM-specific Model - specify_model = Node(SpecifySPMModel( - concatenate_runs = True, input_units = 'secs', output_units = 'secs', - time_repetition = self.tr, high_pass_filter_cutoff = 128), - name='specify_model') + specify_model = Node(SpecifySPMModel(), name = 'specify_model') + specify_model.inputs.concatenate_runs = True + specify_model.inputs.input_units = 'secs' + specify_model.inputs.output_units = 'secs' + specify_model.inputs.time_repetition = TaskInformation()['RepetitionTime'] + specify_model.inputs.high_pass_filter_cutoff = 128 # Level1Design - Generates an SPM design matrix - l1_design = Node(Level1Design( - bases = {'hrf': {'derivs': [0, 0]}}, timing_units = 'secs', - interscan_interval = self.tr), name='l1_design') + model_design = Node(Level1Design(), name = 'model_design') + model_design.inputs.bases = {'hrf': {'derivs': [0, 0]}} + model_design.inputs.timing_units = 'secs' + model_design.inputs.interscan_interval = TaskInformation()['RepetitionTime'] # EstimateModel - estimate the parameters of the model - l1_estimate = Node(EstimateModel( - estimation_method={'Classical': 1}), - name='l1_estimate') - - # Function node get_contrasts - get the contrasts - contrasts = Node(Function( - function = self.get_contrasts, - input_names = [], - output_names = ['contrasts']), - name = 'contrasts') + model_estimate = Node(EstimateModel(), name = 'model_estimate') + model_estimate.inputs.estimation_method = {'Classical': 1} # Function node get_parameters_file - get parameters files parameters = Node(Function( @@ -313,55 +241,57 @@ def get_subject_level_analysis(self): parameters.inputs.working_dir = self.directories.working_dir # EstimateContrast - estimates contrasts - contrast_estimate = Node(EstimateContrast(), - name = 'contrast_estimate') + contrast_estimate = Node(EstimateContrast(), name = 'contrast_estimate') + contrast_estimate.inputs.contrasts = self.subject_level_contrasts # Function node remove_gunzip_files - remove output of the gunzip node - remove_gunzip_files = Node(Function( - function = self.remove_gunzip_files, - input_names = ['_', 'subject_id', 'working_dir'], + remove_gunzip_files = MapNode(Function( + function = remove_file, + input_names = ['_', 'file_name'], output_names = []), - name = 'remove_gunzip_files') - remove_gunzip_files.inputs.working_dir = self.directories.working_dir + name = 'remove_gunzip_files', iterfield = 'file_name') - # Function node remove_smoothed_files - remove output of the smooth node - remove_smoothed_files = Node(Function( - function = self.remove_smoothed_files, - input_names = ['_', 'subject_id', 'working_dir'], + # Function node remove_smoothed_files - remove output of the smoothing node + remove_smoothed_files = MapNode(Function( + function = remove_file, + input_names = ['_', 'file_name'], output_names = []), - name = 'remove_smoothed_files') - remove_smoothed_files.inputs.working_dir = self.directories.working_dir + name = 'remove_smoothed_files', iterfield = 'file_name') # Create l1 analysis workflow and connect its nodes - l1_analysis = Workflow(base_dir = self.directories.working_dir, name = 'l1_analysis') - l1_analysis.connect([ - (infosource, selectfiles, [('subject_id', 'subject_id')]), - (infosource, remove_gunzip_files, [('subject_id', 'subject_id')]), - (infosource, remove_smoothed_files, [('subject_id', 'subject_id')]), - (subject_infos, specify_model, [('subject_info', 'subject_info')]), - (contrasts, contrast_estimate, [('contrasts', 'contrasts')]), - (selectfiles, parameters, [('param', 'filepaths')]), - (selectfiles, subject_infos, [('event', 'event_files')]), - (infosource, parameters, [('subject_id', 'subject_id')]), - (selectfiles, gunzip_func, [('func', 'in_file')]), - (gunzip_func, smooth, [('out_file', 'in_files')]), - (smooth, remove_gunzip_files, [('smoothed_files', '_')]), - (smooth, specify_model, [('smoothed_files', 'functional_runs')]), + subject_level_analysis = Workflow( + base_dir = self.directories.working_dir, + name = 'subject_level_analysis' + ) + subject_level_analysis.connect([ + (information_source, select_files, [('subject_id', 'subject_id')]), + (information_source, remove_gunzip_files, [('subject_id', 'subject_id')]), + (information_source, remove_smoothed_files, [('subject_id', 'subject_id')]), + (subject_information, specify_model, [('subject_info', 'subject_info')]), + (select_files, parameters, [('param', 'filepaths')]), + (select_files, subject_information, [('event', 'event_files')]), + (information_source, parameters, [('subject_id', 'subject_id')]), + (select_files, gunzip_func, [('func', 'in_file')]), + (gunzip_func, smoothing, [('out_file', 'in_files')]), + (gunzip_func, remove_gunzip_files, [('out_file', 'file_name')]), + (smoothing, remove_gunzip_files, [('smoothed_files', '_')]), + (smoothing, remove_smoothed_files, [('smoothed_files', 'file_name')]), + (smoothing, specify_model, [('smoothed_files', 'functional_runs')]), (parameters, specify_model, [('parameters_file', 'realignment_parameters')]), - (specify_model, l1_design, [('session_info', 'session_info')]), - (l1_design, l1_estimate, [('spm_mat_file', 'spm_mat_file')]), - (l1_estimate, contrast_estimate, [ + (specify_model, model_design, [('session_info', 'session_info')]), + (model_design, model_estimate, [('spm_mat_file', 'spm_mat_file')]), + (model_estimate, contrast_estimate, [ ('spm_mat_file', 'spm_mat_file'), ('beta_images', 'beta_images'), ('residual_image', 'residual_image')]), - (contrast_estimate, datasink, [ - ('con_images', 'l1_analysis.@con_images'), - ('spmT_images', 'l1_analysis.@spmT_images'), - ('spm_mat_file', 'l1_analysis.@spm_mat_file')]), + (contrast_estimate, data_sink, [ + ('con_images', 'subject_level_analysis.@con_images'), + ('spmT_images', 'subject_level_analysis.@spmT_images'), + ('spm_mat_file', 'subject_level_analysis.@spm_mat_file')]), (contrast_estimate, remove_smoothed_files, [('spmT_images', '_')]) ]) - return l1_analysis + return subject_level_analysis def get_subject_level_outputs(self): """ Return the names of the files the subject level analysis is supposed to generate. """ @@ -369,18 +299,18 @@ def get_subject_level_outputs(self): # Contrat maps templates = [join( self.directories.output_dir, - 'l1_analysis', '_subject_id_{subject_id}', f'con_{contrast_id}.nii')\ + 'subject_level_analysis', '_subject_id_{subject_id}', f'con_{contrast_id}.nii')\ for contrast_id in self.contrast_list] # SPM.mat file templates += [join( self.directories.output_dir, - 'l1_analysis', '_subject_id_{subject_id}', 'SPM.mat')] + 'subject_level_analysis', '_subject_id_{subject_id}', 'SPM.mat')] # spmT maps templates += [join( self.directories.output_dir, - 'l1_analysis', '_subject_id_{subject_id}', f'spmT_{contrast_id}.nii')\ + 'subject_level_analysis', '_subject_id_{subject_id}', f'spmT_{contrast_id}.nii')\ for contrast_id in self.contrast_list] # Format with subject_ids @@ -390,40 +320,6 @@ def get_subject_level_outputs(self): return return_list - def get_subset_contrasts(file_list, method, subject_list, participants_file): - """ - Parameters : - - file_list : original file list selected by selectfiles node - - subject_list : list of subject IDs that are in the wanted group for the analysis - - participants_file: str, file containing participants characteristics - - method: str, one of 'equalRange', 'equalIndifference' or 'groupComp' - - This function return the file list containing only the files belonging - to subject in the wanted group. - """ - equal_indifference_id = [] - equal_range_id = [] - equal_indifference_files = [] - equal_range_files = [] - - with open(participants_file, 'rt') as file: - next(file) # skip the header - for line in file: - info = line.strip().split() - if info[0][-3:] in subject_list and info[1] == 'equalIndifference': - equal_indifference_id.append(info[0][-3:]) - elif info[0][-3:] in subject_list and info[1] == 'equalRange': - equal_range_id.append(info[0][-3:]) - - for file in file_list: - sub_id = file.split('/') - if sub_id[-1][-7:-4] in equal_indifference_id: - equal_indifference_files.append(file) - elif sub_id[-1][-7:-4] in equal_range_id: - equal_range_files.append(file) - - return equal_indifference_id, equal_range_id, equal_indifference_files, equal_range_files - def get_group_level_analysis(self): """ Return all workflows for the group level analysis. @@ -443,79 +339,94 @@ def get_group_level_analysis_sub_workflow(self, method): - method: one of 'equalRange', 'equalIndifference' or 'groupComp' Returns: - - l2_analysis: nipype.WorkFlow + - group_level_analysis: nipype.WorkFlow """ # Compute the number of participants used to do the analysis nb_subjects = len(self.subject_list) # Infosource - a function free node to iterate over the list of subject names - infosource_groupanalysis = Node( + information_source = Node( IdentityInterface( - fields=['contrast_id', 'subjects'], - subjects = self.subject_list), - name='infosource_groupanalysis') - infosource_groupanalysis.iterables = [('contrast_id', self.contrast_list)] + fields=['contrast_id']), + name='information_source') + information_source.iterables = [('contrast_id', self.contrast_list)] # SelectFiles templates = { # Contrasts for all participants 'contrast' : join(self.directories.output_dir, - 'l1_analysis', '_subject_id_*', 'con_{contrast_id}.nii'), - # Participants file - 'participants' : join(self.directories.dataset_dir, 'participants.tsv') + 'subject_level_analysis', '_subject_id_*', 'con_{contrast_id}.nii') } - selectfiles_groupanalysis = Node(SelectFiles( - templates, - base_directory = self.directories.results_dir, - force_list = True), - name="selectfiles_groupanalysis") + select_files = Node(SelectFiles(templates), name = 'select_files') + select_files.inputs.base_directory = self.directories.results_dir + select_files.inputs.force_list = True # Datasink - save important files - datasink_groupanalysis = Node(DataSink( - base_directory = str(self.directories.output_dir) + data_sink = Node(DataSink(), name = 'data_sink') + data_sink.inputs.base_directory = self.directories.output_dir + + # Function Node get_equal_range_subjects + # Get subjects in the equalRange group and in the subject_list + get_equal_range_subjects = Node(Function( + function = list_intersection, + input_names = ['list_1', 'list_2'], + output_names = ['out_list'] + ), + name = 'get_equal_range_subjects' + ) + get_equal_range_subjects.inputs.list_1 = get_group('equalRange') + get_equal_range_subjects.inputs.list_2 = self.subject_list + + # Function Node get_equal_indifference_subjects + # Get subjects in the equalIndifference group and in the subject_list + get_equal_indifference_subjects = Node(Function( + function = list_intersection, + input_names = ['list_1', 'list_2'], + output_names = ['out_list'] + ), + name = 'get_equal_indifference_subjects' + ) + get_equal_indifference_subjects.inputs.list_1 = get_group('equalIndifference') + get_equal_indifference_subjects.inputs.list_2 = self.subject_list + + # Function Node elements_in_string + # Get contrast files for required subjects + # Note : using a MapNode with elements_in_string requires using clean_list to remove + # None values from the out_list + get_contrasts = MapNode(Function( + function = elements_in_string, + input_names = ['input_str', 'elements'], + output_names = ['out_list'] ), - name = 'datasink_groupanalysis') - - # Node to select subset of contrasts - sub_contrasts = Node(Function( - function = self.get_subset_contrasts, - input_names = ['file_list', 'method', 'subject_list', 'participants_file'], - output_names = [ - 'equalIndifference_id', - 'equalRange_id', - 'equalIndifference_files', - 'equalRange_files']), - name = 'sub_contrasts') - sub_contrasts.inputs.method = method + name = 'get_contrasts', iterfield = 'input_str' + ) # Estimate model - estimate_model = Node(EstimateModel( - estimation_method = {'Classical':1}), - name = 'estimate_model') + estimate_model = Node(EstimateModel(), name = 'estimate_model') + estimate_model.inputs.estimation_method = {'Classical':1} # Estimate contrasts - estimate_contrast = Node(EstimateContrast( - group_contrast = True), - name = 'estimate_contrast') + estimate_contrast = Node(EstimateContrast(), name = 'estimate_contrast') + estimate_contrast.inputs.group_contrast = True ## Create thresholded maps - threshold = MapNode(Threshold( - use_fwe_correction=False, height_threshold = 0.001, extent_fdr_p_threshold = 0.05, - use_topo_fdr = False, force_activation = True), - name = 'threshold', iterfield = ['stat_image', 'contrast_index']) - - l2_analysis = Workflow( + threshold = MapNode(Threshold(), + name = 'threshold', + iterfield = ['stat_image', 'contrast_index']) + threshold.inputs.use_fwe_correction = False + threshold.inputs.height_threshold = 0.001 + threshold.inputs.extent_fdr_p_threshold = 0.05 + threshold.inputs.use_topo_fdr = False + threshold.inputs.force_activation = True + threshold.synchronize = True + + group_level_analysis = Workflow( base_dir = self.directories.working_dir, - name = f'l2_analysis_{method}_nsub_{nb_subjects}') - l2_analysis.connect([ - (infosource_groupanalysis, selectfiles_groupanalysis, [ - ('contrast_id', 'contrast_id')]), - (infosource_groupanalysis, sub_contrasts, [ - ('subjects', 'subject_list')]), - (selectfiles_groupanalysis, sub_contrasts, [ - ('contrast', 'file_list'), - ('participants', 'participants_file')]), + name = f'group_level_analysis_{method}_nsub_{nb_subjects}') + group_level_analysis.connect([ + (information_source, select_files, [('contrast_id', 'contrast_id')]), + (select_files, get_contrasts, [('contrasts', 'input_str')]), (estimate_model, estimate_contrast, [ ('spm_mat_file', 'spm_mat_file'), ('residual_image', 'residual_image'), @@ -523,52 +434,80 @@ def get_group_level_analysis_sub_workflow(self, method): (estimate_contrast, threshold, [ ('spm_mat_file', 'spm_mat_file'), ('spmT_images', 'stat_image')]), - (estimate_model, datasink_groupanalysis, [ - ('mask_image', f'l2_analysis_{method}_nsub_{nb_subjects}.@mask')]), - (estimate_contrast, datasink_groupanalysis, [ - ('spm_mat_file', f'l2_analysis_{method}_nsub_{nb_subjects}.@spm_mat'), - ('spmT_images', f'l2_analysis_{method}_nsub_{nb_subjects}.@T'), - ('con_images', f'l2_analysis_{method}_nsub_{nb_subjects}.@con')]), - (threshold, datasink_groupanalysis, [ - ('thresholded_map', f'l2_analysis_{method}_nsub_{nb_subjects}.@thresh')])]) + (estimate_model, data_sink, [ + ('mask_image', f'group_level_analysis_{method}_nsub_{nb_subjects}.@mask')]), + (estimate_contrast, data_sink, [ + ('spm_mat_file', f'group_level_analysis_{method}_nsub_{nb_subjects}.@spm_mat'), + ('spmT_images', f'group_level_analysis_{method}_nsub_{nb_subjects}.@T'), + ('con_images', f'group_level_analysis_{method}_nsub_{nb_subjects}.@con')]), + (threshold, data_sink, [ + ('thresholded_map', f'group_level_analysis_{method}_nsub_{nb_subjects}.@thresh')])]) if method in ('equalRange', 'equalIndifference'): - contrasts = [('Group', 'T', ['mean'], [1]), ('Group', 'T', ['mean'], [-1])] + estimate_contrast.inputs.contrasts = [ + ('Group', 'T', ['mean'], [1]), + ('Group', 'T', ['mean'], [-1]) + ] + + threshold.inputs.contrast_index = [1, 2] # Specify design matrix one_sample_t_test_design = Node(OneSampleTTestDesign(), name = 'one_sample_t_test_design') + group_level_analysis.connect([ + (get_contrasts, one_sample_t_test_design, [ + (('out_list', clean_list), 'in_files') + ]), + (one_sample_t_test_design, estimate_model, [('spm_mat_file', 'spm_mat_file')]) + ]) - threshold.inputs.contrast_index = [1, 2] - threshold.synchronize = True + if method == 'equalRange': + group_level_analysis.connect([ + (get_equal_range_subjects, get_contrasts, [('out_list', 'elements')]) + ]) - l2_analysis.connect([ - (sub_contrasts, one_sample_t_test_design, [(f'{method}_files', 'in_files')]), - (one_sample_t_test_design, estimate_model, [('spm_mat_file', 'spm_mat_file')])]) + elif method == 'equalIndifference': + group_level_analysis.connect([ + (get_equal_indifference_subjects, get_contrasts, [('out_list', 'elements')]) + ]) elif method == 'groupComp': - contrasts = [ - ('Eq range vs Eq indiff in loss', 'T', ['Group_{1}', 'Group_{2}'], [1, -1])] + estimate_contrast.inputs.contrasts = [ + ('Eq range vs Eq indiff in loss', 'T', ['Group_{1}', 'Group_{2}'], [1, -1]) + ] threshold.inputs.contrast_index = [1] - threshold.synchronize = True - # Node for the design matrix - two_sample_t_test_design = Node(TwoSampleTTestDesign( - unequal_variance=True), - name = 'two_sample_t_test_design') + # Function Node elements_in_string + # Get contrast files for required subjects + # Note : using a MapNode with elements_in_string requires using clean_list to remove + # None values from the out_list + get_contrasts_2 = MapNode(Function( + function = elements_in_string, + input_names = ['input_str', 'elements'], + output_names = ['out_list'] + ), + name = 'get_contrasts_2', iterfield = 'input_str' + ) - l2_analysis.connect([ - (sub_contrasts, two_sample_t_test_design, [ - ('equalRange_files', 'group1_files'), - ('equalIndifference_files', 'group2_files')]), - (two_sample_t_test_design, estimate_model, [ - ('spm_mat_file', 'spm_mat_file')]) + # Specify design matrix + two_sample_t_test_design = Node(TwoSampleTTestDesign(), + name = 'two_sample_t_test_design') + two_sample_t_test_design.inputs.unequal_variance = True + + group_level_analysis.connect([ + (get_equal_range_subjects, get_contrasts, [('out_list', 'elements')]), + (get_equal_indifference_subjects, get_contrasts_2, [('out_list', 'elements')]), + (get_contrasts, two_sample_t_test_design, [ + (('out_list', clean_list), 'group1_files') + ]), + (get_contrasts_2, two_sample_t_test_design, [ + (('out_list', clean_list), 'group2_files') + ]), + (two_sample_t_test_design, estimate_model, [('spm_mat_file', 'spm_mat_file')]) ]) - estimate_contrast.inputs.contrasts = contrasts - - return l2_analysis + return group_level_analysis def get_group_level_outputs(self): """ Return all names for the files the group level analysis is supposed to generate. """ @@ -587,7 +526,7 @@ def get_group_level_outputs(self): parameter_sets = product(*parameters.values()) template = join( self.directories.output_dir, - 'l2_analysis_{method}_nsub_{nb_subjects}', + 'group_level_analysis_{method}_nsub_{nb_subjects}', '_contrast_id_{contrast_id}', '{file}' ) @@ -608,7 +547,7 @@ def get_group_level_outputs(self): parameter_sets = product(*parameters.values()) template = join( self.directories.output_dir, - 'l2_analysis_{method}_nsub_{nb_subjects}', + 'group_level_analysis_{method}_nsub_{nb_subjects}', '_contrast_id_{contrast_id}', '{file}' ) @@ -625,23 +564,41 @@ def get_hypotheses_outputs(self): """ nb_sub = len(self.subject_list) files = [ - join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), - join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0002', 'spmT_0001.nii'), - join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), - join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0002', 'spmT_0001.nii'), - join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), - join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0002', 'spmT_0001.nii'), - join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), - join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0002', 'spmT_0001.nii'), - join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0003', '_threshold1', 'spmT_0002_thr.nii'), - join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0002.nii'), - join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0003', '_threshold1', 'spmT_0001_thr.nii'), - join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0001.nii'), - join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0003', '_threshold0', 'spmT_0001_thr.nii'), - join(f'l2_analysis_equalIndifference_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0001.nii'), - join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0003', '_threshold0', 'spmT_0002_thr.nii'), - join(f'l2_analysis_equalRange_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0002.nii'), - join(f'l2_analysis_groupComp_nsub_{nb_sub}', '_contrast_id_0003', '_threshold0', 'spmT_0001_thr.nii'), - join(f'l2_analysis_groupComp_nsub_{nb_sub}', '_contrast_id_0003', 'spmT_0001.nii') + join(f'group_level_analysis_equalIndifference_nsub_{nb_sub}', + '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), + join(f'group_level_analysis_equalIndifference_nsub_{nb_sub}', + '_contrast_id_0002', 'spmT_0001.nii'), + join(f'group_level_analysis_equalRange_nsub_{nb_sub}', + '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), + join(f'group_level_analysis_equalRange_nsub_{nb_sub}', + '_contrast_id_0002', 'spmT_0001.nii'), + join(f'group_level_analysis_equalIndifference_nsub_{nb_sub}', + '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), + join(f'group_level_analysis_equalIndifference_nsub_{nb_sub}', + '_contrast_id_0002', 'spmT_0001.nii'), + join(f'group_level_analysis_equalRange_nsub_{nb_sub}', + '_contrast_id_0002', '_threshold0', 'spmT_0001_thr.nii'), + join(f'group_level_analysis_equalRange_nsub_{nb_sub}', + '_contrast_id_0002', 'spmT_0001.nii'), + join(f'group_level_analysis_equalIndifference_nsub_{nb_sub}', + '_contrast_id_0003', '_threshold1', 'spmT_0002_thr.nii'), + join(f'group_level_analysis_equalIndifference_nsub_{nb_sub}', + '_contrast_id_0003', 'spmT_0002.nii'), + join(f'group_level_analysis_equalRange_nsub_{nb_sub}', + '_contrast_id_0003', '_threshold1', 'spmT_0001_thr.nii'), + join(f'group_level_analysis_equalRange_nsub_{nb_sub}', + '_contrast_id_0003', 'spmT_0001.nii'), + join(f'group_level_analysis_equalIndifference_nsub_{nb_sub}', + '_contrast_id_0003', '_threshold0', 'spmT_0001_thr.nii'), + join(f'group_level_analysis_equalIndifference_nsub_{nb_sub}', + '_contrast_id_0003', 'spmT_0001.nii'), + join(f'group_level_analysis_equalRange_nsub_{nb_sub}', + '_contrast_id_0003', '_threshold0', 'spmT_0002_thr.nii'), + join(f'group_level_analysis_equalRange_nsub_{nb_sub}', + '_contrast_id_0003', 'spmT_0002.nii'), + join(f'group_level_analysis_groupComp_nsub_{nb_sub}', + '_contrast_id_0003', '_threshold0', 'spmT_0001_thr.nii'), + join(f'group_level_analysis_groupComp_nsub_{nb_sub}', + '_contrast_id_0003', 'spmT_0001.nii') ] return [join(self.directories.output_dir, f) for f in files] diff --git a/narps_open/utils/configuration/testing_config.toml b/narps_open/utils/configuration/testing_config.toml index b1fb28ba..40733c5a 100644 --- a/narps_open/utils/configuration/testing_config.toml +++ b/narps_open/utils/configuration/testing_config.toml @@ -3,9 +3,9 @@ title = "Testing configuration for the NARPS open pipelines project" config_type = "testing" [directories] -dataset = "run/data/ds001734/" +dataset = "data/original/ds001734/" reproduced_results = "run/data/reproduced/" -narps_results = "run/data/results/" +narps_results = "data/results/" test_data = "tests/test_data/" test_runs = "run/" diff --git a/tests/pipelines/test_team_J7F9.py b/tests/pipelines/test_team_J7F9.py index da9023c4..2c2f0704 100644 --- a/tests/pipelines/test_team_J7F9.py +++ b/tests/pipelines/test_team_J7F9.py @@ -10,10 +10,14 @@ pytest -q test_team_J7F9.py pytest -q test_team_J7F9.py -k """ +from os.path import join from pytest import helpers, mark +from numpy import isclose from nipype import Workflow +from nipype.interfaces.base import Bunch +from narps_open.utils.configuration import Configuration from narps_open.pipelines.team_J7F9 import PipelineTeamJ7F9 class TestPipelinesTeamJ7F9: @@ -61,6 +65,50 @@ def test_outputs(): assert len(pipeline.get_group_level_outputs()) == 63 assert len(pipeline.get_hypotheses_outputs()) == 18 + @staticmethod + @mark.unit_test + def test_subject_information(): + """ Test the get_subject_information method """ + + test_event_file = join(Configuration()['directories']['test_data'], 'pipelines', 'events.tsv') + information = PipelineTeamJ7F9.get_subject_information([test_event_file, test_event_file]) + + for run_id in [0, 1]: + bunch = information [run_id] + + assert isinstance(bunch, Bunch) + assert bunch.conditions == ['trial', 'missed'] + + reference_durations = [ + [0.0, 0.0, 0.0, 0.0, 0.0], + [0.0, 0.0, 0.0, 0.0, 0.0] + ] + assert len(reference_durations) == len(bunch.durations) + for reference_array, test_array in zip(reference_durations, bunch.durations): + assert isclose(reference_array, test_array).all() + + reference_onsets = [ + [4.071, 11.834, 19.535, 27.535, 36.435], + [19.535] + ] + assert len(reference_onsets) == len(bunch.onsets) + for reference_array, test_array in zip(reference_onsets, bunch.onsets): + assert isclose(reference_array, test_array).all() + + paramateric_modulation = bunch.pmod[0] + + assert isinstance(paramateric_modulation, Bunch) + assert paramateric_modulation.name == ['gain', 'loss'] + assert paramateric_modulation.poly == [1, 1] + + reference_param = [ + [-8.4, 11.6, 15.6, -12.4, -6.4], + [-8.2, -0.2, 4.8, 0.8, 2.8] + ] + assert len(reference_param) == len(paramateric_modulation.param) + for reference_array, test_array in zip(reference_param, paramateric_modulation.param): + assert isclose(reference_array, test_array).all() + @staticmethod @mark.pipeline_test def test_execution(): diff --git a/tests/test_data/pipelines/events.tsv b/tests/test_data/pipelines/events.tsv new file mode 100644 index 00000000..4b8f04e6 --- /dev/null +++ b/tests/test_data/pipelines/events.tsv @@ -0,0 +1,6 @@ +onset duration gain loss RT participant_response +4.071 4 14 6 2.388 weakly_accept +11.834 4 34 14 2.289 strongly_accept +19.535 4 38 19 0 NoResp +27.535 4 10 15 2.08 strongly_reject +36.435 4 16 17 2.288 weakly_reject \ No newline at end of file From 81cfbaec8519fbc9c5dc6fbbf1635b0ad143780f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 3 Jan 2024 15:20:14 +0100 Subject: [PATCH 12/25] [TEST] test get_confounds method [skip ci] --- narps_open/pipelines/team_J7F9.py | 73 +++++++++---------- tests/pipelines/test_team_J7F9.py | 39 +++++++++- tests/test_data/pipelines/confounds.tsv | 4 + .../pipelines/team_J7F9/confounds.tsv | 3 + 4 files changed, 78 insertions(+), 41 deletions(-) create mode 100644 tests/test_data/pipelines/confounds.tsv create mode 100644 tests/test_data/pipelines/team_J7F9/confounds.tsv diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index b3e85996..0156f733 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -113,57 +113,51 @@ def get_subject_information(event_files): ) ], regressor_names = None, - regressors=None + regressors = None ) ) return subject_information - def get_parameters_file(filepaths, subject_id, working_dir): + def get_confounds_file(filepath, subject_id, run_id, working_dir): """ - Create new tsv files with only desired parameters per subject per run. + Create a new tsv files with only desired confounds per subject per run. Parameters : - - filepaths : paths to subject parameters file (i.e. one per run) - - subject_id : subject for whom the 1st level analysis is made + - filepath : path to the subject confounds file + - subject_id : related subject id + - run_id : related run id - working_dir: str, name of the directory for intermediate results Return : - - parameters_file : paths to new files containing only desired parameters. + - confounds_file : path to new file containing only desired confounds """ from os import makedirs from os.path import join, isdir + from pandas import DataFrame, read_csv from numpy import array, transpose - if not isinstance(filepaths, list): - filepaths = [filepaths] - - parameters_file = [] - for file_id, file in enumerate(filepaths): - data_frame = read_csv(file, sep = '\t', header=0) - - # Parameters we want to use for the model - temp_list = array([ - data_frame['X'], data_frame['Y'], data_frame['Z'], - data_frame['RotX'], data_frame['RotY'], data_frame['RotZ'], - data_frame['CSF'], data_frame['WhiteMatter'], data_frame['GlobalSignal']]) - retained_parameters = DataFrame(transpose(temp_list)) + # Open original confounds file + data_frame = read_csv(filepath, sep = '\t', header=0) - # Write parameters to a parameters file - # TODO : warning !!! filepaths must be ordered (1,2,3,4) for the following code to work - new_path = join(working_dir, 'parameters_file', - f'parameters_file_sub-{subject_id}_run-{str(file_id + 1).zfill(2)}.tsv') + # Extract confounds we want to use for the model + retained_parameters = DataFrame(transpose(array([ + data_frame['X'], data_frame['Y'], data_frame['Z'], + data_frame['RotX'], data_frame['RotY'], data_frame['RotZ'], + data_frame['CSF'], data_frame['WhiteMatter'], data_frame['GlobalSignal']]))) - makedirs(join(working_dir, 'parameters_file'), exist_ok = True) + # Write confounds to a file + confounds_file = join(working_dir, 'confounds_files', + f'confounds_file_sub-{subject_id}_run-{run_id}.tsv') - with open(new_path, 'w') as writer: - writer.write(retained_parameters.to_csv( - sep = '\t', index = False, header = False, na_rep = '0.0')) + makedirs(join(working_dir, 'confounds_files'), exist_ok = True) - parameters_file.append(new_path) + with open(confounds_file, 'w', encoding = 'utf-8') as writer: + writer.write(retained_parameters.to_csv( + sep = '\t', index = False, header = False, na_rep = '0.0')) - return parameters_file + return confounds_file def get_subject_level_analysis(self): """ @@ -232,13 +226,14 @@ def get_subject_level_analysis(self): model_estimate = Node(EstimateModel(), name = 'model_estimate') model_estimate.inputs.estimation_method = {'Classical': 1} - # Function node get_parameters_file - get parameters files - parameters = Node(Function( - function = self.get_parameters_file, - input_names = ['filepaths', 'subject_id', 'working_dir'], - output_names = ['parameters_file']), - name = 'parameters') - parameters.inputs.working_dir = self.directories.working_dir + # Function node get_confounds_file - get confounds files + confounds = MapNode(Function( + function = self.get_confounds_file, + input_names = ['filepath', 'subject_id', 'run_id', 'working_dir'], + output_names = ['confounds_file']), + name = 'confounds', iterfield = 'run_id') + confounds.inputs.working_dir = self.directories.working_dir + confounds.inputs.run_id = self.run_list # EstimateContrast - estimates contrasts contrast_estimate = Node(EstimateContrast(), name = 'contrast_estimate') @@ -268,16 +263,16 @@ def get_subject_level_analysis(self): (information_source, remove_gunzip_files, [('subject_id', 'subject_id')]), (information_source, remove_smoothed_files, [('subject_id', 'subject_id')]), (subject_information, specify_model, [('subject_info', 'subject_info')]), - (select_files, parameters, [('param', 'filepaths')]), + (select_files, confounds, [('param', 'filepaths')]), (select_files, subject_information, [('event', 'event_files')]), - (information_source, parameters, [('subject_id', 'subject_id')]), + (information_source, confounds, [('subject_id', 'subject_id')]), (select_files, gunzip_func, [('func', 'in_file')]), (gunzip_func, smoothing, [('out_file', 'in_files')]), (gunzip_func, remove_gunzip_files, [('out_file', 'file_name')]), (smoothing, remove_gunzip_files, [('smoothed_files', '_')]), (smoothing, remove_smoothed_files, [('smoothed_files', 'file_name')]), (smoothing, specify_model, [('smoothed_files', 'functional_runs')]), - (parameters, specify_model, [('parameters_file', 'realignment_parameters')]), + (confounds, specify_model, [('parameters_file', 'realignment_parameters')]), (specify_model, model_design, [('session_info', 'session_info')]), (model_design, model_estimate, [('spm_mat_file', 'spm_mat_file')]), (model_estimate, contrast_estimate, [ diff --git a/tests/pipelines/test_team_J7F9.py b/tests/pipelines/test_team_J7F9.py index 2c2f0704..a7629dbb 100644 --- a/tests/pipelines/test_team_J7F9.py +++ b/tests/pipelines/test_team_J7F9.py @@ -10,9 +10,12 @@ pytest -q test_team_J7F9.py pytest -q test_team_J7F9.py -k """ -from os.path import join +from os import mkdir +from os.path import join, exists +from shutil import rmtree +from filecmp import cmp -from pytest import helpers, mark +from pytest import helpers, mark, fixture from numpy import isclose from nipype import Workflow from nipype.interfaces.base import Bunch @@ -20,6 +23,17 @@ from narps_open.utils.configuration import Configuration from narps_open.pipelines.team_J7F9 import PipelineTeamJ7F9 +TEMPORARY_DIR = join(Configuration()['directories']['test_runs'], 'test_J7F9') + +@fixture +def remove_test_dir(): + """ A fixture to remove temporary directory created by tests """ + + rmtree(TEMPORARY_DIR, ignore_errors = True) + mkdir(TEMPORARY_DIR) + yield # test runs here + #rmtree(TEMPORARY_DIR, ignore_errors = True) + class TestPipelinesTeamJ7F9: """ A class that contains all the unit tests for the PipelineTeamJ7F9 class.""" @@ -109,6 +123,27 @@ def test_subject_information(): for reference_array, test_array in zip(reference_param, paramateric_modulation.param): assert isclose(reference_array, test_array).all() + @staticmethod + @mark.unit_test + def test_confounds_file(remove_test_dir): + """ Test the get_confounds_file method """ + + confounds_file = join( + Configuration()['directories']['test_data'], 'pipelines', 'confounds.tsv') + reference_file = join( + Configuration()['directories']['test_data'], 'pipelines', 'team_J7F9', 'confounds.tsv') + + # Get new confounds file + PipelineTeamJ7F9.get_confounds_file(confounds_file, 'sid', 'rid', TEMPORARY_DIR) + + # Check confounds file was created + created_confounds_file = join( + TEMPORARY_DIR, 'confounds_files', 'confounds_file_sub-sid_run-rid.tsv') + assert exists(created_confounds_file) + + # Check contents + assert cmp(reference_file, created_confounds_file) + @staticmethod @mark.pipeline_test def test_execution(): diff --git a/tests/test_data/pipelines/confounds.tsv b/tests/test_data/pipelines/confounds.tsv new file mode 100644 index 00000000..f49d4fea --- /dev/null +++ b/tests/test_data/pipelines/confounds.tsv @@ -0,0 +1,4 @@ +CSF WhiteMatter GlobalSignal stdDVARS non-stdDVARS vx-wisestdDVARS FramewiseDisplacement tCompCor00 tCompCor01 tCompCor02 tCompCor03 tCompCor04 tCompCor05 aCompCor00 aCompCor01 aCompCor02 aCompCor03 aCompCor04 aCompCor05 Cosine00 Cosine01 Cosine02 Cosine03 Cosine04 Cosine05 NonSteadyStateOutlier00 X Y Z RotX RotY RotZ +6551.281999999999 6476.4653 9874.576 n/a n/a n/a n/a 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 -0.0 0.0 +6484.7285 6473.4890000000005 9830.212 1.09046686 52.78273392 1.05943739 0.13527900930999998 0.0263099209 -0.0673065879 0.0934882554 -0.0079328884 0.0338007737 -0.011491083999999999 -0.042411347099999996 0.027736422900000002 0.0453303087 -0.07022609490000001 0.0963618709 -0.0200867957 0.0665186088 0.0665174038 0.0665153954 0.0665125838 0.0665089688 0.06650455059999999 0.0 -0.00996895 -0.0313444 -3.00931e-06 0.00132687 -0.000384193 -0.00016819 +6441.5337 6485.7256 9821.212 1.07520139 52.04382706 1.03821933 0.12437666391 -0.0404820317 0.034150583 0.13661184210000002 0.0745358691 -0.0054829985999999995 -0.0217322686 0.046214115199999996 0.005774624 -0.043909359800000006 -0.075619539 0.17546891539999998 -0.0345256763 0.0665153954 0.06650455059999999 0.06648647719999999 0.0664611772 0.0664286533 0.0663889091 0.0 -2.56954e-05 -0.00923735 0.0549667 0.000997278 -0.00019745 -0.000398988 diff --git a/tests/test_data/pipelines/team_J7F9/confounds.tsv b/tests/test_data/pipelines/team_J7F9/confounds.tsv new file mode 100644 index 00000000..81d3ecd9 --- /dev/null +++ b/tests/test_data/pipelines/team_J7F9/confounds.tsv @@ -0,0 +1,3 @@ +0.0 0.0 0.0 0.0 -0.0 0.0 6551.281999999999 6476.4653 9874.576 +-0.00996895 -0.0313444 -3.00931e-06 0.00132687 -0.000384193 -0.00016819 6484.7285 6473.4890000000005 9830.212 +-2.56954e-05 -0.00923735 0.0549667 0.000997278 -0.00019745 -0.000398988 6441.5337 6485.7256 9821.212 From 3f7441c41f6856acc32db400e6041d1871b66bb8 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 3 Jan 2024 15:42:18 +0100 Subject: [PATCH 13/25] [TEST] test get_confounds method [skip ci] --- tests/pipelines/test_team_J7F9.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/pipelines/test_team_J7F9.py b/tests/pipelines/test_team_J7F9.py index a7629dbb..a80eb2cd 100644 --- a/tests/pipelines/test_team_J7F9.py +++ b/tests/pipelines/test_team_J7F9.py @@ -32,7 +32,7 @@ def remove_test_dir(): rmtree(TEMPORARY_DIR, ignore_errors = True) mkdir(TEMPORARY_DIR) yield # test runs here - #rmtree(TEMPORARY_DIR, ignore_errors = True) + rmtree(TEMPORARY_DIR, ignore_errors = True) class TestPipelinesTeamJ7F9: """ A class that contains all the unit tests for the PipelineTeamJ7F9 class.""" From 9b759dff8e20ae5a780fcccb3681216f84ca9f56 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 3 Jan 2024 15:46:32 +0100 Subject: [PATCH 14/25] [PEP8] [skip ci] --- narps_open/pipelines/team_J7F9.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 0156f733..4eff5213 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -133,7 +133,7 @@ def get_confounds_file(filepath, subject_id, run_id, working_dir): - confounds_file : path to new file containing only desired confounds """ from os import makedirs - from os.path import join, isdir + from os.path import join from pandas import DataFrame, read_csv from numpy import array, transpose From 4ec1cf6a511a1834b5da4cf4ec0158f4caee4614 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 3 Jan 2024 17:33:07 +0100 Subject: [PATCH 15/25] [BUG] runs inputs of get_subject_information no longer needed [skip ci] --- narps_open/pipelines/team_J7F9.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 4eff5213..0d5c508f 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -203,10 +203,9 @@ def get_subject_level_analysis(self): # Function node get_subject_information - get subject specific condition information subject_information = Node(Function( function = self.get_subject_information, - input_names = ['event_files', 'runs'], + input_names = ['event_files'], output_names = ['subject_info']), name = 'subject_information') - subject_information.inputs.runs = self.run_list # SpecifyModel - generates SPM-specific Model specify_model = Node(SpecifySPMModel(), name = 'specify_model') From 2dd6525ca6b53723af14e25de4b0fa3764f044a3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 3 Jan 2024 17:36:35 +0100 Subject: [PATCH 16/25] [BUG] get_confounds input [skip ci] --- narps_open/pipelines/team_J7F9.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 0d5c508f..65bac979 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -174,8 +174,7 @@ def get_subject_level_analysis(self): # Templates to select files node template = { - # Parameter file - 'param' : join('derivatives', 'fmriprep', 'sub-{subject_id}', 'func', + 'confounds' : join('derivatives', 'fmriprep', 'sub-{subject_id}', 'func', 'sub-{subject_id}_task-MGT_run-*_bold_confounds.tsv'), 'func' : join('derivatives', 'fmriprep', 'sub-{subject_id}', 'func', 'sub-{subject_id}_task-MGT_run-*_bold_space-MNI152NLin2009cAsym_preproc.nii.gz'), @@ -262,7 +261,7 @@ def get_subject_level_analysis(self): (information_source, remove_gunzip_files, [('subject_id', 'subject_id')]), (information_source, remove_smoothed_files, [('subject_id', 'subject_id')]), (subject_information, specify_model, [('subject_info', 'subject_info')]), - (select_files, confounds, [('param', 'filepaths')]), + (select_files, confounds, [('confounds', 'filepath')]), (select_files, subject_information, [('event', 'event_files')]), (information_source, confounds, [('subject_id', 'subject_id')]), (select_files, gunzip_func, [('func', 'in_file')]), From ccb6cd7d2f61ea8e1725566e5aa9927d817b85c1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Wed, 3 Jan 2024 17:43:32 +0100 Subject: [PATCH 17/25] [BUG] get_confounds input [skip ci] --- narps_open/pipelines/team_J7F9.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 65bac979..fc1a3cf6 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -229,7 +229,7 @@ def get_subject_level_analysis(self): function = self.get_confounds_file, input_names = ['filepath', 'subject_id', 'run_id', 'working_dir'], output_names = ['confounds_file']), - name = 'confounds', iterfield = 'run_id') + name = 'confounds', iterfield = ['filepath', 'run_id']) confounds.inputs.working_dir = self.directories.working_dir confounds.inputs.run_id = self.run_list @@ -270,7 +270,7 @@ def get_subject_level_analysis(self): (smoothing, remove_gunzip_files, [('smoothed_files', '_')]), (smoothing, remove_smoothed_files, [('smoothed_files', 'file_name')]), (smoothing, specify_model, [('smoothed_files', 'functional_runs')]), - (confounds, specify_model, [('parameters_file', 'realignment_parameters')]), + (confounds, specify_model, [('confounds_file', 'realignment_parameters')]), (specify_model, model_design, [('session_info', 'session_info')]), (model_design, model_estimate, [('spm_mat_file', 'spm_mat_file')]), (model_estimate, contrast_estimate, [ From 015e9c724f21ad9ced17997459c14f2d53191b13 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Thu, 4 Jan 2024 11:20:30 +0100 Subject: [PATCH 18/25] Cleaning code [skip ci] --- narps_open/pipelines/team_J7F9.py | 8 +++----- 1 file changed, 3 insertions(+), 5 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index fc1a3cf6..7594f46c 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -85,7 +85,7 @@ def get_subject_information(event_files): weights_loss.append(float(info[3])) # Missed trials - if float(info[4]) < 0.1 or str(info[5]) == 'NoResp': + if float(info[4]) < 0.1 or 'NoResp' in info[5]: onsets_missed.append(float(info[0])) durations_missed.append(0.0) @@ -258,8 +258,6 @@ def get_subject_level_analysis(self): ) subject_level_analysis.connect([ (information_source, select_files, [('subject_id', 'subject_id')]), - (information_source, remove_gunzip_files, [('subject_id', 'subject_id')]), - (information_source, remove_smoothed_files, [('subject_id', 'subject_id')]), (subject_information, specify_model, [('subject_info', 'subject_info')]), (select_files, confounds, [('confounds', 'filepath')]), (select_files, subject_information, [('event', 'event_files')]), @@ -267,7 +265,6 @@ def get_subject_level_analysis(self): (select_files, gunzip_func, [('func', 'in_file')]), (gunzip_func, smoothing, [('out_file', 'in_files')]), (gunzip_func, remove_gunzip_files, [('out_file', 'file_name')]), - (smoothing, remove_gunzip_files, [('smoothed_files', '_')]), (smoothing, remove_smoothed_files, [('smoothed_files', 'file_name')]), (smoothing, specify_model, [('smoothed_files', 'functional_runs')]), (confounds, specify_model, [('confounds_file', 'realignment_parameters')]), @@ -281,7 +278,8 @@ def get_subject_level_analysis(self): ('con_images', 'subject_level_analysis.@con_images'), ('spmT_images', 'subject_level_analysis.@spmT_images'), ('spm_mat_file', 'subject_level_analysis.@spm_mat_file')]), - (contrast_estimate, remove_smoothed_files, [('spmT_images', '_')]) + (data_sink, remove_smoothed_files, [('out_file', '_')]), + (data_sink, remove_gunzip_files, [('out_file', '_')]) ]) return subject_level_analysis From e2bdee471a90cd4c660dac141d443e5566746834 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Thu, 4 Jan 2024 14:58:56 +0100 Subject: [PATCH 19/25] [BUG] real subject ids for search patterns [skip ci] --- narps_open/pipelines/team_J7F9.py | 27 ++++++++++++++++++++------- 1 file changed, 20 insertions(+), 7 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 7594f46c..0d208d4b 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -381,6 +381,11 @@ def get_group_level_analysis_sub_workflow(self, method): get_equal_indifference_subjects.inputs.list_1 = get_group('equalIndifference') get_equal_indifference_subjects.inputs.list_2 = self.subject_list + # Create a function to complete the subject ids out from the get_equal_*_subjects nodes + # If not complete, subject id '001' in search patterns + # would match all contrast files with 'con_0001.nii'. + complete_subject_ids = lambda l : [f'_subject_id_{a}' for a in l] + # Function Node elements_in_string # Get contrast files for required subjects # Note : using a MapNode with elements_in_string requires using clean_list to remove @@ -454,13 +459,17 @@ def get_group_level_analysis_sub_workflow(self, method): if method == 'equalRange': group_level_analysis.connect([ - (get_equal_range_subjects, get_contrasts, [('out_list', 'elements')]) + (get_equal_range_subjects, get_contrasts, [ + (('out_list', complete_subject_ids), 'elements') ]) + ]) elif method == 'equalIndifference': group_level_analysis.connect([ - (get_equal_indifference_subjects, get_contrasts, [('out_list', 'elements')]) + (get_equal_indifference_subjects, get_contrasts, [ + (('out_list', complete_subject_ids), 'elements') ]) + ]) elif method == 'groupComp': estimate_contrast.inputs.contrasts = [ @@ -487,16 +496,20 @@ def get_group_level_analysis_sub_workflow(self, method): two_sample_t_test_design.inputs.unequal_variance = True group_level_analysis.connect([ - (get_equal_range_subjects, get_contrasts, [('out_list', 'elements')]), - (get_equal_indifference_subjects, get_contrasts_2, [('out_list', 'elements')]), + (get_equal_range_subjects, get_contrasts, [ + (('out_list', complete_subject_ids), 'elements') + ]), + (get_equal_indifference_subjects, get_contrasts_2, [ + (('out_list', complete_subject_ids), 'elements') + ]), (get_contrasts, two_sample_t_test_design, [ (('out_list', clean_list), 'group1_files') - ]), + ]), (get_contrasts_2, two_sample_t_test_design, [ (('out_list', clean_list), 'group2_files') - ]), + ]), (two_sample_t_test_design, estimate_model, [('spm_mat_file', 'spm_mat_file')]) - ]) + ]) return group_level_analysis From cfa72791870f22901486f23c8bdc72089dadcb6c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Thu, 4 Jan 2024 16:22:49 +0100 Subject: [PATCH 20/25] Consistent Bunches along runs as output of subject_information [skip ci] --- narps_open/pipelines/team_J7F9.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 0d208d4b..289966fa 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -100,9 +100,9 @@ def get_subject_information(event_files): # Fill Bunch subject_information.append( Bunch( - conditions = ['trial'] if not onsets_missed else ['trial', 'missed'], - onsets = [onsets] if not onsets_missed else [onsets, onsets_missed], - durations = [durations] if not onsets_missed else [durations, durations_missed], + conditions = ['trial', 'missed'], + onsets = [onsets, onsets_missed], + durations = [durations, durations_missed], amplitudes = None, tmod = None, pmod = [ From 4f5ff0a699a81998537973987748c354eace97b1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Thu, 4 Jan 2024 17:43:44 +0100 Subject: [PATCH 21/25] [BUG] with subject_information Bunch [skip ci] --- narps_open/pipelines/team_J7F9.py | 74 ++++++++++++++++++++----------- 1 file changed, 47 insertions(+), 27 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 289966fa..760f8d89 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -60,18 +60,28 @@ def get_subject_information(event_files): subject_information = [] - # Create one Bunch per run - for event_file in event_files: - - # Create empty lists - onsets = [] - durations = [] - weights_gain = [] - weights_loss = [] - onsets_missed = [] - durations_missed = [] + # Create empty dicts + onsets = {} + durations = {} + weights_gain = {} + weights_loss = {} + onsets_missed = {} + durations_missed = {} + + # Run list + run_list = [str(r).zfill(2) for r in range(1, len(event_files) + 1)] + + # Parse event file + for run_id, event_file in zip(run_list, event_files): + + # Init empty lists inside directiries + onsets[run_id] = [] + durations[run_id] = [] + weights_gain[run_id] = [] + weights_loss[run_id] = [] + onsets_missed[run_id] = [] + durations_missed[run_id] = [] - # Parse event file with open(event_file, 'rt') as file: next(file) # skip the header @@ -79,37 +89,47 @@ def get_subject_information(event_files): info = line.strip().split() # Trials - onsets.append(float(info[0])) - durations.append(0.0) - weights_gain.append(float(info[2])) - weights_loss.append(float(info[3])) + onsets[run_id].append(float(info[0])) + durations[run_id].append(0.0) + weights_gain[run_id].append(float(info[2])) + weights_loss[run_id].append(float(info[3])) # Missed trials if float(info[4]) < 0.1 or 'NoResp' in info[5]: - onsets_missed.append(float(info[0])) - durations_missed.append(0.0) + onsets_missed[run_id].append(float(info[0])) + durations_missed[run_id].append(0.0) + + # Compute mean weight values across all runs + mean_gain_weight = mean(list(weights_gain.values())) + mean_loss_weight = mean(list(weights_loss.values())) + + # Check if there are any missed trials across all runs + missed_trials = any(t for t in onsets_missed.values()) + + # Create one Bunch per run + for run_id in run_list: # Mean center gain and loss weights - mean_weight = mean(weights_gain) - for element_id, element in enumerate(weights_gain): - weights_gain[element_id] = element - mean_weight - mean_weight = mean(weights_loss) - for element_id, element in enumerate(weights_loss): - weights_loss[element_id] = element - mean_weight + for element_id, element in enumerate(weights_gain[run_id]): + weights_gain[run_id][element_id] = element - mean_gain_weight + for element_id, element in enumerate(weights_loss[run_id]): + weights_loss[run_id][element_id] = element - mean_loss_weight # Fill Bunch subject_information.append( Bunch( - conditions = ['trial', 'missed'], - onsets = [onsets, onsets_missed], - durations = [durations, durations_missed], + conditions = ['trial', 'missed'] if missed_trials else ['trial'], + onsets = [onsets[run_id], onsets_missed[run_id]] if missed_trials\ + else [onsets[run_id]], + durations = [durations[run_id], durations_missed[run_id]] if missed_trials\ + else [durations[run_id]], amplitudes = None, tmod = None, pmod = [ Bunch( name = ['gain', 'loss'], poly = [1, 1], - param = [weights_gain, weights_loss] + param = [weights_gain[run_id], weights_loss[run_id]] ) ], regressor_names = None, From a9ab13a4163093d1e68cd3661907c1bb0504371e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Thu, 4 Jan 2024 18:09:52 +0100 Subject: [PATCH 22/25] Bug with group level [skip ci] --- narps_open/pipelines/team_J7F9.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 760f8d89..f43fa911 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -365,13 +365,13 @@ def get_group_level_analysis_sub_workflow(self, method): # SelectFiles templates = { # Contrasts for all participants - 'contrast' : join(self.directories.output_dir, + 'contrasts' : join(self.directories.output_dir, 'subject_level_analysis', '_subject_id_*', 'con_{contrast_id}.nii') } select_files = Node(SelectFiles(templates), name = 'select_files') select_files.inputs.base_directory = self.directories.results_dir - select_files.inputs.force_list = True + select_files.inputs.force_lists = True # Datasink - save important files data_sink = Node(DataSink(), name = 'data_sink') From fc37fda87af0db7f10a7d952721d2e3bc9ef031b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Fri, 5 Jan 2024 09:41:11 +0100 Subject: [PATCH 23/25] Remove large files earlier in the process [skip ci] --- narps_open/pipelines/team_J7F9.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index f43fa911..6d8274d1 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -285,8 +285,10 @@ def get_subject_level_analysis(self): (select_files, gunzip_func, [('func', 'in_file')]), (gunzip_func, smoothing, [('out_file', 'in_files')]), (gunzip_func, remove_gunzip_files, [('out_file', 'file_name')]), + (smoothing, remove_gunzip_files, [('smoothed_files', '_')]), (smoothing, remove_smoothed_files, [('smoothed_files', 'file_name')]), (smoothing, specify_model, [('smoothed_files', 'functional_runs')]), + (specify_model, remove_smoothed_files, [('session_info', '_')]), (confounds, specify_model, [('confounds_file', 'realignment_parameters')]), (specify_model, model_design, [('session_info', 'session_info')]), (model_design, model_estimate, [('spm_mat_file', 'spm_mat_file')]), @@ -297,10 +299,9 @@ def get_subject_level_analysis(self): (contrast_estimate, data_sink, [ ('con_images', 'subject_level_analysis.@con_images'), ('spmT_images', 'subject_level_analysis.@spmT_images'), - ('spm_mat_file', 'subject_level_analysis.@spm_mat_file')]), - (data_sink, remove_smoothed_files, [('out_file', '_')]), - (data_sink, remove_gunzip_files, [('out_file', '_')]) + ('spm_mat_file', 'subject_level_analysis.@spm_mat_file') ]) + ]) return subject_level_analysis From 7a03d5269e2b67bbe9c68e9cbcf1c00a4c112f19 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Fri, 5 Jan 2024 12:15:08 +0100 Subject: [PATCH 24/25] [TEST] exhaustive testing of get_subject_information for J7F9 [skip ci] --- narps_open/pipelines/team_J7F9.py | 13 +- tests/pipelines/test_team_J7F9.py | 134 +++++++++++++----- .../pipelines/team_J7F9/events_resp.tsv | 5 + 3 files changed, 111 insertions(+), 41 deletions(-) create mode 100644 tests/test_data/pipelines/team_J7F9/events_resp.tsv diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index 6d8274d1..e478ba16 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -55,7 +55,7 @@ def get_subject_information(event_files): Returns : - subject_information : list of Bunch for subject level analysis. """ - from numpy import mean + from numpy import mean, ravel from nipype.interfaces.base import Bunch subject_information = [] @@ -100,8 +100,15 @@ def get_subject_information(event_files): durations_missed[run_id].append(0.0) # Compute mean weight values across all runs - mean_gain_weight = mean(list(weights_gain.values())) - mean_loss_weight = mean(list(weights_loss.values())) + all_weights_gain = [] + for value in weights_gain.values(): + all_weights_gain += value + mean_gain_weight = mean(all_weights_gain) + + all_weights_loss = [] + for value in weights_loss.values(): + all_weights_loss += value + mean_loss_weight = mean(all_weights_loss) # Check if there are any missed trials across all runs missed_trials = any(t for t in onsets_missed.values()) diff --git a/tests/pipelines/test_team_J7F9.py b/tests/pipelines/test_team_J7F9.py index a80eb2cd..706c0b56 100644 --- a/tests/pipelines/test_team_J7F9.py +++ b/tests/pipelines/test_team_J7F9.py @@ -34,6 +34,14 @@ def remove_test_dir(): yield # test runs here rmtree(TEMPORARY_DIR, ignore_errors = True) +def compare_float_2d_arrays(array_1, array_2): + """ Assert array_1 and array_2 are close enough """ + + assert len(array_1) == len(array_2) + for reference_array, test_array in zip(array_1, array_2): + assert len(reference_array) == len(test_array) + assert isclose(reference_array, test_array).all() + class TestPipelinesTeamJ7F9: """ A class that contains all the unit tests for the PipelineTeamJ7F9 class.""" @@ -84,44 +92,94 @@ def test_outputs(): def test_subject_information(): """ Test the get_subject_information method """ - test_event_file = join(Configuration()['directories']['test_data'], 'pipelines', 'events.tsv') - information = PipelineTeamJ7F9.get_subject_information([test_event_file, test_event_file]) - - for run_id in [0, 1]: - bunch = information [run_id] - - assert isinstance(bunch, Bunch) - assert bunch.conditions == ['trial', 'missed'] - - reference_durations = [ - [0.0, 0.0, 0.0, 0.0, 0.0], - [0.0, 0.0, 0.0, 0.0, 0.0] - ] - assert len(reference_durations) == len(bunch.durations) - for reference_array, test_array in zip(reference_durations, bunch.durations): - assert isclose(reference_array, test_array).all() - - reference_onsets = [ - [4.071, 11.834, 19.535, 27.535, 36.435], - [19.535] - ] - assert len(reference_onsets) == len(bunch.onsets) - for reference_array, test_array in zip(reference_onsets, bunch.onsets): - assert isclose(reference_array, test_array).all() - - paramateric_modulation = bunch.pmod[0] - - assert isinstance(paramateric_modulation, Bunch) - assert paramateric_modulation.name == ['gain', 'loss'] - assert paramateric_modulation.poly == [1, 1] - - reference_param = [ - [-8.4, 11.6, 15.6, -12.4, -6.4], - [-8.2, -0.2, 4.8, 0.8, 2.8] - ] - assert len(reference_param) == len(paramateric_modulation.param) - for reference_array, test_array in zip(reference_param, paramateric_modulation.param): - assert isclose(reference_array, test_array).all() + # Get test files + test_file = join(Configuration()['directories']['test_data'], 'pipelines', 'events.tsv') + test_file_2 = join(Configuration()['directories']['test_data'], + 'pipelines', 'team_J7F9', 'events_resp.tsv') + + # Prepare several scenarii + info_missed = PipelineTeamJ7F9.get_subject_information([test_file, test_file]) + info_ok = PipelineTeamJ7F9.get_subject_information([test_file_2, test_file_2]) + info_half = PipelineTeamJ7F9.get_subject_information([test_file_2, test_file]) + + # Compare bunches to expected + bunch = info_missed[0] + assert isinstance(bunch, Bunch) + assert bunch.conditions == ['trial', 'missed'] + compare_float_2d_arrays(bunch.onsets, [[4.071, 11.834, 19.535, 27.535, 36.435], [19.535]]) + compare_float_2d_arrays(bunch.durations, [[0.0, 0.0, 0.0, 0.0, 0.0], [0.0]]) + assert bunch.amplitudes == None + assert bunch.tmod == None + assert bunch.pmod[0].name == ['gain', 'loss'] + assert bunch.pmod[0].poly == [1, 1] + compare_float_2d_arrays(bunch.pmod[0].param, [[-8.4, 11.6, 15.6, -12.4, -6.4], [-8.2, -0.2, 4.8, 0.8, 2.8]]) + assert bunch.regressor_names == None + assert bunch.regressors == None + + bunch = info_missed[1] + assert isinstance(bunch, Bunch) + assert bunch.conditions == ['trial', 'missed'] + compare_float_2d_arrays(bunch.onsets, [[4.071, 11.834, 19.535, 27.535, 36.435], [19.535]]) + compare_float_2d_arrays(bunch.durations, [[0.0, 0.0, 0.0, 0.0, 0.0], [0.0]]) + assert bunch.amplitudes == None + assert bunch.tmod == None + assert bunch.pmod[0].name == ['gain', 'loss'] + assert bunch.pmod[0].poly == [1, 1] + compare_float_2d_arrays(bunch.pmod[0].param, [[-8.4, 11.6, 15.6, -12.4, -6.4], [-8.2, -0.2, 4.8, 0.8, 2.8]]) + assert bunch.regressor_names == None + assert bunch.regressors == None + + bunch = info_ok[0] + assert isinstance(bunch, Bunch) + assert bunch.conditions == ['trial'] + compare_float_2d_arrays(bunch.onsets, [[4.071, 11.834, 27.535, 36.435]]) + compare_float_2d_arrays(bunch.durations, [[0.0, 0.0, 0.0, 0.0]]) + assert bunch.amplitudes == None + assert bunch.tmod == None + assert bunch.pmod[0].name == ['gain', 'loss'] + assert bunch.pmod[0].poly == [1, 1] + compare_float_2d_arrays(bunch.pmod[0].param, [[-4.5, 15.5, -8.5, -2.5], [-7.0, 1.0, 2.0, 4.0]]) + assert bunch.regressor_names == None + assert bunch.regressors == None + + bunch = info_ok[1] + assert isinstance(bunch, Bunch) + assert bunch.conditions == ['trial'] + compare_float_2d_arrays(bunch.onsets, [[4.071, 11.834, 27.535, 36.435]]) + compare_float_2d_arrays(bunch.durations, [[0.0, 0.0, 0.0, 0.0]]) + assert bunch.amplitudes == None + assert bunch.tmod == None + assert bunch.pmod[0].name == ['gain', 'loss'] + assert bunch.pmod[0].poly == [1, 1] + compare_float_2d_arrays(bunch.pmod[0].param, [[-4.5, 15.5, -8.5, -2.5], [-7.0, 1.0, 2.0, 4.0]]) + assert bunch.regressor_names == None + assert bunch.regressors == None + + bunch = info_half[0] + assert isinstance(bunch, Bunch) + assert bunch.conditions == ['trial', 'missed'] + compare_float_2d_arrays(bunch.onsets, [[4.071, 11.834, 27.535, 36.435], []]) + compare_float_2d_arrays(bunch.durations, [[0.0, 0.0, 0.0, 0.0], []]) + assert bunch.amplitudes == None + assert bunch.tmod == None + assert bunch.pmod[0].name == ['gain', 'loss'] + assert bunch.pmod[0].poly == [1, 1] + compare_float_2d_arrays(bunch.pmod[0].param, [[-6.666666666666668, 13.333333333333332, -10.666666666666668, -4.666666666666668], [-7.666666666666666, 0.3333333333333339, 1.333333333333334, 3.333333333333334]]) + assert bunch.regressor_names == None + assert bunch.regressors == None + + bunch = info_half[1] + assert isinstance(bunch, Bunch) + assert bunch.conditions == ['trial', 'missed'] + compare_float_2d_arrays(bunch.onsets, [[4.071, 11.834, 19.535, 27.535, 36.435], [19.535]]) + compare_float_2d_arrays(bunch.durations, [[0.0, 0.0, 0.0, 0.0, 0.0], [0.0]]) + assert bunch.amplitudes == None + assert bunch.tmod == None + assert bunch.pmod[0].name == ['gain', 'loss'] + assert bunch.pmod[0].poly == [1, 1] + compare_float_2d_arrays(bunch.pmod[0].param, [[-6.666666666666668, 13.333333333333332, 17.333333333333332, -10.666666666666668, -4.666666666666668], [-7.666666666666666, 0.3333333333333339, 5.333333333333334, 1.333333333333334, 3.333333333333334]]) + assert bunch.regressor_names == None + assert bunch.regressors == None @staticmethod @mark.unit_test diff --git a/tests/test_data/pipelines/team_J7F9/events_resp.tsv b/tests/test_data/pipelines/team_J7F9/events_resp.tsv new file mode 100644 index 00000000..dd5ea1a5 --- /dev/null +++ b/tests/test_data/pipelines/team_J7F9/events_resp.tsv @@ -0,0 +1,5 @@ +onset duration gain loss RT participant_response +4.071 4 14 6 2.388 weakly_accept +11.834 4 34 14 2.289 strongly_accept +27.535 4 10 15 2.08 strongly_reject +36.435 4 16 17 2.288 weakly_reject \ No newline at end of file From b2dd3815b4babace941f455834958ee90aec9d84 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Boris=20Cl=C3=A9net?= Date: Fri, 5 Jan 2024 15:28:51 +0100 Subject: [PATCH 25/25] [BUG] get_contrasts_2 iterfield [skip ci] --- narps_open/pipelines/team_J7F9.py | 1 + 1 file changed, 1 insertion(+) diff --git a/narps_open/pipelines/team_J7F9.py b/narps_open/pipelines/team_J7F9.py index e478ba16..d04e66e2 100644 --- a/narps_open/pipelines/team_J7F9.py +++ b/narps_open/pipelines/team_J7F9.py @@ -524,6 +524,7 @@ def get_group_level_analysis_sub_workflow(self, method): two_sample_t_test_design.inputs.unequal_variance = True group_level_analysis.connect([ + (select_files, get_contrasts_2, [('contrasts', 'input_str')]), (get_equal_range_subjects, get_contrasts, [ (('out_list', complete_subject_ids), 'elements') ]),