diff --git a/env.yml b/env.yml index a1b1185..1d65dce 100644 --- a/env.yml +++ b/env.yml @@ -29,6 +29,9 @@ dependencies: - fsl-melodic=2111.3 - fsl-miscmaths=2203.2 - fsl-susan=2111.0 + # Workflow dependencies: ANTs + - ants=2.5 + # Other dependencies - pip - pip: - templateflow==24.0 diff --git a/scripts/fetch_templates.py b/scripts/fetch_templates.py index 04881d6..d622d8a 100755 --- a/scripts/fetch_templates.py +++ b/scripts/fetch_templates.py @@ -25,11 +25,11 @@ def fetch_MNI2009(): """ template = 'MNI152NLin2009cAsym' - tf.get(template, resolution=(1, 2), desc=None, suffix='T1w') - tf.get(template, resolution=(1, 2), desc='brain', suffix='mask') + # tf.get(template, resolution=(1, 2), desc=None, suffix='T1w') + # tf.get(template, resolution=(1, 2), desc='brain', suffix='mask') tf.get(template, resolution=1, atlas=None, desc='carpet', suffix='dseg') - tf.get(template, resolution=2, desc='fMRIPrep', suffix='boldref') - tf.get(template, resolution=1, label='brain', suffix='probseg') + # tf.get(template, resolution=2, desc='fMRIPrep', suffix='boldref') + # tf.get(template, resolution=1, label='brain', suffix='probseg') def fetch_MNI6(): @@ -48,6 +48,10 @@ def fetch_MNI6(): tf.get(template, resolution=(1, 2), desc='brain', suffix='mask') # CIFTI tf.get(template, resolution=2, atlas='HCP', suffix='dseg') + # Transform from MNI152NLin2009cAsym to MNI152NLin6Asym + tf.get( + template, mode='image', suffix='xfm', extension='.h5', **{'from': 'MNI152NLin2009cAsym'} + ) def fetch_OASIS(): @@ -108,7 +112,7 @@ def fetch_fsLR(): def fetch_all(): - # fetch_MNI2009() + fetch_MNI2009() fetch_MNI6() # fetch_OASIS() # fetch_fsaverage() diff --git a/src/fmripost_aroma/cli/workflow.py b/src/fmripost_aroma/cli/workflow.py index 61a5a71..386d020 100644 --- a/src/fmripost_aroma/cli/workflow.py +++ b/src/fmripost_aroma/cli/workflow.py @@ -90,11 +90,12 @@ def build_workflow(config_file, retval): # Called with reports only if config.execution.reports_only: build_log.log(25, 'Running --reports-only on participants %s', ', '.join(subject_list)) - retval['return_code'] = generate_reports( + failed_reports = generate_reports( subject_list=config.execution.participant_label, output_dir=config.execution.output_dir, run_uuid=config.execution.run_uuid, ) + retval['return_code'] = len(failed_reports) return retval # Build main workflow diff --git a/src/fmripost_aroma/data/reports-spec-func.yml b/src/fmripost_aroma/data/reports-spec-func.yml index 97ff86c..29e8513 100644 --- a/src/fmripost_aroma/data/reports-spec-func.yml +++ b/src/fmripost_aroma/data/reports-spec-func.yml @@ -21,11 +21,14 @@ sections: anatomical white matter mask, which appears as a red contour. static: false subtitle: Alignment of functional and anatomical MRI data (coregistration) - - bids: - datatype: figures - desc: '[aggr,nonaggr,orthaggr]carpetplot' - extension: [.html] - suffix: bold + - bids: {datatype: figures, desc: preprocCarpetplot, suffix: bold} + title: Preprocessed BOLD + - bids: {datatype: figures, desc: nonaggrCarpetplot, suffix: bold} + title: Non-Aggressively Denoised BOLD + - bids: {datatype: figures, desc: aggrCarpetplot, suffix: bold} + title: Aggressively Denoised BOLD + - bids: {datatype: figures, desc: orthaggrCarpetplot, suffix: bold} + title: Othogonalized Aggressively Denoised BOLD - name: About nested: true reportlets: diff --git a/src/fmripost_aroma/data/reports-spec.yml b/src/fmripost_aroma/data/reports-spec.yml index 7e0a637..eb346a3 100644 --- a/src/fmripost_aroma/data/reports-spec.yml +++ b/src/fmripost_aroma/data/reports-spec.yml @@ -21,11 +21,14 @@ sections: anatomical white matter mask, which appears as a red contour. static: false subtitle: Alignment of functional and anatomical MRI data (coregistration) - - bids: - datatype: figures - desc: '[aggr,nonaggr,orthaggr]carpetplot' - extension: [.html] - suffix: bold + - bids: {datatype: figures, desc: preprocCarpetplot, suffix: bold} + title: Preprocessed BOLD + - bids: {datatype: figures, desc: nonaggrCarpetplot, suffix: bold} + title: Non-Aggressively Denoised BOLD + - bids: {datatype: figures, desc: aggrCarpetplot, suffix: bold} + title: Aggressively Denoised BOLD + - bids: {datatype: figures, desc: orthaggrCarpetplot, suffix: bold} + title: Othogonalized Aggressively Denoised BOLD - name: About nested: true reportlets: diff --git a/src/fmripost_aroma/interfaces/reportlets.py b/src/fmripost_aroma/interfaces/reportlets.py index 74615b0..7fa726e 100644 --- a/src/fmripost_aroma/interfaces/reportlets.py +++ b/src/fmripost_aroma/interfaces/reportlets.py @@ -67,7 +67,7 @@ ABOUT_TEMPLATE = """\t """ diff --git a/src/fmripost_aroma/tests/test_base.py b/src/fmripost_aroma/tests/test_base.py index 0534018..2044cea 100644 --- a/src/fmripost_aroma/tests/test_base.py +++ b/src/fmripost_aroma/tests/test_base.py @@ -33,5 +33,8 @@ def test_init_denoise_wf(tmp_path_factory): config.execution.output_dir = tempdir / 'out' config.execution.work_dir = tempdir / 'work' - wf = init_denoise_wf(bold_file='sub-01_task-rest_bold.nii.gz') + wf = init_denoise_wf( + bold_file='sub-01_task-rest_bold.nii.gz', + metadata={'RepetitionTime': 2.0}, + ) assert wf.name == 'denoise_task_rest_wf' diff --git a/src/fmripost_aroma/workflows/aroma.py b/src/fmripost_aroma/workflows/aroma.py index 3a903c8..ea04d80 100644 --- a/src/fmripost_aroma/workflows/aroma.py +++ b/src/fmripost_aroma/workflows/aroma.py @@ -388,7 +388,7 @@ def init_ica_aroma_wf( return workflow -def init_denoise_wf(bold_file): +def init_denoise_wf(bold_file, metadata): """Build a workflow that denoises a BOLD series using AROMA confounds. This workflow performs the denoising in the requested output space(s). @@ -397,6 +397,7 @@ def init_denoise_wf(bold_file): from niworkflows.engine.workflows import LiterateWorkflow as Workflow from fmripost_aroma.interfaces.confounds import ICADenoise + from fmripost_aroma.workflows.confounds import init_carpetplot_wf workflow = Workflow(name=_get_wf_name(bold_file, 'denoise')) @@ -405,6 +406,7 @@ def init_denoise_wf(bold_file): fields=[ 'bold_file', 'bold_mask', + 'confounds_file', 'mixing', 'classifications', 'skip_vols', @@ -427,6 +429,22 @@ def init_denoise_wf(bold_file): ]), ]) # fmt:skip + preproc_carpetplot_wf = init_carpetplot_wf( + mem_gb=config.DEFAULT_MEMORY_MIN_GB, + metadata=metadata, + cifti_output=False, + name='preproc_carpet_wf', + ) + preproc_carpetplot_wf.inputs.inputnode.desc = 'preprocCarpetplot' + workflow.connect([ + (inputnode, preproc_carpetplot_wf, [ + ('bold_mask', 'inputnode.bold_mask'), + ('confounds_file', 'inputnode.confounds_file'), + ('skip_vols', 'inputnode.dummy_scans'), + ]), + (inputnode, preproc_carpetplot_wf, [('bold_file', 'inputnode.bold')]), + ]) # fmt:skip + for denoise_method in config.workflow.denoise_method: denoise = pe.Node( ICADenoise(method=denoise_method), @@ -476,6 +494,22 @@ def init_denoise_wf(bold_file): (add_non_steady_state, ds_denoised, [('bold_add', 'in_file')]), ]) # fmt:skip + carpetplot_wf = init_carpetplot_wf( + mem_gb=config.DEFAULT_MEMORY_MIN_GB, + metadata=metadata, + cifti_output=False, + name=f'{denoise_method}_carpet_wf', + ) + carpetplot_wf.inputs.inputnode.desc = f'{denoise_method}Carpetplot' + workflow.connect([ + (inputnode, carpetplot_wf, [ + ('bold_mask', 'inputnode.bold_mask'), + ('confounds_file', 'inputnode.confounds_file'), + ('skip_vols', 'inputnode.dummy_scans'), + ]), + (ds_denoised, carpetplot_wf, [('out_file', 'inputnode.bold')]), + ]) # fmt:skip + return workflow diff --git a/src/fmripost_aroma/workflows/base.py b/src/fmripost_aroma/workflows/base.py index 2ef7825..a0acb82 100644 --- a/src/fmripost_aroma/workflows/base.py +++ b/src/fmripost_aroma/workflows/base.py @@ -201,15 +201,6 @@ def init_single_subject_wf(subject_id: str): # Patch standard-space BOLD files into 'bold' key subject_data['bold'] = listify(subject_data['bold_mni152nlin6asym']) - if not subject_data['bold_mni152nlin6asym']: - task_id = config.execution.task_id - raise RuntimeError( - f"No MNI152NLin6Asym:res-2 BOLD images found for participant {subject_id} and " - f"task {task_id if task_id else ''}. " - "All workflows require MNI152NLin6Asym:res-2 BOLD images. " - f"Please check your BIDS filters: {config.execution.bids_filters}." - ) - # Make sure we always go through these two checks if not subject_data['bold']: task_id = config.execution.task_id @@ -406,10 +397,11 @@ def init_single_run_wf(bold_file): if config.workflow.denoise_method: # Now denoise the output-space BOLD data using ICA-AROMA - denoise_wf = init_denoise_wf(bold_file=bold_file) + denoise_wf = init_denoise_wf(bold_file=bold_file, metadata=bold_metadata) denoise_wf.inputs.inputnode.skip_vols = skip_vols denoise_wf.inputs.inputnode.space = 'MNI152NLin6Asym' denoise_wf.inputs.inputnode.res = '2' + denoise_wf.inputs.inputnode.confounds_file = functional_cache['confounds'] workflow.connect([ (mni6_buffer, denoise_wf, [ @@ -422,6 +414,12 @@ def init_single_run_wf(bold_file): ]), ]) # fmt:skip + # Fill-in datasinks seen so far + for node in workflow.list_node_names(): + if node.split('.')[-1].startswith('ds_'): + workflow.get_node(node).inputs.base_directory = config.execution.output_dir + workflow.get_node(node).inputs.source_file = bold_file + return workflow diff --git a/src/fmripost_aroma/workflows/confounds.py b/src/fmripost_aroma/workflows/confounds.py index c504ce8..d288efa 100644 --- a/src/fmripost_aroma/workflows/confounds.py +++ b/src/fmripost_aroma/workflows/confounds.py @@ -54,17 +54,11 @@ def init_carpetplot_wf( Inputs ------ bold - BOLD image, after the prescribed corrections (STC, HMC and SDC) - when available. + BOLD image, in MNI152NLin6Asym space + 2mm resolution. bold_mask - BOLD series mask + BOLD series mask in same space as ``bold``. confounds_file TSV of all aggregated confounds - boldref2anat_xfm - Affine matrix that maps the BOLD reference space into alignment with - the anatomical (T1w) space - std2anat_xfm - ANTs-compatible affine-and-warp transform file cifti_bold BOLD image in CIFTI format, to be used in place of volumetric BOLD crown_mask @@ -79,6 +73,7 @@ def init_carpetplot_wf( out_carpetplot Path of the generated SVG file """ + from fmriprep.interfaces.confounds import FMRISummary from nipype.interfaces import utility as niu from nipype.pipeline import engine as pe from niworkflows.engine.workflows import LiterateWorkflow as Workflow @@ -86,8 +81,7 @@ def init_carpetplot_wf( from templateflow.api import get as get_template from fmripost_aroma.config import DEFAULT_MEMORY_MIN_GB - from fmripost_aroma.interfaces import DerivativesDataSink - from fmripost_aroma.interfaces.confounds import FMRISummary + from fmripost_aroma.interfaces.bids import DerivativesDataSink inputnode = pe.Node( niu.IdentityInterface( @@ -95,13 +89,10 @@ def init_carpetplot_wf( 'bold', 'bold_mask', 'confounds_file', - 'boldref2anat_xfm', - 'std2anat_xfm', 'cifti_bold', - 'crown_mask', - 'acompcor_mask', 'dummy_scans', - ] + 'desc', + ], ), name='inputnode', ) @@ -113,10 +104,12 @@ def init_carpetplot_wf( FMRISummary( tr=metadata['RepetitionTime'], confounds_list=[ - ('global_signal', None, 'GS'), - ('csf', None, 'CSF'), - ('white_matter', None, 'WM'), - ('std_dvars', None, 'DVARS'), + ('trans_x', 'mm', 'x'), + ('trans_y', 'mm', 'y'), + ('trans_z', 'mm', 'z'), + ('rot_x', 'deg', 'pitch'), + ('rot_y', 'deg', 'roll'), + ('rot_z', 'deg', 'yaw'), ('framewise_displacement', 'mm', 'FD'), ], ), @@ -125,7 +118,6 @@ def init_carpetplot_wf( ) ds_report_bold_conf = pe.Node( DerivativesDataSink( - desc='carpetplot', datatype='figures', extension='svg', dismiss_entities=('echo', 'den', 'res'), @@ -137,10 +129,8 @@ def init_carpetplot_wf( parcels = pe.Node(niu.Function(function=_carpet_parcellation), name='parcels') parcels.inputs.nifti = not cifti_output - # List transforms - mrg_xfms = pe.Node(niu.Merge(2), name='mrg_xfms') - # Warp segmentation into EPI space + # Warp segmentation into MNI152NLin6Asym space resample_parc = pe.Node( ApplyTransforms( dimension=3, @@ -153,7 +143,17 @@ def init_carpetplot_wf( extension=['.nii', '.nii.gz'], ) ), - invert_transform_flags=[True, False], + transforms=[ + str( + get_template( + template='MNI152NLin6Asym', + mode='image', + suffix='xfm', + extension='.h5', + **{'from': 'MNI152NLin2009cAsym'}, + ), + ), + ], interpolation='MultiLabel', args='-u int', ), @@ -165,28 +165,22 @@ def init_carpetplot_wf( workflow.connect(inputnode, 'cifti_bold', conf_plot, 'in_cifti') workflow.connect([ - (inputnode, mrg_xfms, [ - ('boldref2anat_xfm', 'in1'), - ('std2anat_xfm', 'in2'), - ]), (inputnode, resample_parc, [('bold_mask', 'reference_image')]), - (inputnode, parcels, [('crown_mask', 'crown_mask')]), - (inputnode, parcels, [('acompcor_mask', 'acompcor_mask')]), (inputnode, conf_plot, [ ('bold', 'in_nifti'), ('confounds_file', 'confounds_file'), ('dummy_scans', 'drop_trs'), ]), - (mrg_xfms, resample_parc, [('out', 'transforms')]), (resample_parc, parcels, [('output_image', 'segmentation')]), (parcels, conf_plot, [('out', 'in_segm')]), + (inputnode, ds_report_bold_conf, [('desc', 'desc')]), (conf_plot, ds_report_bold_conf, [('out_file', 'in_file')]), (conf_plot, outputnode, [('out_file', 'out_carpetplot')]), ]) # fmt:skip return workflow -def _carpet_parcellation(segmentation, crown_mask, acompcor_mask, nifti=False): +def _carpet_parcellation(segmentation, nifti=False): """Generate the union of two masks.""" from pathlib import Path @@ -202,9 +196,9 @@ def _carpet_parcellation(segmentation, crown_mask, acompcor_mask, nifti=False): lut[255] = 5 if nifti else 0 # Cerebellum # Apply lookup table seg = lut[np.uint16(img.dataobj)] - seg[np.bool_(nb.load(crown_mask).dataobj)] = 6 if nifti else 2 + # seg[np.bool_(nb.load(crown_mask).dataobj)] = 6 if nifti else 2 # Separate deep from shallow WM+CSF - seg[np.bool_(nb.load(acompcor_mask).dataobj)] = 4 if nifti else 1 + # seg[np.bool_(nb.load(acompcor_mask).dataobj)] = 4 if nifti else 1 outimg = img.__class__(seg.astype('uint8'), img.affine, img.header) outimg.set_data_dtype('uint8')