From 963e7927b0608055458b70b86eb72b55305a70ad Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 1 Feb 2024 07:53:26 +0100 Subject: [PATCH 01/47] update: add source space info to Warp data for datagrabbers --- junifer/datagrabber/aomic/id1000.py | 3 +++ junifer/datagrabber/aomic/piop1.py | 3 +++ junifer/datagrabber/aomic/piop2.py | 3 +++ junifer/datagrabber/dmcc13_benchmark.py | 3 +++ junifer/datagrabber/hcp1200/hcp1200.py | 3 +++ 5 files changed, 15 insertions(+) diff --git a/junifer/datagrabber/aomic/id1000.py b/junifer/datagrabber/aomic/id1000.py index e8b4af1cb..0458751a2 100644 --- a/junifer/datagrabber/aomic/id1000.py +++ b/junifer/datagrabber/aomic/id1000.py @@ -154,4 +154,7 @@ class DataladAOMICID1000(PatternDataladDataGrabber): out["T1w"].update({"space": "native"}) else: out["T1w"].update({"space": "MNI152NLin2009cAsym"}) + if out.get("Warp"): + # Add source space information + out["Warp"].update({"src": "MNI152NLin2009cAsym"}) return out diff --git a/junifer/datagrabber/aomic/piop1.py b/junifer/datagrabber/aomic/piop1.py index 87190f466..e65485874 100644 --- a/junifer/datagrabber/aomic/piop1.py +++ b/junifer/datagrabber/aomic/piop1.py @@ -204,6 +204,9 @@ class DataladAOMICPIOP1(PatternDataladDataGrabber): out["T1w"].update({"space": "native"}) else: out["T1w"].update({"space": "MNI152NLin2009cAsym"}) + if out.get("Warp"): + # Add source space information + out["Warp"].update({"src": "MNI152NLin2009cAsym"}) return out def get_elements(self) -> List: diff --git a/junifer/datagrabber/aomic/piop2.py b/junifer/datagrabber/aomic/piop2.py index 9ec991ca8..e5ce18899 100644 --- a/junifer/datagrabber/aomic/piop2.py +++ b/junifer/datagrabber/aomic/piop2.py @@ -204,4 +204,7 @@ class DataladAOMICPIOP2(PatternDataladDataGrabber): out["T1w"].update({"space": "native"}) else: out["T1w"].update({"space": "MNI152NLin2009cAsym"}) + if out.get("Warp"): + # Add source space information + out["Warp"].update({"src": "MNI152NLin2009cAsym"}) return out diff --git a/junifer/datagrabber/dmcc13_benchmark.py b/junifer/datagrabber/dmcc13_benchmark.py index 8ab710afe..1a5c3d29e 100644 --- a/junifer/datagrabber/dmcc13_benchmark.py +++ b/junifer/datagrabber/dmcc13_benchmark.py @@ -271,6 +271,9 @@ class DMCC13Benchmark(PatternDataladDataGrabber): out["T1w"].update({"space": "native"}) else: out["T1w"].update({"space": "MNI152NLin2009cAsym"}) + if out.get("Warp"): + # Add source space information + out["Warp"].update({"src": "MNI152NLin2009cAsym"}) return out def get_elements(self) -> List: diff --git a/junifer/datagrabber/hcp1200/hcp1200.py b/junifer/datagrabber/hcp1200/hcp1200.py index 7c1bc7f57..9ca9e48f1 100644 --- a/junifer/datagrabber/hcp1200/hcp1200.py +++ b/junifer/datagrabber/hcp1200/hcp1200.py @@ -159,6 +159,9 @@ class HCP1200(PatternDataGrabber): # Add space for T1w data type if "T1w" in out: out["T1w"].update({"space": "native"}) + # Add source space for Warp data type + if "Warp" in out: + out["Warp"].update({"src": "MNI152NLin6Asym"}) return out def get_elements(self) -> List: -- 2.52.0 From 8e2b9339a67bfdca8914505cd51222a880ba64e7 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 1 Feb 2024 07:54:48 +0100 Subject: [PATCH 02/47] feat: add support for parcellation space transformation in get_parcellation() --- junifer/data/parcellations.py | 96 ++++++++++++++++++++++++++++------- 1 file changed, 78 insertions(+), 18 deletions(-) diff --git a/junifer/data/parcellations.py b/junifer/data/parcellations.py index 9c2b5d62b..272a535e9 100644 --- a/junifer/data/parcellations.py +++ b/junifer/data/parcellations.py @@ -22,6 +22,7 @@ from nilearn import datasets, image from ..pipeline import WorkDirManager from ..utils import logger, raise_error, run_ext_cmd, warn_with_log +from .template_spaces import get_template, get_xfm from .utils import closest_resolution @@ -242,10 +243,29 @@ def get_parcellation( If ``extra_input`` is None when ``target_data``'s space is native. """ + # Check pre-requirements for space manipulation + target_space = target_data["space"] + # Set target standard space to target space + target_std_space = target_space + # Extra data type requirement check if target space is native + if target_space == "native": + # Check for extra inputs + if extra_input is None: + raise_error( + "No extra input provided, requires `Warp` and `T1w` " + "data types in particular for transformation to " + f"{target_data['space']} space for further computation." + ) + # Set target standard space to warp file space source + target_std_space = extra_input["Warp"]["src"] + # Get the min of the voxels sizes and use it as the resolution target_img = target_data["data"] resolution = np.min(target_img.header.get_zooms()[:3]) + # Create tempdir + tempdir = WorkDirManager().get_tempdir(prefix="parcellations") + # Load the parcellations all_parcellations = [] all_labels = [] @@ -255,14 +275,65 @@ def get_parcellation( name=name, resolution=resolution, ) - # Resample all of them to the image - resampled_img = image.resample_to_img( + # Resample loaded parcellation to target image + img_to_merge = image.resample_to_img( source_img=img, target_img=target_img, interpolation="nearest", copy=True, ) - all_parcellations.append(resampled_img) + # Convert parcellation spaces if required + if space != target_std_space: + # Get xfm file + xfm_file_path = get_xfm(src=space, dst=target_std_space) + # Get target standard space template + target_std_space_template_img = get_template( + space=target_std_space, + target_data=target_data, + extra_input=extra_input, + ) + + # Save resampled parcellation image to a component-scoped tempfile + resampled_parcellation_path = ( + tempdir / f"{name}_resampled_to_{resolution}.nii.gz" + ) + nib.save(img_to_merge, resampled_parcellation_path) + + # Save template + target_std_space_template_path = ( + tempdir / f"{target_std_space}_T1w_{resolution}.nii.gz" + ) + nib.save( + target_std_space_template_img, target_std_space_template_path + ) + + # Set warped parcellation path + warped_resampled_parcellation_path = tempdir / ( + f"{name}_resampled_to_{resolution}_warped_to_" + f"{target_std_space}.nii.gz" + ) + + logger.debug( + f"Using ANTs to warp {name} (resolution={resolution}) " + f"from {space} to {target_std_space}" + ) + # Set antsApplyTransforms command + apply_transforms_cmd = [ + "antsApplyTransforms", + "-d 3", + "-e 2", + "-n 'GenericLabel[NearestNeighbor]'", + f"-i {resampled_parcellation_path.resolve()}", + f"-r {target_std_space_template_path.resolve()}", + f"-t {xfm_file_path}", + f"-o {warped_resampled_parcellation_path}", + ] + # Call antsApplyTransforms + run_ext_cmd(name="antsApplyTransforms", cmd=apply_transforms_cmd) + + img_to_merge = nib.load(warped_resampled_parcellation_path) + + all_parcellations.append(img_to_merge) all_labels.append(labels) all_spaces.append(space) @@ -285,19 +356,8 @@ def get_parcellation( klass=RuntimeError, ) - # Warp parcellation if target data is native - if target_data["space"] == "native": - # Check for extra inputs - if extra_input is None: - raise_error( - "No extra input provided, requires `Warp` and `T1w` " - "data types in particular for transformation to " - f"{target_data['space']} space for further computation." - ) - - # Create component-scoped tempdir - tempdir = WorkDirManager().get_tempdir(prefix="parcellations") - + # Warp parcellation if target space is native + if target_space == "native": # Save parcellation image to a component-scoped tempfile prewarp_parcellation_path = tempdir / "prewarp_parcellation.nii.gz" nib.save(resampled_parcellation_img, prewarp_parcellation_path) @@ -359,8 +419,8 @@ def get_parcellation( # Load nifti resampled_parcellation_img = nib.load(warped_parcellation_path) - # Delete tempdir - WorkDirManager().delete_tempdir(tempdir) + # Delete tempdir + WorkDirManager().delete_tempdir(tempdir) return resampled_parcellation_img, labels # type: ignore -- 2.52.0 From c3dfb08fb8b91380e9dbdc8b1c9048a5ce317029 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 1 Feb 2024 08:14:04 +0100 Subject: [PATCH 03/47] fix: remove same space check before merging in get_parcellation() --- junifer/data/parcellations.py | 24 +++++++----------------- 1 file changed, 7 insertions(+), 17 deletions(-) diff --git a/junifer/data/parcellations.py b/junifer/data/parcellations.py index 272a535e9..f46b24c5b 100644 --- a/junifer/data/parcellations.py +++ b/junifer/data/parcellations.py @@ -237,8 +237,7 @@ def get_parcellation( Raises ------ RuntimeError - If parcellations are in different spaces and they need to be merged or - if warp / transformation file extension is not ".mat" or ".h5". + If warp / transformation file extension is not ".mat" or ".h5". ValueError If ``extra_input`` is None when ``target_data``'s space is native. @@ -269,7 +268,6 @@ def get_parcellation( # Load the parcellations all_parcellations = [] all_labels = [] - all_spaces = [] for name in parcellation: img, labels, _, space = load_parcellation( name=name, @@ -335,26 +333,18 @@ def get_parcellation( all_parcellations.append(img_to_merge) all_labels.append(labels) - all_spaces.append(space) # Avoid merging if there is only one parcellation if len(all_parcellations) == 1: resampled_parcellation_img = all_parcellations[0] labels = all_labels[0] + # Parcellations are already transformed to target standard space else: - # Merge the parcellations only if all parcellations are in the same - # space - if len(set(all_spaces)) == 1: - resampled_parcellation_img, labels = merge_parcellations( - parcellations_list=all_parcellations, - parcellations_names=parcellation, - labels_lists=all_labels, - ) - else: - raise_error( - msg="Parcellations are in different spaces, unable to merge.", - klass=RuntimeError, - ) + resampled_parcellation_img, labels = merge_parcellations( + parcellations_list=all_parcellations, + parcellations_names=parcellation, + labels_lists=all_labels, + ) # Warp parcellation if target space is native if target_space == "native": -- 2.52.0 From b5aa0d0c28cba11b4e3795ed38a3bfc4b02db70a Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 1 Feb 2024 09:33:46 +0100 Subject: [PATCH 04/47] chore: improve commentary for get_parcellation() --- junifer/data/parcellations.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/junifer/data/parcellations.py b/junifer/data/parcellations.py index f46b24c5b..f19bbc029 100644 --- a/junifer/data/parcellations.py +++ b/junifer/data/parcellations.py @@ -262,7 +262,7 @@ def get_parcellation( target_img = target_data["data"] resolution = np.min(target_img.header.get_zooms()[:3]) - # Create tempdir + # Create component-scoped tempdir tempdir = WorkDirManager().get_tempdir(prefix="parcellations") # Load the parcellations -- 2.52.0 From 8631c7202466a16482fae4988f95920bf7cdd1a0 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 1 Feb 2024 09:35:49 +0100 Subject: [PATCH 05/47] feat: add support for mask space transformation in get_mask() --- junifer/data/masks.py | 134 ++++++++++++++++++++++++++++-------------- 1 file changed, 91 insertions(+), 43 deletions(-) diff --git a/junifer/data/masks.py b/junifer/data/masks.py index 908f4c8e6..c12fda253 100644 --- a/junifer/data/masks.py +++ b/junifer/data/masks.py @@ -30,6 +30,7 @@ from nilearn.masking import ( from ..pipeline import WorkDirManager from ..utils import logger, raise_error, run_ext_cmd +from .template_spaces import get_template, get_xfm from .utils import closest_resolution @@ -198,30 +199,43 @@ def get_mask( # noqa: C901 Raises ------ RuntimeError - If masks are in different spaces and they need to be intersected / - unionized or - if warp / transformation file extension is not ".mat" or ".h5". + If warp / transformation file extension is not ".mat" or ".h5". ValueError If extra key is provided in addition to mask name in ``masks`` or if no mask is provided or if ``masks = "inherit"`` but ``extra_input`` is None or ``mask_item`` is None or ``mask_items``'s value is not in ``extra_input`` or if callable parameters are passed to non-callable mask or - if multiple masks are provided and their spaces do not match or if parameters are passed to :func:`nilearn.masking.intersect_masks` when there is only one mask or if ``extra_input`` is None when ``target_data``'s space is native. """ + # Check pre-requirements for space manipulation + target_space = target_data["space"] + # Set target standard space to target space + target_std_space = target_space + # Extra data type requirement check if target space is native + if target_space == "native": + # Check for extra inputs + if extra_input is None: + raise_error( + "No extra input provided, requires `Warp` and `T1w` " + "data types in particular for transformation to " + f"{target_data['space']} space for further computation." + ) + # Set target standard space to warp file space source + target_std_space = extra_input["Warp"]["src"] + # Get the min of the voxels sizes and use it as the resolution target_img = target_data["data"] - inherited_mask_item = target_data.get("mask_item", None) resolution = np.min(target_img.header.get_zooms()[:3]) + # Convert masks to list if not already if not isinstance(masks, list): masks = [masks] - # Check that dicts have only one key + # Check that masks passed as dicts have only one key invalid_elements = [ x for x in masks if isinstance(x, dict) and len(x) != 1 ] @@ -248,9 +262,15 @@ def get_mask( # noqa: C901 if len(true_masks) == 0: raise_error("No mask was passed. At least one mask is required.") + + # Get the data type for the input data type's mask + inherited_mask_item = target_data.get("mask_item", None) + + # Create component-scoped tempdir + tempdir = WorkDirManager().get_tempdir(prefix="masks") + # Get all the masks all_masks = [] - all_spaces = [] for t_mask in true_masks: if isinstance(t_mask, dict): mask_name = next(iter(t_mask.keys())) @@ -281,7 +301,6 @@ def get_mask( # noqa: C901 f"because the item ({inherited_mask_item}) does not exist." ) mask_img = extra_input[inherited_mask_item]["data"] - mask_space = target_data["space"] # Starting with new mask else: # Load mask @@ -290,7 +309,7 @@ def get_mask( # noqa: C901 ) # Replace mask space with target space if mask's space is inherit if mask_space == "inherit": - mask_space = target_data["space"] + mask_space = target_std_space # If mask is callable like from nilearn if callable(mask_object): if mask_params is None: @@ -306,31 +325,72 @@ def get_mask( # noqa: C901 ) # Resample mask to target image mask_img = resample_to_img( - mask_object, - target_img, + source_img=mask_object, + target_img=target_img, interpolation="nearest", copy=True, ) - all_spaces.append(mask_space) + # Convert mask space if required + if mask_space != target_std_space: + # Get xfm file + xfm_file_path = get_xfm(src=mask_space, dst=target_std_space) + # Get target standard space template + target_std_space_template_img = get_template( + space=target_std_space, + target_data=target_data, + extra_input=extra_input, + ) + + # Save resampled mask image to a component-scoped tempfile + resampled_mask_path = ( + tempdir / f"{mask_name}_resampled_to_{resolution}.nii.gz" + ) + nib.save(mask_img, resampled_mask_path) + + # Save template + target_std_space_template_path = ( + tempdir / f"{target_std_space}_T1w_{resolution}.nii.gz" + ) + nib.save( + target_std_space_template_img, + target_std_space_template_path, + ) + + # Set warped mask path + warped_resampled_mask_path = tempdir / ( + f"{mask_name}_resampled_to_{resolution}_warped_to_" + f"{target_std_space}.nii.gz" + ) + + logger.debug( + f"Using ANTs to warp {mask_name} " + f"(resolution={resolution}) from {mask_space} to " + f"{target_std_space}" + ) + # Set antsApplyTransforms command + apply_transforms_cmd = [ + "antsApplyTransforms", + "-d 3", + "-e 2", + "-n 'GenericLabel[NearestNeighbor]'", + f"-i {resampled_mask_path.resolve()}", + f"-r {target_std_space_template_path.resolve()}", + f"-t {xfm_file_path}", + f"-o {warped_resampled_mask_path}", + ] + # Call antsApplyTransforms + run_ext_cmd( + name="antsApplyTransforms", cmd=apply_transforms_cmd + ) + + mask_img = nib.load(warped_resampled_mask_path) + all_masks.append(mask_img) # Multiple masks, need intersection / union if len(all_masks) > 1: - # Make a set of unique spaces - unique_spaces = set(all_spaces) - # Intersect / union of masks only if all masks are in the same space - if len(unique_spaces) == 1: - mask_img = intersect_masks(all_masks, **intersect_params) - # Store the mask space for further checks - mask_space = next(iter(unique_spaces)) - else: - raise_error( - msg=( - f"Masks are in different spaces: {unique_spaces}, " - "unable to merge." - ), - klass=RuntimeError, - ) + # Intersect / union of masks + mask_img = intersect_masks(all_masks, **intersect_params) # Single mask else: if len(intersect_params) > 0: @@ -340,21 +400,9 @@ def get_mask( # noqa: C901 "when there is only one mask." ) mask_img = all_masks[0] - mask_space = all_spaces[0] - - # Warp mask if target data is native and mask space is not native - if target_data["space"] == "native" and target_data["space"] != mask_space: - # Check for extra inputs - if extra_input is None: - raise_error( - "No extra input provided, requires `Warp` and `T1w` " - "data types in particular for transformation to " - f"{target_data['space']} space for further computation." - ) - - # Create component-scoped tempdir - tempdir = WorkDirManager().get_tempdir(prefix="masks") + # Warp mask if target data is native + if target_space == "native": # Save mask image to a component-scoped tempfile prewarp_mask_path = tempdir / "prewarp_mask.nii.gz" nib.save(mask_img, prewarp_mask_path) @@ -413,8 +461,8 @@ def get_mask( # noqa: C901 # Load nifti mask_img = nib.load(warped_mask_path) - # Delete tempdir - WorkDirManager().delete_tempdir(tempdir) + # Delete tempdir + WorkDirManager().delete_tempdir(tempdir) return mask_img # type: ignore -- 2.52.0 From 08ecafe2cb9490a7090ff6980382434a140b6c52 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 2 Feb 2024 07:46:56 +0100 Subject: [PATCH 06/47] update: store resampled + warped individual parcellations in element scoped tempfile --- junifer/data/parcellations.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/junifer/data/parcellations.py b/junifer/data/parcellations.py index f19bbc029..6f7c90c2d 100644 --- a/junifer/data/parcellations.py +++ b/junifer/data/parcellations.py @@ -264,6 +264,12 @@ def get_parcellation( # Create component-scoped tempdir tempdir = WorkDirManager().get_tempdir(prefix="parcellations") + # Create element-scoped tempdir so that warped parcellation is + # available later as nibabel stores file path reference for + # loading on computation + element_tempdir = WorkDirManager().get_element_tempdir( + prefix="parcellations" + ) # Load the parcellations all_parcellations = [] @@ -306,7 +312,7 @@ def get_parcellation( ) # Set warped parcellation path - warped_resampled_parcellation_path = tempdir / ( + warped_resampled_parcellation_path = element_tempdir / ( f"{name}_resampled_to_{resolution}_warped_to_" f"{target_std_space}.nii.gz" ) @@ -352,12 +358,6 @@ def get_parcellation( prewarp_parcellation_path = tempdir / "prewarp_parcellation.nii.gz" nib.save(resampled_parcellation_img, prewarp_parcellation_path) - # Create element-scoped tempdir so that warped parcellation is - # available later as nibabel stores file path reference for - # loading on computation - element_tempdir = WorkDirManager().get_element_tempdir( - prefix="parcellations" - ) # Create an element-scoped tempfile for warped output warped_parcellation_path = ( element_tempdir / "parcellation_warped.nii.gz" -- 2.52.0 From b718e358db6c12546a6c0c416a316de3b1b98e13 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 2 Feb 2024 07:47:58 +0100 Subject: [PATCH 07/47] update: store resampled + warped individual masks in element scoped tempfile --- junifer/data/masks.py | 11 +++++------ 1 file changed, 5 insertions(+), 6 deletions(-) diff --git a/junifer/data/masks.py b/junifer/data/masks.py index c12fda253..e3459fd42 100644 --- a/junifer/data/masks.py +++ b/junifer/data/masks.py @@ -268,6 +268,10 @@ def get_mask( # noqa: C901 # Create component-scoped tempdir tempdir = WorkDirManager().get_tempdir(prefix="masks") + # Create element-scoped tempdir so that warped mask is + # available later as nibabel stores file path reference for + # loading on computation + element_tempdir = WorkDirManager().get_element_tempdir(prefix="masks") # Get all the masks all_masks = [] @@ -357,7 +361,7 @@ def get_mask( # noqa: C901 ) # Set warped mask path - warped_resampled_mask_path = tempdir / ( + warped_resampled_mask_path = element_tempdir / ( f"{mask_name}_resampled_to_{resolution}_warped_to_" f"{target_std_space}.nii.gz" ) @@ -407,11 +411,6 @@ def get_mask( # noqa: C901 prewarp_mask_path = tempdir / "prewarp_mask.nii.gz" nib.save(mask_img, prewarp_mask_path) - # Create element-scoped tempdir so that warped mask is - # available later as nibabel stores file path reference for - # loading on computation - element_tempdir = WorkDirManager().get_element_tempdir(prefix="masks") - # Create an element-scoped tempfile for warped output warped_mask_path = element_tempdir / "mask_warped.nii.gz" -- 2.52.0 From 2c7f4fedf1c759bd863ae7a19495af81eb51fd65 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 2 Feb 2024 07:48:23 +0100 Subject: [PATCH 08/47] chore: improve path handling when running external command for get_mask and get_parcellation --- junifer/data/masks.py | 4 ++-- junifer/data/parcellations.py | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/junifer/data/masks.py b/junifer/data/masks.py index e3459fd42..ff1a7ea4d 100644 --- a/junifer/data/masks.py +++ b/junifer/data/masks.py @@ -379,8 +379,8 @@ def get_mask( # noqa: C901 "-n 'GenericLabel[NearestNeighbor]'", f"-i {resampled_mask_path.resolve()}", f"-r {target_std_space_template_path.resolve()}", - f"-t {xfm_file_path}", - f"-o {warped_resampled_mask_path}", + f"-t {xfm_file_path.resolve()}", + f"-o {warped_resampled_mask_path.resolve()}", ] # Call antsApplyTransforms run_ext_cmd( diff --git a/junifer/data/parcellations.py b/junifer/data/parcellations.py index 6f7c90c2d..8023dd513 100644 --- a/junifer/data/parcellations.py +++ b/junifer/data/parcellations.py @@ -329,8 +329,8 @@ def get_parcellation( "-n 'GenericLabel[NearestNeighbor]'", f"-i {resampled_parcellation_path.resolve()}", f"-r {target_std_space_template_path.resolve()}", - f"-t {xfm_file_path}", - f"-o {warped_resampled_parcellation_path}", + f"-t {xfm_file_path.resolve()}", + f"-o {warped_resampled_parcellation_path.resolve()}", ] # Call antsApplyTransforms run_ext_cmd(name="antsApplyTransforms", cmd=apply_transforms_cmd) -- 2.52.0 From 3fbdaaf4508deaef1199a198452252f7120f78c0 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 2 Feb 2024 09:14:39 +0100 Subject: [PATCH 09/47] fix: resample after warping in get_parcellation and correct cmd for warping --- junifer/data/parcellations.py | 39 +++++++++++++++++------------------ 1 file changed, 19 insertions(+), 20 deletions(-) diff --git a/junifer/data/parcellations.py b/junifer/data/parcellations.py index 8023dd513..605c270b3 100644 --- a/junifer/data/parcellations.py +++ b/junifer/data/parcellations.py @@ -279,13 +279,7 @@ def get_parcellation( name=name, resolution=resolution, ) - # Resample loaded parcellation to target image - img_to_merge = image.resample_to_img( - source_img=img, - target_img=target_img, - interpolation="nearest", - copy=True, - ) + # Convert parcellation spaces if required if space != target_std_space: # Get xfm file @@ -297,11 +291,9 @@ def get_parcellation( extra_input=extra_input, ) - # Save resampled parcellation image to a component-scoped tempfile - resampled_parcellation_path = ( - tempdir / f"{name}_resampled_to_{resolution}.nii.gz" - ) - nib.save(img_to_merge, resampled_parcellation_path) + # Save parcellation image to a component-scoped tempfile + parcellation_path = tempdir / f"{name}.nii.gz" + nib.save(img, parcellation_path) # Save template target_std_space_template_path = ( @@ -312,30 +304,37 @@ def get_parcellation( ) # Set warped parcellation path - warped_resampled_parcellation_path = element_tempdir / ( - f"{name}_resampled_to_{resolution}_warped_to_" - f"{target_std_space}.nii.gz" + warped_parcellation_path = element_tempdir / ( + f"{name}_warped_from_{space}_to_" f"{target_std_space}.nii.gz" ) logger.debug( - f"Using ANTs to warp {name} (resolution={resolution}) " + f"Using ANTs to warp {name} " f"from {space} to {target_std_space}" ) # Set antsApplyTransforms command apply_transforms_cmd = [ "antsApplyTransforms", "-d 3", - "-e 2", + "-e 3", "-n 'GenericLabel[NearestNeighbor]'", - f"-i {resampled_parcellation_path.resolve()}", + f"-i {parcellation_path.resolve()}", f"-r {target_std_space_template_path.resolve()}", f"-t {xfm_file_path.resolve()}", - f"-o {warped_resampled_parcellation_path.resolve()}", + f"-o {warped_parcellation_path.resolve()}", ] # Call antsApplyTransforms run_ext_cmd(name="antsApplyTransforms", cmd=apply_transforms_cmd) - img_to_merge = nib.load(warped_resampled_parcellation_path) + img = nib.load(warped_parcellation_path) + + # Resample parcellation to target image + img_to_merge = image.resample_to_img( + source_img=img, + target_img=target_img, + interpolation="nearest", + copy=True, + ) all_parcellations.append(img_to_merge) all_labels.append(labels) -- 2.52.0 From f32669d2ad671190904921522397d56797b61090 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 2 Feb 2024 09:16:21 +0100 Subject: [PATCH 10/47] fix: correct cmd for warping in get_mask --- junifer/data/masks.py | 23 ++++++++++------------- 1 file changed, 10 insertions(+), 13 deletions(-) diff --git a/junifer/data/masks.py b/junifer/data/masks.py index ff1a7ea4d..99ef6ddc9 100644 --- a/junifer/data/masks.py +++ b/junifer/data/masks.py @@ -345,11 +345,9 @@ def get_mask( # noqa: C901 extra_input=extra_input, ) - # Save resampled mask image to a component-scoped tempfile - resampled_mask_path = ( - tempdir / f"{mask_name}_resampled_to_{resolution}.nii.gz" - ) - nib.save(mask_img, resampled_mask_path) + # Save mask image to a component-scoped tempfile + mask_path = tempdir / f"{mask_name}.nii.gz" + nib.save(mask_img, mask_path) # Save template target_std_space_template_path = ( @@ -361,33 +359,32 @@ def get_mask( # noqa: C901 ) # Set warped mask path - warped_resampled_mask_path = element_tempdir / ( - f"{mask_name}_resampled_to_{resolution}_warped_to_" + warped_mask_path = element_tempdir / ( + f"{mask_name}_warped_from_{mask_space}_to_" f"{target_std_space}.nii.gz" ) logger.debug( f"Using ANTs to warp {mask_name} " - f"(resolution={resolution}) from {mask_space} to " - f"{target_std_space}" + f"from {mask_space} to {target_std_space}" ) # Set antsApplyTransforms command apply_transforms_cmd = [ "antsApplyTransforms", "-d 3", - "-e 2", + "-e 3", "-n 'GenericLabel[NearestNeighbor]'", - f"-i {resampled_mask_path.resolve()}", + f"-i {mask_path.resolve()}", f"-r {target_std_space_template_path.resolve()}", f"-t {xfm_file_path.resolve()}", - f"-o {warped_resampled_mask_path.resolve()}", + f"-o {warped_mask_path.resolve()}", ] # Call antsApplyTransforms run_ext_cmd( name="antsApplyTransforms", cmd=apply_transforms_cmd ) - mask_img = nib.load(warped_resampled_mask_path) + mask_img = nib.load(warped_mask_path) all_masks.append(mask_img) -- 2.52.0 From 9e2215624d5f122dd5e2031b8ff8c1462b543a03 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 2 Feb 2024 11:30:54 +0100 Subject: [PATCH 11/47] chore: add changelog 299.feature --- docs/changes/newsfragments/299.feature | 1 + 1 file changed, 1 insertion(+) create mode 100644 docs/changes/newsfragments/299.feature diff --git a/docs/changes/newsfragments/299.feature b/docs/changes/newsfragments/299.feature new file mode 100644 index 000000000..13219b002 --- /dev/null +++ b/docs/changes/newsfragments/299.feature @@ -0,0 +1 @@ +Add support for on-the-fly template space transformation in :func:`.get_parcellation` and :func:`.get_mask` to allow parcellation and mask in different template spaces to work with a ``DataGrabber``'s data in a specified template space. -- 2.52.0 From 55bf8a2a8176c5d2e6937022d2984b1b90daf659 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 2 Feb 2024 12:03:39 +0100 Subject: [PATCH 12/47] update: run jobs in docs-preview.yml inside ANTs container --- .github/workflows/docs-preview.yml | 1 + 1 file changed, 1 insertion(+) diff --git a/.github/workflows/docs-preview.yml b/.github/workflows/docs-preview.yml index d70e97ae5..154ae47da 100644 --- a/.github/workflows/docs-preview.yml +++ b/.github/workflows/docs-preview.yml @@ -18,6 +18,7 @@ concurrency: preview-${{ github.ref }} jobs: build-docs: runs-on: ubuntu-latest + container: antsx/ants:latest steps: - name: Checkout source uses: actions/checkout@v4 -- 2.52.0 From 7e14371bd23642c9fcefb4cec1749f2a06871811 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Mon, 5 Feb 2024 14:14:25 +0100 Subject: [PATCH 13/47] Revert "update: run jobs in docs-preview.yml inside ANTs container" This reverts commit e0b7e3aac76f25eba26cfcc120e83ca611383b12. --- .github/workflows/docs-preview.yml | 1 - 1 file changed, 1 deletion(-) diff --git a/.github/workflows/docs-preview.yml b/.github/workflows/docs-preview.yml index 154ae47da..d70e97ae5 100644 --- a/.github/workflows/docs-preview.yml +++ b/.github/workflows/docs-preview.yml @@ -18,7 +18,6 @@ concurrency: preview-${{ github.ref }} jobs: build-docs: runs-on: ubuntu-latest - container: antsx/ants:latest steps: - name: Checkout source uses: actions/checkout@v4 -- 2.52.0 From a1d5327c5d599dae30cc16353099c3cb3f33db26 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Mon, 5 Feb 2024 14:23:07 +0100 Subject: [PATCH 14/47] chore: add ANTs installation in docs-preview.yml --- .github/workflows/docs-preview.yml | 14 ++++++++++++++ 1 file changed, 14 insertions(+) diff --git a/.github/workflows/docs-preview.yml b/.github/workflows/docs-preview.yml index d70e97ae5..c7f4d1233 100644 --- a/.github/workflows/docs-preview.yml +++ b/.github/workflows/docs-preview.yml @@ -44,6 +44,20 @@ jobs: run: | git config --global user.email "runner@github.com" git config --global user.name "GITHUB CI Runner" + - name: Install ANTs + run: | + echo "++ Add universe repo" + sudo add-apt-repository -y universe + + echo "++ Update package manager info" + sudo apt-get update -qq + + echo "++ Downloading ANTs" + curl -fsSL -o ants.zip https://github.com/ANTsX/ANTs/releases/download/v2.5.1/ants-2.5.1-ubuntu-22.04-X64-gcc.zip + unzip ants.zip -d /opt + mv /opt/ants-2.5.1/bin/* /opt/ants-2.5.1 + rm ants.zip + echo "/opt/ants-2.5.1" >> $GITHUB_PATH - name: Test build docs run: | BUILDDIR=_build/main make -C docs/ local -- 2.52.0 From 3ed2315b70b8d5dcec09376d6b555c124f116714 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Wed, 14 Feb 2024 18:06:33 +0100 Subject: [PATCH 15/47] update: add yaml for no space conversion example and test data --- examples/yamls/partly_cloudy_agg_mean_tian.yml | 16 ++++++++++++++++ .../tests/data/partly_cloudy_agg_mean_tian.yml | 16 ++++++++++++++++ 2 files changed, 32 insertions(+) create mode 100644 examples/yamls/partly_cloudy_agg_mean_tian.yml create mode 100644 junifer/api/tests/data/partly_cloudy_agg_mean_tian.yml diff --git a/examples/yamls/partly_cloudy_agg_mean_tian.yml b/examples/yamls/partly_cloudy_agg_mean_tian.yml new file mode 100644 index 000000000..da4ee0448 --- /dev/null +++ b/examples/yamls/partly_cloudy_agg_mean_tian.yml @@ -0,0 +1,16 @@ +workdir: /tmp + +with: junifer.testing.registry + +datagrabber: + kind: PartlyCloudyTestingDataGrabber + +markers: + - kind: ParcelAggregation + parcellation: TianxS1x3TxMNInonlinear2009cAsym + method: mean + name: tian-s1-3T_mean + +storage: + kind: HDF5FeatureStorage + uri: /tmp/partly_cloudy_agg_mean_tian.hdf5 diff --git a/junifer/api/tests/data/partly_cloudy_agg_mean_tian.yml b/junifer/api/tests/data/partly_cloudy_agg_mean_tian.yml new file mode 100644 index 000000000..da4ee0448 --- /dev/null +++ b/junifer/api/tests/data/partly_cloudy_agg_mean_tian.yml @@ -0,0 +1,16 @@ +workdir: /tmp + +with: junifer.testing.registry + +datagrabber: + kind: PartlyCloudyTestingDataGrabber + +markers: + - kind: ParcelAggregation + parcellation: TianxS1x3TxMNInonlinear2009cAsym + method: mean + name: tian-s1-3T_mean + +storage: + kind: HDF5FeatureStorage + uri: /tmp/partly_cloudy_agg_mean_tian.hdf5 -- 2.52.0 From 8b7a1358736529f9a4f9593b5744b76458648ab8 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Wed, 14 Feb 2024 18:07:23 +0100 Subject: [PATCH 16/47] fix: update tests in test_cli.py --- junifer/api/tests/test_cli.py | 18 ++++++------------ 1 file changed, 6 insertions(+), 12 deletions(-) diff --git a/junifer/api/tests/test_cli.py b/junifer/api/tests/test_cli.py index c4167db4b..36eb49f9d 100644 --- a/junifer/api/tests/test_cli.py +++ b/junifer/api/tests/test_cli.py @@ -54,16 +54,13 @@ def test_run_and_collect_commands( """ # Get test config - infile = Path(__file__).parent / "data" / "gmd_mean.yaml" + infile = Path(__file__).parent / "data" / "partly_cloudy_agg_mean_tian.yml" # Read test config contents = yaml.load(infile) # Working directory - workdir = tmp_path / "workdir" - contents["workdir"] = str(workdir.resolve()) - # Output directory - outdir = tmp_path / "outdir" + contents["workdir"] = str(tmp_path.resolve()) # Storage - contents["storage"]["uri"] = str(outdir.resolve()) + contents["storage"]["uri"] = str((tmp_path / "out.hdf5").resolve()) # Write new test config outfile = tmp_path / "in.yaml" yaml.dump(contents, stream=outfile) @@ -117,16 +114,13 @@ def test_run_using_element_file(tmp_path: Path, elements: str) -> None: f.write(elements) # Get test config - infile = Path(__file__).parent / "data" / "gmd_mean.yaml" + infile = Path(__file__).parent / "data" / "partly_cloudy_agg_mean_tian.yml" # Read test config contents = yaml.load(infile) # Working directory - workdir = tmp_path / "workdir" - contents["workdir"] = str(workdir.resolve()) - # Output directory - outdir = tmp_path / "outdir" + contents["workdir"] = str(tmp_path.resolve()) # Storage - contents["storage"]["uri"] = str(outdir.resolve()) + contents["storage"]["uri"] = str((tmp_path / "out.hdf5").resolve()) # Write new test config outfile = tmp_path / "in.yaml" yaml.dump(contents, stream=outfile) -- 2.52.0 From d95658323ae9b5cb0ba880f3ee1aae42781ea9b8 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 17:59:19 +0100 Subject: [PATCH 17/47] update: improve tests for test_functions.py --- junifer/api/tests/test_functions.py | 212 +++++++++++++++------------- 1 file changed, 116 insertions(+), 96 deletions(-) diff --git a/junifer/api/tests/test_functions.py b/junifer/api/tests/test_functions.py index 687ffcded..728a2a6b8 100644 --- a/junifer/api/tests/test_functions.py +++ b/junifer/api/tests/test_functions.py @@ -7,7 +7,7 @@ import logging from pathlib import Path -from typing import List, Tuple, Union +from typing import Dict, List, Tuple, Union import pytest from ruamel.yaml import YAML @@ -24,97 +24,104 @@ yaml.default_flow_style = False yaml.allow_unicode = True yaml.indent(mapping=2, sequence=4, offset=2) -# Define datagrabber -datagrabber = { - "kind": "OasisVBMTestingDataGrabber", -} -# Define markers -markers = [ - { - "name": "Schaefer1000x7_Mean", - "kind": "ParcelAggregation", - "parcellation": "Schaefer1000x7", - "method": "mean", - }, - { - "name": "Schaefer1000x7_Std", - "kind": "ParcelAggregation", - "parcellation": "Schaefer1000x7", - "method": "std", - }, -] - -# Define storage -storage = { - "kind": "SQLiteFeatureStorage", -} +@pytest.fixture +def datagrabber() -> Dict[str, str]: + """Return a datagrabber as a dictionary.""" + return { + "kind": "PartlyCloudyTestingDataGrabber", + } -def test_run_single_element(tmp_path: Path) -> None: +@pytest.fixture +def markers() -> List[Dict[str, str]]: + """Return markers as a list of dictionary.""" + return [ + { + "name": "tian-s1-3T_mean", + "kind": "ParcelAggregation", + "parcellation": "TianxS1x3TxMNInonlinear2009cAsym", + "method": "mean", + }, + { + "name": "tian-s1-3T_std", + "kind": "ParcelAggregation", + "parcellation": "TianxS1x3TxMNInonlinear2009cAsym", + "method": "std", + }, + ] + + +@pytest.fixture +def storage() -> Dict[str, str]: + """Return a storage as a dictionary.""" + return { + "kind": "SQLiteFeatureStorage", + } + + +def test_run_single_element( + tmp_path: Path, + datagrabber: Dict[str, str], + markers: List[Dict[str, str]], + storage: Dict[str, str], +) -> None: """Test run function with single element. Parameters ---------- tmp_path : pathlib.Path The path to the test directory. + datagrabber : dict + Testing datagrabber as dictionary. + markers : list of dict + Testing markers as list of dictionary. + storage : dict + Testing storage as dictionary. """ - # Create working directory - workdir = tmp_path / "workdir_single" - workdir.mkdir() - # Create output directory - outdir = workdir / "out" - outdir.mkdir() - # Create storage - uri = outdir / "test.sqlite" - storage["uri"] = uri # type: ignore + # Set storage + storage["uri"] = str((tmp_path / "out.sqlite").resolve()) # Run operations run( - workdir=workdir, + workdir=tmp_path, datagrabber=datagrabber, markers=markers, storage=storage, elements=["sub-01"], ) # Check files - files = list(outdir.glob("*.sqlite")) + files = list(tmp_path.glob("*.sqlite")) assert len(files) == 1 -def test_run_single_element_with_preprocessing(tmp_path: Path) -> None: +def test_run_single_element_with_preprocessing( + tmp_path: Path, + markers: List[Dict[str, str]], + storage: Dict[str, str], +) -> None: """Test run function with single element and pre-processing. Parameters ---------- tmp_path : pathlib.Path The path to the test directory. + markers : list of dict + Testing markers as list of dictionary. + storage : dict + Testing storage as dictionary. """ - # Create working directory - workdir = tmp_path / "workdir_single_with_preprocess" - workdir.mkdir() - # Create output directory - outdir = workdir / "out" - outdir.mkdir() - # Create storage - uri = outdir / "test.sqlite" - storage["uri"] = uri # type: ignore + # Set storage + storage["uri"] = str((tmp_path / "out.sqlite").resolve()) # Run operations run( - workdir=workdir, + workdir=tmp_path, datagrabber={ "kind": "PartlyCloudyTestingDataGrabber", "reduce_confounds": False, }, - markers=[ - { - "name": "Schaefer100x17_mean_FC", - "kind": "FunctionalConnectivityParcels", - "parcellation": "Schaefer100x17", - "agg_method": "mean", - } - ], + markers=markers, storage=storage, preprocessors=[ { @@ -124,97 +131,110 @@ def test_run_single_element_with_preprocessing(tmp_path: Path) -> None: elements=["sub-01"], ) # Check files - files = list(outdir.glob("*.sqlite")) + files = list(tmp_path.glob("*.sqlite")) assert len(files) == 1 -def test_run_multi_element(tmp_path: Path) -> None: - """Test run function with multi element. +def test_run_multi_element_multi_output( + tmp_path: Path, + datagrabber: Dict[str, str], + markers: List[Dict[str, str]], + storage: Dict[str, str], +) -> None: + """Test run function with multi element and multi output. Parameters ---------- tmp_path : pathlib.Path The path to the test directory. + datagrabber : dict + Testing datagrabber as dictionary. + markers : list of dict + Testing markers as list of dictionary. + storage : dict + Testing storage as dictionary. """ - # Create working directory - workdir = tmp_path / "workdir_multi" - workdir.mkdir() - # Create output directory - outdir = workdir / "out" - outdir.mkdir() - # Create storage - uri = outdir / "test.sqlite" - storage["uri"] = uri # type: ignore + # Set storage + storage["uri"] = str((tmp_path / "out.sqlite").resolve()) storage["single_output"] = False # type: ignore # Run operations run( - workdir=workdir, + workdir=tmp_path, datagrabber=datagrabber, markers=markers, storage=storage, elements=["sub-01", "sub-03"], ) # Check files - files = list(outdir.glob("*.sqlite")) + files = list(tmp_path.glob("*.sqlite")) assert len(files) == 2 -def test_run_multi_element_single_output(tmp_path: Path) -> None: - """Test run function with multi element. +def test_run_multi_element_single_output( + tmp_path: Path, + datagrabber: Dict[str, str], + markers: List[Dict[str, str]], + storage: Dict[str, str], +) -> None: + """Test run function with multi element and single output. Parameters ---------- tmp_path : pathlib.Path The path to the test directory. + datagrabber : dict + Testing datagrabber as dictionary. + markers : list of dict + Testing markers as list of dictionary. + storage : dict + Testing storage as dictionary. """ - # Create working directory - workdir = tmp_path / "workdir_multi" - workdir.mkdir() - # Create output directory - outdir = workdir / "out" - outdir.mkdir() - # Create storage - uri = outdir / "test.sqlite" - storage["uri"] = uri # type: ignore + # Set storage + storage["uri"] = str((tmp_path / "out.sqlite").resolve()) storage["single_output"] = True # type: ignore # Run operations run( - workdir=workdir, + workdir=tmp_path, datagrabber=datagrabber, markers=markers, storage=storage, elements=["sub-01", "sub-03"], ) # Check files - files = list(outdir.glob("*.sqlite")) + files = list(tmp_path.glob("*.sqlite")) assert len(files) == 1 - assert files[0].name == "test.sqlite" + assert files[0].name == "out.sqlite" -def test_run_and_collect(tmp_path: Path) -> None: +def test_run_and_collect( + tmp_path: Path, + datagrabber: Dict[str, str], + markers: List[Dict[str, str]], + storage: Dict[str, str], +) -> None: """Test run and collect functions. Parameters ---------- tmp_path : pathlib.Path The path to the test directory. + datagrabber : dict + Testing datagrabber as dictionary. + markers : list of dict + Testing markers as list of dictionary. + storage : dict + Testing storage as dictionary. """ - # Create working directory - workdir = tmp_path / "workdir" - workdir.mkdir() - # Create output directory - outdir = workdir / "out" - outdir.mkdir() - # Create storage - uri = outdir / "test.sqlite" - storage["uri"] = uri # type: ignore + # Set storage + uri = tmp_path / "out.sqlite" + storage["uri"] = str(uri.resolve()) storage["single_output"] = False # type: ignore # Run operations run( - workdir=workdir, + workdir=tmp_path, datagrabber=datagrabber, markers=markers, storage=storage, @@ -225,7 +245,7 @@ def test_run_and_collect(tmp_path: Path) -> None: ) elements = dg.get_elements() # type: ignore # This should create 10 files - files = list(outdir.glob("*.sqlite")) + files = list(tmp_path.glob("*.sqlite")) assert len(files) == len(elements) # But the test.sqlite file should not exist assert not uri.exists() -- 2.52.0 From bdbca3a68ef8e4c4d8c7cc7a108fc4411104f8f5 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:01:59 +0100 Subject: [PATCH 18/47] update: improve tests for test_masks.py --- junifer/data/tests/test_masks.py | 129 ++++++++++++++++++------------- 1 file changed, 74 insertions(+), 55 deletions(-) diff --git a/junifer/data/tests/test_masks.py b/junifer/data/tests/test_masks.py index 4bc23f2da..ff5394ce9 100644 --- a/junifer/data/tests/test_masks.py +++ b/junifer/data/tests/test_masks.py @@ -29,6 +29,7 @@ from junifer.data.masks import ( register_mask, ) from junifer.datareader import DefaultDataReader +from junifer.pipeline.utils import _check_ants from junifer.testing.datagrabbers import ( OasisVBMTestingDataGrabber, SPMAuditoryTestingDataGrabber, @@ -215,18 +216,19 @@ def test_vickery_patil_error() -> None: def test_get_mask() -> None: """Test the get_mask function.""" - reader = DefaultDataReader() with OasisVBMTestingDataGrabber() as dg: - input = dg["sub-01"] - input = reader.fit_transform(input) - vbm_gm = input["VBM_GM"] + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + vbm_gm = element_data["VBM_GM"] vbm_gm_img = vbm_gm["data"] - mask = get_mask(masks="GM_prob0.2", target_data=vbm_gm) + mask = get_mask(masks="compute_brain_mask", target_data=vbm_gm) assert mask.shape == vbm_gm_img.shape assert_array_equal(mask.affine, vbm_gm_img.affine) - raw_mask_img, _, _ = load_mask("GM_prob0.2", resolution=1.5) + raw_mask_callable, _, _ = load_mask( + "compute_brain_mask", resolution=1.5 + ) + raw_mask_img = raw_mask_callable(vbm_gm_img) # type: ignore res_mask_img = resample_to_img( raw_mask_img, vbm_gm_img, @@ -245,13 +247,11 @@ def test_mask_callable() -> None: _available_masks["identity"] = { "family": "Callable", "func": ident, - "space": "MNI", + "space": "MNI152Lin", } - reader = DefaultDataReader() with OasisVBMTestingDataGrabber() as dg: - input = dg["sub-01"] - input = reader.fit_transform(input) - vbm_gm = input["VBM_GM"] + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + vbm_gm = element_data["VBM_GM"] vbm_gm_img = vbm_gm["data"] mask = get_mask(masks="identity", target_data=vbm_gm) @@ -262,11 +262,9 @@ def test_mask_callable() -> None: def test_get_mask_errors() -> None: """Test passing wrong parameters to get_mask.""" - reader = DefaultDataReader() with OasisVBMTestingDataGrabber() as dg: - input = dg["sub-01"] - input = reader.fit_transform(input) - vbm_gm = input["VBM_GM"] + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + vbm_gm = element_data["VBM_GM"] # Test wrong masks definitions (more than one key per dict) with pytest.raises(ValueError, match=r"only one key"): get_mask(masks={"GM_prob0.2": {}, "Other": {}}, target_data=vbm_gm) @@ -286,7 +284,8 @@ def test_get_mask_errors() -> None: ValueError, match=r"parameters to the intersection" ): get_mask( - masks=["GM_prob0.2", {"threshold": 1}], target_data=vbm_gm + masks=["compute_brain_mask", {"threshold": 1}], + target_data=vbm_gm, ) # Test "inherited" masks errors @@ -317,12 +316,6 @@ def test_get_mask_errors() -> None: ("compute_brain_mask", compute_brain_mask, {"threshold": 0.2}, False), ("compute_background_mask", compute_background_mask, None, False), ("compute_epi_mask", compute_epi_mask, None, False), - ( - "fetch_icbm152_brain_gm_mask", - fetch_icbm152_brain_gm_mask, - None, - True, - ), ], ) def test_nilearn_compute_masks( @@ -345,11 +338,9 @@ def test_nilearn_compute_masks( Whether to resample the mask to the target data. """ - reader = DefaultDataReader() with SPMAuditoryTestingDataGrabber() as dg: - input = dg["sub001"] - input = reader.fit_transform(input) - bold = input["BOLD"] + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + bold = element_data["BOLD"] bold_img = bold["data"] if params is None: @@ -376,29 +367,59 @@ def test_nilearn_compute_masks( assert_array_equal(mask.get_fdata(), ni_mask.get_fdata()) +@pytest.mark.skipif( + _check_ants() is False, reason="requires ANTs to be in PATH" +) +def test_nilearn_compute_masks_with_space_transformation() -> None: + """Test nilearn compute mask function with space transformation.""" + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + bold = element_data["BOLD"] + bold_img = bold["data"] + + # Get mask + mask = get_mask(masks="fetch_icbm152_brain_gm_mask", target_data=bold) + # Check affine matrices are almost equal + assert_array_equal(mask.affine, bold_img.affine) + + # Mask needs resampling + ni_mask = resample_to_img( + fetch_icbm152_brain_gm_mask(), + bold_img, + interpolation="nearest", + copy=True, + ) + assert_array_equal(mask.get_fdata(), ni_mask.get_fdata()) + + def test_get_mask_inherit() -> None: """Test using the inherit mask functionality.""" - reader = DefaultDataReader() with SPMAuditoryTestingDataGrabber() as dg: - input = dg["sub001"] - input = reader.fit_transform(input) + element_data = DefaultDataReader().fit_transform(dg["sub001"]) # Compute brain mask using nilearn - gm_mask = compute_brain_mask(input["BOLD"]["data"], threshold=0.2) + gm_mask = compute_brain_mask( + element_data["BOLD"]["data"], threshold=0.2 + ) # Get mask using the compute_brain_mask function mask1 = get_mask( masks={"compute_brain_mask": {"threshold": 0.2}}, - target_data=input["BOLD"], + target_data=element_data["BOLD"], ) # Now get the mask using the inherit functionality, passing the # computed mask as extra data extra_input = { - "BOLD_MASK": {"data": gm_mask, "space": input["BOLD"]["space"]} + "BOLD_MASK": { + "data": gm_mask, + "space": element_data["BOLD"]["space"], + } } - input["BOLD"]["mask_item"] = "BOLD_MASK" + element_data["BOLD"]["mask_item"] = "BOLD_MASK" mask2 = get_mask( - masks="inherit", target_data=input["BOLD"], extra_input=extra_input + masks="inherit", + target_data=element_data["BOLD"], + extra_input=extra_input, ) # Both masks should be equal @@ -408,7 +429,6 @@ def test_get_mask_inherit() -> None: @pytest.mark.parametrize( "masks,params", [ - (["GM_prob0.2", "GM_prob0.2_cortex"], {}), (["compute_brain_mask", "compute_background_mask"], {}), (["compute_brain_mask", "compute_epi_mask"], {}), ], @@ -426,10 +446,8 @@ def test_get_mask_multiple( Parameters to pass to the intersect_masks function. """ - reader = DefaultDataReader() with SPMAuditoryTestingDataGrabber() as dg: - input = dg["sub001"] - input = reader.fit_transform(input) + element_data = DefaultDataReader().fit_transform(dg["sub001"]) if not isinstance(masks, list): junifer_masks = [masks] else: @@ -438,10 +456,12 @@ def test_get_mask_multiple( # Convert params to junifer style (one dict per param) junifer_params = [{k: params[k]} for k in params.keys()] junifer_masks.extend(junifer_params) - target_img = input["BOLD"]["data"] + target_img = element_data["BOLD"]["data"] resolution = np.min(target_img.header.get_zooms()[:3]) - computed = get_mask(masks=junifer_masks, target_data=input["BOLD"]) + computed = get_mask( + masks=junifer_masks, target_data=element_data["BOLD"] + ) masks_names = [ next(iter(x.keys())) if isinstance(x, dict) else x for x in masks @@ -480,19 +500,18 @@ def test_get_mask_multiple( assert_array_equal(computed.get_fdata(), expected.get_fdata()) -def test_get_mask_multiple_incorrect_space() -> None: - """Test incorrect space error for getting multiple masks.""" - reader = DefaultDataReader() +@pytest.mark.skipif( + _check_ants() is False, reason="requires ANTs to be in PATH" +) +def test_get_mask_multiple_different_space() -> None: + """Test tailored multiple mask fetch in different space.""" with SPMAuditoryTestingDataGrabber() as dg: - input = dg["sub001"] - input = reader.fit_transform(input) - - with pytest.raises(RuntimeError, match="unable to merge."): - get_mask( - masks=[ - "GM_prob0.2", - "compute_brain_mask", - "fetch_icbm152_brain_gm_mask", - ], - target_data=input["BOLD"], - ) + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + # Get tailored mask + get_mask( + masks=[ + "compute_brain_mask", + "fetch_icbm152_brain_gm_mask", + ], + target_data=element_data["BOLD"], + ) -- 2.52.0 From 6f47e4b5c8fcd099a63537b23b80cd28db81341b Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:02:20 +0100 Subject: [PATCH 19/47] update: improve tests for test_parcellations.py --- junifer/data/tests/test_parcellations.py | 115 +++++++++++++---------- 1 file changed, 66 insertions(+), 49 deletions(-) diff --git a/junifer/data/tests/test_parcellations.py b/junifer/data/tests/test_parcellations.py index a19cac8f6..e730be664 100644 --- a/junifer/data/tests/test_parcellations.py +++ b/junifer/data/tests/test_parcellations.py @@ -30,7 +30,11 @@ from junifer.data.parcellations import ( register_parcellation, ) from junifer.datareader import DefaultDataReader -from junifer.testing.datagrabbers import OasisVBMTestingDataGrabber +from junifer.pipeline.utils import _check_ants +from junifer.testing.datagrabbers import ( + OasisVBMTestingDataGrabber, + PartlyCloudyTestingDataGrabber, +) def test_register_parcellation_built_in_check() -> None: @@ -58,7 +62,7 @@ def test_register_parcellation_already_registered() -> None: name="testparc", parcellation_path="testparc.nii.gz", parcels_labels=["1", "2", "3"], - space="MNI", + space="MNI152Lin", ) assert ( load_parcellation("testparc", path_only=True)[2].name @@ -71,13 +75,13 @@ def test_register_parcellation_already_registered() -> None: name="testparc", parcellation_path="testparc.nii.gz", parcels_labels=["1", "2", "3"], - space="MNI", + space="MNI152Lin", ) register_parcellation( name="testparc", parcellation_path="testparc2.nii.gz", parcels_labels=["1", "2", "3"], - space="MNI", + space="MNI152Lin", overwrite=True, ) @@ -100,14 +104,16 @@ def test_parcellation_wrong_labels_values(tmp_path: Path) -> None: assert schaefer is not None # Test wrong number of labels - register_parcellation("WrongLabels", schaefer_path, labels[:10], "MNI") + register_parcellation( + "WrongLabels", schaefer_path, labels[:10], "MNI152Lin" + ) with pytest.raises(ValueError, match=r"has 100 parcels but 10"): load_parcellation("WrongLabels") # Test wrong number of labels register_parcellation( - "WrongLabels2", schaefer_path, [*labels, "wrong"], "MNI" + "WrongLabels2", schaefer_path, [*labels, "wrong"], "MNI152Lin" ) with pytest.raises(ValueError, match=r"has 100 parcels but 101"): @@ -119,7 +125,9 @@ def test_parcellation_wrong_labels_values(tmp_path: Path) -> None: new_schaefer_img = new_img_like(schaefer, schaefer_data) nib.save(new_schaefer_img, new_schaefer_path) - register_parcellation("WrongValues", new_schaefer_path, labels[:-1], "MNI") + register_parcellation( + "WrongValues", new_schaefer_path, labels[:-1], "MNI152Lin" + ) with pytest.raises(ValueError, match=r"the range [0, 99]"): load_parcellation("WrongValues") @@ -129,7 +137,9 @@ def test_parcellation_wrong_labels_values(tmp_path: Path) -> None: new_schaefer_img = new_img_like(schaefer, schaefer_data) nib.save(new_schaefer_img, new_schaefer_path) - register_parcellation("WrongValues2", new_schaefer_path, labels, "MNI") + register_parcellation( + "WrongValues2", new_schaefer_path, labels, "MNI152Lin" + ) with pytest.raises(ValueError, match=r"the range [0, 100]"): load_parcellation("WrongValues2") @@ -137,13 +147,25 @@ def test_parcellation_wrong_labels_values(tmp_path: Path) -> None: @pytest.mark.parametrize( "name, parcellation_path, parcels_labels, space, overwrite", [ - ("testparc_1", "testparc_1.nii.gz", ["1", "2", "3"], "MNI", True), - ("testparc_2", "testparc_2.nii.gz", ["1", "2", "6"], "MNI", True), + ( + "testparc_1", + "testparc_1.nii.gz", + ["1", "2", "3"], + "MNI152Lin", + True, + ), + ( + "testparc_2", + "testparc_2.nii.gz", + ["1", "2", "6"], + "MNI152Lin", + True, + ), ( "testparc_3", Path("testparc_3.nii.gz"), ["1", "2", "6"], - "MNI", + "MNI152Lin", True, ), ], @@ -1172,28 +1194,26 @@ def test_merge_parcellations_3D_multiple_duplicated_labels() -> None: def test_get_parcellation_single() -> None: """Test tailored single parcellation fetch.""" - reader = DefaultDataReader() - with OasisVBMTestingDataGrabber() as dg: - element = dg["sub-01"] - element_data = reader.fit_transform(element) - vbm_gm = element_data["VBM_GM"] - vbm_gm_img = vbm_gm["data"] + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + bold = element_data["BOLD"] + bold_img = bold["data"] # Get tailored parcellation tailored_parcellation, tailored_labels = get_parcellation( - parcellation=["Schaefer100x7"], - target_data=vbm_gm, + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=bold, ) # Check shape and affine with original element data - assert tailored_parcellation.shape == vbm_gm_img.shape - assert_array_equal(tailored_parcellation.affine, vbm_gm_img.affine) + assert tailored_parcellation.shape == bold_img.shape[:3] + assert_array_equal(tailored_parcellation.affine, bold_img.affine) # Get raw parcellation raw_parcellation, raw_labels, _, _ = load_parcellation( - "Schaefer100x7", + "TianxS1x3TxMNInonlinear2009cAsym", resolution=1.5, ) resampled_raw_parcellation = resample_to_img( source_img=raw_parcellation, - target_img=vbm_gm_img, + target_img=bold_img, interpolation="nearest", copy=True, ) @@ -1207,36 +1227,34 @@ def test_get_parcellation_single() -> None: def test_get_parcellation_multi_same_space() -> None: """Test tailored multi parcellation fetch in same space.""" - reader = DefaultDataReader() - with OasisVBMTestingDataGrabber() as dg: - element = dg["sub-01"] - element_data = reader.fit_transform(element) - vbm_gm = element_data["VBM_GM"] - vbm_gm_img = vbm_gm["data"] + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + bold = element_data["BOLD"] + bold_img = bold["data"] # Get tailored parcellation tailored_parcellation, tailored_labels = get_parcellation( parcellation=[ - "Schaefer100x7", - "TianxS2x3TxMNI6thgeneration", + "Shen_2015_268", + "TianxS1x3TxMNInonlinear2009cAsym", ], - target_data=vbm_gm, + target_data=bold, ) # Check shape and affine with original element data - assert tailored_parcellation.shape == vbm_gm_img.shape - assert_array_equal(tailored_parcellation.affine, vbm_gm_img.affine) + assert tailored_parcellation.shape == bold_img.shape[:3] + assert_array_equal(tailored_parcellation.affine, bold_img.affine) # Get raw parcellations raw_parcellations = [] raw_labels = [] parcellations_names = [ - "Schaefer100x7", - "TianxS2x3TxMNI6thgeneration", + "Shen_2015_268", + "TianxS1x3TxMNInonlinear2009cAsym", ] for name in parcellations_names: img, labels, _, _ = load_parcellation(name=name, resolution=1.5) # Resample raw parcellations resampled_img = resample_to_img( source_img=img, - target_img=vbm_gm_img, + target_img=bold_img, interpolation="nearest", copy=True, ) @@ -1256,19 +1274,18 @@ def test_get_parcellation_multi_same_space() -> None: assert tailored_labels == merged_labels +@pytest.mark.skipif( + _check_ants() is False, reason="requires ANTs to be in PATH" +) def test_get_parcellation_multi_different_space() -> None: """Test tailored multi parcellation fetch in different space.""" - reader = DefaultDataReader() with OasisVBMTestingDataGrabber() as dg: - element = dg["sub-01"] - element_data = reader.fit_transform(element) - vbm_gm = element_data["VBM_GM"] + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) # Get tailored parcellation - with pytest.raises(RuntimeError, match="unable to merge."): - get_parcellation( - parcellation=[ - "Schaefer100x7", - "SUITxSUIT", - ], - target_data=vbm_gm, - ) + get_parcellation( + parcellation=[ + "Schaefer100x7", + "TianxS1x3TxMNInonlinear2009cAsym", + ], + target_data=element_data["VBM_GM"], + ) -- 2.52.0 From 5130ef59da5e3da6dedf746f28d7ed4a709e73ff Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:02:47 +0100 Subject: [PATCH 20/47] update: improve tests for test_falff_parcels.py --- junifer/markers/falff/tests/test_falff_parcels.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/junifer/markers/falff/tests/test_falff_parcels.py b/junifer/markers/falff/tests/test_falff_parcels.py index 4356c4954..0d847c8f4 100644 --- a/junifer/markers/falff/tests/test_falff_parcels.py +++ b/junifer/markers/falff/tests/test_falff_parcels.py @@ -18,7 +18,7 @@ from junifer.storage import SQLiteFeatureStorage from junifer.testing.datagrabbers import PartlyCloudyTestingDataGrabber -PARCELLATION = "Schaefer100x7" +PARCELLATION = "TianxS1x3TxMNInonlinear2009cAsym" def test_ALFFParcels(caplog: pytest.LogCaptureFixture, tmp_path: Path) -> None: @@ -59,7 +59,7 @@ def test_ALFFParcels(caplog: pytest.LogCaptureFixture, tmp_path: Path) -> None: output_bold_data = output_bold["data"] # Assert BOLD output data dimension assert output_bold_data.ndim == 2 - assert output_bold_data.shape == (1, 100) + assert output_bold_data.shape == (1, 16) # Reset log capture caplog.clear() -- 2.52.0 From 771b406d01600d40b899b5c6704eed0afe1fafd1 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:03:11 +0100 Subject: [PATCH 21/47] update: improve tests for test_falff_spheres.py --- junifer/markers/falff/tests/test_falff_spheres.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/junifer/markers/falff/tests/test_falff_spheres.py b/junifer/markers/falff/tests/test_falff_spheres.py index 115c6c205..d231d7f56 100644 --- a/junifer/markers/falff/tests/test_falff_spheres.py +++ b/junifer/markers/falff/tests/test_falff_spheres.py @@ -82,14 +82,14 @@ def test_ALFFSpheres(caplog: pytest.LogCaptureFixture, tmp_path: Path) -> None: "fractional", [True, False], ids=["fractional", "non-fractional"] ) def test_ALFFSpheres_comparison(tmp_path: Path, fractional: bool) -> None: - """Test ALFFSpheres using afni. + """Test ALFFSpheres implementation comparison. Parameters ---------- tmp_path : pathlib.Path - The path to the test directory. + The path to the test directory. fractional : bool - Whether to compute fractional ALFF or not. + Whether to compute fractional ALFF or not. """ with PartlyCloudyTestingDataGrabber() as dg: -- 2.52.0 From 250b95a6570da3e43343f405d8a98cab689e1b7c Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:03:32 +0100 Subject: [PATCH 22/47] update: improve tests for test_crossparcellation_functional_connectivity.py --- ...ossparcellation_functional_connectivity.py | 91 ++++++++++--------- 1 file changed, 46 insertions(+), 45 deletions(-) diff --git a/junifer/markers/functional_connectivity/tests/test_crossparcellation_functional_connectivity.py b/junifer/markers/functional_connectivity/tests/test_crossparcellation_functional_connectivity.py index 720b817fa..8fde0adf3 100644 --- a/junifer/markers/functional_connectivity/tests/test_crossparcellation_functional_connectivity.py +++ b/junifer/markers/functional_connectivity/tests/test_crossparcellation_functional_connectivity.py @@ -2,14 +2,17 @@ # Authors: Leonard Sasse # Kaustubh R. Patil +# Synchon Mandal # License: AGPL from pathlib import Path import pytest -from nilearn import image +from junifer.datareader import DefaultDataReader from junifer.markers.functional_connectivity import CrossParcellationFC +from junifer.pipeline import WorkDirManager +from junifer.pipeline.utils import _check_ants from junifer.storage import SQLiteFeatureStorage from junifer.testing.datagrabbers import SPMAuditoryTestingDataGrabber @@ -18,32 +21,53 @@ parcellation_one = "Schaefer100x17" parcellation_two = "Schaefer200x17" -def test_compute() -> None: - """Test CrossParcellationFC compute().""" +def test_init() -> None: + """Test CrossParcellationFC init().""" + with pytest.raises(ValueError, match="must be different"): + CrossParcellationFC( + parcellation_one="a", + parcellation_two="a", + correlation_method="pearson", + ) + +def test_get_output_type() -> None: + """Test CrossParcellationFC get_output_type().""" + crossparcellation = CrossParcellationFC( + parcellation_one=parcellation_one, parcellation_two=parcellation_two + ) + assert "matrix" == crossparcellation.get_output_type("BOLD") + + +@pytest.mark.skipif( + _check_ants() is False, reason="requires ANTs to be in PATH" +) +def test_compute(tmp_path: Path) -> None: + """Test CrossParcellationFC compute(). + + Parameters + ---------- + tmp_path : pathlib.Path + The path to the test directory. + + """ with SPMAuditoryTestingDataGrabber() as dg: - out = dg["sub001"] - niimg = image.load_img(str(out["BOLD"]["path"].absolute())) - input_dict = { - "BOLD": { - "data": niimg, - "path": out["BOLD"]["path"], - "meta": {"element": "sub001"}, - "space": "MNI", - } - } - + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + WorkDirManager().workdir = tmp_path crossparcellation = CrossParcellationFC( parcellation_one=parcellation_one, parcellation_two=parcellation_two, correlation_method="spearman", ) - out = crossparcellation.compute(input_dict["BOLD"]) + out = crossparcellation.compute(element_data["BOLD"]) assert out["data"].shape == (200, 100) assert len(out["col_names"]) == 100 assert len(out["row_names"]) == 200 +@pytest.mark.skipif( + _check_ants() is False, reason="requires ANTs to be in PATH" +) def test_store(tmp_path: Path) -> None: """Test CrossParcellationFC store(). @@ -53,43 +77,20 @@ def test_store(tmp_path: Path) -> None: The path to the test directory. """ - with SPMAuditoryTestingDataGrabber() as dg: - input_dict = dg["sub001"] - niimg = image.load_img(str(input_dict["BOLD"]["path"].absolute())) - - input_dict["BOLD"]["data"] = niimg - + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + WorkDirManager().workdir = tmp_path crossparcellation = CrossParcellationFC( parcellation_one=parcellation_one, parcellation_two=parcellation_two, correlation_method="spearman", ) - uri = tmp_path / "test_crossparcellation.sqlite" - storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") - crossparcellation.fit_transform(input_dict, storage=storage) + storage = SQLiteFeatureStorage( + uri=tmp_path / "test_crossparcellation.sqlite", upsert="ignore" + ) + # Fit transform marker on data with storage + crossparcellation.fit_transform(input=element_data, storage=storage) features = storage.list_features() assert any( x["name"] == "BOLD_CrossParcellationFC" for x in features.values() ) - - -def test_get_output_type() -> None: - """Test CrossParcellationFC get_output_type().""" - - crossparcellation = CrossParcellationFC( - parcellation_one=parcellation_one, parcellation_two=parcellation_two - ) - input_ = "BOLD" - output = crossparcellation.get_output_type(input_) - assert output == "matrix" - - -def test_init_() -> None: - """Test CrossParcellationFC init().""" - with pytest.raises(ValueError, match="must be different"): - CrossParcellationFC( - parcellation_one="a", - parcellation_two="a", - correlation_method="pearson", - ) -- 2.52.0 From 3960fb8fc6a641add1f704fd208c85f451cde0e1 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:04:11 +0100 Subject: [PATCH 23/47] update: improve tests for test_edge_functional_connectivity_parcels.py --- ...st_edge_functional_connectivity_parcels.py | 69 +++++++++---------- 1 file changed, 31 insertions(+), 38 deletions(-) diff --git a/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_parcels.py b/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_parcels.py index 7822652fc..2aabae185 100644 --- a/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_parcels.py +++ b/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_parcels.py @@ -6,10 +6,10 @@ from pathlib import Path -from nilearn import datasets, image - +from junifer.datareader import DefaultDataReader from junifer.markers.functional_connectivity import EdgeCentricFCParcels from junifer.storage import SQLiteFeatureStorage +from junifer.testing.datagrabbers import PartlyCloudyTestingDataGrabber def test_EdgeCentricFCParcels(tmp_path: Path) -> None: @@ -21,42 +21,35 @@ def test_EdgeCentricFCParcels(tmp_path: Path) -> None: The path to the test directory. """ - # get a dataset - ni_data = datasets.fetch_spm_auditory(subject_id="sub001") - fmri_img = image.concat_imgs(ni_data.func) # type: ignore + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + marker = EdgeCentricFCParcels( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + cor_method_params={"empirical": True}, + ) + # Check correct output + assert marker.get_output_type("BOLD") == "matrix" - # Check empirical correlation method parameters - efc = EdgeCentricFCParcels( - parcellation="TianxS1x3TxMNInonlinear2009cAsym", - cor_method_params={"empirical": True}, - ) - all_out = efc.fit_transform( - {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - ) + # Fit-transform the data + edge_fc = marker.fit_transform(element_data) + edge_fc_bold = edge_fc["BOLD"] - out = all_out["BOLD"] + # For 16 ROIs we should get (16 * (16 -1) / 2) edges in the ETS + n_edges = int(16 * (16 - 1) / 2) + assert "data" in edge_fc_bold + assert "row_names" in edge_fc_bold + assert "col_names" in edge_fc_bold + assert edge_fc_bold["data"].shape[0] == n_edges + assert edge_fc_bold["data"].shape[1] == n_edges + assert len(set(edge_fc_bold["row_names"])) == n_edges + assert len(set(edge_fc_bold["col_names"])) == n_edges - # for 16 ROIs we should get (16 * (16 -1) / 2) edges in the ETS - n_edges = int(16 * (16 - 1) / 2) - assert "data" in out - assert "row_names" in out - assert "col_names" in out - assert out["data"].shape[0] == n_edges - assert out["data"].shape[1] == n_edges - assert len(set(out["row_names"])) == n_edges - assert len(set(out["col_names"])) == n_edges - - # check correct output - assert efc.get_output_type("BOLD") == "matrix" - - uri = tmp_path / "test_fc_parcellation.sqlite" - # Single storage, must be the uri - storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") - meta = {"element": {"subject": "test"}, "dependencies": {"numpy"}} - input = {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - all_out = efc.fit_transform(input, storage=storage) - - features = storage.list_features() - assert any( - x["name"] == "BOLD_EdgeCentricFCParcels" for x in features.values() - ) + # Store + storage = SQLiteFeatureStorage( + uri=tmp_path / "test_edge_fc_parcels.sqlite", upsert="ignore" + ) + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "BOLD_EdgeCentricFCParcels" for x in features.values() + ) -- 2.52.0 From b39eb0db302086c403f0e1f2ee671970a0b25616 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:04:40 +0100 Subject: [PATCH 24/47] update: improve tests for test_edge_functional_connectivity_spheres.py --- ...st_edge_functional_connectivity_spheres.py | 90 ++++++++----------- 1 file changed, 37 insertions(+), 53 deletions(-) diff --git a/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_spheres.py b/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_spheres.py index 8695e55f5..38b871ebf 100644 --- a/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_spheres.py +++ b/junifer/markers/functional_connectivity/tests/test_edge_functional_connectivity_spheres.py @@ -6,10 +6,10 @@ from pathlib import Path -from nilearn import datasets, image - +from junifer.datareader import DefaultDataReader from junifer.markers.functional_connectivity import EdgeCentricFCSpheres from junifer.storage import SQLiteFeatureStorage +from junifer.testing.datagrabbers import SPMAuditoryTestingDataGrabber def test_EdgeCentricFCSpheres(tmp_path: Path) -> None: @@ -21,57 +21,41 @@ def test_EdgeCentricFCSpheres(tmp_path: Path) -> None: The path to the test directory. """ - # get a dataset - ni_data = datasets.fetch_spm_auditory(subject_id="sub001") - fmri_img = image.concat_imgs(ni_data.func) # type: ignore + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + marker = EdgeCentricFCSpheres( + coords="DMNBuckner", radius=5.0, cor_method="correlation" + ) + # Check correct output + assert marker.get_output_type("BOLD") == "matrix" - efc = EdgeCentricFCSpheres( - coords="DMNBuckner", radius=5.0, cor_method="correlation" - ) - all_out = efc.fit_transform( - {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - ) + # Fit-transform the data + edge_fc = marker.fit_transform(element_data) + edge_fc_bold = edge_fc["BOLD"] - out = all_out["BOLD"] + # There are six DMNBuckner coordinates, so + # for 6 ROIs we should get (6 * (6 -1) / 2) edges in the ETS + n_edges = int(6 * (6 - 1) / 2) + assert "data" in edge_fc_bold + assert "row_names" in edge_fc_bold + assert "col_names" in edge_fc_bold + assert edge_fc_bold["data"].shape == (n_edges, n_edges) + assert len(set(edge_fc_bold["row_names"])) == n_edges + assert len(set(edge_fc_bold["col_names"])) == n_edges - # There are six DMNBuckner coordinates, so - # for 6 ROIs we should get (6 * (6 -1) / 2) edges in the ETS - n_edges = int(6 * (6 - 1) / 2) - assert "data" in out - assert "row_names" in out - assert "col_names" in out - assert out["data"].shape[0] == n_edges - assert out["data"].shape[1] == n_edges - assert len(set(out["row_names"])) == n_edges - assert len(set(out["col_names"])) == n_edges - - # check correct output - assert efc.get_output_type("BOLD") == "matrix" - - # Check empirical correlation method parameters - efc = EdgeCentricFCSpheres( - coords="DMNBuckner", - radius=5.0, - cor_method="correlation", - cor_method_params={"empirical": True}, - ) - - meta = { - "element": {"subject": "sub001"}, - "dependencies": {"nilearn"}, - } - all_out = efc.fit_transform( - {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - ) - - uri = tmp_path / "test_fc_parcellation.sqlite" - # Single storage, must be the uri - storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") - meta = {"element": {"subject": "test"}, "dependencies": {"numpy"}} - input = {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - all_out = efc.fit_transform(input, storage=storage) - - features = storage.list_features() - assert any( - x["name"] == "BOLD_EdgeCentricFCSpheres" for x in features.values() - ) + # Check empirical correlation method parameters + marker = EdgeCentricFCSpheres( + coords="DMNBuckner", + radius=5.0, + cor_method="correlation", + cor_method_params={"empirical": True}, + ) + # Store + storage = SQLiteFeatureStorage( + uri=tmp_path / "test_edge_fc_spheres.sqlite", upsert="ignore" + ) + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "BOLD_EdgeCentricFCSpheres" for x in features.values() + ) -- 2.52.0 From 959ab616d1b8f5140e10723ca04fc31bb9dd44a9 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:05:11 +0100 Subject: [PATCH 25/47] update: improve tests for test_functional_connectivity_parcels.py --- .../test_functional_connectivity_parcels.py | 126 ++++++++---------- 1 file changed, 57 insertions(+), 69 deletions(-) diff --git a/junifer/markers/functional_connectivity/tests/test_functional_connectivity_parcels.py b/junifer/markers/functional_connectivity/tests/test_functional_connectivity_parcels.py index 560635f26..01036507b 100644 --- a/junifer/markers/functional_connectivity/tests/test_functional_connectivity_parcels.py +++ b/junifer/markers/functional_connectivity/tests/test_functional_connectivity_parcels.py @@ -2,20 +2,22 @@ # Authors: Amir Omidvarnia # Kaustubh R. Patil +# Synchon Mandal # License: AGPL from pathlib import Path -from nilearn import datasets, image from nilearn.connectome import ConnectivityMeasure from nilearn.maskers import NiftiLabelsMasker -from numpy.testing import assert_array_almost_equal, assert_array_equal +from numpy.testing import assert_array_almost_equal +from junifer.data import get_parcellation +from junifer.datareader import DefaultDataReader from junifer.markers.functional_connectivity import ( FunctionalConnectivityParcels, ) -from junifer.markers.parcel_aggregation import ParcelAggregation from junifer.storage import SQLiteFeatureStorage +from junifer.testing.datagrabbers import PartlyCloudyTestingDataGrabber def test_FunctionalConnectivityParcels(tmp_path: Path) -> None: @@ -27,74 +29,60 @@ def test_FunctionalConnectivityParcels(tmp_path: Path) -> None: The path to the test directory. """ - # get a dataset - ni_data = datasets.fetch_spm_auditory(subject_id="sub001") - fmri_img = image.concat_imgs(ni_data.func) # type: ignore + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + marker = FunctionalConnectivityParcels( + parcellation="TianxS1x3TxMNInonlinear2009cAsym" + ) + # Check correct output + assert marker.get_output_type("BOLD") == "matrix" - fc = FunctionalConnectivityParcels(parcellation="Schaefer100x7") - all_out = fc.fit_transform( - {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - ) + # Fit-transform the data + fc = marker.fit_transform(element_data) + fc_bold = fc["BOLD"] - out = all_out["BOLD"] + assert "data" in fc_bold + assert "row_names" in fc_bold + assert "col_names" in fc_bold + assert fc_bold["data"].shape == (16, 16) + assert len(set(fc_bold["row_names"])) == 16 + assert len(set(fc_bold["col_names"])) == 16 - assert "data" in out - assert "row_names" in out - assert "col_names" in out - assert out["data"].shape[0] == 100 - assert out["data"].shape[1] == 100 - assert len(set(out["row_names"])) == 100 - assert len(set(out["col_names"])) == 100 + # Compare with nilearn + # Load testing parcellation for the target data + testing_parcellation, _ = get_parcellation( + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=element_data["BOLD"], + ) + # Extract timeseries + nifti_labels_masker = NiftiLabelsMasker( + labels_img=testing_parcellation, standardize=False + ) + extracted_timeseries = nifti_labels_masker.fit_transform( + element_data["BOLD"]["data"] + ) + # Compute the connectivity measure + connectivity_measure = ConnectivityMeasure( + kind="covariance" + ).fit_transform([extracted_timeseries])[0] - # get the timeseries using pa - pa = ParcelAggregation( - parcellation="Schaefer100x7", method="mean", on="BOLD" - ) - meta = { - "element": {"subject": "sub001"}, - "dependencies": {"nilearn"}, - } - ts = pa.compute({"data": fmri_img, "meta": meta, "space": "MNI"}) + # Check that FC are almost equal + assert_array_almost_equal( + connectivity_measure, fc_bold["data"], decimal=3 + ) - # compare with nilearn - # Get the testing parcellation (for nilearn) - parcellation = datasets.fetch_atlas_schaefer_2018( - n_rois=100, yeo_networks=7, resolution_mm=2 - ) - masker = NiftiLabelsMasker( - labels_img=parcellation["maps"], standardize=False - ) - ts_ni = masker.fit_transform(fmri_img) - - # check the TS are almost equal - assert_array_equal(ts_ni, ts["data"]) - - # Check that FC are almost equal - cm = ConnectivityMeasure(kind="covariance") - out_ni = cm.fit_transform([ts_ni])[0] - assert_array_almost_equal(out_ni, out["data"], decimal=3) - - # check correct output - assert fc.get_output_type("BOLD") == "matrix" - - # Check empirical correlation method parameters - fc = FunctionalConnectivityParcels( - parcellation="Schaefer100x7", cor_method_params={"empirical": True} - ) - - all_out = fc.fit_transform( - {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - ) - - uri = tmp_path / "test_fc_parcellation.sqlite" - # Single storage, must be the uri - storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") - meta = {"element": {"subject": "test"}, "dependencies": {"numpy"}} - input = {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - all_out = fc.fit_transform(input, storage=storage) - - features = storage.list_features() - assert any( - x["name"] == "BOLD_FunctionalConnectivityParcels" - for x in features.values() - ) + # Check empirical correlation method parameters + marker = FunctionalConnectivityParcels( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + cor_method_params={"empirical": True}, + ) + # Store + storage = SQLiteFeatureStorage( + uri=tmp_path / "test_fc_parcels.sqlite", upsert="ignore" + ) + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "BOLD_FunctionalConnectivityParcels" + for x in features.values() + ) -- 2.52.0 From 02237def5b9a4c98e56722135c06d5d867d4f12d Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:05:44 +0100 Subject: [PATCH 26/47] update: improve tests for test_functional_connectivity_spheres.py --- .../test_functional_connectivity_spheres.py | 170 +++++++++--------- 1 file changed, 90 insertions(+), 80 deletions(-) diff --git a/junifer/markers/functional_connectivity/tests/test_functional_connectivity_spheres.py b/junifer/markers/functional_connectivity/tests/test_functional_connectivity_spheres.py index 6f695c787..fe75c0503 100644 --- a/junifer/markers/functional_connectivity/tests/test_functional_connectivity_spheres.py +++ b/junifer/markers/functional_connectivity/tests/test_functional_connectivity_spheres.py @@ -3,21 +3,24 @@ # Authors: Amir Omidvarnia # Kaustubh R. Patil # Federico Raimondo +# Synchon Mandal # License: AGPL from pathlib import Path import pytest -from nilearn import datasets, image from nilearn.connectome import ConnectivityMeasure +from nilearn.maskers import NiftiSpheresMasker from numpy.testing import assert_array_almost_equal from sklearn.covariance import EmpiricalCovariance +from junifer.data import get_coordinates +from junifer.datareader import DefaultDataReader from junifer.markers.functional_connectivity import ( FunctionalConnectivitySpheres, ) -from junifer.markers.sphere_aggregation import SphereAggregation from junifer.storage import SQLiteFeatureStorage +from junifer.testing.datagrabbers import SPMAuditoryTestingDataGrabber def test_FunctionalConnectivitySpheres(tmp_path: Path) -> None: @@ -29,56 +32,57 @@ def test_FunctionalConnectivitySpheres(tmp_path: Path) -> None: The path to the test directory. """ - # get a dataset - ni_data = datasets.fetch_spm_auditory(subject_id="sub001") - fmri_img = image.concat_imgs(ni_data.func) # type: ignore + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + marker = FunctionalConnectivitySpheres( + coords="DMNBuckner", radius=5.0, cor_method="correlation" + ) + # Check correct output + assert marker.get_output_type("BOLD") == "matrix" - fc = FunctionalConnectivitySpheres( - coords="DMNBuckner", radius=5.0, cor_method="correlation" - ) - all_out = fc.fit_transform( - {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - ) + # Fit-transform the data + fc = marker.fit_transform(element_data) + fc_bold = fc["BOLD"] - out = all_out["BOLD"] + assert "data" in fc_bold + assert "row_names" in fc_bold + assert "col_names" in fc_bold + assert fc_bold["data"].shape == (6, 6) + assert len(set(fc_bold["row_names"])) == 6 + assert len(set(fc_bold["col_names"])) == 6 - assert "data" in out - assert "row_names" in out - assert "col_names" in out - assert out["data"].shape[0] == 6 - assert out["data"].shape[1] == 6 - assert len(set(out["row_names"])) == 6 - assert len(set(out["col_names"])) == 6 + # Compare with nilearn + # Load testing coordinates for the target data + testing_coords, _ = get_coordinates( + coords="DMNBuckner", target_data=element_data["BOLD"] + ) + # Extract timeseries + nifti_spheres_masker = NiftiSpheresMasker( + seeds=testing_coords, radius=5.0 + ) + extracted_timeseries = nifti_spheres_masker.fit_transform( + element_data["BOLD"]["data"] + ) + # Compute the connectivity measure + connectivity_measure = ConnectivityMeasure( + kind="correlation" + ).fit_transform([extracted_timeseries])[0] - # get the timeseries using sa - sa = SphereAggregation( - coords="DMNBuckner", radius=5.0, method="mean", on="BOLD" - ) - ts = sa.compute({"data": fmri_img, "meta": {}, "space": "MNI"}) + # Check that FC are almost equal + assert_array_almost_equal( + connectivity_measure, fc_bold["data"], decimal=3 + ) - # Check that FC are almost equal when using nileran - cm = ConnectivityMeasure(kind="correlation") - out_ni = cm.fit_transform([ts["data"]])[0] - assert_array_almost_equal(out_ni, out["data"], decimal=3) - - # check correct output - assert fc.get_output_type("BOLD") == "matrix" - - uri = tmp_path / "test_fc_parcel.sqlite" - # Single storage, must be the uri - storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") - meta = { - "element": {"subject": "test"}, - "dependencies": {"numpy", "nilearn"}, - } - input = {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - all_out = fc.fit_transform(input, storage=storage) - - features = storage.list_features() - assert any( - x["name"] == "BOLD_FunctionalConnectivitySpheres" - for x in features.values() - ) + # Store + storage = SQLiteFeatureStorage( + uri=tmp_path / "test_fc_spheres.sqlite", upsert="ignore" + ) + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "BOLD_FunctionalConnectivitySpheres" + for x in features.values() + ) def test_FunctionalConnectivitySpheres_empirical(tmp_path: Path) -> None: @@ -90,43 +94,49 @@ def test_FunctionalConnectivitySpheres_empirical(tmp_path: Path) -> None: The path to the test directory. """ + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + marker = FunctionalConnectivitySpheres( + coords="DMNBuckner", + radius=5.0, + cor_method="correlation", + cor_method_params={"empirical": True}, + ) + # Check correct output + assert marker.get_output_type("BOLD") == "matrix" - # get a dataset - ni_data = datasets.fetch_spm_auditory(subject_id="sub001") - fmri_img = image.concat_imgs(ni_data.func) # type: ignore + # Fit-transform the data + fc = marker.fit_transform(element_data) + fc_bold = fc["BOLD"] - fc = FunctionalConnectivitySpheres( - coords="DMNBuckner", - radius=5.0, - cor_method="correlation", - cor_method_params={"empirical": True}, - ) - all_out = fc.fit_transform( - {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - ) + assert "data" in fc_bold + assert "row_names" in fc_bold + assert "col_names" in fc_bold + assert fc_bold["data"].shape == (6, 6) + assert len(set(fc_bold["row_names"])) == 6 + assert len(set(fc_bold["col_names"])) == 6 - out = all_out["BOLD"] + # Compare with nilearn + # Load testing coordinates for the target data + testing_coords, _ = get_coordinates( + coords="DMNBuckner", target_data=element_data["BOLD"] + ) + # Extract timeseries + nifti_spheres_masker = NiftiSpheresMasker( + seeds=testing_coords, radius=5.0 + ) + extracted_timeseries = nifti_spheres_masker.fit_transform( + element_data["BOLD"]["data"] + ) + # Compute the connectivity measure + connectivity_measure = ConnectivityMeasure( + cov_estimator=EmpiricalCovariance(), kind="correlation" # type: ignore + ).fit_transform([extracted_timeseries])[0] - assert "data" in out - assert "row_names" in out - assert "col_names" in out - assert out["data"].shape[0] == 6 - assert out["data"].shape[1] == 6 - assert len(set(out["row_names"])) == 6 - assert len(set(out["col_names"])) == 6 - - # get the timeseries using sa - sa = SphereAggregation( - coords="DMNBuckner", radius=5.0, method="mean", on="BOLD" - ) - ts = sa.compute({"data": fmri_img, "space": "MNI"}) - - # Check that FC are almost equal when using nileran - cm = ConnectivityMeasure( - cov_estimator=EmpiricalCovariance(), kind="correlation" # type: ignore - ) - out_ni = cm.fit_transform([ts["data"]])[0] - assert_array_almost_equal(out_ni, out["data"], decimal=3) + # Check that FC are almost equal + assert_array_almost_equal( + connectivity_measure, fc_bold["data"], decimal=3 + ) def test_FunctionalConnectivitySpheres_error() -> None: -- 2.52.0 From b18b5daae60b4f0dfdf5d4d9cf42f5fcf1be06ab Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:06:15 +0100 Subject: [PATCH 27/47] update: improve tests for test_reho_parcels.py --- .../markers/reho/tests/test_reho_parcels.py | 26 ++++++++++++------- 1 file changed, 16 insertions(+), 10 deletions(-) diff --git a/junifer/markers/reho/tests/test_reho_parcels.py b/junifer/markers/reho/tests/test_reho_parcels.py index 531fe4fc8..94c854d49 100644 --- a/junifer/markers/reho/tests/test_reho_parcels.py +++ b/junifer/markers/reho/tests/test_reho_parcels.py @@ -12,12 +12,12 @@ import scipy as sp from junifer.datareader import DefaultDataReader from junifer.markers import ReHoParcels from junifer.pipeline import WorkDirManager -from junifer.pipeline.utils import _check_afni +from junifer.pipeline.utils import _check_afni, _check_ants from junifer.storage import SQLiteFeatureStorage -from junifer.testing.datagrabbers import SPMAuditoryTestingDataGrabber - - -PARCELLATION = "Schaefer100x7" +from junifer.testing.datagrabbers import ( + PartlyCloudyTestingDataGrabber, + SPMAuditoryTestingDataGrabber, +) def test_ReHoParcels(caplog: pytest.LogCaptureFixture, tmp_path: Path) -> None: @@ -32,13 +32,16 @@ def test_ReHoParcels(caplog: pytest.LogCaptureFixture, tmp_path: Path) -> None: """ with caplog.at_level(logging.DEBUG): - with SPMAuditoryTestingDataGrabber() as dg: - element_data = DefaultDataReader().fit_transform(dg["sub001"]) + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) # Update workdir to current test's tmp_path WorkDirManager().workdir = tmp_path # Initialize marker - marker = ReHoParcels(parcellation=PARCELLATION, using="junifer") + marker = ReHoParcels( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + using="junifer", + ) # Fit transform marker on data output = marker.fit_transform(element_data) @@ -72,6 +75,9 @@ def test_ReHoParcels(caplog: pytest.LogCaptureFixture, tmp_path: Path) -> None: assert "Creating cache" not in caplog.text +@pytest.mark.skipif( + _check_ants() is False, reason="requires ANTs to be in PATH" +) @pytest.mark.skipif( _check_afni() is False, reason="requires AFNI to be in PATH" ) @@ -91,7 +97,7 @@ def test_ReHoParcels_comparison(tmp_path: Path) -> None: # Initialize marker junifer_marker = ReHoParcels( - parcellation=PARCELLATION, using="junifer" + parcellation="Schaefer100x7", using="junifer" ) # Fit transform marker on data junifer_output = junifer_marker.fit_transform(element_data) @@ -99,7 +105,7 @@ def test_ReHoParcels_comparison(tmp_path: Path) -> None: junifer_output_bold = junifer_output["BOLD"] # Initialize marker - afni_marker = ReHoParcels(parcellation=PARCELLATION, using="afni") + afni_marker = ReHoParcels(parcellation="Schaefer100x7", using="afni") # Fit transform marker on data afni_output = afni_marker.fit_transform(element_data) # Get BOLD output -- 2.52.0 From 031d2026764ed1023a36518fc04744d976fde6ab Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:07:02 +0100 Subject: [PATCH 28/47] update: improve tests for test_temporal_snr_parcels.py --- .../tests/test_temporal_snr_parcels.py | 75 ++++++++++--------- 1 file changed, 38 insertions(+), 37 deletions(-) diff --git a/junifer/markers/temporal_snr/tests/test_temporal_snr_parcels.py b/junifer/markers/temporal_snr/tests/test_temporal_snr_parcels.py index e888a74ff..e696e7220 100644 --- a/junifer/markers/temporal_snr/tests/test_temporal_snr_parcels.py +++ b/junifer/markers/temporal_snr/tests/test_temporal_snr_parcels.py @@ -1,18 +1,39 @@ """Provide tests for temporal signal-to-noise ratio using parcellation.""" # Authors: Leonard Sasse +# Synchon Mandal # License: AGPL from pathlib import Path -from nilearn import datasets, image - +from junifer.datareader import DefaultDataReader from junifer.markers.temporal_snr import TemporalSNRParcels -from junifer.storage import SQLiteFeatureStorage +from junifer.storage import HDF5FeatureStorage +from junifer.testing.datagrabbers import PartlyCloudyTestingDataGrabber -def test_TemporalSNRParcels(tmp_path: Path) -> None: - """Test TemporalSNRParcels. +def test_TemporalSNRParcels_computation() -> None: + """Test TemporalSNRParcels fit-transform.""" + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + marker = TemporalSNRParcels( + parcellation="TianxS1x3TxMNInonlinear2009cAsym" + ) + # Check correct output + assert marker.get_output_type("BOLD") == "vector" + + # Fit-transform the data + tsnr_parcels = marker.fit_transform(element_data) + tsnr_parcels_bold = tsnr_parcels["BOLD"] + + assert "data" in tsnr_parcels_bold + assert "col_names" in tsnr_parcels_bold + assert tsnr_parcels_bold["data"].shape == (1, 16) + assert len(set(tsnr_parcels_bold["col_names"])) == 16 + + +def test_TemporalSNRParcels_storage(tmp_path: Path) -> None: + """Test TemporalSNRParcels store. Parameters ---------- @@ -20,35 +41,15 @@ def test_TemporalSNRParcels(tmp_path: Path) -> None: The path to the test directory. """ - # get a dataset - ni_data = datasets.fetch_spm_auditory(subject_id="sub001") - fmri_img = image.concat_imgs(ni_data.func) # type: ignore - - tsnr_parcels = TemporalSNRParcels(parcellation="Schaefer100x7") - all_out = tsnr_parcels.fit_transform( - {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - ) - - out = all_out["BOLD"] - - assert "data" in out - assert "col_names" in out - - assert out["data"].shape[0] == 1 - assert out["data"].shape[1] == 100 - assert len(set(out["col_names"])) == 100 - - # check correct output - assert tsnr_parcels.get_output_type("BOLD") == "vector" - - uri = tmp_path / "test_tsnr_parcellation.sqlite" - # Single storage, must be the uri - storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") - meta = {"element": {"subject": "test"}, "dependencies": {"numpy"}} - input = {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - all_out = tsnr_parcels.fit_transform(input, storage=storage) - - features = storage.list_features() - assert any( - x["name"] == "BOLD_TemporalSNRParcels" for x in features.values() - ) + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + marker = TemporalSNRParcels( + parcellation="TianxS1x3TxMNInonlinear2009cAsym" + ) + # Store + storage = HDF5FeatureStorage(tmp_path / "test_tsnr_parcels.hdf5") + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "BOLD_TemporalSNRParcels" for x in features.values() + ) -- 2.52.0 From 678b38901b46abbadcd0da940539d8e3fe47a814 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:07:30 +0100 Subject: [PATCH 29/47] update: improve tests for test_temporal_snr_spheres.py --- .../tests/test_temporal_snr_spheres.py | 72 +++++++++---------- 1 file changed, 34 insertions(+), 38 deletions(-) diff --git a/junifer/markers/temporal_snr/tests/test_temporal_snr_spheres.py b/junifer/markers/temporal_snr/tests/test_temporal_snr_spheres.py index 70ad1752e..c25d2a0a5 100644 --- a/junifer/markers/temporal_snr/tests/test_temporal_snr_spheres.py +++ b/junifer/markers/temporal_snr/tests/test_temporal_snr_spheres.py @@ -1,19 +1,39 @@ """Provide tests for temporal signal-to-noise spheres.""" # Authors: Leonard Sasse +# Synchon Mandal # License: AGPL from pathlib import Path import pytest -from nilearn import datasets, image +from junifer.datareader import DefaultDataReader from junifer.markers.temporal_snr import TemporalSNRSpheres -from junifer.storage import SQLiteFeatureStorage +from junifer.storage import HDF5FeatureStorage +from junifer.testing.datagrabbers import SPMAuditoryTestingDataGrabber -def test_TemporalSNRSpheres(tmp_path: Path) -> None: - """Test TemporalSNRSpheres. +def test_TemporalSNRSpheres_computation() -> None: + """Test TemporalSNRSpheres fit-transform.""" + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + marker = TemporalSNRSpheres(coords="DMNBuckner", radius=5.0) + # Check correct output + assert marker.get_output_type("BOLD") == "vector" + + # Fit-transform the data + tsnr_spheres = marker.fit_transform(element_data) + tsnr_spheres_bold = tsnr_spheres["BOLD"] + + assert "data" in tsnr_spheres_bold + assert "col_names" in tsnr_spheres_bold + assert tsnr_spheres_bold["data"].shape == (1, 6) + assert len(set(tsnr_spheres_bold["col_names"])) == 6 + + +def test_TemporalSNRSpheres_storage(tmp_path: Path) -> None: + """Test TemporalSNRSpheres store. Parameters ---------- @@ -21,40 +41,16 @@ def test_TemporalSNRSpheres(tmp_path: Path) -> None: The path to the test directory. """ - # get a dataset - ni_data = datasets.fetch_spm_auditory(subject_id="sub001") - fmri_img = image.concat_imgs(ni_data.func) # type: ignore - - tsnr_spheres = TemporalSNRSpheres(coords="DMNBuckner", radius=5.0) - all_out = tsnr_spheres.fit_transform( - {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - ) - - out = all_out["BOLD"] - - assert "data" in out - assert "col_names" in out - assert out["data"].shape[0] == 1 - assert out["data"].shape[1] == 6 - assert len(set(out["col_names"])) == 6 - - # check correct output - assert tsnr_spheres.get_output_type("BOLD") == "vector" - - uri = tmp_path / "test_tsnr_coords.sqlite" - # Single storage, must be the uri - storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") - meta = { - "element": {"subject": "test"}, - "dependencies": {"numpy", "nilearn"}, - } - input = {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - all_out = tsnr_spheres.fit_transform(input, storage=storage) - - features = storage.list_features() - assert any( - x["name"] == "BOLD_TemporalSNRSpheres" for x in features.values() - ) + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + marker = TemporalSNRSpheres(coords="DMNBuckner", radius=5.0) + # Store + storage = HDF5FeatureStorage(tmp_path / "test_tsnr_spheres.hdf5") + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "BOLD_TemporalSNRSpheres" for x in features.values() + ) def test_TemporalSNRSpheres_error() -> None: -- 2.52.0 From 27ce405488c0b4e651efef3e2265d92ce063aa62 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:07:59 +0100 Subject: [PATCH 30/47] update: improve tests for test_collection.py --- junifer/markers/tests/test_collection.py | 75 ++++++++++++------------ 1 file changed, 38 insertions(+), 37 deletions(-) diff --git a/junifer/markers/tests/test_collection.py b/junifer/markers/tests/test_collection.py index 5f978ee70..c63065f2f 100644 --- a/junifer/markers/tests/test_collection.py +++ b/junifer/markers/tests/test_collection.py @@ -19,7 +19,6 @@ from junifer.pipeline import PipelineStepMixin from junifer.preprocess import fMRIPrepConfoundRemover from junifer.storage import SQLiteFeatureStorage from junifer.testing.datagrabbers import ( - OasisVBMTestingDataGrabber, PartlyCloudyTestingDataGrabber, ) @@ -46,20 +45,20 @@ def test_marker_collection() -> None: """Test MarkerCollection.""" markers = [ ParcelAggregation( - parcellation="Schaefer100x7", + parcellation="TianxS2x3TxMNInonlinear2009cAsym", method="mean", - name="gmd_schaefer100x7_mean", + name="tian_mean", ), ParcelAggregation( - parcellation="Schaefer100x7", + parcellation="TianxS2x3TxMNInonlinear2009cAsym", method="std", - name="gmd_schaefer100x7_std", + name="tian_std", ), ParcelAggregation( - parcellation="Schaefer100x7", + parcellation="TianxS2x3TxMNInonlinear2009cAsym", method="trim_mean", method_params={"proportiontocut": 0.1}, - name="gmd_schaefer100x7_trim_mean90", + name="tian_trim_mean90", ), ] mc = MarkerCollection(markers=markers) # type: ignore @@ -69,7 +68,7 @@ def test_marker_collection() -> None: assert isinstance(mc._datareader, DefaultDataReader) # Create testing datagrabber - dg = OasisVBMTestingDataGrabber() + dg = PartlyCloudyTestingDataGrabber() mc.validate(dg) with dg: @@ -78,17 +77,17 @@ def test_marker_collection() -> None: assert out is not None assert isinstance(out, dict) assert len(out) == 3 - assert "gmd_schaefer100x7_mean" in out - assert "gmd_schaefer100x7_std" in out - assert "gmd_schaefer100x7_trim_mean90" in out + assert "tian_mean" in out + assert "tian_std" in out + assert "tian_trim_mean90" in out for t_marker in markers: t_name = t_marker.name - assert "VBM_GM" in out[t_name] - t_vbm = out[t_name]["VBM_GM"] - assert "data" in t_vbm - assert "col_names" in t_vbm - assert "meta" in t_vbm + assert "BOLD" in out[t_name] + t_bold = out[t_name]["BOLD"] + assert "data" in t_bold + assert "col_names" in t_bold + assert "meta" in t_bold # Test preprocessing class BypassPreprocessing(PipelineStepMixin): @@ -108,7 +107,7 @@ def test_marker_collection() -> None: for t_marker in markers: t_name = t_marker.name assert_array_equal( - out[t_name]["VBM_GM"]["data"], out2[t_name]["VBM_GM"]["data"] + out[t_name]["BOLD"]["data"], out2[t_name]["BOLD"]["data"] ) @@ -151,27 +150,28 @@ def test_marker_collection_storage(tmp_path: Path) -> None: """ markers = [ ParcelAggregation( - parcellation="Schaefer100x7", + parcellation="TianxS2x3TxMNInonlinear2009cAsym", method="mean", - name="gmd_schaefer100x7_mean", + name="tian_mean", ), ParcelAggregation( - parcellation="Schaefer100x7", + parcellation="TianxS2x3TxMNInonlinear2009cAsym", method="std", - name="gmd_schaefer100x7_std", + name="tian_std", ), ParcelAggregation( - parcellation="Schaefer100x7", + parcellation="TianxS2x3TxMNInonlinear2009cAsym", method="trim_mean", method_params={"proportiontocut": 0.1}, - name="gmd_schaefer100x7_trim_mean90", + name="tian_trim_mean90", ), ] - # Test storage - dg = OasisVBMTestingDataGrabber() - - uri = tmp_path / "test_marker_collection_storage.sqlite" - storage = SQLiteFeatureStorage(uri=uri) + # Setup datagrabber + dg = PartlyCloudyTestingDataGrabber() + # Setup storage + storage = SQLiteFeatureStorage( + tmp_path / "test_marker_collection_storage.sqlite" + ) mc = MarkerCollection( markers=markers, # type: ignore storage=storage, @@ -197,23 +197,24 @@ def test_marker_collection_storage(tmp_path: Path) -> None: features = storage.list_features() assert len(features) == 3 + feature_md5 = next(iter(features.keys())) t_feature = storage.read_df(feature_md5=feature_md5) - fname = "gmd_schaefer100x7_mean" - t_data = out[fname]["VBM_GM"]["data"] # type: ignore - cols = out[fname]["VBM_GM"]["col_names"] # type: ignore + fname = "tian_mean" + t_data = out[fname]["BOLD"]["data"] # type: ignore + cols = out[fname]["BOLD"]["col_names"] # type: ignore assert_array_equal(t_feature[cols].values, t_data) # type: ignore feature_md5 = list(features.keys())[1] t_feature = storage.read_df(feature_md5=feature_md5) - fname = "gmd_schaefer100x7_std" - t_data = out[fname]["VBM_GM"]["data"] # type: ignore - cols = out[fname]["VBM_GM"]["col_names"] # type: ignore + fname = "tian_std" + t_data = out[fname]["BOLD"]["data"] # type: ignore + cols = out[fname]["BOLD"]["col_names"] # type: ignore assert_array_equal(t_feature[cols].values, t_data) # type: ignore feature_md5 = list(features.keys())[2] t_feature = storage.read_df(feature_md5=feature_md5) - fname = "gmd_schaefer100x7_trim_mean90" - t_data = out[fname]["VBM_GM"]["data"] # type: ignore - cols = out[fname]["VBM_GM"]["col_names"] # type: ignore + fname = "tian_trim_mean90" + t_data = out[fname]["BOLD"]["data"] # type: ignore + cols = out[fname]["BOLD"]["col_names"] # type: ignore assert_array_equal(t_feature[cols].values, t_data) # type: ignore -- 2.52.0 From 9d68b444bd5f8d523e4510555b1bdfa9bf9abe6f Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:08:22 +0100 Subject: [PATCH 31/47] update: improve tests for test_ets_rss.py --- junifer/markers/tests/test_ets_rss.py | 68 +++++++++++---------------- 1 file changed, 28 insertions(+), 40 deletions(-) diff --git a/junifer/markers/tests/test_ets_rss.py b/junifer/markers/tests/test_ets_rss.py index 04c102893..10d9f8de9 100644 --- a/junifer/markers/tests/test_ets_rss.py +++ b/junifer/markers/tests/test_ets_rss.py @@ -8,52 +8,47 @@ from pathlib import Path -from nilearn import image from nilearn.maskers import NiftiLabelsMasker -from junifer.data import load_parcellation +from junifer.data import get_parcellation +from junifer.datareader import DefaultDataReader from junifer.markers.ets_rss import RSSETSMarker from junifer.storage import SQLiteFeatureStorage -from junifer.testing.datagrabbers import SPMAuditoryTestingDataGrabber +from junifer.testing.datagrabbers import PartlyCloudyTestingDataGrabber # Set parcellation -PARCELLATION = "Schaefer100x17" +PARCELLATION = "TianxS1x3TxMNInonlinear2009cAsym" def test_compute() -> None: """Test RSS ETS compute().""" - with SPMAuditoryTestingDataGrabber() as dg: - # Fetch element - out = dg["sub001"] - # Load BOLD image - niimg = image.load_img(str(out["BOLD"]["path"].absolute())) - # Create input data - input_dict = { - "data": niimg, - "path": out["BOLD"]["path"], - "space": "MNI", - } + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) # Compute the RSSETSMarker - ets_rss_marker = RSSETSMarker(parcellation=PARCELLATION) - new_out = ets_rss_marker.compute(input_dict) + marker = RSSETSMarker(parcellation=PARCELLATION) + rss_ets = marker.compute(element_data["BOLD"]) - # Load parcellation - test_parcellation, _, _, _ = load_parcellation(PARCELLATION) - # Compute the NiftiLabelsMasker - test_masker = NiftiLabelsMasker(test_parcellation) - test_ts = test_masker.fit_transform(niimg) + # Compare with nilearn + # Load testing parcellation + test_parcellation, _ = get_parcellation( + parcellation=[PARCELLATION], + target_data=element_data["BOLD"], + ) + # Extract timeseries + nifti_labels_masker = NiftiLabelsMasker(labels_img=test_parcellation) + extacted_timeseries = nifti_labels_masker.fit_transform( + element_data["BOLD"]["data"] + ) # Assert the dimension of timeseries - n_time, _ = test_ts.shape - assert n_time == len(new_out["data"]) + assert extacted_timeseries.shape[0] == len(rss_ets["data"]) def test_get_output_type() -> None: """Test RSS ETS get_output_type().""" - ets_rss_marker = RSSETSMarker(parcellation=PARCELLATION) - input_ = "BOLD" - output = ets_rss_marker.get_output_type(input_) - assert output == "timeseries" + assert "timeseries" == RSSETSMarker( + parcellation=PARCELLATION + ).get_output_type("BOLD") def test_store(tmp_path: Path) -> None: @@ -65,20 +60,13 @@ def test_store(tmp_path: Path) -> None: The path to the test directory. """ - with SPMAuditoryTestingDataGrabber() as dg: - # Fetch element - elem = dg["sub001"] - # Load BOLD image - niimg = image.load_img(str(elem["BOLD"]["path"].absolute())) - elem["BOLD"]["data"] = niimg + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) # Compute the RSSETSMarker - ets_rss_marker = RSSETSMarker(parcellation=PARCELLATION) + marker = RSSETSMarker(parcellation=PARCELLATION) # Create storage - storage = SQLiteFeatureStorage( - uri=str((tmp_path / "test.sqlite").absolute()) - ) + storage = SQLiteFeatureStorage(tmp_path / "test_rss_ets.sqlite") # Store - ets_rss_marker.fit_transform(input=elem, storage=storage) - + marker.fit_transform(input=element_data, storage=storage) features = storage.list_features() assert any(x["name"] == "BOLD_RSSETSMarker" for x in features.values()) -- 2.52.0 From efe53eae26c4190c026b64ce310d217ac79b0183 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 18:08:43 +0100 Subject: [PATCH 32/47] update: improve tests for test_sphere_aggregation.py --- .../markers/tests/test_sphere_aggregation.py | 352 ++++++++++-------- 1 file changed, 199 insertions(+), 153 deletions(-) diff --git a/junifer/markers/tests/test_sphere_aggregation.py b/junifer/markers/tests/test_sphere_aggregation.py index 04b9de474..7416764dc 100644 --- a/junifer/markers/tests/test_sphere_aggregation.py +++ b/junifer/markers/tests/test_sphere_aggregation.py @@ -1,22 +1,23 @@ """Provide tests for sphere aggregation.""" # Authors: Federico Raimondo +# Synchon Mandal # License: AGPL -import typing from pathlib import Path -from typing import Dict -import nibabel as nib import pytest -from nilearn import datasets -from nilearn.image import concat_imgs from nilearn.maskers import NiftiSpheresMasker from numpy.testing import assert_array_equal -from junifer.data import load_coordinates, load_mask +from junifer.data import get_coordinates, get_mask +from junifer.datareader import DefaultDataReader from junifer.markers.sphere_aggregation import SphereAggregation from junifer.storage import SQLiteFeatureStorage +from junifer.testing.datagrabbers import ( + OasisVBMTestingDataGrabber, + SPMAuditoryTestingDataGrabber, +) # Define common variables @@ -36,51 +37,66 @@ def test_SphereAggregation_input_output() -> None: def test_SphereAggregation_3D() -> None: """Test SphereAggregation object on 3D images.""" - # Get the testing coordinates (for nilearn) - coordinates, _, _ = load_coordinates(COORDS) + with OasisVBMTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + # Create SphereAggregation object + marker = SphereAggregation( + coords=COORDS, method="mean", radius=RADIUS, on="VBM_GM" + ) + sphere_agg_vbm_gm_data = marker.fit_transform(element_data)["VBM_GM"][ + "data" + ] - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) + # Compare with nilearn + # Load testing coordinates + testing_coords, _ = get_coordinates( + coords=COORDS, target_data=element_data["VBM_GM"] + ) + # Extract data + nifti_spheres_masker = NiftiSpheresMasker( + seeds=testing_coords, radius=RADIUS + ) + nifti_spheres_masked_vbm_gm = nifti_spheres_masker.fit_transform( + element_data["VBM_GM"]["data"] + ) - # Create NiftSpheresMasker - nifti_masker = NiftiSpheresMasker(seeds=coordinates, radius=RADIUS) - auto4d = nifti_masker.fit_transform(img) - - # Create SphereAggregation object - marker = SphereAggregation( - coords=COORDS, method="mean", radius=RADIUS, on="VBM_GM" - ) - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - jun_values4d = marker.fit_transform(input)["VBM_GM"]["data"] - - assert jun_values4d.ndim == 2 - assert_array_equal(auto4d.shape, jun_values4d.shape) - assert_array_equal(auto4d, jun_values4d) + assert sphere_agg_vbm_gm_data.ndim == 2 + assert_array_equal( + nifti_spheres_masked_vbm_gm.shape, sphere_agg_vbm_gm_data.shape + ) + assert_array_equal(nifti_spheres_masked_vbm_gm, sphere_agg_vbm_gm_data) def test_SphereAggregation_4D() -> None: """Test SphereAggregation object on 4D images.""" - # Get the testing coordinates (for nilearn) - coordinates, _, _ = load_coordinates(COORDS) + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + # Create SphereAggregation object + marker = SphereAggregation( + coords=COORDS, method="mean", radius=RADIUS, on="BOLD" + ) + sphere_agg_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] - # Get the SPM auditory data - subject_data = datasets.fetch_spm_auditory() - fmri_img = concat_imgs(subject_data.func) # type: ignore + # Compare with nilearn + # Load testing coordinates + testing_coords, _ = get_coordinates( + coords=COORDS, target_data=element_data["BOLD"] + ) + # Extract data + nifti_spheres_masker = NiftiSpheresMasker( + seeds=testing_coords, radius=RADIUS + ) + nifti_spheres_masked_bold = nifti_spheres_masker.fit_transform( + element_data["BOLD"]["data"] + ) - # Create NiftSpheresMasker - nifti_masker = NiftiSpheresMasker(seeds=coordinates, radius=RADIUS) - auto4d = nifti_masker.fit_transform(fmri_img) - - # Create SphereAggregation object - marker = SphereAggregation(coords=COORDS, method="mean", radius=RADIUS) - input = {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - jun_values4d = marker.fit_transform(input)["BOLD"]["data"] - - assert jun_values4d.ndim == 2 - assert_array_equal(auto4d.shape, jun_values4d.shape) - assert_array_equal(auto4d, jun_values4d) + assert sphere_agg_bold_data.ndim == 2 + assert_array_equal( + nifti_spheres_masked_bold.shape, sphere_agg_bold_data.shape + ) + assert_array_equal(nifti_spheres_masked_bold, sphere_agg_bold_data) def test_SphereAggregation_storage(tmp_path: Path) -> None: @@ -92,124 +108,143 @@ def test_SphereAggregation_storage(tmp_path: Path) -> None: The path to the test directory. """ - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) - uri = tmp_path / "test_sphere_storage_3D.sqlite" + # Store 3D + with OasisVBMTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + storage = SQLiteFeatureStorage( + uri=tmp_path / "test_sphere_storage_3D.sqlite", upsert="ignore" + ) + marker = SphereAggregation( + coords=COORDS, method="mean", radius=RADIUS, on="VBM_GM" + ) + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "VBM_GM_SphereAggregation" for x in features.values() + ) - storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") - meta = { - "element": {"subject": "sub-01", "session": "ses-01"}, - "dependencies": {"nilearn", "nibabel"}, - } - input = {"VBM_GM": {"data": img, "meta": meta, "space": "MNI"}} - marker = SphereAggregation( - coords=COORDS, method="mean", radius=RADIUS, on="VBM_GM" - ) - - marker.fit_transform(input, storage=storage) - - features: Dict = typing.cast(Dict, storage.list_features()) - assert any( - x["name"] == "VBM_GM_SphereAggregation" for x in features.values() - ) - - meta = { - "element": {"subject": "sub-01", "session": "ses-01"}, - "dependencies": {"nilearn", "nibabel"}, - } - # Get the SPM auditory data - subject_data = datasets.fetch_spm_auditory() - fmri_img = concat_imgs(subject_data.func) # type: ignore - input = {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - marker = SphereAggregation( - coords=COORDS, method="mean", radius=RADIUS, on="BOLD" - ) - - marker.fit_transform(input, storage=storage) - features: Dict = typing.cast(Dict, storage.list_features()) - assert any( - x["name"] == "BOLD_SphereAggregation" for x in features.values() - ) + # Store 4D + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + storage = SQLiteFeatureStorage( + uri=tmp_path / "test_sphere_storage_4D.sqlite", upsert="ignore" + ) + marker = SphereAggregation( + coords=COORDS, method="mean", radius=RADIUS, on="BOLD" + ) + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "BOLD_SphereAggregation" for x in features.values() + ) def test_SphereAggregation_3D_mask() -> None: """Test SphereAggregation object on 3D images using mask.""" - # Get the testing coordinates (for nilearn) - coordinates, _, _ = load_coordinates(COORDS) + with OasisVBMTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + # Create SphereAggregation object + marker = SphereAggregation( + coords=COORDS, + method="mean", + radius=RADIUS, + on="VBM_GM", + masks="compute_brain_mask", + ) + sphere_agg_vbm_gm_data = marker.fit_transform(element_data)["VBM_GM"][ + "data" + ] - # Get one mask - mask_img, _, _ = load_mask("GM_prob0.2") + # Compare with nilearn + # Load testing coordinates + testing_coords, _ = get_coordinates( + coords=COORDS, target_data=element_data["VBM_GM"] + ) + # Load mask + mask_img = get_mask( + "compute_brain_mask", target_data=element_data["VBM_GM"] + ) + # Extract data + nifti_spheres_masker = NiftiSpheresMasker( + seeds=testing_coords, radius=RADIUS, mask_img=mask_img + ) + nifti_spheres_masked_vbm_agg = nifti_spheres_masker.fit_transform( + element_data["VBM_GM"]["data"] + ) - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) - - # Create NiftSpheresMasker - nifti_masker = NiftiSpheresMasker( - seeds=coordinates, radius=RADIUS, mask_img=mask_img - ) - auto4d = nifti_masker.fit_transform(img) - - # Create SphereAggregation object - marker = SphereAggregation( - coords=COORDS, - method="mean", - radius=RADIUS, - on="VBM_GM", - masks="GM_prob0.2", - ) - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - jun_values4d = marker.fit_transform(input)["VBM_GM"]["data"] - - assert jun_values4d.ndim == 2 - assert_array_equal(auto4d.shape, jun_values4d.shape) - assert_array_equal(auto4d, jun_values4d) + assert sphere_agg_vbm_gm_data.ndim == 2 + assert_array_equal( + nifti_spheres_masked_vbm_agg.shape, + nifti_spheres_masked_vbm_agg.shape, + ) + assert_array_equal( + nifti_spheres_masked_vbm_agg, nifti_spheres_masked_vbm_agg + ) def test_SphereAggregation_4D_agg_time() -> None: """Test SphereAggregation object on 4D images, aggregating time.""" - # Get the testing coordinates (for nilearn) - coordinates, _, _ = load_coordinates(COORDS) + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + # Create SphereAggregation object + marker = SphereAggregation( + coords=COORDS, + method="mean", + radius=RADIUS, + time_method="mean", + on="BOLD", + ) + sphere_agg_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] - # Get the SPM auditory data - subject_data = datasets.fetch_spm_auditory() - fmri_img = concat_imgs(subject_data.func) # type: ignore + # Compare with nilearn + # Load testing coordinates + testing_coords, _ = get_coordinates( + coords=COORDS, target_data=element_data["BOLD"] + ) + # Extract data + nifti_spheres_masker = NiftiSpheresMasker( + seeds=testing_coords, radius=RADIUS + ) + nifti_spheres_masked_bold = nifti_spheres_masker.fit_transform( + element_data["BOLD"]["data"] + ) + nifti_spheres_masked_bold_mean = nifti_spheres_masked_bold.mean(axis=0) - # Create NiftSpheresMasker - nifti_masker = NiftiSpheresMasker(seeds=coordinates, radius=RADIUS) - auto4d = nifti_masker.fit_transform(fmri_img) - auto_mean = auto4d.mean(axis=0) + assert sphere_agg_bold_data.ndim == 1 + assert_array_equal( + nifti_spheres_masked_bold_mean.shape, sphere_agg_bold_data.shape + ) + assert_array_equal( + nifti_spheres_masked_bold_mean, sphere_agg_bold_data + ) - # Create SphereAggregation object - marker = SphereAggregation( - coords=COORDS, method="mean", radius=RADIUS, time_method="mean" - ) - input = {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - jun_values4d = marker.fit_transform(input)["BOLD"]["data"] + # Test picking first time point + nifti_spheres_masked_bold_pick_0 = nifti_spheres_masked_bold[:1, :] + marker = SphereAggregation( + coords=COORDS, + method="mean", + radius=RADIUS, + time_method="select", + time_method_params={"pick": [0]}, + on="BOLD", + ) + sphere_agg_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] - assert jun_values4d.ndim == 1 - assert_array_equal(auto_mean.shape, jun_values4d.shape) - assert_array_equal(auto_mean, jun_values4d) + assert sphere_agg_bold_data.ndim == 2 + assert_array_equal( + nifti_spheres_masked_bold_pick_0.shape, sphere_agg_bold_data.shape + ) + assert_array_equal( + nifti_spheres_masked_bold_pick_0, sphere_agg_bold_data + ) - auto_pick_0 = auto4d[:1, :] - marker = SphereAggregation( - coords=COORDS, - method="mean", - radius=RADIUS, - time_method="select", - time_method_params={"pick": [0]}, - ) - - input = {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - jun_values4d = marker.fit_transform(input)["BOLD"]["data"] - - assert jun_values4d.ndim == 2 - assert_array_equal(auto_pick_0.shape, jun_values4d.shape) - assert_array_equal(auto_pick_0, jun_values4d) +def test_SphereAggregation_errors() -> None: + """Test errors for SphereAggregation.""" with pytest.raises(ValueError, match="can only be used with BOLD data"): SphereAggregation( coords=COORDS, @@ -231,12 +266,23 @@ def test_SphereAggregation_4D_agg_time() -> None: on="VBM_GM", ) - with pytest.warns(RuntimeWarning, match="No time dimension to aggregate"): - input = { - "BOLD": { - "data": fmri_img.slicer[..., 0:1], - "meta": {}, - "space": "MNI", - } - } - marker.fit_transform(input) + +def test_SphereAggregation_warning() -> None: + """Test warning for SphereAggregation.""" + with SPMAuditoryTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub001"]) + with pytest.warns( + RuntimeWarning, match="No time dimension to aggregate" + ): + marker = SphereAggregation( + coords=COORDS, + method="mean", + radius=RADIUS, + time_method="select", + time_method_params={"pick": [0]}, + on="BOLD", + ) + element_data["BOLD"]["data"] = element_data["BOLD"]["data"].slicer[ + ..., 0:1 + ] + marker.fit_transform(element_data) -- 2.52.0 From 7b1241d363c380c2fbd7c3fc125c81a9bbef3414 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 15 Feb 2024 19:17:48 +0100 Subject: [PATCH 33/47] update: improve tests for test_parcel_aggregation.py --- .../markers/tests/test_parcel_aggregation.py | 1043 +++++++++-------- 1 file changed, 566 insertions(+), 477 deletions(-) diff --git a/junifer/markers/tests/test_parcel_aggregation.py b/junifer/markers/tests/test_parcel_aggregation.py index 12bd01c3f..d7532f181 100644 --- a/junifer/markers/tests/test_parcel_aggregation.py +++ b/junifer/markers/tests/test_parcel_aggregation.py @@ -10,16 +10,17 @@ from pathlib import Path import nibabel as nib import numpy as np import pytest -from nilearn import datasets -from nilearn.image import concat_imgs, math_img, new_img_like, resample_to_img +from nilearn.image import math_img, new_img_like from nilearn.maskers import NiftiLabelsMasker, NiftiMasker from nilearn.masking import compute_brain_mask from numpy.testing import assert_array_almost_equal, assert_array_equal from scipy.stats import trim_mean -from junifer.data import load_mask, load_parcellation, register_parcellation +from junifer.data import get_mask, get_parcellation, register_parcellation +from junifer.datareader import DefaultDataReader from junifer.markers.parcel_aggregation import ParcelAggregation from junifer.storage import SQLiteFeatureStorage +from junifer.testing.datagrabbers import PartlyCloudyTestingDataGrabber def test_ParcelAggregation_input_output() -> None: @@ -36,123 +37,145 @@ def test_ParcelAggregation_input_output() -> None: def test_ParcelAggregation_3D() -> None: """Test ParcelAggregation object on 3D images.""" - # Get the testing parcellation (for nilearn) - parcellation = datasets.fetch_atlas_schaefer_2018(n_rois=100) - - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) - - # Mask parcellation manually - parcellation_img_res = resample_to_img( - parcellation.maps, - img, - interpolation="nearest", - ) - parcellation_bin = math_img( - "img != 0", - img=parcellation_img_res, - ) - - # Create NiftiMasker - masker = NiftiMasker(parcellation_bin, target_affine=img.affine) - data = masker.fit_transform(img) - parcellation_values = masker.transform(parcellation_img_res) - parcellation_values = np.squeeze(parcellation_values).astype(int) - - # Compute the mean manually - manual = [] - for t_v in sorted(np.unique(parcellation_values)): - t_values = np.mean(data[:, parcellation_values == t_v]) - manual.append(t_values) - manual = np.array(manual)[np.newaxis, :] - - # Create NiftiLabelsMasker - nifti_masker = NiftiLabelsMasker(labels_img=parcellation.maps) - auto = nifti_masker.fit_transform(img) - - # Check that arrays are almost equal - assert_array_almost_equal(auto, manual) - - # Use the ParcelAggregation object - marker = ParcelAggregation( - parcellation="Schaefer100x7", - method="mean", - name="gmd_schaefer100x7_mean", - on="VBM_GM", - ) # Test passing "on" as a keyword argument - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - jun_values3d_mean = marker.fit_transform(input)["VBM_GM"]["data"] - - assert jun_values3d_mean.ndim == 2 - assert jun_values3d_mean.shape[0] == 1 - assert_array_equal(manual, jun_values3d_mean) - - # Test using another function (std) - manual = [] - for t_v in sorted(np.unique(parcellation_values)): - t_values = np.std(data[:, parcellation_values == t_v]) - manual.append(t_values) - manual = np.array(manual)[np.newaxis, :] - - # Use the ParcelAggregation object - marker = ParcelAggregation(parcellation="Schaefer100x7", method="std") - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - jun_values3d_std = marker.fit_transform(input)["VBM_GM"]["data"] - - assert jun_values3d_std.ndim == 2 - assert jun_values3d_std.shape[0] == 1 - assert_array_equal(manual, jun_values3d_std) - - # Test using another function with parameters - manual = [] - for t_v in sorted(np.unique(parcellation_values)): - t_values = trim_mean( - data[:, parcellation_values == t_v], - proportiontocut=0.1, - axis=None, # type: ignore + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + # Create ParcelAggregation object + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + on="BOLD", ) - manual.append(t_values) - manual = np.array(manual)[np.newaxis, :] + element_data["BOLD"]["data"] = element_data["BOLD"]["data"].slicer[ + ..., 0:1 + ] - # Use the ParcelAggregation object - marker = ParcelAggregation( - parcellation="Schaefer100x7", - method="trim_mean", - method_params={"proportiontocut": 0.1}, - ) - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - jun_values3d_tm = marker.fit_transform(input)["VBM_GM"]["data"] + # Compare with nilearn + # Load testing parcellation + testing_parcellation, _ = get_parcellation( + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=element_data["BOLD"], + ) + # Binarize parcellation + testing_parcellation_bin = math_img( + "img != 0", + img=testing_parcellation, + ) + # Create NiftiMasker + masker = NiftiMasker( + testing_parcellation_bin, + target_affine=element_data["BOLD"]["data"].affine, + ) + data = masker.fit_transform(element_data["BOLD"]["data"]) + parcellation_values = np.squeeze( + masker.transform(testing_parcellation) + ).astype(int) + # Compute the mean manually + manual = [] + for t_v in sorted(np.unique(parcellation_values)): + t_values = np.mean(data[:, parcellation_values == t_v]) + manual.append(t_values) + manual = np.array(manual)[np.newaxis, :] - assert jun_values3d_tm.ndim == 2 - assert jun_values3d_tm.shape[0] == 1 - assert_array_equal(manual, jun_values3d_tm) + # Create NiftiLabelsMasker + nifti_labels_masker = NiftiLabelsMasker( + labels_img=testing_parcellation + ) + nifti_labels_masked_bold = nifti_labels_masker.fit_transform( + element_data["BOLD"]["data"].slicer[..., 0:1] + ) + + parcel_agg_mean_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] + # Check that arrays are almost equal + assert_array_equal(parcel_agg_mean_bold_data, manual) + assert_array_almost_equal(nifti_labels_masked_bold, manual) + + # Check further + assert parcel_agg_mean_bold_data.ndim == 2 + assert parcel_agg_mean_bold_data.shape[0] == 1 + assert_array_equal( + nifti_labels_masked_bold.shape, parcel_agg_mean_bold_data.shape + ) + assert_array_equal(nifti_labels_masked_bold, parcel_agg_mean_bold_data) + + # Compute std manually + manual = [] + for t_v in sorted(np.unique(parcellation_values)): + t_values = np.std(data[:, parcellation_values == t_v]) + manual.append(t_values) + manual = np.array(manual)[np.newaxis, :] + + # Create ParcelAggregation object + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="std", + on="BOLD", + ) + parcel_agg_std_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] + assert parcel_agg_std_bold_data.ndim == 2 + assert parcel_agg_std_bold_data.shape[0] == 1 + assert_array_equal(parcel_agg_std_bold_data, manual) + + # Test using another function with parameters + manual = [] + for t_v in sorted(np.unique(parcellation_values)): + t_values = trim_mean( + data[:, parcellation_values == t_v], + proportiontocut=0.1, + axis=None, # type: ignore + ) + manual.append(t_values) + manual = np.array(manual)[np.newaxis, :] + + # Create ParcelAggregation object + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="trim_mean", + method_params={"proportiontocut": 0.1}, + on="BOLD", + ) + parcel_agg_trim_mean_bold_data = marker.fit_transform(element_data)[ + "BOLD" + ]["data"] + assert parcel_agg_trim_mean_bold_data.ndim == 2 + assert parcel_agg_trim_mean_bold_data.shape[0] == 1 + assert_array_equal(parcel_agg_trim_mean_bold_data, manual) def test_ParcelAggregation_4D(): """Test ParcelAggregation object on 4D images.""" - # Get the testing parcellation (for nilearn) - parcellation = datasets.fetch_atlas_schaefer_2018( - n_rois=100, yeo_networks=7, resolution_mm=2 - ) + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + # Create ParcelAggregation object + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", method="mean" + ) + parcel_agg_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] - # Get the SPM auditory data: - subject_data = datasets.fetch_spm_auditory() - fmri_img = concat_imgs(subject_data.func) # type: ignore + # Compare with nilearn + # Load testing parcellation + testing_parcellation, _ = get_parcellation( + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=element_data["BOLD"], + ) + # Extract data + nifti_labels_masker = NiftiLabelsMasker( + labels_img=testing_parcellation + ) + nifti_labels_masked_bold = nifti_labels_masker.fit_transform( + element_data["BOLD"]["data"] + ) - # Create NiftiLabelsMasker - nifti_masker = NiftiLabelsMasker(labels_img=parcellation.maps) - auto4d = nifti_masker.fit_transform(fmri_img) - - # Create ParcelAggregation object - marker = ParcelAggregation(parcellation="Schaefer100x7", method="mean") - input = {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - jun_values4d = marker.fit_transform(input)["BOLD"]["data"] - - assert jun_values4d.ndim == 2 - assert_array_equal(auto4d.shape, jun_values4d.shape) - assert_array_equal(auto4d, jun_values4d) + assert parcel_agg_bold_data.ndim == 2 + assert_array_equal( + nifti_labels_masked_bold.shape, parcel_agg_bold_data.shape + ) + assert_array_equal(nifti_labels_masked_bold, parcel_agg_bold_data) def test_ParcelAggregation_storage(tmp_path: Path) -> None: @@ -164,130 +187,148 @@ def test_ParcelAggregation_storage(tmp_path: Path) -> None: The path to the test directory. """ - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) - uri = tmp_path / "test_sphere_storage_3D.sqlite" + # Store 3D + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + storage = SQLiteFeatureStorage( + uri=tmp_path / "test_parcel_storage_3D.sqlite", upsert="ignore" + ) + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + on="BOLD", + ) + element_data["BOLD"]["data"] = element_data["BOLD"]["data"].slicer[ + ..., 0:1 + ] + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "BOLD_ParcelAggregation" for x in features.values() + ) - storage = SQLiteFeatureStorage(uri=uri, upsert="ignore") - meta = { - "element": {"subject": "sub-01", "session": "ses-01"}, - "dependencies": {"nilearn", "nibabel"}, - } - input = {"VBM_GM": {"data": img, "meta": meta, "space": "MNI"}} - marker = ParcelAggregation( - parcellation="Schaefer100x7", method="mean", on="VBM_GM" - ) - - marker.fit_transform(input, storage=storage) - - features = storage.list_features() - assert any( - x["name"] == "VBM_GM_ParcelAggregation" for x in features.values() - ) - - meta = { - "element": {"subject": "sub-01", "session": "ses-01"}, - "dependencies": {"nilearn", "nibabel"}, - } - # Get the SPM auditory data - subject_data = datasets.fetch_spm_auditory() - fmri_img = concat_imgs(subject_data.func) # type: ignore - input = {"BOLD": {"data": fmri_img, "meta": meta, "space": "MNI"}} - marker = ParcelAggregation( - parcellation="Schaefer100x7", method="mean", on="BOLD" - ) - - marker.fit_transform(input, storage=storage) - features = storage.list_features() - assert any( - x["name"] == "BOLD_ParcelAggregation" for x in features.values() - ) + # Store 4D + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + storage = SQLiteFeatureStorage( + uri=tmp_path / "test_parcel_storage_4D.sqlite", upsert="ignore" + ) + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + on="BOLD", + ) + marker.fit_transform(input=element_data, storage=storage) + features = storage.list_features() + assert any( + x["name"] == "BOLD_ParcelAggregation" for x in features.values() + ) def test_ParcelAggregation_3D_mask() -> None: """Test ParcelAggregation object on 3D images with mask.""" + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + # Create ParcelAggregation object + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + name="tian_mean", + on="BOLD", + masks="compute_brain_mask", + ) + element_data["BOLD"]["data"] = element_data["BOLD"]["data"].slicer[ + ..., 0:1 + ] + parcel_agg_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] - # Get the testing parcellation (for nilearn) - parcellation = datasets.fetch_atlas_schaefer_2018(n_rois=100) + # Compare with nilearn + # Load testing parcellation + testing_parcellation, _ = get_parcellation( + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=element_data["BOLD"], + ) + # Load mask + mask_img = get_mask( + "compute_brain_mask", target_data=element_data["BOLD"] + ) + # Extract data + nifti_labels_masker = NiftiLabelsMasker( + labels_img=testing_parcellation, mask_img=mask_img + ) + nifti_labels_masked_bold = nifti_labels_masker.fit_transform( + element_data["BOLD"]["data"].slicer[..., 0:1] + ) - # Get one mask - mask_img, _, _ = load_mask("GM_prob0.2") - - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) - - # Create NiftiLabelsMasker - nifti_masker = NiftiLabelsMasker( - labels_img=parcellation.maps, mask_img=mask_img - ) - auto = nifti_masker.fit_transform(img) - - # Use the ParcelAggregation object - marker = ParcelAggregation( - parcellation="Schaefer100x7", - method="mean", - masks="GM_prob0.2", - name="gmd_schaefer100x7_mean", - on="VBM_GM", - ) # Test passing "on" as a keyword argument - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - jun_values3d_mean = marker.fit_transform(input)["VBM_GM"]["data"] - - assert jun_values3d_mean.ndim == 2 - assert jun_values3d_mean.shape[0] == 1 - assert_array_almost_equal(auto, jun_values3d_mean) + assert parcel_agg_bold_data.ndim == 2 + assert_array_equal( + nifti_labels_masked_bold.shape, parcel_agg_bold_data.shape + ) + assert_array_equal(nifti_labels_masked_bold, parcel_agg_bold_data) def test_ParcelAggregation_3D_mask_computed() -> None: """Test ParcelAggregation object on 3D images with computed masks.""" + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + element_data["BOLD"]["data"] = element_data["BOLD"]["data"].slicer[ + ..., 0:1 + ] - # Get the testing parcellation (for nilearn) - parcellation = datasets.fetch_atlas_schaefer_2018(n_rois=100) + # Compare with nilearn + # Load testing parcellation + testing_parcellation, _ = get_parcellation( + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=element_data["BOLD"], + ) + # Get a mask + mask_img = compute_brain_mask( + element_data["BOLD"]["data"], threshold=0.2 + ) + # Create NiftiLabelsMasker + nifti_labels_masker = NiftiLabelsMasker( + labels_img=testing_parcellation, mask_img=mask_img + ) + nifti_labels_masked_bold_good = nifti_labels_masker.fit_transform( + element_data["BOLD"]["data"] + ) - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) + # Get another mask + mask_img = compute_brain_mask( + element_data["BOLD"]["data"], threshold=0.5 + ) + # Create NiftiLabelsMasker + nifti_labels_masker = NiftiLabelsMasker( + labels_img=testing_parcellation, mask_img=mask_img + ) + nifti_labels_masked_bold_bad = nifti_labels_masker.fit_transform( + mask_img + ) - # Get one mask - mask_img = compute_brain_mask(img, threshold=0.2) + # Use the ParcelAggregation object + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + masks={"compute_brain_mask": {"threshold": 0.2}}, + name="tian_mean", + on="BOLD", + ) + parcel_agg_mean_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] - # Create NiftiLabelsMasker - nifti_masker = NiftiLabelsMasker( - labels_img=parcellation.maps, mask_img=mask_img - ) - auto = nifti_masker.fit_transform(img) + assert parcel_agg_mean_bold_data.ndim == 2 + assert parcel_agg_mean_bold_data.shape[0] == 1 + assert_array_almost_equal( + nifti_labels_masked_bold_good, parcel_agg_mean_bold_data + ) - # Get one mask - mask_img = compute_brain_mask(img, threshold=0.5) - - # Create NiftiLabelsMasker - nifti_masker = NiftiLabelsMasker( - labels_img=parcellation.maps, mask_img=mask_img - ) - auto_bad = nifti_masker.fit_transform(img) - - # Use the ParcelAggregation object - marker = ParcelAggregation( - parcellation="Schaefer100x7", - method="mean", - masks={"compute_brain_mask": {"threshold": 0.2}}, - name="gmd_schaefer100x7_mean", - on="VBM_GM", - ) # Test passing "on" as a keyword argument - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - jun_values3d_mean = marker.fit_transform(input)["VBM_GM"]["data"] - - assert jun_values3d_mean.ndim == 2 - assert jun_values3d_mean.shape[0] == 1 - assert_array_almost_equal(auto, jun_values3d_mean) - - with pytest.raises(AssertionError): - assert_array_almost_equal(jun_values3d_mean, auto_bad) + with pytest.raises(AssertionError): + assert_array_almost_equal( + parcel_agg_mean_bold_data, nifti_labels_masked_bold_bad + ) def test_ParcelAggregation_3D_multiple_non_overlapping(tmp_path: Path) -> None: @@ -299,89 +340,93 @@ def test_ParcelAggregation_3D_multiple_non_overlapping(tmp_path: Path) -> None: The path to the test directory. """ + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + element_data["BOLD"]["data"] = element_data["BOLD"]["data"].slicer[ + ..., 0:1 + ] - # Get the testing parcellation - parcellation, labels, _, _ = load_parcellation("Schaefer100x7") + # Load testing parcellation + testing_parcellation, labels = get_parcellation( + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=element_data["BOLD"], + ) - assert parcellation is not None + # Create two parcellations from it + parcellation_data = testing_parcellation.get_fdata() + parcellation1_data = parcellation_data.copy() + parcellation1_data[parcellation1_data > 8] = 0 + parcellation2_data = parcellation_data.copy() + parcellation2_data[parcellation2_data <= 8] = 0 + parcellation2_data[parcellation2_data > 0] -= 8 + labels1 = labels[:8] + labels2 = labels[8:] - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) + parcellation1_img = new_img_like( + testing_parcellation, parcellation1_data + ) + parcellation2_img = new_img_like( + testing_parcellation, parcellation2_data + ) - # Create two parcellations from it - parcellation_data = parcellation.get_fdata() - parcellation1_data = parcellation_data.copy() - parcellation1_data[parcellation1_data > 50] = 0 - parcellation2_data = parcellation_data.copy() - parcellation2_data[parcellation2_data <= 50] = 0 - parcellation2_data[parcellation2_data > 0] -= 50 - labels1 = labels[:50] - labels2 = labels[50:] + parcellation1_path = tmp_path / "parcellation1.nii.gz" + parcellation2_path = tmp_path / "parcellation2.nii.gz" - parcellation1_img = new_img_like(parcellation, parcellation1_data) - parcellation2_img = new_img_like(parcellation, parcellation2_data) + nib.save(parcellation1_img, parcellation1_path) + nib.save(parcellation2_img, parcellation2_path) - parcellation1_path = tmp_path / "parcellation1.nii.gz" - parcellation2_path = tmp_path / "parcellation2.nii.gz" + register_parcellation( + name="TianxS1x3TxMNInonlinear2009cAsym_low", + parcellation_path=parcellation1_path, + parcels_labels=labels1, + space="MNI152NLin2009cAsym", + overwrite=True, + ) + register_parcellation( + name="TianxS1x3TxMNInonlinear2009cAsym_high", + parcellation_path=parcellation2_path, + parcels_labels=labels2, + space="MNI152NLin2009cAsym", + overwrite=True, + ) - nib.save(parcellation1_img, parcellation1_path) - nib.save(parcellation2_img, parcellation2_path) + # Use the ParcelAggregation object on the original parcellation + marker_original = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + name="tian_mean", + on="BOLD", + ) + orig_mean = marker_original.fit_transform(element_data)["BOLD"] - register_parcellation( - name="Schaefer100x7_low", - parcellation_path=parcellation1_path, - parcels_labels=labels1, - space="MNI", - overwrite=True, - ) - register_parcellation( - name="Schaefer100x7_high", - parcellation_path=parcellation2_path, - parcels_labels=labels2, - space="MNI", - overwrite=True, - ) + orig_mean_data = orig_mean["data"] + assert orig_mean_data.ndim == 2 + assert orig_mean_data.shape == (1, 16) - # Use the ParcelAggregation object on the original parcellation - marker_original = ParcelAggregation( - parcellation="Schaefer100x7", - method="mean", - name="gmd_schaefer100x7_mean", - on="VBM_GM", - ) # Test passing "on" as a keyword argument - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - orig_mean = marker_original.fit_transform(input)["VBM_GM"] + # Use the ParcelAggregation object on the two parcellations + marker_split = ParcelAggregation( + parcellation=[ + "TianxS1x3TxMNInonlinear2009cAsym_low", + "TianxS1x3TxMNInonlinear2009cAsym_high", + ], + method="mean", + name="tian_mean", + on="BOLD", + ) - orig_mean_data = orig_mean["data"] - assert orig_mean_data.ndim == 2 - assert orig_mean_data.shape[0] == 1 - assert orig_mean_data.shape[1] == 100 - # assert_array_almost_equal(auto, jun_values3d_mean) + # No warnings should be raised + with warnings.catch_warnings(): + warnings.simplefilter("error", category=UserWarning) + split_mean = marker_split.fit_transform(element_data)["BOLD"] - # Use the ParcelAggregation object on the two parcellations - marker_split = ParcelAggregation( - parcellation=["Schaefer100x7_low", "Schaefer100x7_high"], - method="mean", - name="gmd_schaefer100x7_mean", - on="VBM_GM", - ) # Test passing "on" as a keyword argument - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} + split_mean_data = split_mean["data"] - # No warnings should be raised - with warnings.catch_warnings(): - warnings.simplefilter("error", category=UserWarning) - split_mean = marker_split.fit_transform(input)["VBM_GM"] - split_mean_data = split_mean["data"] + assert split_mean_data.ndim == 2 + assert split_mean_data.shape == (1, 16) - assert split_mean_data.ndim == 2 - assert split_mean_data.shape[0] == 1 - assert split_mean_data.shape[1] == 100 - - # Data and labels should be the same - assert_array_equal(orig_mean_data, split_mean_data) - assert orig_mean["col_names"] == split_mean["col_names"] + # Data and labels should be the same + assert_array_equal(orig_mean_data, split_mean_data) + assert orig_mean["col_names"] == split_mean["col_names"] def test_ParcelAggregation_3D_multiple_overlapping(tmp_path: Path) -> None: @@ -393,95 +438,100 @@ def test_ParcelAggregation_3D_multiple_overlapping(tmp_path: Path) -> None: The path to the test directory. """ + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + element_data["BOLD"]["data"] = element_data["BOLD"]["data"].slicer[ + ..., 0:1 + ] - # Get the testing parcellation - parcellation, labels, _, _ = load_parcellation("Schaefer100x7") + # Load testing parcellation + testing_parcellation, labels = get_parcellation( + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=element_data["BOLD"], + ) - assert parcellation is not None + # Create two parcellations from it + parcellation_data = testing_parcellation.get_fdata() + parcellation1_data = parcellation_data.copy() + parcellation1_data[parcellation1_data > 8] = 0 + parcellation2_data = parcellation_data.copy() - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) + # Make the second parcellation overlap with the first + parcellation2_data[parcellation2_data <= 6] = 0 + parcellation2_data[parcellation2_data > 0] -= 6 + labels1 = [f"low_{x}" for x in labels[:8]] # Change the labels + labels2 = [f"high_{x}" for x in labels[6:]] # Change the labels - # Create two parcellations from it - parcellation_data = parcellation.get_fdata() - parcellation1_data = parcellation_data.copy() - parcellation1_data[parcellation1_data > 50] = 0 - parcellation2_data = parcellation_data.copy() + parcellation1_img = new_img_like( + testing_parcellation, parcellation1_data + ) + parcellation2_img = new_img_like( + testing_parcellation, parcellation2_data + ) - # Make the second parcellation overlap with the first - parcellation2_data[parcellation2_data <= 45] = 0 - parcellation2_data[parcellation2_data > 0] -= 45 - labels1 = [f"low_{x}" for x in labels[:50]] # Change the labels - labels2 = [f"high_{x}" for x in labels[45:]] # Change the labels + parcellation1_path = tmp_path / "parcellation1.nii.gz" + parcellation2_path = tmp_path / "parcellation2.nii.gz" - parcellation1_img = new_img_like(parcellation, parcellation1_data) - parcellation2_img = new_img_like(parcellation, parcellation2_data) + nib.save(parcellation1_img, parcellation1_path) + nib.save(parcellation2_img, parcellation2_path) - parcellation1_path = tmp_path / "parcellation1.nii.gz" - parcellation2_path = tmp_path / "parcellation2.nii.gz" + register_parcellation( + name="TianxS1x3TxMNInonlinear2009cAsym_low", + parcellation_path=parcellation1_path, + parcels_labels=labels1, + space="MNI152NLin2009cAsym", + overwrite=True, + ) + register_parcellation( + name="TianxS1x3TxMNInonlinear2009cAsym_high", + parcellation_path=parcellation2_path, + parcels_labels=labels2, + space="MNI152NLin2009cAsym", + overwrite=True, + ) - nib.save(parcellation1_img, parcellation1_path) - nib.save(parcellation2_img, parcellation2_path) + # Use the ParcelAggregation object on the original parcellation + marker_original = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + name="tian_mean", + on="BOLD", + ) + orig_mean = marker_original.fit_transform(element_data)["BOLD"] - register_parcellation( - name="Schaefer100x7_low2", - parcellation_path=parcellation1_path, - parcels_labels=labels1, - space="MNI", - overwrite=True, - ) - register_parcellation( - name="Schaefer100x7_high2", - parcellation_path=parcellation2_path, - parcels_labels=labels2, - space="MNI", - overwrite=True, - ) + orig_mean_data = orig_mean["data"] + assert orig_mean_data.ndim == 2 + assert orig_mean_data.shape == (1, 16) - # Use the ParcelAggregation object on the original parcellation - marker_original = ParcelAggregation( - parcellation="Schaefer100x7", - method="mean", - name="gmd_schaefer100x7_mean", - on="VBM_GM", - ) # Test passing "on" as a keyword argument - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - orig_mean = marker_original.fit_transform(input)["VBM_GM"] + # Use the ParcelAggregation object on the two parcellations + marker_split = ParcelAggregation( + parcellation=[ + "TianxS1x3TxMNInonlinear2009cAsym_low", + "TianxS1x3TxMNInonlinear2009cAsym_high", + ], + method="mean", + name="tian_mean", + on="BOLD", + ) + # Warning should be raised + with pytest.warns(RuntimeWarning, match="overlapping voxels"): + split_mean = marker_split.fit_transform(element_data)["BOLD"] - orig_mean_data = orig_mean["data"] - assert orig_mean_data.ndim == 2 - assert orig_mean_data.shape[0] == 1 - assert orig_mean_data.shape[1] == 100 - # assert_array_almost_equal(auto, jun_values3d_mean) + split_mean_data = split_mean["data"] - # Use the ParcelAggregation object on the two parcellations - marker_split = ParcelAggregation( - parcellation=["Schaefer100x7_low2", "Schaefer100x7_high2"], - method="mean", - name="gmd_schaefer100x7_mean", - on="VBM_GM", - ) # Test passing "on" as a keyword argument - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - with pytest.warns(RuntimeWarning, match="overlapping voxels"): - split_mean = marker_split.fit_transform(input)["VBM_GM"] - split_mean_data = split_mean["data"] + assert split_mean_data.ndim == 2 + assert split_mean_data.shape == (1, 18) - assert split_mean_data.ndim == 2 - assert split_mean_data.shape[0] == 1 - assert split_mean_data.shape[1] == 105 + # Overlapping voxels should be NaN + assert np.isnan(split_mean_data[:, 8:10]).all() - # Overlapping voxels should be NaN - assert np.isnan(split_mean_data[:, 50:55]).all() + non_nan = split_mean_data[~np.isnan(split_mean_data)] + # Data should be the same + assert_array_equal(orig_mean_data, non_nan[None, :]) - non_nan = split_mean_data[~np.isnan(split_mean_data)] - # Data should be the same - assert_array_equal(orig_mean_data, non_nan[None, :]) - - # Labels should be "low" for the first 50 and "high" for the second 50 - assert all(x.startswith("low") for x in split_mean["col_names"][:50]) - assert all(x.startswith("high") for x in split_mean["col_names"][50:]) + # Labels should be "low" for the first 8 and "high" for the second 8 + assert all(x.startswith("low") for x in split_mean["col_names"][:8]) + assert all(x.startswith("high") for x in split_mean["col_names"][8:]) def test_ParcelAggregation_3D_multiple_duplicated_labels( @@ -495,135 +545,164 @@ def test_ParcelAggregation_3D_multiple_duplicated_labels( The path to the test directory. """ + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + element_data["BOLD"]["data"] = element_data["BOLD"]["data"].slicer[ + ..., 0:1 + ] - # Get the testing parcellation - parcellation, labels, _, _ = load_parcellation("Schaefer100x7") + # Load testing parcellation + testing_parcellation, labels = get_parcellation( + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=element_data["BOLD"], + ) - assert parcellation is not None + # Create two parcellations from it + parcellation_data = testing_parcellation.get_fdata() + parcellation1_data = parcellation_data.copy() + parcellation1_data[parcellation1_data > 8] = 0 + parcellation2_data = parcellation_data.copy() + parcellation2_data[parcellation2_data <= 8] = 0 + parcellation2_data[parcellation2_data > 0] -= 8 + labels1 = labels[:8] + labels2 = labels[7:-1] # One label is duplicated - # Get the oasis VBM data - oasis_dataset = datasets.fetch_oasis_vbm(n_subjects=1) - vbm = oasis_dataset.gray_matter_maps[0] - img = nib.load(vbm) + parcellation1_img = new_img_like( + testing_parcellation, parcellation1_data + ) + parcellation2_img = new_img_like( + testing_parcellation, parcellation2_data + ) - # Create two parcellations from it - parcellation_data = parcellation.get_fdata() - parcellation1_data = parcellation_data.copy() - parcellation1_data[parcellation1_data > 50] = 0 - parcellation2_data = parcellation_data.copy() - parcellation2_data[parcellation2_data <= 50] = 0 - parcellation2_data[parcellation2_data > 0] -= 50 - labels1 = labels[:50] - labels2 = labels[49:-1] # One label is duplicated + parcellation1_path = tmp_path / "parcellation1.nii.gz" + parcellation2_path = tmp_path / "parcellation2.nii.gz" - parcellation1_img = new_img_like(parcellation, parcellation1_data) - parcellation2_img = new_img_like(parcellation, parcellation2_data) + nib.save(parcellation1_img, parcellation1_path) + nib.save(parcellation2_img, parcellation2_path) - parcellation1_path = tmp_path / "parcellation1.nii.gz" - parcellation2_path = tmp_path / "parcellation2.nii.gz" + register_parcellation( + name="TianxS1x3TxMNInonlinear2009cAsym_low", + parcellation_path=parcellation1_path, + parcels_labels=labels1, + space="MNI152NLin2009cAsym", + overwrite=True, + ) + register_parcellation( + name="TianxS1x3TxMNInonlinear2009cAsym_high", + parcellation_path=parcellation2_path, + parcels_labels=labels2, + space="MNI152NLin2009cAsym", + overwrite=True, + ) - nib.save(parcellation1_img, parcellation1_path) - nib.save(parcellation2_img, parcellation2_path) + # Use the ParcelAggregation object on the original parcellation + marker_original = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + name="tian_mean", + on="BOLD", + ) + orig_mean = marker_original.fit_transform(element_data)["BOLD"] - register_parcellation( - name="Schaefer100x7_low", - parcellation_path=parcellation1_path, - parcels_labels=labels1, - space="MNI", - overwrite=True, - ) - register_parcellation( - name="Schaefer100x7_high", - parcellation_path=parcellation2_path, - parcels_labels=labels2, - space="MNI", - overwrite=True, - ) + orig_mean_data = orig_mean["data"] + assert orig_mean_data.ndim == 2 + assert orig_mean_data.shape == (1, 16) - # Use the ParcelAggregation object on the original parcellation - marker_original = ParcelAggregation( - parcellation="Schaefer100x7", - method="mean", - name="gmd_schaefer100x7_mean", - on="VBM_GM", - ) # Test passing "on" as a keyword argument - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} - orig_mean = marker_original.fit_transform(input)["VBM_GM"] + # Use the ParcelAggregation object on the two parcellations + marker_split = ParcelAggregation( + parcellation=[ + "TianxS1x3TxMNInonlinear2009cAsym_low", + "TianxS1x3TxMNInonlinear2009cAsym_high", + ], + method="mean", + name="tian_mean", + on="BOLD", + ) - orig_mean_data = orig_mean["data"] - assert orig_mean_data.ndim == 2 - assert orig_mean_data.shape[0] == 1 - assert orig_mean_data.shape[1] == 100 - # assert_array_almost_equal(auto, jun_values3d_mean) + # Warning should be raised + with pytest.warns(RuntimeWarning, match="duplicated labels."): + split_mean = marker_split.fit_transform(element_data)["BOLD"] - # Use the ParcelAggregation object on the two parcellations - marker_split = ParcelAggregation( - parcellation=["Schaefer100x7_low", "Schaefer100x7_high"], - method="mean", - name="gmd_schaefer100x7_mean", - on="VBM_GM", - ) # Test passing "on" as a keyword argument - input = {"VBM_GM": {"data": img, "meta": {}, "space": "MNI"}} + split_mean_data = split_mean["data"] - with pytest.warns(RuntimeWarning, match="duplicated labels."): - split_mean = marker_split.fit_transform(input)["VBM_GM"] - split_mean_data = split_mean["data"] + assert split_mean_data.ndim == 2 + assert split_mean_data.shape == (1, 16) - assert split_mean_data.ndim == 2 - assert split_mean_data.shape[0] == 1 - assert split_mean_data.shape[1] == 100 + # Data should be the same + assert_array_equal(orig_mean_data, split_mean_data) - # Data should be the same - assert_array_equal(orig_mean_data, split_mean_data) - - # Labels should be prefixed with the parcellation name - col_names = [f"Schaefer100x7_low_{x}" for x in labels1] - col_names += [f"Schaefer100x7_high_{x}" for x in labels2] - assert col_names == split_mean["col_names"] + # Labels should be prefixed with the parcellation name + col_names = [ + f"TianxS1x3TxMNInonlinear2009cAsym_low_{x}" for x in labels1 + ] + col_names += [ + f"TianxS1x3TxMNInonlinear2009cAsym_high_{x}" for x in labels2 + ] + assert col_names == split_mean["col_names"] def test_ParcelAggregation_4D_agg_time(): """Test ParcelAggregation object on 4D images, aggregating time.""" - # Get the testing parcellation (for nilearn) - parcellation = datasets.fetch_atlas_schaefer_2018( - n_rois=100, yeo_networks=7, resolution_mm=2 - ) + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + # Create ParcelAggregation object + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + time_method="mean", + on="BOLD", + ) + parcel_agg_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] - # Get the SPM auditory data: - subject_data = datasets.fetch_spm_auditory() - fmri_img = concat_imgs(subject_data.func) # type: ignore + # Compare with nilearn + # Loading testing parcellation + testing_parcellation, _ = get_parcellation( + parcellation=["TianxS1x3TxMNInonlinear2009cAsym"], + target_data=element_data["BOLD"], + ) + # Extract data + nifti_labels_masker = NiftiLabelsMasker( + labels_img=testing_parcellation + ) + nifti_labels_masked_bold = nifti_labels_masker.fit_transform( + element_data["BOLD"]["data"] + ) + nifti_labels_masked_bold_mean = nifti_labels_masked_bold.mean(axis=0) - # Create NiftiLabelsMasker - nifti_masker = NiftiLabelsMasker(labels_img=parcellation.maps) - auto4d = nifti_masker.fit_transform(fmri_img) - auto_mean = auto4d.mean(axis=0) + assert parcel_agg_bold_data.ndim == 1 + assert_array_equal( + nifti_labels_masked_bold_mean.shape, parcel_agg_bold_data.shape + ) + assert_array_almost_equal( + nifti_labels_masked_bold_mean, parcel_agg_bold_data, decimal=2 + ) - # Create ParcelAggregation object - marker = ParcelAggregation( - parcellation="Schaefer100x7", method="mean", time_method="mean" - ) - input = {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - jun_values4d = marker.fit_transform(input)["BOLD"]["data"] + # Test picking first time point + nifti_labels_masked_bold_pick_0 = nifti_labels_masked_bold[:1, :] + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + time_method="select", + time_method_params={"pick": [0]}, + on="BOLD", + ) + parcel_agg_bold_data = marker.fit_transform(element_data)["BOLD"][ + "data" + ] - assert jun_values4d.ndim == 1 - assert_array_equal(auto_mean.shape, jun_values4d.shape) - assert_array_almost_equal(auto_mean, jun_values4d, decimal=2) + assert parcel_agg_bold_data.ndim == 2 + assert_array_equal( + nifti_labels_masked_bold_pick_0.shape, parcel_agg_bold_data.shape + ) + assert_array_equal( + nifti_labels_masked_bold_pick_0, parcel_agg_bold_data + ) - auto_pick_0 = auto4d[:1, :] - marker = ParcelAggregation( - parcellation="Schaefer100x7", - method="mean", - time_method="select", - time_method_params={"pick": [0]}, - ) - - input = {"BOLD": {"data": fmri_img, "meta": {}, "space": "MNI"}} - jun_values4d = marker.fit_transform(input)["BOLD"]["data"] - - assert jun_values4d.ndim == 2 - assert_array_equal(auto_pick_0.shape, jun_values4d.shape) - assert_array_equal(auto_pick_0, jun_values4d) +def test_ParcelAggregation_errors() -> None: + """Test errors for ParcelAggregation.""" with pytest.raises(ValueError, match="can only be used with BOLD data"): ParcelAggregation( parcellation="Schaefer100x7", @@ -643,12 +722,22 @@ def test_ParcelAggregation_4D_agg_time(): on="VBM_GM", ) - with pytest.warns(RuntimeWarning, match="No time dimension to aggregate"): - input = { - "BOLD": { - "data": fmri_img.slicer[..., 0:1], - "meta": {}, - "space": "MNI", - } - } - marker.fit_transform(input) + +def test_ParcelAggregation_warning() -> None: + """Test warning for ParcelAggregation.""" + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + with pytest.warns( + RuntimeWarning, match="No time dimension to aggregate" + ): + marker = ParcelAggregation( + parcellation="TianxS1x3TxMNInonlinear2009cAsym", + method="mean", + time_method="select", + time_method_params={"pick": [0]}, + on="BOLD", + ) + element_data["BOLD"]["data"] = element_data["BOLD"]["data"].slicer[ + ..., 0:1 + ] + marker.fit_transform(element_data) -- 2.52.0 From e9601f56bdb3ddc1eabc97c492b4d8b13b9a5e35 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 16 Feb 2024 10:26:43 +0100 Subject: [PATCH 34/47] update: improve tests for test_functions.py --- junifer/api/tests/test_functions.py | 14 +++++++++++++- 1 file changed, 13 insertions(+), 1 deletion(-) diff --git a/junifer/api/tests/test_functions.py b/junifer/api/tests/test_functions.py index 728a2a6b8..2d8909862 100644 --- a/junifer/api/tests/test_functions.py +++ b/junifer/api/tests/test_functions.py @@ -259,6 +259,9 @@ def test_queue_correct_yaml_config( tmp_path: Path, monkeypatch: pytest.MonkeyPatch, caplog: pytest.LogCaptureFixture, + datagrabber: Dict[str, str], + markers: List[Dict[str, str]], + storage: Dict[str, str], ) -> None: """Test proper YAML config generation for queueing. @@ -270,6 +273,12 @@ def test_queue_correct_yaml_config( The pytest.MonkeyPatch object. caplog : pytest.LogCaptureFixture The pytest.LogCaptureFixture object. + datagrabber : dict + Testing datagrabber as dictionary. + markers : list of dict + Testing markers as list of dictionary. + storage : dict + Testing storage as dictionary. """ with monkeypatch.context() as m: @@ -281,7 +290,7 @@ def test_queue_correct_yaml_config( "workdir": str(tmp_path.resolve()), "datagrabber": datagrabber, "markers": markers, - "storage": {"kind": "SQLiteFeatureStorage"}, + "storage": storage, "env": { "kind": "conda", "name": "junifer", @@ -499,6 +508,7 @@ def test_queue_without_elements( tmp_path: Path, monkeypatch: pytest.MonkeyPatch, caplog: pytest.LogCaptureFixture, + datagrabber: Dict[str, str], ) -> None: """Test queue without elements. @@ -510,6 +520,8 @@ def test_queue_without_elements( The pytest.MonkeyPatch object. caplog : pytest.LogCaptureFixture The pytest.LogCaptureFixture object. + datagrabber : dict + Testing datagrabber as dictionary. """ with monkeypatch.context() as m: -- 2.52.0 From 39e0eec9ba7f2319a051e3e29e760014649381d4 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 16 Feb 2024 10:27:46 +0100 Subject: [PATCH 35/47] chore: fix ReHo comparison test to a different Pearson's r tolerance --- junifer/markers/reho/tests/test_reho_parcels.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/junifer/markers/reho/tests/test_reho_parcels.py b/junifer/markers/reho/tests/test_reho_parcels.py index 94c854d49..746aff615 100644 --- a/junifer/markers/reho/tests/test_reho_parcels.py +++ b/junifer/markers/reho/tests/test_reho_parcels.py @@ -116,4 +116,4 @@ def test_ReHoParcels_comparison(tmp_path: Path) -> None: junifer_output_bold["data"].flatten(), afni_output_bold["data"].flatten(), ) - assert r >= 0.3 # this is very bad, but they differ... + assert r >= 0.2 # this is very bad, but they differ... -- 2.52.0 From 127b6a520629bfd4ff0f3f407d9aa352c0f188ec Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 16 Feb 2024 14:33:13 +0100 Subject: [PATCH 36/47] update: add template_type parameter to get_template() --- junifer/data/template_spaces.py | 39 ++++++++++++++--- junifer/data/tests/test_template_spaces.py | 49 ++++++++++++++++++---- 2 files changed, 75 insertions(+), 13 deletions(-) diff --git a/junifer/data/template_spaces.py b/junifer/data/template_spaces.py index fb8b18122..eb4a3e6da 100644 --- a/junifer/data/template_spaces.py +++ b/junifer/data/template_spaces.py @@ -99,6 +99,7 @@ def get_template( space: str, target_data: Dict[str, Any], extra_input: Optional[Dict[str, Any]] = None, + template_type: str = "T1w", ) -> nib.Nifti1Image: """Get template for the space, tailored for the target image. @@ -112,6 +113,8 @@ def get_template( extra_input : dict, optional The other fields in the data object. Useful for accessing other data types (default None). + template_type : {"T1w", "whole", "gm", "wm", "csf"}, optional + The template type to retrieve (default "T1w"). Returns ------- @@ -121,15 +124,19 @@ def get_template( Raises ------ ValueError - If ``space`` is invalid. + If ``space`` or ``template_type`` is invalid. RuntimeError - If template in the required resolution is not found. + If required template is not found. """ # Check for invalid space; early check to raise proper error if space not in tflow.templates(): raise_error(f"Unknown template space: {space}") + # Check for template type + if template_type not in ["T1w", "whole", "gm", "wm", "csf"]: + raise_error(f"Unknown template type: {template_type}") + # Get the min of the voxels sizes and use it as the resolution target_img = target_data["data"] resolution = np.min(target_img.header.get_zooms()[:3]).astype(int) @@ -145,18 +152,38 @@ def get_template( logger.info(f"Downloading template {space} in resolution {resolution}") # Retrieve template try: + suffix = None + desc = None + label = None + if template_type == "T1w": + suffix = template_type + desc = None + label = None + elif template_type == "whole": + suffix = "mask" + desc = "brain" + label = None + elif template_type in ["gm", "wm", "csf"]: + suffix = "probseg" + desc = None + label = template_type.upper() + # Set kwargs for fetching + kwargs = { + "suffix": suffix, + "desc": desc, + "label": label, + } template_path = tflow.get( space, raise_empty=True, resolution=resolution, - suffix="T1w", - desc=None, extension="nii.gz", + **kwargs, ) except Exception: # noqa: BLE001 raise_error( - f"Template {space} not found in the required resolution " - f"{resolution}", + f"Template {space} ({template_type}) with resolution {resolution} " + "not found", klass=RuntimeError, ) else: diff --git a/junifer/data/tests/test_template_spaces.py b/junifer/data/tests/test_template_spaces.py index cae6cf1be..7470823bf 100644 --- a/junifer/data/tests/test_template_spaces.py +++ b/junifer/data/tests/test_template_spaces.py @@ -11,7 +11,10 @@ import pytest from junifer.data import get_template, get_xfm from junifer.datareader import DefaultDataReader -from junifer.testing.datagrabbers import OasisVBMTestingDataGrabber +from junifer.testing.datagrabbers import ( + OasisVBMTestingDataGrabber, + PartlyCloudyTestingDataGrabber, +) @pytest.mark.skipif( @@ -33,15 +36,32 @@ def test_get_xfm(tmp_path: Path) -> None: assert isinstance(xfm_path, Path) -def test_get_template() -> None: - """Test tailored template image fetch.""" - with OasisVBMTestingDataGrabber() as dg: +@pytest.mark.parametrize( + "template_type", + [ + "T1w", + "whole", + "gm", + "wm", + "csf", + ], +) +def test_get_template(template_type: str) -> None: + """Test tailored template image fetch. + + Parameters + ---------- + template_type : str + The parametrized template type. + + """ + with PartlyCloudyTestingDataGrabber() as dg: element = dg["sub-01"] element_data = DefaultDataReader().fit_transform(element) - vbm_gm = element_data["VBM_GM"] + bold = element_data["BOLD"] # Get tailored parcellation tailored_template = get_template( - space=vbm_gm["space"], target_data=vbm_gm + space=bold["space"], target_data=bold, template_type=template_type ) assert isinstance(tailored_template, nib.Nifti1Image) @@ -54,7 +74,22 @@ def test_get_template_invalid_space() -> None: vbm_gm = element_data["VBM_GM"] # Get tailored parcellation with pytest.raises(ValueError, match="Unknown template space:"): - _ = get_template(space="andromeda", target_data=vbm_gm) + get_template(space="andromeda", target_data=vbm_gm) + + +def test_get_template_invalid_template_type() -> None: + """Test invalid template type check for template fetch.""" + with OasisVBMTestingDataGrabber() as dg: + element = dg["sub-01"] + element_data = DefaultDataReader().fit_transform(element) + vbm_gm = element_data["VBM_GM"] + # Get tailored parcellation + with pytest.raises(ValueError, match="Unknown template type:"): + get_template( + space=vbm_gm["space"], + target_data=vbm_gm, + template_type="xenon", + ) def test_get_template_closest_resolution() -> None: -- 2.52.0 From 55eb12db2f312d9a228887d8a7fd22895a9c94af Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 16 Feb 2024 15:05:25 +0100 Subject: [PATCH 37/47] update: rewrite compute_brain_mask() to allow variable template fetching via templateflow --- junifer/data/masks.py | 91 +++++++++++++++++++++++++++++- junifer/data/tests/test_masks.py | 96 +++++++++++++++++++++++++++++--- 2 files changed, 177 insertions(+), 10 deletions(-) diff --git a/junifer/data/masks.py b/junifer/data/masks.py index 99ef6ddc9..bff06a813 100644 --- a/junifer/data/masks.py +++ b/junifer/data/masks.py @@ -20,10 +20,9 @@ from typing import ( import nibabel as nib import numpy as np from nilearn.datasets import fetch_icbm152_brain_gm_mask -from nilearn.image import resample_to_img +from nilearn.image import get_data, new_img_like, resample_to_img from nilearn.masking import ( compute_background_mask, - compute_brain_mask, compute_epi_mask, intersect_masks, ) @@ -41,6 +40,87 @@ if TYPE_CHECKING: _masks_path = Path(__file__).parent / "masks" +def compute_brain_mask( + target_data: Dict[str, Any], + extra_input: Optional[Dict[str, Any]] = None, + mask_type: str = "whole", + threshold: float = 0.5, +) -> "Nifti1Image": + """Compute the whole-brain, grey-matter or white-matter mask. + + This mask is calculated using the template space and resolution as found + in the ``target_data``. + + Parameters + ---------- + target_data : dict + The corresponding item of the data object for which mask will be + loaded. + extra_input : dict, optional + The other fields in the data object. Useful for accessing other data + types (default None). + mask_type : {"whole", "gm", "wm"}, optional + Type of mask to be computed: + + * "whole" : whole-brain mask + * "gm" : grey-matter mask + * "wm" : white-matter mask + + (default "whole"). + threshold : float, optional + The value under which the template is cut off (default 0.5). + + Returns + ------- + Nifti1Image + The mask (3D image). + + Raises + ------ + ValueError + If ``mask_type`` is invalid or + if ``extra_input`` is None when ``target_data``'s space is native. + + """ + logger.debug(f"Computing {mask_type} mask") + + if mask_type not in ["whole", "gm", "wm"]: + raise_error(f"Unknown mask type: {mask_type}") + + # Check pre-requirements for space manipulation + target_space = target_data["space"] + # Set target standard space to target space + target_std_space = target_space + # Extra data type requirement check if target space is native + if target_space == "native": + # Check for extra inputs + if extra_input is None: + raise_error( + "No extra input provided, requires `Warp` " + "data type to infer target template space." + ) + # Set target standard space to warp file space source + target_std_space = extra_input["Warp"]["src"] + + # Fetch template in closest resolution + template = get_template( + space=target_std_space, + target_data=target_data, + extra_input=extra_input, + template_type=mask_type if mask_type in ["gm", "wm"] else "T1w", + ) + # Resample template to target image + target_img = target_data["data"] + resampled_template = resample_to_img( + source_img=template, target_img=target_img + ) + + # Threshold and get mask + mask = (get_data(resampled_template) >= threshold).astype("int8") + + return new_img_like(target_img, mask) # type: ignore + + def _fetch_icbm152_brain_gm_mask( target_img: "Nifti1Image", **kwargs, @@ -318,7 +398,12 @@ def get_mask( # noqa: C901 if callable(mask_object): if mask_params is None: mask_params = {} - mask_img = mask_object(target_img, **mask_params) + # From nilearn + if mask_name != "compute_brain_mask": + mask_img = mask_object(target_img, **mask_params) + # Not from nilearn + else: + mask_img = mask_object(target_data, **mask_params) # Mask is a Nifti1Image else: # Mask params provided diff --git a/junifer/data/tests/test_masks.py b/junifer/data/tests/test_masks.py index ff5394ce9..bd31f2b40 100644 --- a/junifer/data/tests/test_masks.py +++ b/junifer/data/tests/test_masks.py @@ -5,16 +5,17 @@ # Synchon Mandal # License: AGPL +import socket from pathlib import Path from typing import Callable, Dict, List, Optional, Union +import nibabel as nib import numpy as np import pytest from nilearn.datasets import fetch_icbm152_brain_gm_mask from nilearn.image import resample_to_img from nilearn.masking import ( compute_background_mask, - compute_brain_mask, compute_epi_mask, intersect_masks, ) @@ -23,19 +24,97 @@ from numpy.testing import assert_array_almost_equal, assert_array_equal from junifer.data.masks import ( _available_masks, _load_vickery_patil_mask, + compute_brain_mask, get_mask, list_masks, load_mask, register_mask, ) +from junifer.datagrabber import DMCC13Benchmark from junifer.datareader import DefaultDataReader from junifer.pipeline.utils import _check_ants from junifer.testing.datagrabbers import ( OasisVBMTestingDataGrabber, + PartlyCloudyTestingDataGrabber, SPMAuditoryTestingDataGrabber, ) +@pytest.mark.parametrize( + "mask_type, threshold", + [ + ("whole", 0.2), + ("whole", 0.5), + ("whole", 0.8), + ("gm", 0.2), + ("gm", 0.5), + ("gm", 0.8), + ("wm", 0.2), + ("wm", 0.5), + ("wm", 0.8), + ], +) +def test_compute_brain_mask(mask_type: str, threshold: float) -> None: + """Test compute_brain_mask(). + + Parameters + ---------- + mask_type : str + The parametrized mask type. + threshold : float + The parametrized threshold. + + """ + with PartlyCloudyTestingDataGrabber() as dg: + element_data = DefaultDataReader().fit_transform(dg["sub-01"]) + mask = compute_brain_mask( + target_data=element_data["BOLD"], + extra_input=None, + mask_type=mask_type, + ) + assert isinstance(mask, nib.Nifti1Image) + + +@pytest.mark.skipif( + socket.gethostname() != "juseless", + reason="only for juseless", +) +@pytest.mark.parametrize( + "mask_type", + [ + "whole", + "gm", + "wm", + ], +) +def test_compute_brain_mask_for_native(mask_type: str) -> None: + """Test compute_brain_mask(). + + Parameters + ---------- + mask_type : str + The parametrized mask type. + + """ + with DMCC13Benchmark( + types=["BOLD"], + sessions=["wave1bas"], + tasks=["Rest"], + phase_encodings=["AP"], + runs=["1"], + native_t1w=True, + ) as dg: + element_data = DefaultDataReader().fit_transform( + dg[("f1031ax", "wave1bas", "Rest", "AP", "1")] + ) + mask = compute_brain_mask( + target_data=element_data["BOLD"], + extra_input=None, + mask_type=mask_type, + ) + assert isinstance(mask, nib.Nifti1Image) + + def test_register_mask_built_in_check() -> None: """Test mask registration check for built-in masks.""" with pytest.raises(ValueError, match=r"built-in mask"): @@ -228,7 +307,7 @@ def test_get_mask() -> None: raw_mask_callable, _, _ = load_mask( "compute_brain_mask", resolution=1.5 ) - raw_mask_img = raw_mask_callable(vbm_gm_img) # type: ignore + raw_mask_img = raw_mask_callable(vbm_gm) # type: ignore res_mask_img = resample_to_img( raw_mask_img, vbm_gm_img, @@ -313,7 +392,6 @@ def test_get_mask_errors() -> None: @pytest.mark.parametrize( "mask_name,function,params,resample", [ - ("compute_brain_mask", compute_brain_mask, {"threshold": 0.2}, False), ("compute_background_mask", compute_background_mask, None, False), ("compute_epi_mask", compute_epi_mask, None, False), ], @@ -397,9 +475,7 @@ def test_get_mask_inherit() -> None: with SPMAuditoryTestingDataGrabber() as dg: element_data = DefaultDataReader().fit_transform(dg["sub001"]) # Compute brain mask using nilearn - gm_mask = compute_brain_mask( - element_data["BOLD"]["data"], threshold=0.2 - ) + gm_mask = compute_brain_mask(element_data["BOLD"], threshold=0.2) # Get mask using the compute_brain_mask function mask1 = get_mask( @@ -484,7 +560,13 @@ def test_get_mask_multiple( ] for t_func in mask_funcs: - mask_imgs.append(_available_masks[t_func]["func"](target_img)) + # Bypass for custom mask + if t_func == "compute_brain_mask": + mask_imgs.append( + _available_masks[t_func]["func"](element_data["BOLD"]) + ) + else: + mask_imgs.append(_available_masks[t_func]["func"](target_img)) mask_imgs = [ resample_to_img( -- 2.52.0 From 65219c709bd55f0f456f51e9f27267949078383b Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 16 Feb 2024 15:06:33 +0100 Subject: [PATCH 38/47] chore: deprecate fetch_icbm152_brain_gm_mask() --- junifer/data/masks.py | 36 ++++++++++++++++++++-- junifer/data/tests/test_masks.py | 52 +++++--------------------------- 2 files changed, 41 insertions(+), 47 deletions(-) diff --git a/junifer/data/masks.py b/junifer/data/masks.py index bff06a813..0791c97a8 100644 --- a/junifer/data/masks.py +++ b/junifer/data/masks.py @@ -28,7 +28,7 @@ from nilearn.masking import ( ) from ..pipeline import WorkDirManager -from ..utils import logger, raise_error, run_ext_cmd +from ..utils import logger, raise_error, run_ext_cmd, warn_with_log from .template_spaces import get_template, get_xfm from .utils import closest_resolution @@ -124,7 +124,7 @@ def compute_brain_mask( def _fetch_icbm152_brain_gm_mask( target_img: "Nifti1Image", **kwargs, -): +) -> "Nifti1Image": """Fetch ICBM152 brain mask and resample. Parameters @@ -140,7 +140,20 @@ def _fetch_icbm152_brain_gm_mask( nibabel.Nifti1Image The resampled mask. + Warns + ----- + DeprecationWarning + If this function is used. + """ + warn_with_log( + msg=( + "It is recommended to use ``compute_brain_mask`` with " + "``mask_type='gm'``. This function will be removed in the next " + "release. For now, it's available for backward compatibility." + ), + category=DeprecationWarning, + ) mask = fetch_icbm152_brain_gm_mask(**kwargs) mask = resample_to_img( mask, target_img, interpolation="nearest", copy=True @@ -279,7 +292,9 @@ def get_mask( # noqa: C901 Raises ------ RuntimeError - If warp / transformation file extension is not ".mat" or ".h5". + If warp / transformation file extension is not ".mat" or ".h5" or + if fetch_icbm152_brain_gm_mask is used and requires warping to + other template space. ValueError If extra key is provided in addition to mask name in ``masks`` or if no mask is provided or @@ -387,6 +402,21 @@ def get_mask( # noqa: C901 mask_img = extra_input[inherited_mask_item]["data"] # Starting with new mask else: + # Restrict fetch_icbm152_brain_gm_mask if target std space doesn't + # match + if ( + mask_name == "fetch_icbm152_brain_gm_mask" + and target_std_space != "MNI152NLin2009aAsym" + ): + raise_error( + ( + "``fetch_icbm152_brain_gm_mask`` is deprecated and " + "space transformation to any other template space is " + "prohibited as it will lead to unforeseen errors. " + "``compute_brain_mask`` is a better alternative." + ), + klass=RuntimeError, + ) # Load mask mask_object, _, mask_space = load_mask( mask_name, path_only=False, resolution=resolution diff --git a/junifer/data/tests/test_masks.py b/junifer/data/tests/test_masks.py index bd31f2b40..9386d9f05 100644 --- a/junifer/data/tests/test_masks.py +++ b/junifer/data/tests/test_masks.py @@ -12,7 +12,6 @@ from typing import Callable, Dict, List, Optional, Union import nibabel as nib import numpy as np import pytest -from nilearn.datasets import fetch_icbm152_brain_gm_mask from nilearn.image import resample_to_img from nilearn.masking import ( compute_background_mask, @@ -32,7 +31,6 @@ from junifer.data.masks import ( ) from junifer.datagrabber import DMCC13Benchmark from junifer.datareader import DefaultDataReader -from junifer.pipeline.utils import _check_ants from junifer.testing.datagrabbers import ( OasisVBMTestingDataGrabber, PartlyCloudyTestingDataGrabber, @@ -388,6 +386,14 @@ def test_get_mask_errors() -> None: masks="inherit", target_data=vbm_gm, extra_input=extra_input ) + # Block fetch_icbm152_brain_gm_mask space transformation + with pytest.raises(RuntimeError, match="prohibited"): + get_mask( + masks="fetch_icbm152_brain_gm_mask", + target_data=vbm_gm, + extra_input=extra_input, + ) + @pytest.mark.parametrize( "mask_name,function,params,resample", @@ -445,31 +451,6 @@ def test_nilearn_compute_masks( assert_array_equal(mask.get_fdata(), ni_mask.get_fdata()) -@pytest.mark.skipif( - _check_ants() is False, reason="requires ANTs to be in PATH" -) -def test_nilearn_compute_masks_with_space_transformation() -> None: - """Test nilearn compute mask function with space transformation.""" - with SPMAuditoryTestingDataGrabber() as dg: - element_data = DefaultDataReader().fit_transform(dg["sub001"]) - bold = element_data["BOLD"] - bold_img = bold["data"] - - # Get mask - mask = get_mask(masks="fetch_icbm152_brain_gm_mask", target_data=bold) - # Check affine matrices are almost equal - assert_array_equal(mask.affine, bold_img.affine) - - # Mask needs resampling - ni_mask = resample_to_img( - fetch_icbm152_brain_gm_mask(), - bold_img, - interpolation="nearest", - copy=True, - ) - assert_array_equal(mask.get_fdata(), ni_mask.get_fdata()) - - def test_get_mask_inherit() -> None: """Test using the inherit mask functionality.""" with SPMAuditoryTestingDataGrabber() as dg: @@ -580,20 +561,3 @@ def test_get_mask_multiple( expected = intersect_masks(mask_imgs, **params) assert_array_equal(computed.get_fdata(), expected.get_fdata()) - - -@pytest.mark.skipif( - _check_ants() is False, reason="requires ANTs to be in PATH" -) -def test_get_mask_multiple_different_space() -> None: - """Test tailored multiple mask fetch in different space.""" - with SPMAuditoryTestingDataGrabber() as dg: - element_data = DefaultDataReader().fit_transform(dg["sub001"]) - # Get tailored mask - get_mask( - masks=[ - "compute_brain_mask", - "fetch_icbm152_brain_gm_mask", - ], - target_data=element_data["BOLD"], - ) -- 2.52.0 From cca884aeeb7c9be7a6e99836c04e2a6ce6848063 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 16 Feb 2024 15:07:24 +0100 Subject: [PATCH 39/47] chore: add changelog 299.{enh,change} --- docs/changes/newsfragments/299.change | 1 + docs/changes/newsfragments/299.enh | 1 + 2 files changed, 2 insertions(+) create mode 100644 docs/changes/newsfragments/299.change create mode 100644 docs/changes/newsfragments/299.enh diff --git a/docs/changes/newsfragments/299.change b/docs/changes/newsfragments/299.change new file mode 100644 index 000000000..c42f194f3 --- /dev/null +++ b/docs/changes/newsfragments/299.change @@ -0,0 +1 @@ +Add ``template_type`` parameter to :func:`.get_template` by `Synchon Mandal`_ diff --git a/docs/changes/newsfragments/299.enh b/docs/changes/newsfragments/299.enh new file mode 100644 index 000000000..eadf2622a --- /dev/null +++ b/docs/changes/newsfragments/299.enh @@ -0,0 +1 @@ +Rewrite :func:`.compute_brain_mask` to allow variable template fetching via templateflow, according to target data by `Synchon Mandal`_ -- 2.52.0 From b667c80ea25ad67c373608b1f88a438778be8e40 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Wed, 21 Feb 2024 17:31:54 +0100 Subject: [PATCH 40/47] chore: improve docstring for register_parcellation() --- junifer/data/parcellations.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/junifer/data/parcellations.py b/junifer/data/parcellations.py index 605c270b3..4331d41be 100644 --- a/junifer/data/parcellations.py +++ b/junifer/data/parcellations.py @@ -155,7 +155,7 @@ def register_parcellation( parcels_labels : list of str The list of labels for the parcellation. space : str - The space of the parcellation. + The template space of the parcellation, for e.g., "MNI152NLin6Asym". overwrite : bool, optional If True, overwrite an existing parcellation with the same name. Does not apply to built-in parcellations (default False). -- 2.52.0 From 1692454870d80d7c4e75d1360e358c8db20bdd30 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Wed, 21 Feb 2024 17:32:16 +0100 Subject: [PATCH 41/47] chore: improve docstring for register_mask() --- junifer/data/masks.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/junifer/data/masks.py b/junifer/data/masks.py index 0791c97a8..9595a1c53 100644 --- a/junifer/data/masks.py +++ b/junifer/data/masks.py @@ -217,7 +217,7 @@ def register_mask( mask_path : str or pathlib.Path The path to the mask file. space : str - The space of the mask. + The space of the mask, for e.g., "MNI152NLin6Asym". overwrite : bool, optional If True, overwrite an existing mask with the same name. Does not apply to built-in mask (default False). -- 2.52.0 From 7518411e23a6c7aba17bb318e977ea3e301e1716 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Wed, 21 Feb 2024 17:32:33 +0100 Subject: [PATCH 42/47] chore: improve docstring for register_coordinates() --- junifer/data/coordinates.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/junifer/data/coordinates.py b/junifer/data/coordinates.py index 0f7cacaba..b07b42ce4 100644 --- a/junifer/data/coordinates.py +++ b/junifer/data/coordinates.py @@ -138,7 +138,7 @@ def register_coordinates( voi_names : list of str The names of the VOIs. space : str - The space of the coordinates. + The space of the coordinates, for e.g., "MNI". overwrite : bool, optional If True, overwrite an existing list of coordinates with the same name. Does not apply to built-in coordinates (default False). -- 2.52.0 From a0782b5113bde418d7c9826f43b4cd95caf94261 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Tue, 19 Mar 2024 10:46:56 +0100 Subject: [PATCH 43/47] fix: update tests of api.tests.test_functions.py --- junifer/api/tests/test_functions.py | 26 ++++++++++++++++++++++++-- 1 file changed, 24 insertions(+), 2 deletions(-) diff --git a/junifer/api/tests/test_functions.py b/junifer/api/tests/test_functions.py index 2d8909862..728dc510b 100644 --- a/junifer/api/tests/test_functions.py +++ b/junifer/api/tests/test_functions.py @@ -534,13 +534,24 @@ def test_queue_without_elements( assert "Queue done" in caplog.text -def test_reset_run(tmp_path: Path) -> None: +def test_reset_run( + tmp_path: Path, + datagrabber: Dict[str, str], + markers: List[Dict[str, str]], + storage: Dict[str, str], +) -> None: """Test reset function for run. Parameters ---------- tmp_path : pathlib.Path The path to the test directory. + datagrabber : dict + Testing datagrabber as dictionary. + markers : list of dict + Testing markers as list of dictionary. + storage : dict + Testing storage as dictionary. """ # Create storage @@ -567,7 +578,12 @@ def test_reset_run(tmp_path: Path) -> None: ), ) def test_reset_queue( - tmp_path: Path, monkeypatch: pytest.MonkeyPatch, job_name: str + tmp_path: Path, + monkeypatch: pytest.MonkeyPatch, + datagrabber: Dict[str, str], + markers: List[Dict[str, str]], + storage: Dict[str, str], + job_name: str, ) -> None: """Test reset function for queue. @@ -577,6 +593,12 @@ def test_reset_queue( The path to the test directory. monkeypatch : pytest.MonkeyPatch The pytest.MonkeyPatch object. + datagrabber : dict + Testing datagrabber as dictionary. + markers : list of dict + Testing markers as list of dictionary. + storage : dict + Testing storage as dictionary. job_name : str The parametrized job name. -- 2.52.0 From bbd6e4a4c6cf00057cd45fad8e344c4650c3846c Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Fri, 22 Mar 2024 13:47:16 +0100 Subject: [PATCH 44/47] fix: correct file for queue test --- junifer/api/tests/test_cli.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/junifer/api/tests/test_cli.py b/junifer/api/tests/test_cli.py index 36eb49f9d..4efd8c912 100644 --- a/junifer/api/tests/test_cli.py +++ b/junifer/api/tests/test_cli.py @@ -222,7 +222,7 @@ def test_queue(tmp_path: Path, monkeypatch: pytest.MonkeyPatch) -> None: @pytest.mark.parametrize( "action, action_file", [ - (run, "gmd_mean.yaml"), + (run, "partly_cloudy_agg_mean_tian.yml"), (queue, "gmd_mean_htcondor.yaml"), ], ) -- 2.52.0 From e07305fcbb8b7c997bc7081d5fe3b097e6925780 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 28 Mar 2024 10:37:12 +0100 Subject: [PATCH 45/47] update: rename whole to brain for get_template() and compute_brain_mask() --- junifer/data/masks.py | 10 +++++----- junifer/data/template_spaces.py | 6 +++--- junifer/data/tests/test_masks.py | 8 ++++---- 3 files changed, 12 insertions(+), 12 deletions(-) diff --git a/junifer/data/masks.py b/junifer/data/masks.py index 9595a1c53..a7be41ff7 100644 --- a/junifer/data/masks.py +++ b/junifer/data/masks.py @@ -43,7 +43,7 @@ _masks_path = Path(__file__).parent / "masks" def compute_brain_mask( target_data: Dict[str, Any], extra_input: Optional[Dict[str, Any]] = None, - mask_type: str = "whole", + mask_type: str = "brain", threshold: float = 0.5, ) -> "Nifti1Image": """Compute the whole-brain, grey-matter or white-matter mask. @@ -59,14 +59,14 @@ def compute_brain_mask( extra_input : dict, optional The other fields in the data object. Useful for accessing other data types (default None). - mask_type : {"whole", "gm", "wm"}, optional + mask_type : {"brain", "gm", "wm"}, optional Type of mask to be computed: - * "whole" : whole-brain mask + * "brain" : whole-brain mask * "gm" : grey-matter mask * "wm" : white-matter mask - (default "whole"). + (default "brain"). threshold : float, optional The value under which the template is cut off (default 0.5). @@ -84,7 +84,7 @@ def compute_brain_mask( """ logger.debug(f"Computing {mask_type} mask") - if mask_type not in ["whole", "gm", "wm"]: + if mask_type not in ["brain", "gm", "wm"]: raise_error(f"Unknown mask type: {mask_type}") # Check pre-requirements for space manipulation diff --git a/junifer/data/template_spaces.py b/junifer/data/template_spaces.py index eb4a3e6da..251854773 100644 --- a/junifer/data/template_spaces.py +++ b/junifer/data/template_spaces.py @@ -113,7 +113,7 @@ def get_template( extra_input : dict, optional The other fields in the data object. Useful for accessing other data types (default None). - template_type : {"T1w", "whole", "gm", "wm", "csf"}, optional + template_type : {"T1w", "brain", "gm", "wm", "csf"}, optional The template type to retrieve (default "T1w"). Returns @@ -134,7 +134,7 @@ def get_template( raise_error(f"Unknown template space: {space}") # Check for template type - if template_type not in ["T1w", "whole", "gm", "wm", "csf"]: + if template_type not in ["T1w", "brain", "gm", "wm", "csf"]: raise_error(f"Unknown template type: {template_type}") # Get the min of the voxels sizes and use it as the resolution @@ -159,7 +159,7 @@ def get_template( suffix = template_type desc = None label = None - elif template_type == "whole": + elif template_type == "brain": suffix = "mask" desc = "brain" label = None diff --git a/junifer/data/tests/test_masks.py b/junifer/data/tests/test_masks.py index 9386d9f05..ddf448618 100644 --- a/junifer/data/tests/test_masks.py +++ b/junifer/data/tests/test_masks.py @@ -41,9 +41,9 @@ from junifer.testing.datagrabbers import ( @pytest.mark.parametrize( "mask_type, threshold", [ - ("whole", 0.2), - ("whole", 0.5), - ("whole", 0.8), + ("brain", 0.2), + ("brain", 0.5), + ("brain", 0.8), ("gm", 0.2), ("gm", 0.5), ("gm", 0.8), @@ -80,7 +80,7 @@ def test_compute_brain_mask(mask_type: str, threshold: float) -> None: @pytest.mark.parametrize( "mask_type", [ - "whole", + "brain", "gm", "wm", ], -- 2.52.0 From 5723247e9461730f53e686a3c0548dddc78e309e Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Thu, 28 Mar 2024 10:37:12 +0100 Subject: [PATCH 46/47] update: rename whole to brain for get_template() and compute_brain_mask() --- junifer/data/tests/test_template_spaces.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/junifer/data/tests/test_template_spaces.py b/junifer/data/tests/test_template_spaces.py index 7470823bf..d52964259 100644 --- a/junifer/data/tests/test_template_spaces.py +++ b/junifer/data/tests/test_template_spaces.py @@ -40,7 +40,7 @@ def test_get_xfm(tmp_path: Path) -> None: "template_type", [ "T1w", - "whole", + "brain", "gm", "wm", "csf", -- 2.52.0 From 78d63d677aa19941494486debf4b91692fd621a4 Mon Sep 17 00:00:00 2001 From: Synchon Mandal Date: Tue, 2 Apr 2024 10:02:15 +0200 Subject: [PATCH 47/47] chore: update test_ALFFParcels_comparison in test_falff_parcels.py --- junifer/markers/falff/tests/test_falff_parcels.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/junifer/markers/falff/tests/test_falff_parcels.py b/junifer/markers/falff/tests/test_falff_parcels.py index 0d847c8f4..a4a1119ee 100644 --- a/junifer/markers/falff/tests/test_falff_parcels.py +++ b/junifer/markers/falff/tests/test_falff_parcels.py @@ -123,4 +123,4 @@ def test_ALFFParcels_comparison(tmp_path: Path, fractional: bool) -> None: junifer_output_bold["data"][0], afni_output_bold["data"][0], ) - assert r > 0.99 + assert r > 0.97 -- 2.52.0