Update evaluation code and add script for resizing annotations

Author: constantinpape

flamingo_tools/validation.py

@@ -176,13 +176,21 @@ def compute_matches_for_annotated_slice(
     segmentation_ids = np.unique(segmentation)[1:]
 
     # Crop to the minimal enclosing bounding box of points and segmented objects.
-    bb_seg = np.where(segmentation != 0)
-    bb_seg = tuple(slice(int(bb.min()), int(bb.max())) for bb in bb_seg)
-    bb_points = tuple(
-        slice(int(np.floor(annotations[coords].min())), int(np.ceil(annotations[coords].max())) + 1)
-        for coords in coordinates
-    )
-    bbox = tuple(slice(min(bbs.start, bbp.start), max(bbs.stop, bbp.stop)) for bbs, bbp in zip(bb_seg, bb_points))
+    seg_mask = segmentation != 0
+    if seg_mask.sum() > 0:
+        bb_seg = np.where(seg_mask)
+        bb_seg = tuple(slice(int(bb.min()), int(bb.max())) for bb in bb_seg)
+        bb_points = tuple(
+            slice(int(np.floor(annotations[coords].min())), int(np.ceil(annotations[coords].max())) + 1)
+            for coords in coordinates
+        )
+        bbox = tuple(slice(min(bbs.start, bbp.start), max(bbs.stop, bbp.stop)) for bbs, bbp in zip(bb_seg, bb_points))
+    else:
+        print("The segmentation is empty!!!")
+        bbox = tuple(
+            slice(int(np.floor(annotations[coords].min())), int(np.ceil(annotations[coords].max())) + 1)
+            for coords in coordinates
+        )
     segmentation = segmentation[bbox]
 
     annotations = annotations.copy()
@@ -231,6 +239,11 @@ def compute_scores_for_annotated_slice(
     return {"tp": tp, "fp": fp, "fn": fn}
 
 
+# TODO
+def create_consensus_annotations():
+    pass
+
+
 def for_visualization(segmentation, annotations, matches):
     green_red = ["#00FF00", "#FF0000"]
 

scripts/validation/SGNs/analyze.py

@@ -24,9 +24,15 @@ def compute_scores(table, annotator=None):
 def main():
     parser = argparse.ArgumentParser()
     parser.add_argument("result_file")
+    parser.add_argument("--all", action="store_true")
     args = parser.parse_args()
 
     table = pd.read_csv(args.result_file)
+    if args.all:
+        print(table)
+        print()
+        print()
+
     annotators = pd.unique(table.annotator)
 
     results = []

scripts/validation/SGNs/rescale_annotations.py

@@ -0,0 +1,84 @@
+import os
+import shutil
+from glob import glob
+
+import numpy as np
+import pandas as pd
+import tifffile
+import zarr
+
+from flamingo_tools.s3_utils import get_s3_path, BUCKET_NAME, SERVICE_ENDPOINT
+
+
+def get_scale_factor():
+    original_path = "/mnt/ceph-hdd/cold/nim00007/cochlea-lightsheet/M_LR_000169_R/MLR000169R_PV.tif"
+    original_shape = tifffile.memmap(original_path).shape
+
+    cochlea = "M_LR_000169_R"
+    internal_path = os.path.join(cochlea, "images",  "ome-zarr", "SGN_v2.ome.zarr")
+    s3_store, fs = get_s3_path(internal_path, bucket_name=BUCKET_NAME, service_endpoint=SERVICE_ENDPOINT)
+
+    input_key = "s0"
+    with zarr.open(s3_store, mode="r") as f:
+        new_shape = f[input_key].shape
+
+    scale_factor = tuple(
+        float(nsh) / float(osh) for nsh, osh in zip(new_shape, original_shape)
+    )
+    return scale_factor
+
+
+def rescale_annotations(input_path, scale_factor, bkp_folder):
+    annotations = pd.read_csv(input_path)
+
+    annotations_rescaled = annotations.copy()
+    annotations_rescaled["axis-1"] = annotations["axis-1"] * scale_factor[1]
+    annotations_rescaled["axis-2"] = annotations["axis-2"] * scale_factor[2]
+
+    fname = os.path.basename(input_path)
+    name_components = fname.split("_")
+    z = int(name_components[2][1:])
+    new_z = int(np.round(z * scale_factor[0]))
+
+    name_components[2] = f"z{new_z}"
+    name_components = name_components[:-1] + ["rescaled"] + name_components[-1:]
+    new_fname = "_".join(name_components)
+
+    input_folder = os.path.split(input_path)[0]
+    out_path = os.path.join(input_folder, new_fname)
+    bkp_path = os.path.join(bkp_folder, fname)
+
+    # print(input_path)
+    # print(out_path)
+    # print(bkp_path)
+    # print()
+    # return
+
+    shutil.move(input_path, bkp_path)
+    annotations_rescaled.to_csv(out_path, index=False)
+
+
+def main():
+    # scale_factor = get_scale_factor()
+    # print(scale_factor)
+    scale_factor = (2.6314,) * 3
+
+    root = "/mnt/vast-nhr/projects/nim00007/data/moser/cochlea-lightsheet/AnnotatedImageCrops/F1ValidationSGNs"
+    annotation_folders = ["AnnotationsEK", "AnnotationsAMD", "AnnotationsLR"]
+    for folder in annotation_folders:
+        bkp_folder = os.path.join(root, folder, "rescaled_bkp")
+        os.makedirs(bkp_folder, exist_ok=True)
+
+        files = glob(os.path.join(root, folder, "*.csv"))
+        for annotation_file in files:
+            fname = os.path.basename(annotation_file)
+            if not fname.startswith(("MLR169R_PV_z722", "MLR169R_PV_z979")):
+                continue
+            print("Rescaling", annotation_file)
+            rescale_annotations(annotation_file, scale_factor, bkp_folder)
+
+
+# Rescale the point annotations for the cochlea MLR169R, which was
+# annotated at the original scale, but then rescaled for segmentation.
+if __name__ == "__main__":
+    main()

scripts/validation/SGNs/run_evaluation.py

@@ -7,7 +7,7 @@
 )
 
 ROOT = "/mnt/vast-nhr/projects/nim00007/data/moser/cochlea-lightsheet/AnnotatedImageCrops/F1ValidationSGNs"
-ANNOTATION_FOLDERS = ["AnnotationsEK", "AnnotationsAMD", "AnnotationLR"]
+ANNOTATION_FOLDERS = ["AnnotationsEK", "AnnotationsAMD", "AnnotationsLR"]
 
 
 def run_evaluation(root, annotation_folders, result_file, cache_folder):
@@ -27,11 +27,7 @@ def run_evaluation(root, annotation_folders, result_file, cache_folder):
         annotator = folder[len("Annotations"):]
         annotations = sorted(glob(os.path.join(root, folder, "*.csv")))
         for annotation_path in annotations:
-            print(annotation_path)
             cochlea, slice_id = parse_annotation_path(annotation_path)
-            # We don't have this cochlea in MoBIE yet
-            if cochlea == "M_LR_000169_R":
-                continue
 
             print("Run evaluation for", annotator, cochlea, "z=", slice_id)
             segmentation, annotations = fetch_data_for_evaluation(

scripts/validation/SGNs/visualize_validation.py

@@ -30,8 +30,11 @@ def visualize_anotation(annotation_path, cache_folder):
 
     image_path = _match_image_path(annotation_path)
 
+    # For debugging.
+    components = [1]
+    # components = None
     segmentation, annotations = fetch_data_for_evaluation(
-        annotation_path, cache_path=cache_path, components_for_postprocessing=[1],
+        annotation_path, cache_path=cache_path, components_for_postprocessing=components,
     )
 
     image = tifffile.memmap(image_path)