Merge pull request #1270 from effigies/enh/sbref

ENH: Use single-band reference images when available

Author: effigies

docs/workflows.rst

@@ -294,6 +294,9 @@ BOLD reference image estimation
 
 This workflow estimates a reference image for a
 :abbr:`BOLD (blood-oxygen level-dependent)` series.
+If a single-band reference ("sbref") image associated with the BOLD series is
+available, then it is used directly.
+If not, a reference image is estimated from the BOLD series as follows:
 When T1-saturation effects ("dummy scans" or non-steady state volumes) are
 detected, they are averaged and used as reference due to their
 superior tissue contrast.

fmriprep/cli/run.py

@@ -107,7 +107,7 @@ def get_parser():
     g_conf = parser.add_argument_group('Workflow configuration')
     g_conf.add_argument(
         '--ignore', required=False, action='store', nargs="+", default=[],
-        choices=['fieldmaps', 'slicetiming'],
+        choices=['fieldmaps', 'slicetiming', 'sbref'],
         help='ignore selected aspects of the input dataset to disable corresponding '
              'parts of the workflow (a space delimited list)')
     g_conf.add_argument(

fmriprep/workflows/bold/base.py

@@ -242,6 +242,7 @@ def init_func_preproc_wf(bold_file, ignore, freesurfer,
                     'Memory resampled/largemem=%.2f/%.2f GB.'),
                ref_file, mem_gb['filesize'], bold_tlen, mem_gb['resampled'], mem_gb['largemem'])
 
+    sbref_file = None
     # For doc building purposes
     if layout is None or bold_file == 'bold_preprocesing':
         LOGGER.log(25, 'No valid layout: building empty workflow.')
@@ -259,6 +260,28 @@ def init_func_preproc_wf(bold_file, ignore, freesurfer,
         run_stc = True
         multiecho = False
     else:
+        # Find associated sbref, if possible
+        entities = layout.parse_file_entities(ref_file)
+        entities['type'] = 'sbref'
+        files = layout.get(**entities)
+        refbase = os.path.basename(ref_file)
+        if 'sbref' in ignore:
+            LOGGER.info("Single-band reference files ignored.")
+        elif files and multiecho:
+            LOGGER.warning("Single-band reference found, but not supported in "
+                           "multi-echo workflows at this time. Ignoring.")
+        elif files:
+            sbref_file = files[0].filename
+            sbbase = os.path.basename(sbref_file)
+            if len(files) > 1:
+                LOGGER.warning(
+                    "Multiple single-band reference files found for {}; using "
+                    "{}".format(refbase, sbbase))
+            else:
+                LOGGER.log(25, "Using single-band reference file {}".format(sbbase))
+        else:
+            LOGGER.log(25, "No single-band-reference found for {}".format(refbase))
+
         metadata = layout.get_metadata(ref_file)
 
         # Find fieldmaps. Options: (phase1|phase2|phasediff|epi|fieldmap|syn)
@@ -312,13 +335,15 @@ def init_func_preproc_wf(bold_file, ignore, freesurfer,
 """
 
     inputnode = pe.Node(niu.IdentityInterface(
-        fields=['bold_file', 'subjects_dir', 'subject_id',
+        fields=['bold_file', 'sbref_file', 'subjects_dir', 'subject_id',
                 't1_preproc', 't1_brain', 't1_mask', 't1_seg', 't1_tpms',
                 't1_aseg', 't1_aparc',
                 't1_2_mni_forward_transform', 't1_2_mni_reverse_transform',
                 't1_2_fsnative_forward_transform', 't1_2_fsnative_reverse_transform']),
         name='inputnode')
     inputnode.inputs.bold_file = bold_file
+    if sbref_file is not None:
+        inputnode.inputs.sbref_file = sbref_file
 
     outputnode = pe.Node(niu.IdentityInterface(
         fields=['bold_t1', 'bold_mask_t1', 'bold_aseg_t1', 'bold_aparc_t1', 'cifti_variant',
@@ -441,7 +466,8 @@ def init_func_preproc_wf(bold_file, ignore, freesurfer,
     # MAIN WORKFLOW STRUCTURE #######################################################
     workflow.connect([
         # Generate early reference
-        (inputnode, bold_reference_wf, [('bold_file', 'inputnode.bold_file')]),
+        (inputnode, bold_reference_wf, [('bold_file', 'inputnode.bold_file'),
+                                        ('sbref_file', 'inputnode.sbref_file')]),
         # BOLD buffer has slice-time corrected if it was run, original otherwise
         (boldbuffer, bold_split, [('bold_file', 'in_file')]),
         # HMC
@@ -549,7 +575,8 @@ def init_func_preproc_wf(bold_file, ignore, freesurfer,
         # Replace reference with the echo selected with FirstEcho
         workflow.disconnect([
             (inputnode, bold_reference_wf, [
-                ('bold_file', 'inputnode.bold_file')]),
+                ('bold_file', 'inputnode.bold_file'),
+                ('sbref_file', 'inputnode.sbref_file')]),
             (bold_reference_wf, boldbuffer, [
                 ('outputnode.bold_file', 'bold_file')]),
         ])

fmriprep/workflows/bold/util.py

@@ -87,7 +87,8 @@ def init_bold_reference_wf(omp_nthreads, bold_file=None, name='bold_reference_wf
 First, a reference volume and its skull-stripped version were generated
 using a custom methodology of *fMRIPrep*.
 """
-    inputnode = pe.Node(niu.IdentityInterface(fields=['bold_file']), name='inputnode')
+    inputnode = pe.Node(niu.IdentityInterface(fields=['bold_file', 'sbref_file']),
+                        name='inputnode')
     outputnode = pe.Node(
         niu.IdentityInterface(fields=['bold_file', 'raw_ref_image', 'skip_vols', 'ref_image',
                                       'ref_image_brain', 'bold_mask', 'validation_report',
@@ -102,16 +103,20 @@ def init_bold_reference_wf(omp_nthreads, bold_file=None, name='bold_reference_wf
 
     gen_ref = pe.Node(EstimateReferenceImage(), name="gen_ref",
                       mem_gb=1)  # OE: 128x128x128x50 * 64 / 8 ~ 900MB.
+    # Re-run validation; no effect if no sbref; otherwise apply same validation to sbref as bold
+    validate_ref = pe.Node(ValidateImage(), name='validate_ref', mem_gb=DEFAULT_MEMORY_MIN_GB)
     enhance_and_skullstrip_bold_wf = init_enhance_and_skullstrip_bold_wf(omp_nthreads=omp_nthreads)
 
     workflow.connect([
         (inputnode, validate, [('bold_file', 'in_file')]),
+        (inputnode, gen_ref, [('sbref_file', 'sbref_file')]),
         (validate, gen_ref, [('out_file', 'in_file')]),
-        (gen_ref, enhance_and_skullstrip_bold_wf, [('ref_image', 'inputnode.in_file')]),
+        (gen_ref, validate_ref, [('ref_image', 'in_file')]),
+        (validate_ref, enhance_and_skullstrip_bold_wf, [('out_file', 'inputnode.in_file')]),
         (validate, outputnode, [('out_file', 'bold_file'),
                                 ('out_report', 'validation_report')]),
-        (gen_ref, outputnode, [('ref_image', 'raw_ref_image'),
-                               ('n_volumes_to_discard', 'skip_vols')]),
+        (gen_ref, outputnode, [('n_volumes_to_discard', 'skip_vols')]),
+        (validate_ref, outputnode, [('out_file', 'raw_ref_image')]),
         (enhance_and_skullstrip_bold_wf, outputnode, [
             ('outputnode.bias_corrected_file', 'ref_image'),
             ('outputnode.mask_file', 'bold_mask'),