Merge branch 'develop' into many_pipes

Author: sgiavasis

CHANGELOG.md

@@ -26,6 +26,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Optional positional argument "cfg" to `create_lesion_preproc`.
 - Allow enabling `overwrite_transform` only when the registration method is `ANTS`.
 - `resource_inventory` utility to inventory NodeBlock function inputs and outputs.
+- New switch `mask_sbref` under `func_input_prep` in functional registration and set to default `on`.
+- New resource `desc-head_bold` as non skull-stripped bold from nodeblock `bold_masking`.
+- `censor_file_path` from `offending_timepoints_connector` in the `build_nuisance_regressor` node.
 
 ### Changed
 
@@ -34,17 +37,20 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Made orientation configurable (was hard-coded as "RPI").
 - Resource-not-found errors now include information about where to source those resources.
 - Moved `ref_mask_res_2` and `T1w_template_res-2` fields from registration into surface under `abcd_prefreesurfer_prep`.
+- Moved `find_censors node` inside `create_nuisance_regression_workflow` into its own function/subworkflow as `offending_timepoints_connector`.
 - [FSL-AFNI subworkflow](https://github.com/FCP-INDI/C-PAC/blob/4bdd6c410ef0a9b90f53100ea005af1f7d6e76c0/CPAC/func_preproc/func_preproc.py#L1052C4-L1231C25)
   - Moved `FSL-AFNI subworkflow` from inside a `bold_mask_fsl_afni` nodeblock into a separate function.
   - Renamed `desc-ref_bold` created in this workflow to `desc-unifized_bold`.
   - `coregistration_prep_fmriprep` nodeblock now checks if `desc-unifized_bold` exists in the Resource Pool, if not it runs the `FSL-AFNI subworkflow` to create it.
+- Input `desc-brain_bold` to `desc-preproc_bold` for `sbref` generation nodeblock `coregistration_prep_vol`.
 
 ### Fixed
 
 - A bug in which AWS S3 encryption was looked for in Nipype config instead of pipeline config (only affected uploading logs).
 - Restored `bids-validator` functionality.
 - Fixed empty `shell` variable in cluster run scripts.
 - A bug in which bandpass filters always assumed 1D regressor files have exactly 5 header rows.
+- Removed an erroneous connection to AFNI 3dTProject in nuisance denoising that would unnecessarily send a spike regressor as a censor. This would sometimes cause TRs to unnecessarily be dropped from the timeseries as if scrubbing were being performed.
 
 ### Removed
 
@@ -57,7 +63,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
   - as output from FNIRT registration.
   - as inputs from Nodeblocks requesting it and, replaced with `space-template_desc-brain_mask`.
   - from outputs tsv.
-
+- Inputs `[desc-motion_bold, bold]` from `coregistration_prep_vol` nodeblock.
+- `input` field from `coregistration` in blank and default config.
+- `reg_with_skull` swtich from `func_input_prep` in blank and default config.
 
 ## [1.8.7] - 2024-05-03
 

CPAC/func_preproc/func_preproc.py

@@ -1830,6 +1830,10 @@ def bold_mask_ccs(wf, cfg, strat_pool, pipe_num, opt=None):
             "Description": "The skull-stripped BOLD time-series.",
             "SkullStripped": True,
         },
+        "desc-head_bold": {
+            "Description": "The non skull-stripped BOLD time-series.",
+            "SkullStripped": False,
+        },
     },
 )
 def bold_masking(wf, cfg, strat_pool, pipe_num, opt=None):
@@ -1841,15 +1845,16 @@ def bold_masking(wf, cfg, strat_pool, pipe_num, opt=None):
     func_edge_detect.inputs.expr = "a*b"
     func_edge_detect.inputs.outputtype = "NIFTI_GZ"
 
-    node, out = strat_pool.get_data("desc-preproc_bold")
-    wf.connect(node, out, func_edge_detect, "in_file_a")
+    node_head_bold, out_head_bold = strat_pool.get_data("desc-preproc_bold")
+    wf.connect(node_head_bold, out_head_bold, func_edge_detect, "in_file_a")
 
     node, out = strat_pool.get_data("space-bold_desc-brain_mask")
     wf.connect(node, out, func_edge_detect, "in_file_b")
 
     outputs = {
         "desc-preproc_bold": (func_edge_detect, "out_file"),
         "desc-brain_bold": (func_edge_detect, "out_file"),
+        "desc-head_bold": (node_head_bold, out_head_bold),
     }
 
     return (wf, outputs)

CPAC/nuisance/nuisance.py

@@ -496,6 +496,104 @@ def gather_nuisance(
     return output_file_path, censor_indices
 
 
+def offending_timepoints_connector(
+    nuisance_selectors, name="offending_timepoints_connector"
+):
+    inputspec = pe.Node(
+        util.IdentityInterface(
+            fields=[
+                "fd_j_file_path",
+                "fd_p_file_path",
+                "dvars_file_path",
+            ]
+        ),
+        name="inputspec",
+    )
+
+    wf = pe.Workflow(name=name)
+
+    outputspec = pe.Node(
+        util.IdentityInterface(fields=["out_file"]),
+        name="outputspec",
+    )
+
+    censor_selector = nuisance_selectors.get("Censor")
+
+    find_censors = pe.Node(
+        Function(
+            input_names=[
+                "fd_j_file_path",
+                "fd_j_threshold",
+                "fd_p_file_path",
+                "fd_p_threshold",
+                "dvars_file_path",
+                "dvars_threshold",
+                "number_of_previous_trs_to_censor",
+                "number_of_subsequent_trs_to_censor",
+            ],
+            output_names=["out_file"],
+            function=find_offending_time_points,
+            as_module=True,
+        ),
+        name="find_offending_time_points",
+    )
+
+    if not censor_selector.get("thresholds"):
+        msg = "Censoring requested, but thresh_metric not provided."
+        raise ValueError(msg)
+
+    for threshold in censor_selector["thresholds"]:
+        if "type" not in threshold or threshold["type"] not in [
+            "DVARS",
+            "FD_J",
+            "FD_P",
+        ]:
+            msg = "Censoring requested, but with invalid threshold type."
+            raise ValueError(msg)
+
+        if "value" not in threshold:
+            msg = "Censoring requested, but threshold not provided."
+            raise ValueError(msg)
+
+        if threshold["type"] == "FD_J":
+            find_censors.inputs.fd_j_threshold = threshold["value"]
+            wf.connect(inputspec, "fd_j_file_path", find_censors, "fd_j_file_path")
+
+        if threshold["type"] == "FD_P":
+            find_censors.inputs.fd_p_threshold = threshold["value"]
+            wf.connect(inputspec, "fd_p_file_path", find_censors, "fd_p_file_path")
+
+        if threshold["type"] == "DVARS":
+            find_censors.inputs.dvars_threshold = threshold["value"]
+            wf.connect(inputspec, "dvars_file_path", find_censors, "dvars_file_path")
+
+    if (
+        censor_selector.get("number_of_previous_trs_to_censor")
+        and censor_selector["method"] != "SpikeRegression"
+    ):
+        find_censors.inputs.number_of_previous_trs_to_censor = censor_selector[
+            "number_of_previous_trs_to_censor"
+        ]
+
+    else:
+        find_censors.inputs.number_of_previous_trs_to_censor = 0
+
+    if (
+        censor_selector.get("number_of_subsequent_trs_to_censor")
+        and censor_selector["method"] != "SpikeRegression"
+    ):
+        find_censors.inputs.number_of_subsequent_trs_to_censor = censor_selector[
+            "number_of_subsequent_trs_to_censor"
+        ]
+
+    else:
+        find_censors.inputs.number_of_subsequent_trs_to_censor = 0
+
+    wf.connect(find_censors, "out_file", outputspec, "out_file")
+
+    return wf
+
+
 def create_regressor_workflow(
     nuisance_selectors,
     use_ants,
@@ -1547,6 +1645,38 @@ def create_regressor_workflow(
         "functional_file_path",
     )
 
+    if nuisance_selectors.get("Censor"):
+        if nuisance_selectors["Censor"]["method"] == "SpikeRegression":
+            offending_timepoints_connector_wf = offending_timepoints_connector(
+                nuisance_selectors
+            )
+            nuisance_wf.connect(
+                [
+                    (
+                        inputspec,
+                        offending_timepoints_connector_wf,
+                        [("fd_j_file_path", "inputspec.fd_j_file_path")],
+                    ),
+                    (
+                        inputspec,
+                        offending_timepoints_connector_wf,
+                        [("fd_p_file_path", "inputspec.fd_p_file_path")],
+                    ),
+                    (
+                        inputspec,
+                        offending_timepoints_connector_wf,
+                        [("dvars_file_path", "inputspec.dvars_file_path")],
+                    ),
+                ]
+            )
+
+            nuisance_wf.connect(
+                offending_timepoints_connector_wf,
+                "outputspec.out_file",
+                build_nuisance_regressors,
+                "censor_file_path",
+            )
+
     build_nuisance_regressors.inputs.selector = nuisance_selectors
 
     # Check for any regressors to combine into files
@@ -1656,93 +1786,28 @@ def create_nuisance_regression_workflow(nuisance_selectors, name="nuisance_regre
     nuisance_wf = pe.Workflow(name=name)
 
     if nuisance_selectors.get("Censor"):
-        censor_methods = ["Kill", "Zero", "Interpolate", "SpikeRegression"]
-
-        censor_selector = nuisance_selectors.get("Censor")
-        if censor_selector.get("method") not in censor_methods:
-            msg = (
-                "Improper censoring method specified ({0}), "
-                "should be one of {1}.".format(
-                    censor_selector.get("method"), censor_methods
-                )
-            )
-            raise ValueError(msg)
-
-        find_censors = pe.Node(
-            Function(
-                input_names=[
-                    "fd_j_file_path",
-                    "fd_j_threshold",
-                    "fd_p_file_path",
-                    "fd_p_threshold",
-                    "dvars_file_path",
-                    "dvars_threshold",
-                    "number_of_previous_trs_to_censor",
-                    "number_of_subsequent_trs_to_censor",
-                ],
-                output_names=["out_file"],
-                function=find_offending_time_points,
-                as_module=True,
-            ),
-            name="find_offending_time_points",
+        offending_timepoints_connector_wf = offending_timepoints_connector(
+            nuisance_selectors
         )
-
-        if not censor_selector.get("thresholds"):
-            msg = "Censoring requested, but thresh_metric not provided."
-            raise ValueError(msg)
-
-        for threshold in censor_selector["thresholds"]:
-            if "type" not in threshold or threshold["type"] not in [
-                "DVARS",
-                "FD_J",
-                "FD_P",
-            ]:
-                msg = "Censoring requested, but with invalid threshold type."
-                raise ValueError(msg)
-
-            if "value" not in threshold:
-                msg = "Censoring requested, but threshold not provided."
-                raise ValueError(msg)
-
-            if threshold["type"] == "FD_J":
-                find_censors.inputs.fd_j_threshold = threshold["value"]
-                nuisance_wf.connect(
-                    inputspec, "fd_j_file_path", find_censors, "fd_j_file_path"
-                )
-
-            if threshold["type"] == "FD_P":
-                find_censors.inputs.fd_p_threshold = threshold["value"]
-                nuisance_wf.connect(
-                    inputspec, "fd_p_file_path", find_censors, "fd_p_file_path"
-                )
-
-            if threshold["type"] == "DVARS":
-                find_censors.inputs.dvars_threshold = threshold["value"]
-                nuisance_wf.connect(
-                    inputspec, "dvars_file_path", find_censors, "dvars_file_path"
-                )
-
-        if (
-            censor_selector.get("number_of_previous_trs_to_censor")
-            and censor_selector["method"] != "SpikeRegression"
-        ):
-            find_censors.inputs.number_of_previous_trs_to_censor = censor_selector[
-                "number_of_previous_trs_to_censor"
-            ]
-
-        else:
-            find_censors.inputs.number_of_previous_trs_to_censor = 0
-
-        if (
-            censor_selector.get("number_of_subsequent_trs_to_censor")
-            and censor_selector["method"] != "SpikeRegression"
-        ):
-            find_censors.inputs.number_of_subsequent_trs_to_censor = censor_selector[
-                "number_of_subsequent_trs_to_censor"
+        nuisance_wf.connect(
+            [
+                (
+                    inputspec,
+                    offending_timepoints_connector_wf,
+                    [("fd_j_file_path", "inputspec.fd_j_file_path")],
+                ),
+                (
+                    inputspec,
+                    offending_timepoints_connector_wf,
+                    [("fd_p_file_path", "inputspec.fd_p_file_path")],
+                ),
+                (
+                    inputspec,
+                    offending_timepoints_connector_wf,
+                    [("dvars_file_path", "inputspec.dvars_file_path")],
+                ),
             ]
-
-        else:
-            find_censors.inputs.number_of_subsequent_trs_to_censor = 0
+        )
 
     # Use 3dTproject to perform nuisance variable regression
     nuisance_regression = pe.Node(
@@ -1757,17 +1822,19 @@ def create_nuisance_regression_workflow(nuisance_selectors, name="nuisance_regre
     nuisance_regression.inputs.norm = False
 
     if nuisance_selectors.get("Censor"):
-        if nuisance_selectors["Censor"]["method"] == "SpikeRegression":
-            nuisance_wf.connect(find_censors, "out_file", nuisance_regression, "censor")
-        else:
-            if nuisance_selectors["Censor"]["method"] == "Interpolate":
-                nuisance_regression.inputs.cenmode = "NTRP"
-            else:
-                nuisance_regression.inputs.cenmode = nuisance_selectors["Censor"][
-                    "method"
-                ].upper()
+        if nuisance_selectors["Censor"]["method"] != "SpikeRegression":
+            nuisance_regression.inputs.cenmode = (
+                "NTRP"
+                if nuisance_selectors["Censor"]["method"] == "Interpolate"
+                else nuisance_selectors["Censor"]["method"].upper()
+            )
 
-            nuisance_wf.connect(find_censors, "out_file", nuisance_regression, "censor")
+            nuisance_wf.connect(
+                offending_timepoints_connector_wf,
+                "outputspec.out_file",
+                nuisance_regression,
+                "censor",
+            )
 
     if nuisance_selectors.get("PolyOrt"):
         if not nuisance_selectors["PolyOrt"].get("degree"):

CPAC/nuisance/utils/utils.py

@@ -139,7 +139,7 @@ def find_offending_time_points(
     censor_vector[extended_censors] = 0
 
     out_file_path = os.path.join(os.getcwd(), "censors.tsv")
-    np.savetxt(out_file_path, censor_vector, fmt="%d", header="censor", comments="")
+    np.savetxt(out_file_path, censor_vector, fmt="%d", comments="")
 
     return out_file_path
 

CPAC/pipeline/cpac_pipeline.py

@@ -148,6 +148,7 @@
     coregistration_prep_vol,
     create_func_to_T1template_symmetric_xfm,
     create_func_to_T1template_xfm,
+    mask_sbref,
     overwrite_transform_anat_to_template,
     register_ANTs_anat_to_template,
     register_ANTs_EPI_to_template,
@@ -1288,6 +1289,7 @@ def build_workflow(subject_id, sub_dict, cfg, pipeline_name=None):
                 coregistration_prep_mean,
                 coregistration_prep_fmriprep,
             ],
+            mask_sbref,
         ]
 
         # Distortion/Susceptibility Correction