[ENH]: Allow BOLDWarper to warp between MNI spaces #302

Merged
synchon merged 10 commits from update/multi-mni-bold-warper into main 2024-02-21 14:48:40 +00:00
3 changed files with 203 additions and 45 deletions

View file

@ -0,0 +1 @@
Allow :class:`.BOLDWarper` to warp BOLD data to other MNI spaces by `Synchon Mandal`_

View file

@ -13,8 +13,13 @@ from typing import (
Union, Union,
) )
import nibabel as nib
from templateflow import api as tflow
from ..api.decorators import register_preprocessor from ..api.decorators import register_preprocessor
from ..utils import logger, raise_error from ..data import get_template, get_xfm
from ..pipeline import WorkDirManager
from ..utils import logger, raise_error, run_ext_cmd
from .ants.ants_apply_transforms_warper import _AntsApplyTransformsWarper from .ants.ants_apply_transforms_warper import _AntsApplyTransformsWarper
from .base import BasePreprocessor from .base import BasePreprocessor
from .fsl.apply_warper import _ApplyWarper from .fsl.apply_warper import _ApplyWarper
@ -27,7 +32,13 @@ class BOLDWarper(BasePreprocessor):
Parameters Parameters
---------- ----------
reference : str reference : str
The data type to use as reference for warping. The data type to use as reference for warping, can be either a data
type like "T1w" or a template space like "MNI152NLin2009cAsym".
Raises
------
ValueError
If ``reference`` is invalid.
""" """
@ -49,9 +60,15 @@ class BOLDWarper(BasePreprocessor):
def __init__(self, reference: str) -> None: def __init__(self, reference: str) -> None:
"""Initialize the class.""" """Initialize the class."""
self.ref = reference self.ref = reference
# Initialize superclass based on reference
if self.ref == "T1w":
super().__init__( super().__init__(
on="BOLD", required_data_types=["BOLD", self.ref, "Warp"] on="BOLD", required_data_types=["BOLD", self.ref, "Warp"]
) )
elif self.ref in tflow.templates():
super().__init__(on="BOLD", required_data_types=["BOLD"])
else:
raise_error(f"Unknown reference: {self.ref}")
def get_valid_inputs(self) -> List[str]: def get_valid_inputs(self) -> List[str]:
"""Get valid data types for input. """Get valid data types for input.
@ -97,7 +114,8 @@ class BOLDWarper(BasePreprocessor):
The BOLD input from the Junifer Data object. The BOLD input from the Junifer Data object.
extra_input : dict, optional extra_input : dict, optional
The other fields in the Junifer Data object. Must include the The other fields in the Junifer Data object. Must include the
``Warp`` and ``ref`` value's keys. ``Warp`` and ``ref`` value's keys if native space transformation is
needed.
Returns Returns
------- -------
@ -110,17 +128,23 @@ class BOLDWarper(BasePreprocessor):
Raises Raises
------ ------
ValueError ValueError
If ``extra_input`` is None. If ``extra_input`` is None when transforming to native space
i.e., using "T1w" as reference.
RuntimeError RuntimeError
If warp / transformation file extension is not ".mat" or ".h5". If warp / transformation file extension is not ".mat" or ".h5"
when transforming to native space or
if the BOLD data is in the correct space and does not require
warping.
fraimondo commented 2024-02-14 15:01:56 +00:00 (Migrated from github.com)

This should be a warning if not an error. Why is the BOLDWarper set if you don't need to warp? I prefer an error though.

This should be a warning if not an error. Why is the BOLDWarper set if you don't need to warp? I prefer an error though.
synchon commented 2024-02-14 15:09:04 +00:00 (Migrated from github.com)

I don't understand why it should fail with an error, can make it a warning though. But I prefer a logged message instead of a warning in these cases.

I don't understand why it should fail with an error, can make it a warning though. But I prefer a logged message instead of a warning in these cases.
""" """
logger.debug("Warping BOLD using BOLDWarper") logger.info(f"Warping BOLD to {self.ref} space using BOLDWarper")
# Transform to native space
if self.ref == "T1w":
# Check for extra inputs # Check for extra inputs
if extra_input is None: if extra_input is None:
raise_error( raise_error(
f"No extra input provided, requires `Warp` and `{self.ref}` " "No extra input provided, requires `Warp` and "
"data types in particular." f"`{self.ref}` data types in particular."
) )
# Check for warp file type to use correct tool # Check for warp file type to use correct tool
warp_file_ext = extra_input["Warp"]["path"].suffix warp_file_ext = extra_input["Warp"]["path"].suffix
@ -152,4 +176,69 @@ class BOLDWarper(BasePreprocessor):
), ),
klass=RuntimeError, klass=RuntimeError,
) )
# Transform to template space
else:
# Check pre-requirements for space manipulation
if self.ref == input["space"]:
raise_error(
(
f"Skipped warping as the BOLD data is in {self.ref} "
"space which would mean that you can remove the "
"BOLDWarper from the preprocess step."
),
klass=RuntimeError,
)
else:
# Get xfm file
xfm_file_path = get_xfm(src=input["space"], dst=self.ref)
# Get template space image
template_space_img = get_template(
space=self.ref,
target_data=input,
extra_input=None,
)
# Create component-scoped tempdir
tempdir = WorkDirManager().get_tempdir(prefix="boldwarper")
# Create element-scoped tempdir so that warped BOLD is
# available later as nibabel stores file path reference for
# loading on computation
element_tempdir = WorkDirManager().get_element_tempdir(
prefix="boldwarper"
)
# Save template
template_space_img_path = tempdir / f"{self.ref}_T1w.nii.gz"
nib.save(template_space_img, template_space_img_path)
# Create a tempfile for warped output
warped_bold_path = (
element_tempdir
/ f"bold_warped_from_{input['space']}_to_{self.ref}.nii.gz"
)
logger.debug(
f"Using ANTs to warp BOLD "
f"from {input['space']} to {self.ref}"
)
# Set antsApplyTransforms command
apply_transforms_cmd = [
"antsApplyTransforms",
"-d 3",
"-e 3",
"-n LanczosWindowedSinc",
f"-i {input['path'].resolve()}",
f"-r {template_space_img_path.resolve()}",
f"-t {xfm_file_path.resolve()}",
f"-o {warped_bold_path.resolve()}",
]
# Call antsApplyTransforms
run_ext_cmd(
name="antsApplyTransforms", cmd=apply_transforms_cmd
)
# Modify target data
input["data"] = nib.load(warped_bold_path)
input["space"] = self.ref
return "BOLD", input return "BOLD", input

View file

@ -7,6 +7,7 @@ import socket
from typing import TYPE_CHECKING, List, Tuple from typing import TYPE_CHECKING, List, Tuple
import pytest import pytest
from numpy.testing import assert_array_equal, assert_raises
from junifer.datagrabber import DataladHCP1200, DMCC13Benchmark from junifer.datagrabber import DataladHCP1200, DMCC13Benchmark
from junifer.datareader import DefaultDataReader from junifer.datareader import DefaultDataReader
@ -83,7 +84,7 @@ def test_BOLDWarper_get_output_type(input_: List[str]) -> None:
socket.gethostname() != "juseless", socket.gethostname() != "juseless",
reason="only for juseless", reason="only for juseless",
) )
def test_BOLDWarper_preprocess( def test_BOLDWarper_preprocess_to_native(
datagrabber: "BaseDataGrabber", element: Tuple[str, ...] datagrabber: "BaseDataGrabber", element: Tuple[str, ...]
) -> None: ) -> None:
"""Test BOLDWarper preprocess. """Test BOLDWarper preprocess.
@ -106,3 +107,70 @@ def test_BOLDWarper_preprocess(
) )
assert data_type == "BOLD" assert data_type == "BOLD"
assert isinstance(data, dict) assert isinstance(data, dict)
@pytest.mark.parametrize(
"datagrabber, element, space",
[
[
DMCC13Benchmark(
types=["BOLD"],
sessions=["wave1bas"],
tasks=["Rest"],
phase_encodings=["AP"],
runs=["1"],
native_t1w=False,
),
("f9057kp", "wave1bas", "Rest", "AP", "1"),
"MNI152NLin2009aAsym",
],
[
DMCC13Benchmark(
types=["BOLD"],
sessions=["wave1bas"],
tasks=["Rest"],
phase_encodings=["AP"],
runs=["1"],
native_t1w=False,
),
("f9057kp", "wave1bas", "Rest", "AP", "1"),
"MNI152NLin6Asym",
],
],
)
@pytest.mark.skipif(
_check_ants() is False, reason="requires ANTs to be in PATH"
)
@pytest.mark.skipif(
socket.gethostname() != "juseless",
reason="only for juseless",
)
def test_BOLDWarper_preprocess_to_multi_mni(
datagrabber: "BaseDataGrabber", element: Tuple[str, ...], space: str
) -> None:
"""Test BOLDWarper preprocess.
Parameters
----------
datagrabber : DataGrabber-like object
The parametrized DataGrabber objects.
element : tuple of str
The parametrized elements.
space : str
The parametrized template space to transform to.
"""
with datagrabber as dg:
# Read data
element_data = DefaultDataReader().fit_transform(dg[element])
pre_xfm_data = element_data["BOLD"]["data"].get_fdata().copy()
# Preprocess data
data_type, data = BOLDWarper(reference=space).preprocess(
input=element_data["BOLD"],
extra_input=element_data,
)
assert data_type == "BOLD"
assert isinstance(data, dict)
assert data["space"] == space
with assert_raises(AssertionError):
assert_array_equal(pre_xfm_data, data["data"])