From 0064f01a273375b0b95cbddfc267095bf348133f Mon Sep 17 00:00:00 2001 From: frheault Date: Thu, 1 Aug 2024 14:59:02 -0400 Subject: [PATCH 1/6] lesions_info now use labels as input --- scripts/scil_lesions_info.py | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) diff --git a/scripts/scil_lesions_info.py b/scripts/scil_lesions_info.py index 96adae6f8..c8004c406 100755 --- a/scripts/scil_lesions_info.py +++ b/scripts/scil_lesions_info.py @@ -20,7 +20,6 @@ import nibabel as nib import numpy as np -import scipy.ndimage as ndi from scilpy.image.labels import get_data_as_labels from scilpy.io.image import get_data_as_mask @@ -40,7 +39,7 @@ def _build_arg_parser(): formatter_class=argparse.RawTextHelpFormatter) p.add_argument('in_lesion', - help='Binary mask of the lesion(s) (.nii.gz).') + help='Lesions file as labels (.nii.gz).') p.add_argument('out_json', help='Output file for lesion information (.json).') p1 = p.add_mutually_exclusive_group() @@ -89,7 +88,8 @@ def main(): reference=args.reference) lesion_img = nib.load(args.in_lesion) - lesion_data = get_data_as_mask(lesion_img, dtype=bool) + voxel_sizes = lesion_img.header.get_zooms()[0:3] + lesion_atlas = get_data_as_labels(lesion_img) if args.bundle: bundle_name, _ = split_name_with_nii(os.path.basename(args.bundle)) @@ -98,7 +98,7 @@ def main(): sft.to_corner() streamlines = sft.get_streamlines_copy() map_data = compute_tract_counts_map(streamlines, - lesion_data.shape) + lesion_atlas.shape) map_data[map_data > 0] = 1 elif args.bundle_mask: bundle_name, _ = split_name_with_nii( @@ -112,8 +112,6 @@ def main(): map_data = get_data_as_labels(map_img) is_single_label = args.bundle_labels_map is None - voxel_sizes = lesion_img.header.get_zooms()[0:3] - lesion_atlas, _ = ndi.label(lesion_data) lesion_load_dict = compute_lesion_stats( map_data, lesion_atlas, single_label=is_single_label, From 9e4c7a9b012f62d377a0a755eaaf09a1819413aa Mon Sep 17 00:00:00 2001 From: frheault Date: Thu, 1 Aug 2024 15:12:10 -0400 Subject: [PATCH 2/6] Improved doc --- scripts/scil_lesions_info.py | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/scripts/scil_lesions_info.py b/scripts/scil_lesions_info.py index c8004c406..b866a59ba 100755 --- a/scripts/scil_lesions_info.py +++ b/scripts/scil_lesions_info.py @@ -3,11 +3,15 @@ """ This script will output informations about lesion load in bundle(s). -The input can either be streamlines, binary bundle map, or a bundle voxel +The input can either be streamlines, binary bundle mask, or a bundle voxel label map. +The input lesion file is a labeled volume (.nii.gz) where each lesion is +represented by a unique label. To label a lesion file use +scil_labels_from_mask.py + To be considered a valid lesion, the lesion volume must be at least -min_lesion_vol mm3. This avoid the detection of thousand of single voxel +--min_lesion_vol mm3. This avoid the detection of thousand of single voxel lesions if an automatic lesion segmentation tool is used. Formerly: scil_analyse_lesions_load.py From 8093f95981eb127a1220e7759049dbb26e5890c2 Mon Sep 17 00:00:00 2001 From: frheault Date: Mon, 5 Aug 2024 08:50:32 -0400 Subject: [PATCH 3/6] Transfer arg --- scilpy/image/labels.py | 12 +++++++++++- scripts/scil_labels_from_mask.py | 10 +++++++++- scripts/scil_lesions_info.py | 2 -- 3 files changed, 20 insertions(+), 4 deletions(-) diff --git a/scilpy/image/labels.py b/scilpy/image/labels.py index 56dc8482d..6e9266d28 100644 --- a/scilpy/image/labels.py +++ b/scilpy/image/labels.py @@ -68,7 +68,8 @@ def get_binary_mask_from_labels(atlas, label_list): return mask -def get_labels_from_mask(mask_data, labels=None, background_label=0): +def get_labels_from_mask(mask_data, labels=None, background_label=0, + min_voxel_count=0): """ Get labels from a binary mask which contains multiple blobs. Each blob will be assigned a label, by default starting from 1. Background will @@ -83,6 +84,9 @@ def get_labels_from_mask(mask_data, labels=None, background_label=0): label. background_label: int Label for the background. + min_voxel_count: int, optional + Minimum number of voxels for a blob to be considered. Blobs with fewer + voxels will be ignored. Returns ------- @@ -91,6 +95,12 @@ def get_labels_from_mask(mask_data, labels=None, background_label=0): """ # Get the number of structures and assign labels to each blob label_map, nb_structures = ndi.label(mask_data) + if min_voxel_count: + for label in range(1, nb_structures + 1): + if np.count_nonzero(label_map == label) < min_voxel_count: + label_map[label_map == label] = 0 + label_map, nb_structures = ndi.label(mask_data) + # Assign labels to each blob if provided if labels: # Only keep the first nb_structures labels if the number of labels diff --git a/scripts/scil_labels_from_mask.py b/scripts/scil_labels_from_mask.py index 3f188a829..e91e60b62 100755 --- a/scripts/scil_labels_from_mask.py +++ b/scripts/scil_labels_from_mask.py @@ -34,6 +34,10 @@ def _build_arg_parser(): p.add_argument('--background_label', default=0, type=int, help='Label to assign to the background. [%(default)s]') + p.add_argument('--min_volume', type=float, default=7, + help='Minimum volume in mm3 [%(default)s],' + 'Useful for lesions.') + add_verbose_arg(p) add_overwrite_arg(p) @@ -50,9 +54,13 @@ def main(): # Load mask and get data mask_img = nib.load(args.in_mask) mask_data = get_data_as_mask(mask_img) + voxel_volume = np.prod(np.diag(mask_img.affine)[:3]) + min_voxel_count = args.min_volume // voxel_volume + # Get labels from mask label_map = get_labels_from_mask( - mask_data, args.labels, args.background_label) + mask_data, args.labels, args.background_label, + min_voxel_count=min_voxel_count) # Save result out_img = nib.Nifti1Image(label_map.astype(np.uint16), mask_img.affine) nib.save(out_img, args.out_labels) diff --git a/scripts/scil_lesions_info.py b/scripts/scil_lesions_info.py index b866a59ba..9760c6e68 100755 --- a/scripts/scil_lesions_info.py +++ b/scripts/scil_lesions_info.py @@ -54,8 +54,6 @@ def _build_arg_parser(): p1.add_argument('--bundle_labels_map', help='Path of the bundle labels map (.nii.gz).') - p.add_argument('--min_lesion_vol', type=float, default=7, - help='Minimum lesion volume in mm3 [%(default)s].') p.add_argument('--out_lesion_atlas', metavar='FILE', help='Save the labelized lesion(s) map (.nii.gz).') p.add_argument('--out_lesion_stats', metavar='FILE', From be71344e48103967149f94496c1c5e8646f6c59a Mon Sep 17 00:00:00 2001 From: frheault Date: Tue, 6 Aug 2024 11:48:22 -0400 Subject: [PATCH 4/6] Improved doc --- scripts/scil_lesions_info.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/scripts/scil_lesions_info.py b/scripts/scil_lesions_info.py index 9760c6e68..3bc8112cb 100755 --- a/scripts/scil_lesions_info.py +++ b/scripts/scil_lesions_info.py @@ -3,13 +3,13 @@ """ This script will output informations about lesion load in bundle(s). -The input can either be streamlines, binary bundle mask, or a bundle voxel -label map. - The input lesion file is a labeled volume (.nii.gz) where each lesion is represented by a unique label. To label a lesion file use scil_labels_from_mask.py +Then, the second input can either be streamlines, binary bundle mask, or a +bundle voxel label map. + To be considered a valid lesion, the lesion volume must be at least --min_lesion_vol mm3. This avoid the detection of thousand of single voxel lesions if an automatic lesion segmentation tool is used. From e27c3577544641484c5f80d8934703e2fddf28a8 Mon Sep 17 00:00:00 2001 From: frheault Date: Tue, 6 Aug 2024 12:01:13 -0400 Subject: [PATCH 5/6] Add missing line --- scilpy/image/labels.py | 1 + 1 file changed, 1 insertion(+) diff --git a/scilpy/image/labels.py b/scilpy/image/labels.py index 6e9266d28..9eaf14736 100644 --- a/scilpy/image/labels.py +++ b/scilpy/image/labels.py @@ -99,6 +99,7 @@ def get_labels_from_mask(mask_data, labels=None, background_label=0, for label in range(1, nb_structures + 1): if np.count_nonzero(label_map == label) < min_voxel_count: label_map[label_map == label] = 0 + mask_data = label_map > 0 label_map, nb_structures = ndi.label(mask_data) # Assign labels to each blob if provided From ca306f8567633c032e1c9a20d2f88069a7a35789 Mon Sep 17 00:00:00 2001 From: frheault Date: Tue, 6 Aug 2024 13:44:18 -0400 Subject: [PATCH 6/6] Emmanuel comments --- scilpy/image/labels.py | 10 ++++++++-- scripts/scil_lesions_info.py | 4 ---- scripts/tests/test_labels_from_mask.py | 2 +- 3 files changed, 9 insertions(+), 7 deletions(-) diff --git a/scilpy/image/labels.py b/scilpy/image/labels.py index 9eaf14736..98b97a989 100644 --- a/scilpy/image/labels.py +++ b/scilpy/image/labels.py @@ -96,11 +96,17 @@ def get_labels_from_mask(mask_data, labels=None, background_label=0, # Get the number of structures and assign labels to each blob label_map, nb_structures = ndi.label(mask_data) if min_voxel_count: + new_count = 0 for label in range(1, nb_structures + 1): if np.count_nonzero(label_map == label) < min_voxel_count: label_map[label_map == label] = 0 - mask_data = label_map > 0 - label_map, nb_structures = ndi.label(mask_data) + else: + new_count += 1 + label_map[label_map == label] = new_count + logging.debug( + f"Ignored blob {nb_structures-new_count} with fewer " + "than {min_voxel_count} voxels") + nb_structures = new_count # Assign labels to each blob if provided if labels: diff --git a/scripts/scil_lesions_info.py b/scripts/scil_lesions_info.py index 3bc8112cb..b4d9301c3 100755 --- a/scripts/scil_lesions_info.py +++ b/scripts/scil_lesions_info.py @@ -10,10 +10,6 @@ Then, the second input can either be streamlines, binary bundle mask, or a bundle voxel label map. -To be considered a valid lesion, the lesion volume must be at least ---min_lesion_vol mm3. This avoid the detection of thousand of single voxel -lesions if an automatic lesion segmentation tool is used. - Formerly: scil_analyse_lesions_load.py """ diff --git a/scripts/tests/test_labels_from_mask.py b/scripts/tests/test_labels_from_mask.py index 27aa8d1d2..335abc195 100644 --- a/scripts/tests/test_labels_from_mask.py +++ b/scripts/tests/test_labels_from_mask.py @@ -24,7 +24,7 @@ def test_execution(script_runner, monkeypatch): 'bundle_4_head_tail_offset.nii.gz') ret = script_runner.run('scil_labels_from_mask.py', in_mask, 'labels_from_mask.nii.gz', - '-f') + '--min_volume', '0', '-f') assert ret.success