Fix tonotopic mapping of IHC segmentation

Author: schilling40

flamingo_tools/segmentation/cochlea_mapping.py

@@ -146,7 +146,7 @@ def measure_run_length_sgns(centroids: np.ndarray, scale_factor=10):
     return total_distance, path, path_dict
 
 
-def measure_run_length_ihcs(centroids):
+def measure_run_length_ihcs(centroids, max_edge_distance=50):
     """Measure the run lengths of the IHC segmentation
     by finding the shortest path between the most distant nodes in a Steiner Tree.
 
@@ -159,15 +159,25 @@ def measure_run_length_ihcs(centroids):
         A dictionary containing the position and the length fraction of each point in the path.
     """
     graph = nx.Graph()
-    for num, pos in enumerate(centroids):
+    coords = {}
+    labels = [int(i) for i in range(len(centroids))]
+    for index, element in zip(labels, centroids):
+        coords[index] = element
+
+    for num, pos in coords.items():
         graph.add_node(num, pos=pos)
-    # approximate Steiner tree and find shortest path between the two most distant nodes
-    terminals = set(graph.nodes())  # All nodes are required
-    # Approximate Steiner Tree over all nodes
-    T = steiner_tree(graph, terminals)
-    u, v = find_most_distant_nodes(T)
-    path = nx.shortest_path(T, source=u, target=v)
-    total_distance = nx.path_weight(T, path, weight="weight")
+
+    # create edges between points whose distance is less than threshold max_edge_distance
+    for num_i, pos_i in coords.items():
+        for num_j, pos_j in coords.items():
+            if num_i < num_j:
+                dist = math.dist(pos_i, pos_j)
+                if dist <= max_edge_distance:
+                    graph.add_edge(num_i, num_j, weight=dist)
+
+    u, v = find_most_distant_nodes(graph)
+    path = nx.shortest_path(graph, source=u, target=v)
+    total_distance = nx.path_weight(graph, path, weight="weight")
 
     # assign relative distance to points on path
     path_dict = {}
@@ -180,6 +190,9 @@ def measure_run_length_ihcs(centroids):
         path_dict[num + 1] = {"pos": graph.nodes[p]["pos"], "length_fraction": rel_dist}
     path_dict[len(path)] = {"pos": graph.nodes[path[-1]]["pos"], "length_fraction": 1}
 
+    path_pos = np.array([graph.nodes[p]["pos"] for p in path])
+    path = moving_average_3d(path_pos, window=5)
+
     return total_distance, path, path_dict
 
 

reproducibility/block_extraction/Alphatag_MAMDOTOF1L.json

@@ -0,0 +1,54 @@
+[
+	{
+		"cochlea": "M_AMD_OTOF1_L",
+		"image_channel": [
+			"Apha",
+			"Vglut3",
+			"IHC_v4"
+		],
+		"segmentation_channel": "IHC_v4",
+		"type": "ihc",
+		"n_blocks": 6,
+		"component_list": [
+			5,
+			9
+		],
+		"halo_size": [
+			256,
+			256,
+			50
+		],
+		"crop_centers": [
+			[
+				915,
+				1358,
+				870
+			],
+			[
+				1307,
+				776,
+				764
+			],
+			[
+				687,
+				828,
+				362
+			],
+			[
+				598,
+				1142,
+				1014
+			],
+			[
+				1011,
+				469,
+				1118
+			],
+			[
+				710,
+				35,
+				534
+			]
+		]
+	}
+]

reproducibility/block_extraction/repro_equidistant_centers.py

@@ -32,7 +32,6 @@ def repro_equidistant_centers(
 
     for dic in param_dicts:
         cochlea = dic["cochlea"]
-        img_channel = dic["image_channel"]
         seg_channel = dic["segmentation_channel"]
 
         s3_path = os.path.join(f"{cochlea}", "tables", f"{seg_channel}", "default.tsv")
@@ -50,8 +49,7 @@ def repro_equidistant_centers(
 
         centers = equidistant_centers(table, component_label=component_list, cell_type=cell_type, n_blocks=n_blocks)
         centers = [[int(c) for c in center] for center in centers]
-        ddict = {"cochlea": cochlea}
-        ddict["image_channel"] = img_channel
+        ddict = dic.copy()
         ddict["crop_centers"] = centers
         ddict["halo_size"] = halo_size
         out_dict.append(ddict)