refactor: remove immutable edge notion, upd coord conversions (#51)

Co-authored-by: anna-grim <[email protected]>

src/deep_neurographs/feature_extraction.py

@@ -97,7 +97,8 @@ def generate_img_chunks(
 def generate_img_chunks_via_multithreads(
     neurograph, img_path, labels_path, proposals=None
 ):
-    img = utils.open_tensorstore(img_path, "zarr")
+    driver = "n5" if ".n5" in img_path else "zarr"
+    img = utils.open_tensorstore(img_path, driver)
     labels = utils.open_tensorstore(labels_path, "neuroglancer_precomputed")
     with ThreadPoolExecutor() as executor:
         # Assign Threads
@@ -155,9 +156,9 @@ def generate_img_chunks_via_superchunk(
         # Compute image coordinates
         i, j = tuple(edge)
         xyz = gutils.get_edge_attr(neurograph, edge, "xyz")
-        xyz_i = utils.world_to_img(neurograph, xyz[0])
-        xyz_j = utils.world_to_img(neurograph, xyz[1])
-        chunk, profile = get_img_chunk_features(img, labels, xyz_i, xyz_j)
+        coord_i = utils.to_img(xyz[0])
+        coord_j = utils.to_img(xyz[1])
+        chunk, profile = get_img_chunk_features(img, labels, coord_i, coord_j)
         chunk_features[edge] = chunk
         profile_features[edge] = profile
     return chunk_features, profile_features
@@ -203,19 +204,18 @@ def generate_img_profiles(neurograph, path, proposals=None):
         )
 
 
-def generate_img_profiles_via_multithreads(
-    neurograph, img_path, proposals=None
-):
+def generate_img_profiles_via_multithreads(neurograph, path, proposals=None):
     profile_features = dict()
-    img = utils.open_tensorstore(img_path, "zarr")
+    driver = "n5" if ".n5" in path else "zarr"
+    img = utils.open_tensorstore(path, driver)
     with ThreadPoolExecutor() as executor:
         # Assign threads
         threads = [None] * len(proposals)
         for i, edge in enumerate(proposals):
             xyz_i, xyz_j = gutils.get_edge_attr(neurograph, edge, "xyz")
-            xyz_i = utils.world_to_img(neurograph, xyz_i)
-            xyz_j = utils.world_to_img(neurograph, xyz_j)
-            line = geometry.make_line(xyz_i, xyz_j, N_PROFILE_POINTS)
+            coord_i = utils.to_img(xyz_i)
+            coord_j = utils.to_img(xyz_j)
+            line = geometry.make_line(coord_i, coord_j, N_PROFILE_POINTS)
             threads[i] = executor.submit(geometry.get_profile, img, line, edge)
 
         # Store result
@@ -249,16 +249,16 @@ def generate_img_profiles_via_superchunk(neurograph, path, proposals=None):
 
     """
     features = dict()
-    origin = utils.apply_anisotropy(neurograph.origin, return_int=True)
+    driver = "n5" if ".n5" in path else "zarr"
     img = utils.get_superchunk(
-        path, "zarr", origin, neurograph.shape, from_center=False
+        path, driver, neurograph.origin, neurograph.shape, from_center=False
     )
     img = utils.normalize_img(img)
     for edge in neurograph.mutable_edges:
         xyz_i, xyz_j = neurograph.get_edge_attr(edge, "xyz")
-        xyz_i = utils.world_to_img(neurograph, xyz_i)
-        xyz_j = utils.world_to_img(neurograph, xyz_j)
-        path = geometry.make_line(xyz_i, xyz_j, N_PROFILE_POINTS)
+        coord_i = utils.to_img(xyz_i) - neurograph.origin
+        coord_j = utils.to_img(xyz_j) - neurograph.origin
+        path = geometry.make_line(coord_i, coord_j, N_PROFILE_POINTS)
         features[edge] = geometry.get_profile(img, path, window=WINDOW)
     return features
 
@@ -270,8 +270,8 @@ def generate_skel_features(neurograph, proposals=None):
         features[edge] = np.concatenate(
             (
                 neurograph.compute_length(edge),
-                neurograph.immutable_degree(i),
-                neurograph.immutable_degree(j),
+                neurograph.get_degree_temp(i),
+                neurograph.get_degree_temp(j),
                 get_radii(neurograph, edge),
                 get_avg_radii(neurograph, edge),
                 get_avg_branch_lens(neurograph, edge),
@@ -413,7 +413,7 @@ def get_feature_vectors(neurograph, features, shift=0):
     features = combine_features(features)
     features.keys()
     key = sample(list(features.keys()), 1)[0]
-    n_edges = neurograph.num_mutables()
+    n_edges = neurograph.n_mutables()
     n_features = len(features[key])
 
     # Build
@@ -429,7 +429,7 @@ def get_feature_vectors(neurograph, features, shift=0):
 
 def get_multimodal_features(neurograph, features, shift=0):
     idx_to_edge = dict()
-    n_edges = neurograph.num_mutables()
+    n_edges = neurograph.n_mutables()
     X = np.zeros(((n_edges, 2) + tuple(CHUNK_SIZE)))
     x = np.zeros((n_edges, N_SKEL_FEATURES + N_PROFILE_POINTS))
     y = np.zeros((n_edges))
@@ -445,7 +445,7 @@ def get_multimodal_features(neurograph, features, shift=0):
 
 def get_img_chunks(neurograph, features, shift=0):
     idx_to_edge = dict()
-    n_edges = neurograph.num_mutables()
+    n_edges = neurograph.n_mutables()
     X = np.zeros(((n_edges, 2) + tuple(CHUNK_SIZE)))
     y = np.zeros((n_edges))
     for i, edge in enumerate(features["img_chunks"].keys()):

src/deep_neurographs/geometry.py

@@ -22,7 +22,7 @@ def get_directional(
         elif branch.shape[0] <= d:
             xyz = deepcopy(branch)
         else:
-            xyz = deepcopy(branch[d : window + d, :])
+            xyz = deepcopy(branch[d: window + d, :])
         directionals.append(compute_tangent(xyz))
 
     # Determine best
@@ -38,6 +38,28 @@ def get_directional(
 
 
 def compute_svd(xyz):
+    """
+    Compute singular value decomposition (svd) of an NxD array where N is the
+    number of points and D is the dimension of the space.
+
+    Parameters
+    ----------
+    xyz : numpy.ndarray
+        Array containing data points.
+
+    Returns
+    -------
+    numpy.ndarry
+        Unitary matrix having left singular vectors as columns. Of shape
+        (N, N) or (N, min(N, D)), depending on full_matrices.
+    numpy.ndarray
+        Singular values, sorted in non-increasing order. Of shape (K,), with
+        K = min(N, D).
+    numpy.ndarray
+        Unitary matrix having right singular vectors as rows. Of shape (D, D)
+        or (K, D) depending on full_matrices.
+
+    """
     xyz = xyz - np.mean(xyz, axis=0)
     return svd(xyz)
 
@@ -242,9 +264,9 @@ def align(neurograph, img, branch_1, branch_2, depth):
     best_d2 = None
     best_score = 0
     for d1 in range(min(depth, len(branch_1) - 1)):
-        xyz_1 = neurograph.to_img(branch_1[d1])
+        xyz_1 = neurograph.to_img(branch_1[d1], shift=True)
         for d2 in range(min(depth, len(branch_2) - 1)):
-            xyz_2 = neurograph.to_img(branch_2[d2])
+            xyz_2 = neurograph.to_img(branch_2[d2], shift=True)
             line = make_line(xyz_1, xyz_2, 10)
             score = np.mean(get_profile(img, line, window=[3, 3, 3]))
             if score > best_score:
@@ -481,6 +503,5 @@ def make_line(xyz_1, xyz_2, n_steps):
     return np.array([(1 - t) * xyz_1 + t * xyz_2 for t in t_steps], dtype=int)
 
 
-def normalize(x, norm="l2"):
-    zero_vec = np.zeros((3))
-    return x / abs(dist(zero_vec, x, metric=norm))
+def normalize(vec, norm="l2"):
+    return vec / abs(dist(np.zeros((3)), vec, metric=norm))

src/deep_neurographs/graph_utils.py

@@ -316,7 +316,7 @@ def init_edge_attrs(swc_dict, i):
         Edge attribute dictionary.
 
     """
-    j = swc_dict["idx"][i] 
+    j = swc_dict["idx"][i]
     return {"radius": [swc_dict["radius"][j]], "xyz": [swc_dict["xyz"][j]]}
 
 
@@ -339,7 +339,7 @@ def upd_edge_attrs(swc_dict, attrs, i):
         Edge attribute dictionary.
 
     """
-    j = swc_dict["idx"][i] 
+    j = swc_dict["idx"][i]
     attrs["radius"].append(swc_dict["radius"][j])
     attrs["xyz"].append(swc_dict["xyz"][j])
     return attrs
@@ -406,7 +406,7 @@ def set_node_attrs(swc_dict, nodes):
     """
     attrs = dict()
     for i in nodes:
-        j = swc_dict["idx"][i] 
+        j = swc_dict["idx"][i]
         attrs[i] = {"radius": swc_dict["radius"][j], "xyz": swc_dict["xyz"][j]}
     return attrs
 
@@ -438,10 +438,38 @@ def upd_node_attrs(swc_dict, leafs, junctions, i):
         Updated dictionary if "i" was contained in "junctions.keys()".
 
     """
-    j = swc_dict["idx"][i] 
+    j = swc_dict["idx"][i]
     upd_attrs = {"radius": swc_dict["radius"][j], "xyz": swc_dict["xyz"][j]}
     if i in leafs:
         leafs[i] = upd_attrs
     else:
         junctions[i] = upd_attrs
     return leafs, junctions
+
+
+# -- miscellaneous --
+def creates_cycle(graph, edge):
+    """
+    Checks whether adding "edge" to "graph" creates a cycle.
+
+    Paramaters
+    ----------
+    graph : networkx.Graph
+        Graph to be checked for cycles.
+    edge : tuple
+        Edge to be added to "graph"
+
+    Returns
+    -------
+    bool
+        Indication of whether adding "edge" to graph creates a cycle.
+
+    """
+    graph.add_edges_from([edge])
+    try:
+        nx.find_cycle(graph)
+        graph.remove_edges_from([edge])
+        return True
+    except:
+        graph.remove_edges_from([edge])
+        return False

src/deep_neurographs/intake.py

@@ -48,14 +48,14 @@ def build_neurograph_from_local(
 ):
     # Process swc files
     assert swc_dir or swc_paths, "Provide swc_dir or swc_paths!"
-    bbox = utils.get_bbox(img_patch_origin, img_patch_shape)
+    img_bbox = utils.get_img_bbox(img_patch_origin, img_patch_shape)
     paths = get_paths(swc_dir) if swc_dir else swc_paths
-    swc_dicts = process_local_paths(paths, min_size, bbox=bbox)
+    swc_dicts = process_local_paths(paths, min_size, img_bbox=img_bbox)
 
     # Build neurograph
     neurograph = build_neurograph(
         swc_dicts,
-        bbox=bbox,
+        img_bbox=img_bbox,
         img_path=img_path,
         swc_paths=paths,
         progress_bar=progress_bar,
@@ -207,8 +207,6 @@ def download_gcs_zips(bucket_name, cloud_path, min_size):
             cnt, t1 = report_progress(
                 i, len(zip_paths), chunk_size, cnt, t0, t1
             )
-        if len(swc_dicts) > 2000:
-            stop
     return swc_dicts
 
 
@@ -223,7 +221,7 @@ def count_files_in_zips(bucket, zip_paths):
 # -- Build neurograph ---
 def build_neurograph(
     swc_dicts,
-    bbox=None,
+    img_bbox=None,
     img_path=None,
     swc_paths=None,
     progress_bar=True,
@@ -249,13 +247,16 @@ def build_neurograph(
         print("(2) Combine irreducibles...")
         print("# nodes:", utils.reformat_number(n_nodes))
         print("# edges:", utils.reformat_number(n_edges))
-    neurograph = NeuroGraph(bbox=bbox, img_path=img_path, swc_paths=swc_paths)
+
+    neurograph = NeuroGraph(
+        img_bbox=img_bbox, img_path=img_path, swc_paths=swc_paths
+    )
     t0, t1 = utils.init_timers()
     chunk_size = max(int(n_components * 0.05), 1)
     cnt, i = 1, 0
     while len(irreducibles):
         key, irreducible_set = irreducibles.popitem()
-        neurograph.add_immutables(irreducible_set, key)
+        neurograph.add_component(irreducible_set, key)
         if i > cnt * chunk_size and progress_bar:
             cnt, t1 = report_progress(i, n_components, chunk_size, cnt, t0, t1)
         i += 1