ENH: Add option to exclude projecting high variance voxels to surface #2855

madisoth · 2022-09-15T17:15:53Z

Changes proposed in this pull request

Adds option --project-goodvoxels to address #2822 .

If enabled, voxels whose timeseries have locally high coefficient of variation, or which lie outside the cortex (based on signed distance volumes generated from the WM and pial surfaces), are excluded from volume-to-surface sampling of BOLD timeseries.

Note that unlike HCP / ciftify, no dilation or spatial smoothing is applied to the surface-resampled BOLD timeseries.

Documentation that should be reviewed

…mmits) Squashed commits: [b1597c2] output goodvoxels-masked BOLD surfs in bold_surf_wf

…n a separate branch)

… separate branch)

welcome · 2022-09-15T17:15:55Z

Thanks for opening this pull request! It looks like this is your first time contributing to fMRIPrep. 😄
We invite you to list yourself as an fMRIPrep contributor. To learn more about what that entails and how we credit our contributors, please check out the contributing guidelines. If your name is not already on the list, please insert it, in alphabetical order, into the contributors.json file. Example:

   "name": "Contributor, New FMRIPrep",
   "affiliation": "Department of fMRI prep'ing, Open Science Made-Up University",
   "orcid": "<your id>"
}, ```

Of course, if you want to opt-out this time there is no problem at all with adding your name later. You will be always welcome to add it in the future whenever you feel it should be listed.

…urfaces for test sub)

effigies

Some initial thoughts. This is still marked draft. Are you planning on doing more?

effigies · 2022-09-22T17:34:37Z

docs/workflows.rst

@@ -488,6 +489,13 @@ surface, is sampled at 6 intervals and averaged.

 Surfaces are generated for the "subject native" surface, as well as transformed to the
 ``fsaverage`` template space.
+
+The flag ``--project_goodvoxels``, when enabled, excludes voxels whose timeseries have


Suggested change

The flag ``--project_goodvoxels``, when enabled, excludes voxels whose timeseries have

The flag ``--project-goodvoxels``, when enabled, excludes voxels whose timeseries have

effigies · 2022-09-22T17:35:58Z

fmriprep/interfaces/volume.py

+    surface = File(
+        exists=True,
+        mandatory=True,
+        argstr="%s ",


Suggested change

argstr="%s ",

argstr="%s",

effigies · 2022-09-22T17:39:10Z

fmriprep/interfaces/volume.py

+    surface = File(
+        exists=True,
+        mandatory=True,
+        argstr="%s ",
+        position=0,
+        desc="the input surface",


A couple things:

Nipype convention to use _file to indicate files.

No need for an extra space in argstr.

I think connectome_wb is pretty finicky about the type and extension, so we can try to indicate it in the doc.

Suggested change

surface = File(

exists=True,

mandatory=True,

argstr="%s ",

position=0,

desc="the input surface",

surf_file = File(

exists=True,

mandatory=True,

argstr="%s",

position=0,

desc="Input surface GIFTI file (.surf.gii)",

effigies · 2022-09-22T17:40:09Z

fmriprep/interfaces/volume.py

+    ref_space = File(
+        exists=True,
+        mandatory=True,
+        argstr="%s ",
+        position=1,
+        desc="a volume in the desired output space (dims, spacing, origin)",
+    )


Suggested change

ref_space = File(

exists=True,

mandatory=True,

argstr="%s ",

position=1,

desc="a volume in the desired output space (dims, spacing, origin)",

)

ref_file = File(

exists=True,

mandatory=True,

argstr="%s",

position=1,

desc="NIfTI volume in the desired output space (dims, spacing, origin)",

)

effigies · 2022-09-22T17:40:57Z

fmriprep/interfaces/volume.py

+    out_vol = File(
+        name_source=["surface"],
+        name_template="%s.distvol.nii.gz",
+        argstr="%s ",
+        position=2,
+        desc="output - the output volume",
+    )


I would generally use _ to separate a descriptive tag, rather than a .:

Suggested change

out_vol = File(

name_source=["surface"],

name_template="%s.distvol.nii.gz",

argstr="%s ",

position=2,

desc="output - the output volume",

)

out_file = File(

name_source=["surface"],

name_template="%s_distvol.nii.gz",

argstr="%s",

position=2,

desc="Name of output volume containing signed distances",

)

effigies · 2022-09-22T17:53:13Z

fmriprep/interfaces/volume.py

+        argstr="-winding %s ",
+        usedefault=True,
+        position=8,        
+        desc="winding method for point inside surface test, choices: "
+             "EVEN_ODD (default) "
+             "NEGATIVE "
+             "NONZERO "
+             "NORMALS "


Suggested change

argstr="-winding %s ",

usedefault=True,

position=8,

desc="winding method for point inside surface test, choices: "

"EVEN_ODD (default) "

"NEGATIVE "

"NONZERO "

"NORMALS "

argstr="-winding %s",

usedefault=True,

desc="winding method for point inside surface test"

effigies · 2022-09-22T17:54:37Z

fmriprep/interfaces/volume.py

+    out_vol = File(exists=True, desc="output - the output volume")
+    roi_out = File(desc="output roi volume of where the output has a computed value")


exists=True is asking for trouble: nipy/nipype#3493

Suggested change

out_vol = File(exists=True, desc="output - the output volume")

roi_out = File(desc="output roi volume of where the output has a computed value")

out_file = File(desc="Name of output volume containing signed distances")

out_mask = File(desc="Name of file to store a mask where the ``out_file`` has a computed value")

effigies · 2022-09-22T17:54:57Z

fmriprep/interfaces/volume.py

+    """
+CREATE SIGNED DISTANCE VOLUME FROM SURFACE


Suggested change

"""

CREATE SIGNED DISTANCE VOLUME FROM SURFACE

"""CREATE SIGNED DISTANCE VOLUME FROM SURFACE

effigies · 2022-09-22T17:55:19Z

fmriprep/interfaces/volume.py

+
+    def _list_outputs(self):
+        outputs = super(CreateSignedDistanceVolume, self)._list_outputs()
+        return outputs


Suggested change

def _list_outputs(self):

outputs = super(CreateSignedDistanceVolume, self)._list_outputs()

return outputs

effigies · 2022-09-22T18:13:17Z

fmriprep/workflows/bold/resampling.py

+    # 0, 1 = wm; 2, 3 = pial; 6, 7 = mid
+    # note that order of lh / rh within each surf type is not guaranteed due to use
+    # of unsorted glob by FreeSurferSource prior, but we can do a sort with
+    # _sorted_by_basename to ensure consistent ordering
+    select_wm = pe.Node(
+        niu.Select(index=[0, 1]),
+        name="select_wm",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    select_pial = pe.Node(
+        niu.Select(index=[2, 3]),
+        name="select_pial",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    select_midthick = pe.Node(
+        niu.Select(index=[6, 7]),
+        name="select_midthick",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    create_wm_distvol = pe.MapNode(
+        CreateSignedDistanceVolume(),
+        iterfield=["surface"],


This should really be its own workflow at this point. In this workflow, we could then have a small adaptation:

if project_goodvoxels: ribbon_mask_wf = init_ribbon_mask_wf() apply_ribbon_mask = pe.Node( fsl.ApplyMask(), name_source=['in_file'], name="apply_goodvoxels_ribbon_mask", mem_gb=mem_gb * 3, ) wf.connect([ (inputnode, ribbon_mask_wf, [('anat_giftis', 'anat_giftis'), ('t1w_mask', 't1w_mask')]), (rename_src, apply_ribbon_mask, [('out_file', 'in_file')]), (ribbon_mask_wf, apply_ribbon_mask, [('mask_file', 'mask_file')]), (apply_ribbon_mask, sampler, [('out_file', 'in_file')]), ]) else: wf.connect([ (rename_src, sampler, [('out_file', 'in_file')]), ])

mgxd

My initial review on nibabies overlapped a lot with Chris's comments, so only ported over novel points - There are a lot of styling changes that hopefully will be addressed once #2857 is merged.

mgxd · 2022-09-23T20:22:56Z

fmriprep/workflows/bold/resampling.py

+    select_wm = pe.Node(
+        niu.Select(index=[0, 1]),
+        name="select_wm",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,


Generally, we can run most nipype utility interfaces on the master process since they are not intense operations - no need to potentially submit to a job queue.

Suggested change

mem_gb=DEFAULT_MEMORY_MIN_GB,

mem_gb=DEFAULT_MEMORY_MIN_GB,

run_without_submitting=True,

effigies

It might be worth taking the opportunity to rewrite some of this workflow instead of setting up a long fslmaths chain. Just a quick sketch, below. Happy to discuss at a convenient time.

effigies · 2022-09-22T18:55:58Z

fmriprep/workflows/bold/resampling.py

+    thresh_wm_distvol = pe.MapNode(
+        fsl.maths.MathsCommand(args="-thr 0 -bin -mul 255"),
+        iterfield=["in_file"],
+        name="thresh_wm_distvol",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    uthresh_pial_distvol = pe.MapNode(
+        fsl.maths.MathsCommand(args="-uthr 0 -abs -bin -mul 255"),
+        iterfield=["in_file"],
+        name="uthresh_pial_distvol",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    bin_wm_distvol = pe.MapNode(
+        fsl.maths.UnaryMaths(operation="bin"),
+        iterfield=["in_file"],
+        name="bin_wm_distvol",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    bin_pial_distvol = pe.MapNode(
+        fsl.maths.UnaryMaths(operation="bin"),
+        iterfield=["in_file"],
+        name="bin_pial_distvol",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    split_wm_distvol = pe.Node(
+        niu.Split(splits=[1, 1]),
+        name="split_wm_distvol",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    merge_wm_distvol_no_flatten = pe.Node(
+        niu.Merge(2),
+        no_flatten=True,
+        name="merge_wm_distvol_no_flatten",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    make_ribbon_vol = pe.MapNode(
+        fsl.maths.MultiImageMaths(op_string="-mas %s -mul 255 "),
+        iterfield=["in_file", "operand_files"],
+        name="make_ribbon_vol",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    bin_ribbon_vol = pe.MapNode(
+        fsl.maths.UnaryMaths(operation="bin"),
+        iterfield=["in_file"],
+        name="bin_ribbon_vol",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    split_squeeze_ribbon_vol = pe.Node(
+        niu.Split(splits=[1, 1], squeeze=True),
+        name="split_squeeze_ribbon",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )
+
+    combine_ribbon_vol_hemis = pe.Node(
+        fsl.maths.BinaryMaths(operation="add"),
+        name="combine_ribbon_vol_hemis",
+        mem_gb=DEFAULT_MEMORY_MIN_GB,
+    )


I suspect this workflow could be made more understandable with a couple pure Python interfaces. For instance, I think we could replace all of this with one function:

def make_ribbon(pial_distances, white_distances): dist_lp = nb.load(pial_distances[0]) dist_rp = nb.load(pial_distances[1]) dist_lw = nb.load(white_distances[0]) dist_rw = nb.load(white_distances[1]) lmask = (dist_lw.get_fdata(dtype='f4') > 0) & (dist_lp.get_fdata(dtype='f4') < 0) rmask = (dist_rw.get_fdata(dtype='f4') > 0) & (dist_rp.get_fdata(dtype='f4') < 0) mask = np.uint8(lmask | rmask) mask_img = nb.Nifti1Image(mask, dist_lp.affine, dist_lp.header) mask_img.set_data_dtype(np.uint8) return mask_img

Functions can then be promoted into interfaces with:

class MakeRibbon(SimpleInterface): input_spec = ... output_spec = ... def _run_interface(self, runtime): mask_img = make_ribbon(self.inputs.pial_distances, self.inputs.white_distances) mask_img.to_filename("ribbon.nii.gz") self._outputs["out_file"] = str(Path(runtime.cwd) / "ribbon.nii.gz") return runtime

- Implemented the formatting suggestions by Mathias at the URLs below: - nipreps#235 - nipreps/fmriprep#2855 - Split the 3 workflow-connect blocks in the project-goodvoxels section of resampling.py into 3 different functions: ribbon, outliers, and vol-to-surf - Incomplete comments in resampling.py

effigies · 2023-01-26T14:16:17Z

This needs rebasing/merging with master.

@madisoth Otherwise this is just waiting on templateflow/tpl-fsaverage#2?

FIX: Minor anat wf & workbench.py fixes (+5 squashed commit) Squashed commit: [e995ec4] RF: Moved ribbon wf (from bold to anat) and reformatted [9676f14] Updated Docker version tag [532875e] STY: black and isort [946dc7e] Formatting fixes and resampling.py modularization - Implemented the formatting suggestions by Mathias at the URLs below: - nipreps#235 - nipreps/fmriprep#2855 - Split the 3 workflow-connect blocks in the project-goodvoxels section of resampling.py into 3 different functions: ribbon, outliers, and vol-to-surf - Incomplete comments in resampling.py [114f982] Modularized goodvoxels code in resampling.py Removed duplicate code by abstracting it into functions, clarifying the goodvoxels logic

… (update of #2855) (#2956) Adds option --project-goodvoxels to address #2822. If enabled, voxels whose timeseries have locally high coefficient of variation, or which lie outside the cortex (based on signed distance volumes generated from the WM and pial surfaces), are excluded from volume-to-surface sampling of BOLD timeseries. Not enabled by default. - adds MetricDilate interface for Workbench's wb_command -metric-dilate, which, when using the --project-goodvoxels option, fills zero-value vertices in the BOLD surface timeseries by dilating in data from the nearest "good" (nonzero, nonempty) vertex.

Eric Feczko and others added 8 commits July 28, 2022 16:59

adding goodvoxels freesurfer branch

26eb6e1

added workflow connectors for next stages -- one more needs to be added

0e52bf3

Create init workflow for masking medial wall outliers

e06c412

output "goodvoxels"-masked BOLD surfs in bold_surf_wf (+1 squashed co…

64e61e7

…mmits) Squashed commits: [b1597c2] output goodvoxels-masked BOLD surfs in bold_surf_wf

goodvoxels_gifti work

260c7a6

progress on workbench / gifti surface wf

3041e25

revert work on workbench-based surface sampling wf (to be continued i…

30ae531

…n a separate branch)

revert work on workench-based surface sampling wf (to be continued in…

c555283

… separate branch)

madisoth added 3 commits September 15, 2022 10:43

volume.py docstring import fix

ba7f699

remove CreateSignedDistanceVolume doctest for now (needs anatomical s…

a680447

…urfaces for test sub)

remove unused itk_bold_to_t1 from bold_surf_wf inputnode

17fba92

effigies reviewed Sep 22, 2022

View reviewed changes

mgxd reviewed Sep 23, 2022

View reviewed changes

mgxd assigned madisoth Sep 23, 2022

mgxd changed the title ~~Enh/project goodvoxels~~ ENH: A Sep 23, 2022

mgxd changed the title ~~ENH: A~~ ENH: Add option to exclude projecting high variance voxels to surface Sep 23, 2022

effigies reviewed Sep 25, 2022

View reviewed changes

GregConan mentioned this pull request Dec 5, 2022

Reorganized/cleaned goodvoxels resampling.py madisoth/nibabies#2

Merged

Untested 1st draft integrating madisoth/nibabies:enh/project-goodvoxels

ca9d5d4

effigies added this to the 23.0.0 milestone Dec 7, 2022

madisoth added 2 commits December 15, 2022 12:40

refactor to handle anat_ribbon_wf moving to sMRIPrep

8bf03e2

refactor to handle anat_ribbon_wf moving to sMRIPrep

97e2b79

effigies mentioned this pull request Feb 15, 2023

ENH: Add option to exclude projecting high variance voxels to surface #2953

Closed

1 task

effigies closed this Feb 17, 2023

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

ENH: Add option to exclude projecting high variance voxels to surface #2855

ENH: Add option to exclude projecting high variance voxels to surface #2855

madisoth commented Sep 15, 2022

welcome bot commented Sep 15, 2022

effigies left a comment

effigies Sep 22, 2022

effigies Sep 22, 2022

effigies Sep 22, 2022

effigies Sep 22, 2022

effigies Sep 22, 2022

effigies Sep 22, 2022

effigies Sep 22, 2022

effigies Sep 22, 2022

effigies Sep 22, 2022

effigies Sep 22, 2022

mgxd left a comment

mgxd Sep 23, 2022

effigies left a comment

effigies Sep 22, 2022

effigies commented Jan 26, 2023

	The flag ``--project_goodvoxels``, when enabled, excludes voxels whose timeseries have
	The flag ``--project-goodvoxels``, when enabled, excludes voxels whose timeseries have

		out_vol = File(exists=True, desc="output - the output volume")
		roi_out = File(desc="output roi volume of where the output has a computed value")

	mem_gb=DEFAULT_MEMORY_MIN_GB,
	mem_gb=DEFAULT_MEMORY_MIN_GB,
	run_without_submitting=True,

ENH: Add option to exclude projecting high variance voxels to surface #2855

ENH: Add option to exclude projecting high variance voxels to surface #2855

Conversation

madisoth commented Sep 15, 2022

Changes proposed in this pull request

Documentation that should be reviewed

welcome bot commented Sep 15, 2022

effigies left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mgxd left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

effigies left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

effigies commented Jan 26, 2023