diff --git a/docs/changes/newsfragments/161.feature b/docs/changes/newsfragments/161.feature new file mode 100644 index 0000000000..6c1e3742e5 --- /dev/null +++ b/docs/changes/newsfragments/161.feature @@ -0,0 +1 @@ +Introduce :class:`.Smoothing` for smoothing / blurring images as a preprocessing step by `Synchon Mandal`_ diff --git a/junifer/preprocess/__init__.py b/junifer/preprocess/__init__.py index 15042f0eec..d60c6efddb 100644 --- a/junifer/preprocess/__init__.py +++ b/junifer/preprocess/__init__.py @@ -9,3 +9,4 @@ from .confounds import fMRIPrepConfoundRemover from .bold_warper import BOLDWarper from .warping import SpaceWarper +from .smoothing import Smoothing diff --git a/junifer/preprocess/smoothing/__init__.py b/junifer/preprocess/smoothing/__init__.py new file mode 100644 index 0000000000..0236a63d1e --- /dev/null +++ b/junifer/preprocess/smoothing/__init__.py @@ -0,0 +1,6 @@ +"""Provide imports for smoothing sub-package.""" + +# Authors: Synchon Mandal +# License: AGPL + +from .smoothing import Smoothing diff --git a/junifer/preprocess/smoothing/_afni_smoothing.py b/junifer/preprocess/smoothing/_afni_smoothing.py new file mode 100644 index 0000000000..50fd27161b --- /dev/null +++ b/junifer/preprocess/smoothing/_afni_smoothing.py @@ -0,0 +1,119 @@ +"""Provide class for smoothing via AFNI.""" + +# Authors: Synchon Mandal +# License: AGPL + +from typing import ( + TYPE_CHECKING, + ClassVar, + Dict, + List, + Set, + Union, +) + +import nibabel as nib + +from ...pipeline import WorkDirManager +from ...utils import logger, run_ext_cmd + + +if TYPE_CHECKING: + from nibabel import Nifti1Image + + +__all__ = ["AFNISmoothing"] + + +class AFNISmoothing: + """Class for smoothing via AFNI. + + This class uses AFNI's 3dBlurToFWHM. + + """ + + _EXT_DEPENDENCIES: ClassVar[List[Dict[str, Union[str, List[str]]]]] = [ + { + "name": "afni", + "commands": ["3dBlurToFWHM"], + }, + ] + + _DEPENDENCIES: ClassVar[Set[str]] = {"nibabel"} + + def preprocess( + self, + data: "Nifti1Image", + fwhm: Union[int, float], + ) -> "Nifti1Image": + """Preprocess using AFNI. + + Parameters + ---------- + data : Niimg-like object + Image(s) to preprocess. + fwhm : int or float + Smooth until the value. AFNI estimates the smoothing and then + applies smoothing to reach ``fwhm``. + + Returns + ------- + Niimg-like object + The preprocessed image(s). + + Notes + ----- + For more information on ``3dBlurToFWHM``, check: + https://afni.nimh.nih.gov/pub/dist/doc/program_help/3dBlurToFWHM.html + + As the process also depends on the conversion of AFNI files to NIfTI + via AFNI's ``3dAFNItoNIFTI``, the help for that can be found at: + https://afni.nimh.nih.gov/pub/dist/doc/program_help/3dAFNItoNIFTI.html + + """ + logger.info("Smoothing using AFNI") + + # Create component-scoped tempdir + tempdir = WorkDirManager().get_tempdir(prefix="afni_smoothing") + + # Save target data to a component-scoped tempfile + nifti_in_file_path = tempdir / "input.nii" # needs to be .nii + nib.save(data, nifti_in_file_path) + + # Set 3dBlurToFWHM command + blur_out_path_prefix = tempdir / "blur" + blur_cmd = [ + "3dBlurToFWHM", + f"-input {nifti_in_file_path.resolve()}", + f"-prefix {blur_out_path_prefix.resolve()}", + "-automask", + f"-FWHM {fwhm}", + ] + # Call 3dBlurToFWHM + run_ext_cmd(name="3dBlurToFWHM", cmd=blur_cmd) + + # Create element-scoped tempdir so that the blurred output is + # available later as nibabel stores file path reference for + # loading on computation + element_tempdir = WorkDirManager().get_element_tempdir( + prefix="afni_blur" + ) + # Convert afni to nifti + blur_afni_to_nifti_out_path = ( + element_tempdir / "output.nii" # needs to be .nii + ) + convert_cmd = [ + "3dAFNItoNIFTI", + f"-prefix {blur_afni_to_nifti_out_path.resolve()}", + f"{blur_out_path_prefix}+tlrc.BRIK", + ] + # Call 3dAFNItoNIFTI + run_ext_cmd(name="3dAFNItoNIFTI", cmd=convert_cmd) + + # Load nifti + output_data = nib.load(blur_afni_to_nifti_out_path) + + # Delete tempdir + WorkDirManager().delete_tempdir(tempdir) + + return output_data # type: ignore diff --git a/junifer/preprocess/smoothing/_fsl_smoothing.py b/junifer/preprocess/smoothing/_fsl_smoothing.py new file mode 100644 index 0000000000..2e02ea2226 --- /dev/null +++ b/junifer/preprocess/smoothing/_fsl_smoothing.py @@ -0,0 +1,116 @@ +"""Provide class for smoothing via FSL.""" + +# Authors: Synchon Mandal +# License: AGPL + +from typing import ( + TYPE_CHECKING, + ClassVar, + Dict, + List, + Set, + Union, +) + +import nibabel as nib + +from ...pipeline import WorkDirManager +from ...utils import logger, run_ext_cmd + + +if TYPE_CHECKING: + from nibabel import Nifti1Image + + +__all__ = ["FSLSmoothing"] + + +class FSLSmoothing: + """Class for smoothing via FSL. + + This class uses FSL's susan. + + """ + + _EXT_DEPENDENCIES: ClassVar[List[Dict[str, Union[str, List[str]]]]] = [ + { + "name": "fsl", + "commands": ["susan"], + }, + ] + + _DEPENDENCIES: ClassVar[Set[str]] = {"nibabel"} + + def preprocess( + self, + data: "Nifti1Image", + brightness_threshold: float, + fwhm: float, + ) -> "Nifti1Image": + """Preprocess using FSL. + + Parameters + ---------- + data : Niimg-like object + Image(s) to preprocess. + brightness_threshold : float + Threshold to discriminate between noise and the underlying image. + The value should be set greater than the noise level and less than + the contrast of the underlying image. + fwhm : float + Spatial extent of smoothing. + + Returns + ------- + Niimg-like object + The preprocessed image(s). + + Notes + ----- + For more information on ``SUSAN``, check [1]_ + + References + ---------- + .. [1] Smith, S.M. and Brady, J.M. (1997). + SUSAN - a new approach to low level image processing. + International Journal of Computer Vision, Volume 23(1), + Pages 45-78. + + """ + logger.info("Smoothing using FSL") + + # Create component-scoped tempdir + tempdir = WorkDirManager().get_tempdir(prefix="fsl_smoothing") + + # Save target data to a component-scoped tempfile + nifti_in_file_path = tempdir / "input.nii.gz" + nib.save(data, nifti_in_file_path) + + # Create element-scoped tempdir so that the output is + # available later as nibabel stores file path reference for + # loading on computation + element_tempdir = WorkDirManager().get_element_tempdir( + prefix="fsl_susan" + ) + susan_out_path = element_tempdir / "output.nii.gz" + # Set susan command + susan_cmd = [ + "susan", + f"{nifti_in_file_path.resolve()}", + f"{brightness_threshold}", + f"{fwhm}", + "3", # dimension + "1", # use median when no neighbourhood is found + "0", # use input image to find USAN + f"{susan_out_path.resolve()}", + ] + # Call susan + run_ext_cmd(name="susan", cmd=susan_cmd) + + # Load nifti + output_data = nib.load(susan_out_path) + + # Delete tempdir + WorkDirManager().delete_tempdir(tempdir) + + return output_data # type: ignore diff --git a/junifer/preprocess/smoothing/_nilearn_smoothing.py b/junifer/preprocess/smoothing/_nilearn_smoothing.py new file mode 100644 index 0000000000..1ef0498ee0 --- /dev/null +++ b/junifer/preprocess/smoothing/_nilearn_smoothing.py @@ -0,0 +1,69 @@ +"""Provide class for smoothing via nilearn.""" + +# Authors: Synchon Mandal +# License: AGPL + +from typing import ( + TYPE_CHECKING, + ClassVar, + Literal, + Set, + Union, +) + +from nilearn import image as nimg +from numpy.typing import ArrayLike + +from ...utils import logger + + +if TYPE_CHECKING: + from nibabel import Nifti1Image + + +__all__ = ["NilearnSmoothing"] + + +class NilearnSmoothing: + """Class for smoothing via nilearn. + + This class uses :func:`nilearn.image.smooth_img` to smooth image(s). + + """ + + _DEPENDENCIES: ClassVar[Set[str]] = {"nilearn"} + + def preprocess( + self, + data: "Nifti1Image", + fwhm: Union[int, float, ArrayLike, Literal["fast"], None], + ) -> "Nifti1Image": + """Preprocess using nilearn. + + Parameters + ---------- + data : Niimg-like object + Image(s) to preprocess. + fwhm : scalar, ``numpy.ndarray``, tuple or list of scalar, "fast" or \ + None + Smoothing strength, as a full-width at half maximum, in + millimeters: + + * If nonzero scalar, width is identical in all 3 directions. + * If ``numpy.ndarray``, tuple, or list, it must have 3 elements, + giving the FWHM along each axis. If any of the elements is 0 or + None, smoothing is not performed along that axis. + * If ``"fast"``, a fast smoothing will be performed with a filter + ``[0.2, 1, 0.2]`` in each direction and a normalisation to + preserve the local average value. + * If None, no filtering is performed (useful when just removal of + non-finite values is needed). + + Returns + ------- + Niimg-like object + The preprocessed image(s). + + """ + logger.info("Smoothing using nilearn") + return nimg.smooth_img(imgs=data, fwhm=fwhm) # type: ignore diff --git a/junifer/preprocess/smoothing/smoothing.py b/junifer/preprocess/smoothing/smoothing.py new file mode 100644 index 0000000000..16f71a3943 --- /dev/null +++ b/junifer/preprocess/smoothing/smoothing.py @@ -0,0 +1,174 @@ +"""Provide class for smoothing.""" + +# Authors: Synchon Mandal +# License: AGPL + +from typing import Any, ClassVar, Dict, List, Optional, Tuple, Type, Union + +from ...api.decorators import register_preprocessor +from ...utils import logger, raise_error +from ..base import BasePreprocessor +from ._afni_smoothing import AFNISmoothing +from ._fsl_smoothing import FSLSmoothing +from ._nilearn_smoothing import NilearnSmoothing + + +__all__ = ["Smoothing"] + + +@register_preprocessor +class Smoothing(BasePreprocessor): + """Class for smoothing. + + Parameters + ---------- + using : {"nilearn", "afni", "fsl"} + Implementation to use for smoothing: + + * "nilearn" : Use :func:`nilearn.image.smooth_img` + * "afni" : Use AFNI's ``3dBlurToFWHM`` + * "fsl" : Use FSL SUSAN's ``susan`` + + on : {"T1w", "T2w", "BOLD"} or list of the options + The data type to apply smoothing to. + smoothing_params : dict, optional + Extra parameters for smoothing as a dictionary (default None). + If ``using="nilearn"``, then the valid keys are: + + * ``fmhw`` : scalar, ``numpy.ndarray``, tuple or list of scalar, \ + "fast" or None + Smoothing strength, as a full-width at half maximum, in + millimeters: + + - If nonzero scalar, width is identical in all 3 directions. + - If ``numpy.ndarray``, tuple, or list, it must have 3 elements, + giving the FWHM along each axis. If any of the elements is 0 or + None, smoothing is not performed along that axis. + - If ``"fast"``, a fast smoothing will be performed with a filter + ``[0.2, 1, 0.2]`` in each direction and a normalisation to + preserve the local average value. + - If None, no filtering is performed (useful when just removal of + non-finite values is needed). + + else if ``using="afni"``, then the valid keys are: + + * ``fwhm`` : int or float + Smooth until the value. AFNI estimates the smoothing and then + applies smoothing to reach ``fwhm``. + + else if ``using="fsl"``, then the valid keys are: + + * ``brightness_threshold`` : float + Threshold to discriminate between noise and the underlying image. + The value should be set greater than the noise level and less than + the contrast of the underlying image. + * ``fwhm`` : float + Spatial extent of smoothing. + + """ + + _CONDITIONAL_DEPENDENCIES: ClassVar[List[Dict[str, Union[str, Type]]]] = [ + { + "using": "nilearn", + "depends_on": NilearnSmoothing, + }, + { + "using": "afni", + "depends_on": AFNISmoothing, + }, + { + "using": "fsl", + "depends_on": FSLSmoothing, + }, + ] + + def __init__( + self, + using: str, + on: Union[List[str], str], + smoothing_params: Optional[Dict] = None, + ) -> None: + """Initialize the class.""" + # Validate `using` parameter + valid_using = [dep["using"] for dep in self._CONDITIONAL_DEPENDENCIES] + if using not in valid_using: + raise_error( + f"Invalid value for `using`, should be one of: {valid_using}" + ) + self.using = using + self.smoothing_params = ( + smoothing_params if smoothing_params is not None else {} + ) + super().__init__(on=on) + + def get_valid_inputs(self) -> List[str]: + """Get valid data types for input. + + Returns + ------- + list of str + The list of data types that can be used as input for this + preprocessor. + + """ + return ["T1w", "T2w", "BOLD"] + + def get_output_type(self, input_type: str) -> str: + """Get output type. + + Parameters + ---------- + input_type : str + The data type input to the preprocessor. + + Returns + ------- + str + The data type output by the preprocessor. + + """ + # Does not add any new keys + return input_type + + def preprocess( + self, + input: Dict[str, Any], + extra_input: Optional[Dict[str, Any]] = None, + ) -> Tuple[Dict[str, Any], Optional[Dict[str, Dict[str, Any]]]]: + """Preprocess. + + Parameters + ---------- + input : dict + The input from the Junifer Data object. + extra_input : dict, optional + The other fields in the Junifer Data object. + + Returns + ------- + dict + The computed result as dictionary. + None + Extra "helper" data types as dictionary to add to the Junifer Data + object. + + """ + logger.debug("Smoothing") + + # Conditional preprocessor + if self.using == "nilearn": + preprocessor = NilearnSmoothing() + elif self.using == "afni": + preprocessor = AFNISmoothing() + elif self.using == "fsl": + preprocessor = FSLSmoothing() + # Smooth + output = preprocessor.preprocess( # type: ignore + data=input["data"], + **self.smoothing_params, + ) + + # Modify target data + input["data"] = output + + return input, None diff --git a/junifer/preprocess/smoothing/tests/test_smoothing.py b/junifer/preprocess/smoothing/tests/test_smoothing.py new file mode 100644 index 0000000000..6fe66c55b2 --- /dev/null +++ b/junifer/preprocess/smoothing/tests/test_smoothing.py @@ -0,0 +1,94 @@ +"""Provide tests for Smoothing.""" + +# Authors: Synchon Mandal +# License: AGPL + + +import pytest + +from junifer.datareader import DefaultDataReader +from junifer.pipeline.utils import _check_afni, _check_fsl +from junifer.preprocess import Smoothing +from junifer.testing.datagrabbers import SPMAuditoryTestingDataGrabber + + +@pytest.mark.parametrize( + "data_type", + ["T1w", "BOLD"], +) +def test_Smoothing_nilearn(data_type: str) -> None: + """Test Smoothing using nilearn. + + Parameters + ---------- + data_type : str + The parametrized data type. + + """ + with SPMAuditoryTestingDataGrabber() as dg: + # Read data + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + # Preprocess data + output = Smoothing( + using="nilearn", + on=data_type, + smoothing_params={"fwhm": "fast"}, + ).fit_transform(element_data) + + assert isinstance(output, dict) + + +@pytest.mark.parametrize( + "data_type", + ["T1w", "BOLD"], +) +@pytest.mark.skipif( + _check_afni() is False, reason="requires AFNI to be in PATH" +) +def test_Smoothing_afni(data_type: str) -> None: + """Test Smoothing using AFNI. + + Parameters + ---------- + data_type : str + The parametrized data type. + + """ + with SPMAuditoryTestingDataGrabber() as dg: + # Read data + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + # Preprocess data + output = Smoothing( + using="afni", + on=data_type, + smoothing_params={"fwhm": 3}, + ).fit_transform(element_data) + + assert isinstance(output, dict) + + +@pytest.mark.parametrize( + "data_type", + ["T1w", "BOLD"], +) +@pytest.mark.skipif(_check_fsl() is False, reason="requires FSL to be in PATH") +def test_Smoothing_fsl(data_type: str) -> None: + """Test Smoothing using FSL. + + Parameters + ---------- + data_type : str + The parametrized data type. + + """ + with SPMAuditoryTestingDataGrabber() as dg: + # Read data + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + # Preprocess data + output = Smoothing( + using="fsl", + on=data_type, + smoothing_params={"brightness_threshold": 10.0, "fwhm": 3.0}, + ).fit_transform(element_data) + + assert isinstance(output, dict)