Plot updates

Author: constantinpape

flamingo_tools/extract_block_util.py

@@ -64,6 +64,7 @@ def extract_block(
 
     if output_dir == "":
         output_dir = input_dir
+    os.makedirs(output_dir, exist_ok=True)
 
     if tif:
         if output_key is None:

reproducibility/templates_processing/segment_unet_IHC_template.sbatch

@@ -8,7 +8,8 @@
 #SBATCH --mem 400G
 
 source ~/.bashrc
-micromamba activate micro-sam_gpu
+# micromamba activate micro-sam_gpu
+micromamba activate sam
 
 # Print out some info.
 echo "Submitting job with sbatch from directory: ${SLURM_SUBMIT_DIR}"
@@ -19,7 +20,8 @@ echo "Current node: ${SLURM_NODELIST}"
 # Run the script
 #python myprogram.py $SLURM_ARRAY_TASK_ID
 
-SCRIPT_REPO=/user/schilling40/u15000/flamingo-tools
+# SCRIPT_REPO=/user/schilling40/u15000/flamingo-tools
+SCRIPT_REPO=/user/pape41/u12086/Work/my_projects/flamingo-tools
 cd "$SCRIPT_REPO"/flamingo_tools/segmentation/ || exit
 
 # name of cochlea, as it appears in MoBIE and the NHR

scripts/figures/plot_fig2.py

@@ -66,7 +66,7 @@ def plot_seg_crop(img_path, seg_path, save_path, xlim1, xlim2, ylim1, ylim2, bou
         plt.close()
 
 
-def fig_02b_sgn(save_dir, plot=False):
+def fig_02a_sgn(save_dir, plot=False):
     """Plot crops of SGN segmentation of CochleaNet, Cellpose and micro-sam.
     """
     cochlea_dir = "/mnt/vast-nhr/projects/nim00007/data/moser/cochlea-lightsheet"
@@ -402,9 +402,12 @@ def main():
 
     os.makedirs(args.figure_dir, exist_ok=True)
 
-    # Panel C: Evaluation of the segmentation results:
-    fig_02b_sgn(save_dir=args.figure_dir, plot=args.plot)
+    # Panes A and B: Qualitative comparison of visualization results.
+    fig_02a_sgn(save_dir=args.figure_dir, plot=args.plot)
+    return
     fig_02b_ihc(save_dir=args.figure_dir, plot=args.plot)
+
+    # Panel C: Evaluation of the segmentation results:
     fig_02c(save_path=os.path.join(args.figure_dir, "fig_02c"), plot=args.plot, all_versions=False)
 
     # Panel D: The number of SGNs, IHCs and average number of ribbon synapses per IHC

scripts/figures/segmentation_comparison.py

@@ -0,0 +1,113 @@
+import os
+from glob import glob
+from pathlib import Path
+import json
+
+import imageio.v3 as imageio
+import napari
+import numpy as np
+from skimage.segmentation import find_boundaries
+
+FOR_COMPARISON = ["distance_unet", "micro-sam", "cellpose3"]
+
+
+def _eval_seg(seg, eval_path):
+    with open(eval_path, "r") as f:
+        eval_res = json.load(f)
+
+    correct, wrong = eval_res["tp_objects"], eval_res["fp"]
+    all_ids = correct + wrong
+    seg[~np.isin(seg, all_ids)] = 0
+
+    eva_mask = np.zeros_like(seg)
+
+    eva_mask[np.isin(seg, correct)] = 1
+    eva_mask[np.isin(seg, wrong)] = 2
+
+    bd = find_boundaries(seg)
+    return bd, eva_mask
+
+
+def sgn_comparison():
+    z = 10
+
+    cochlea_dir = "/mnt/vast-nhr/projects/nim00007/data/moser/cochlea-lightsheet"
+    val_sgn_dir = f"{cochlea_dir}/predictions/val_sgn"
+    image_dir = f"{cochlea_dir}/AnnotatedImageCrops/F1ValidationSGNs/for_consensus_annotation"
+
+    image_paths = sorted(glob(os.path.join(image_dir, "*.tif")))
+
+    for path in image_paths:
+        image = imageio.imread(path)[z]
+
+        seg_fname = Path(path).stem + "_seg.tif"
+        eval_fname = Path(path).stem + "_dic.json"
+
+        segmentations, boundaries, eval_im = {}, {}, {}
+        for seg_name in FOR_COMPARISON:
+            seg_path = os.path.join(val_sgn_dir, seg_name, seg_fname)
+            eval_path = os.path.join(val_sgn_dir, seg_name, eval_fname)
+            assert os.path.exists(seg_path), seg_path
+
+            seg = imageio.imread(seg_path)[z]
+
+            bd, eva = _eval_seg(seg, eval_path)
+            segmentations[seg_name] = seg
+            boundaries[seg_name] = bd
+            eval_im[seg_name] = eva
+
+        v = napari.Viewer()
+        v.add_image(image)
+        for seg_name, bd in boundaries.items():
+            v.add_labels(bd, name=seg_name, colormap={1: "cyan"})
+            v.add_labels(eval_im[seg_name], name=f"{seg_name}_eval", colormap={1: "green", 2: "red"})
+        v.title = Path(path).stem
+        napari.run()
+
+
+def ihc_comparison():
+    z = 10
+
+    cochlea_dir = "/mnt/vast-nhr/projects/nim00007/data/moser/cochlea-lightsheet"
+    val_sgn_dir = f"{cochlea_dir}/predictions/val_ihc"
+    image_dir = f"{cochlea_dir}/AnnotatedImageCrops/F1ValidationIHCs"
+
+    image_paths = sorted(glob(os.path.join(image_dir, "*.tif")))
+
+    for path in image_paths:
+        image = imageio.imread(path)[z]
+
+        seg_fname = Path(path).stem + "_seg.tif"
+        eval_fname = Path(path).stem + "_dic.json"
+
+        segmentations, boundaries, eval_im = {}, {}, {}
+        for seg_name in FOR_COMPARISON:
+            # FIXME distance_unet_v4b is missing the eval files
+            seg_name_ = "distance_unet_v3" if seg_name == "distance_unet" else seg_name
+            seg_path = os.path.join(val_sgn_dir, seg_name_, seg_fname)
+            eval_path = os.path.join(val_sgn_dir, seg_name_, eval_fname)
+            assert os.path.exists(seg_path), seg_path
+
+            seg = imageio.imread(seg_path)[z]
+
+            bd, eva = _eval_seg(seg, eval_path)
+            segmentations[seg_name] = seg
+            boundaries[seg_name] = bd
+            eval_im[seg_name] = eva
+
+        v = napari.Viewer()
+        v.add_image(image)
+        for seg_name, bd in boundaries.items():
+            v.add_labels(bd, name=seg_name, colormap={1: "cyan"})
+            v.add_labels(eval_im[seg_name], name=f"{seg_name}_eval", colormap={1: "green", 2: "red"})
+        v.title = Path(path).stem
+        napari.run()
+
+
+def main():
+    # sgn_comparison()
+    ihc_comparison()
+
+
+if __name__ == "__main__":
+    main()

scripts/more-annotations/extract_blocks_da_mlr99l.py

@@ -0,0 +1,26 @@
+from flamingo_tools.extract_block_util import extract_block
+
+
+positions_with_signal = [
+    [1215.074063453085, 912.697256780485, 1036.814204517708],
+    [1030.351117830933, 1262.3358840155736, 1123.2581736686361],
+    [1192.167776682008, 354.058713359485, 767.1544606203263],
+    [916.9294364078347, 754.7061965177552, 923.607923806173],
+]
+positions_without_signal = [
+    [1383.4288658807268, 783.0008672288084, 467.5426478786816],
+]
+
+halo = [256, 256, 64]
+
+
+for pos in positions_with_signal + positions_without_signal:
+    extract_block(
+        input_path="M_LR_000099_L/images/ome-zarr/PV.ome.zarr",
+        coords=pos,
+        output_dir="./MLR99L_for_DA",
+        input_key="s0",
+        roi_halo=halo,
+        tif=True,
+        s3=True
+    )