Add custom censoring options

xcp_d/cli/parser.py

@@ -8,6 +8,7 @@
 import os
 import sys
 
+import xcp_d.cli.parser_utils as types
 from xcp_d import config
 
 
@@ -19,7 +20,6 @@ def _build_parser():
 
     from packaging.version import Version
 
-    from xcp_d.cli.parser_utils import _float_or_auto, _int_or_auto, _restricted_float
     from xcp_d.cli.version import check_latest, is_flagged
     from xcp_d.utils.atlas import select_atlases
 
@@ -249,7 +249,7 @@ def _bids_filter(value, parser):
         "--dummy_scans",
         dest="dummy_scans",
         default=0,
-        type=_int_or_auto,
+        type=types._int_or_auto,
         metavar="{{auto,INT}}",
         help=(
             "Number of volumes to remove from the beginning of each run. "
@@ -391,7 +391,7 @@ def _bids_filter(value, parser):
         "--head-radius",
         "--head_radius",
         default=50,
-        type=_float_or_auto,
+        type=types._float_or_auto,
         help=(
             "Head radius used to calculate framewise displacement, in mm. "
             "The default value is 50 mm, which is recommended for adults. "
@@ -404,15 +404,67 @@ def _bids_filter(value, parser):
         "-f",
         "--fd-thresh",
         "--fd_thresh",
-        default=0.3,
+        default=[0.3, 0],
+        metavar="FLOAT",
         type=float,
+        nargs="+",
         help=(
-            "Framewise displacement threshold for censoring. "
-            "Any volumes with an FD value greater than the threshold will be removed from the "
-            "denoised BOLD data. "
+            "Framewise displacement thresholds for censoring. "
+            "This may be a single value or a pair of values. "
+            "If two values are provided, the first one will determine which volumes are "
+            "removed from the denoised BOLD data, and the second one will determine which "
+            "volumes are removed from the interpolated BOLD data. "
             "A threshold of <=0 will disable censoring completely."
         ),
     )
+    g_censor.add_argument(
+        "--dvars-thresh",
+        "--dvars_thresh",
+        default=[0, 0],
+        metavar="FLOAT",
+        type=float,
+        nargs="+",
+        help=(
+            "DVARS threshold for censoring. "
+            "This may be a single value or a pair of values. "
+            "If two values are provided, the first one will determine which volumes are "
+            "removed from the denoised BOLD data, and the second one will determine which "
+            "volumes are removed from the interpolated BOLD data. "
+            "A threshold of <=0 will disable censoring completely."
+        ),
+    )
+    g_censor.add_argument(
+        "--censor-before",
+        "--censor_before",
+        default=[0, 0],
+        metavar="INT",
+        nargs="+",
+        type=int,
+        help="The number of volumes to remove before any outlier volumes.",
+    )
+    g_censor.add_argument(
+        "--censor-after",
+        "--censor_after",
+        default=[0, 0],
+        metavar="INT",
+        nargs="+",
+        type=int,
+        help="The number of volumes to remove after any outlier volumes.",
+    )
+    g_censor.add_argument(
+        "--censor-between",
+        "--censor_between",
+        default=[0, 0],
+        metavar="INT",
+        nargs="+",
+        type=int,
+        help=(
+            "If any short sets of contiguous non-outliers are found between outliers, "
+            "this parameter will remove them. "
+            "For example, if the value is set to 1, then any cases where only one non-outlier "
+            "volume exists between two outlier volumes will be censored."
+        ),
+    )
     g_censor.add_argument(
         "--min-time",
         "--min_time",
@@ -509,7 +561,7 @@ def _bids_filter(value, parser):
         "--min_coverage",
         required=False,
         default=0.5,
-        type=_restricted_float,
+        type=types._restricted_float,
         help=(
             "Coverage threshold to apply to parcels in each atlas. "
             "Any parcels with lower coverage than the threshold will be replaced with NaNs. "
@@ -875,10 +927,23 @@ def _validate_parameters(opts, build_log, parser):
         build_log.warning("Bandpass filtering is disabled. ALFF outputs will not be generated.")
 
     # Scrubbing parameters
-    if opts.fd_thresh <= 0 and opts.min_time > 0:
+    opts.fd_thresh = types._check_censoring_thresholds(opts.fd_thresh, parser, "--fd-thresh")
+    opts.dvars_thresh = types._check_censoring_thresholds(
+        opts.dvars_thresh, parser, "--dvars-thresh"
+    )
+    opts.censor_before = types._check_censoring_numbers(
+        opts.censor_before, parser, "--censor-before"
+    )
+    opts.censor_after = types._check_censoring_numbers(opts.censor_after, parser, "--censor-after")
+    opts.censor_between = types._check_censoring_numbers(
+        opts.censor_between, parser, "--censor-between"
+    )
+
+    nocensor = any(t <= 0 for t in opts.fd_thresh + opts.dvars_thresh)
+    if nocensor and opts.min_time > 0:
         ignored_params = "\n\t".join(["--min-time"])
         build_log.warning(
-            "Framewise displacement-based scrubbing is disabled. "
+            "Censoring is disabled. "
             f"The following parameters will have no effect:\n\t{ignored_params}"
         )
         opts.min_time = 0

xcp_d/cli/parser_utils.py

@@ -73,6 +73,52 @@ def _float_or_auto(string, is_parser=True):
     return floatarg
 
 
+def _one_or_two_ints(string, is_parser=True):
+    """Check if argument is one or two integers >= 0."""
+    error = argparse.ArgumentTypeError if is_parser else ValueError
+    try:
+        ints = [int(i) for i in string.split()]
+    except ValueError:
+        msg = "Argument must be one or two integers."
+        raise error(msg)
+
+    if len(ints) == 1:
+        ints.append(0)
+    elif len(ints) != 2:
+        raise error("Argument must be one or two integers.")
+
+    if ints[0] < 0 or ints[1] < 0:
+        raise error(
+            "Int arguments must be nonnegative. "
+            "If you wish to disable censoring, set the value to 0."
+        )
+
+    return ints
+
+
+def _one_or_two_floats(string, is_parser=True):
+    """Check if argument is one or two floats >= 0."""
+    error = argparse.ArgumentTypeError if is_parser else ValueError
+    try:
+        floats = [float(i) for i in string.split()]
+    except ValueError:
+        msg = "Argument must be one or two floats."
+        raise error(msg)
+
+    if len(floats) == 1:
+        floats.append(0.0)
+    elif len(floats) != 2:
+        raise error("Argument must be one or two floats.")
+
+    if floats[0] < 0 or floats[1] < 0:
+        raise error(
+            "Float arguments must be nonnegative. "
+            "If you wish to disable censoring, set the value to 0."
+        )
+
+    return floats
+
+
 def _restricted_float(x):
     """From https://stackoverflow.com/a/12117065/2589328."""
     try:
@@ -101,3 +147,25 @@ def __call__(self, parser, namespace, values, option_string=None):  # noqa: U100
             f"{self.__version__}. "
         )
         setattr(namespace, self.dest, values)
+
+
+def _check_censoring_thresholds(values, parser, arg_name):
+    if len(values) == 1:
+        values = [values[0], values[0]]
+    elif len(values) > 2:
+        parser.error(
+            f"Invalid number of values for '{arg_name}': {len(values)}. "
+            "Please provide either one or two values."
+        )
+    return values
+
+
+def _check_censoring_numbers(values, parser, arg_name):
+    if len(values) == 1:
+        values = [values[0], 0]
+    elif len(values) > 2:
+        parser.error(
+            f"Invalid number of values for '{arg_name}': {len(values)}. "
+            "Please provide either one or two values."
+        )
+    return values

xcp_d/config.py

@@ -544,8 +544,16 @@ class workflow(_Config):
     """Order of the filter to apply to the motion regressors."""
     head_radius = 50
     """Radius of the head in mm."""
-    fd_thresh = 0.3
+    fd_thresh = [0.3, 0.0]
     """Framewise displacement threshold for censoring."""
+    dvars_thresh = [0.0, 0.0]
+    """DVARS threshold for censoring."""
+    censor_before = [0, 0]
+    """Number of volumes to remove before FD or DVARS outliers."""
+    censor_after = [0, 0]
+    """Number of volumes to remove after FD or DVARS outliers."""
+    censor_between = [0, 0]
+    """Number of volumes to remove between any FD or DVARS outliers."""
     min_time = 240
     """Post-scrubbing threshold to apply to individual runs in the dataset."""
     bandpass_filter = True