@@ -39,14 +39,71 @@ physical units (nm / nm² / nm³) using the voxel size.
3939
4040---
4141
42+ ## Command-line interface
43+
44+ The same analysis is available headlessly via the ` synapse_net.run_cristae_analysis ` console script
45+ (` cristae_analysis_cli ` in ` synapse_net/tools/cli.py ` ). It takes a crista segmentation and a
46+ mitochondria instance segmentation, computes ` compute_mito_crista_statistics ` , and writes one
47+ ` <stem>_cristae_analysis.csv ` per input pair. It accepts a single file ** or** a directory for each
48+ input (batch); for a directory the input subfolder structure is mirrored under the output, and
49+ already-present tables are skipped unless ` --force ` is given.
50+
51+ Segmentations are read as tif when no key is given, or from an hdf5 dataset when a key is given — so
52+ the crista and mito segmentations can live in separate files or in the * same* hdf5 file under
53+ different keys.
54+
55+ | Flag | Default | Meaning |
56+ | ---| ---| ---|
57+ | ` --crista_path ` / ` -c ` | required | Crista segmentation — file or directory. |
58+ | ` --mito_path ` / ` -m ` | required | Mitochondria instance segmentation — file or directory. |
59+ | ` --output_path ` / ` -o ` | required | Output directory for the result table(s). |
60+ | ` --crista_key ` | — | HDF5 dataset key for the crista segmentation (omit → tif). |
61+ | ` --mito_key ` | — | HDF5 dataset key for the mito segmentation (omit → tif). |
62+ | ` --voxel_size ` | — | Voxel size in ** nm** , applied to all inputs. |
63+ | ` --tomogram_path ` | — | Raw tomogram (mrc/rec), file or directory, used to read the voxel size when ` --voxel_size ` is omitted. |
64+ | ` --membrane_thickness ` | ` 8.0 ` | Membrane shell thickness (nm). |
65+ | ` --border_gap ` | thickness | Distance from volume faces where the membrane is suppressed (nm). |
66+ | ` --method ` | ` skip ` | Orientation anisotropy mode: ` skip ` / ` fast ` / ` exact ` . |
67+ | ` --membrane_mode ` | ` slice_2d ` | Membrane shell construction: ` slice_2d ` / ` shell_3d ` . |
68+ | ` --n_jobs ` | ` -1 ` | Workers for the per-mitochondrion computation (` -1 ` = all cores). |
69+ | ` --force ` | off | Over-write already present result tables. |
70+ | ` --verbose ` / ` -v ` | off | Show a progress bar over the mitochondria of each file. |
71+
72+ One of ` --voxel_size ` or ` --tomogram_path ` must be given (the physical-unit results are wrong without
73+ the correct voxel size — see the units caveat above).
74+
75+ Examples:
76+
77+ ``` bash
78+ # Single pair of tif segmentations, explicit voxel size.
79+ synapse_net.run_cristae_analysis -c crista.tif -m mito.tif -o results/ --voxel_size 0.868 -v
80+
81+ # Both segmentations in one hdf5 file under different keys.
82+ synapse_net.run_cristae_analysis \
83+ -c seg.h5 --crista_key labels/cristae \
84+ -m seg.h5 --mito_key labels/mitochondria \
85+ -o results/ --voxel_size 0.868
86+
87+ # Directory batch, reading the voxel size from the matching raw tomograms.
88+ synapse_net.run_cristae_analysis \
89+ -c crista_dir/ -m mito_dir/ --tomogram_path tomo_dir/ -o results/
90+ ```
91+
92+ The output CSV has one row per mitochondrion with the columns documented in ** Output columns** below.
93+
94+ ---
95+
4296## Computational steps & math
4397
4498### 1. Membrane approximation — ` approximate_membrane `
4599The membrane is approximated as the ** outer boundary shell of the mitochondrion** . For 3D data it is
46100done ** per Z‑slice** with a 2D disk of radius ` mm_thickness / voxel_xy ` :
47101` membrane = mito & ~binary_erosion(mito, disk(r)) ` — i.e. the rim removed by the erosion. Per‑slice
48- (not a 3D ball) prevents the shell bleeding across slices when the shape changes in Z. Voxels within
49- ` border_gap ` of any volume face are removed (so clipped mito edges aren't treated as membrane).
102+ (not a 3D ball) prevents the shell bleeding across slices when the shape changes in Z. A separable
103+ Z‑only erosion (radius ` mm_thickness / voxel_z ` ) then adds ** Z‑caps** where a mito column truly ends in
104+ Z; ends clipped by a volume Z‑face stay uncapped (` border_value=1 ` plus the ` border_gap ` trim below).
105+ Voxels within ` border_gap ` of any volume face are removed (so clipped mito edges aren't treated as
106+ membrane).
50107The eroded interior is the ** lumen** , whose inner‑wall surface is the mesh the junction geodesics run
51108along (see §5 and the ** Membrane Mesh** layer). At mesh time the lumen surface is trimmed to the
52109certain region — the border zone within ` border_gap ` of the volume faces is cropped off, so it ends
@@ -61,8 +118,7 @@ Computes the **structure-tensor eigenvalues** of the binary crista mask via
61118derivative scale (a minimal 1 voxel), ` outer_sigma = neighborhood_size_nm / voxel ` is the integration
62119scale (default 30 nm). The per‑voxel eigenvalues (descending) give
63120** anisotropy = λ_max / (λ_min + ε)** , and the reported ` crista_orientation_anisotropy ` is the mean
64- over crista voxels. When bioimage‑cpp lacks this function, a ** NumPy fallback** builds the tensor with
65- ` np.gradient ` + Gaussian smoothing and a low‑memory chunked ` eigvalsh ` .
121+ over crista voxels.
66122* Interpretation:* high → strongly directional (parallel lamellae), ~ 1 → isotropic/tubular. It is a
67123** magnitude, rotation‑invariant** — it says * how* laminar, not * which* direction (the orientation
68124direction / eigenvectors are not computed).
@@ -94,14 +150,15 @@ nearest‑neighbour distances (`mean_/median_nn_junction_distance_nm`; `NaN` whe
94150junctions) and a ** Clark‑Evans clustering index** ` R = mean_NN / (0.5 · √(A/n)) ` with ` A ` = mito
95151surface area and ` n ` = junction count (` junction_clustering_index ` : <1 clustered, ≈1 random,
96152> 1 dispersed; flat‑surface CSR approximation).
97- > When the bioimage‑cpp geodesic API is unavailable the junction‑distance columns are ` NaN ` and a
98- > single warning is emitted — ** there is no other backend** .
153+ > The bioimage‑cpp geodesic backend (` bioimage_cpp.distance.geodesic_distances_mesh ` ) is a hard
154+ > dependency — ** there is no other backend** . A mito with no usable lumen mesh (empty/degenerate)
155+ > yields ` NaN ` junction‑distance columns.
99156
100157** "Missing" pairs.** Two junctions on ** different connected components** of the lumen surface have no
101158along‑surface path → that pair is ` NaN ` and is ignored by the nearest‑neighbour / clustering
102- summaries. ` membrane_mode="slice_2d" ` (the default) has no Z‑caps and fragments across slices more
103- readily than ` "shell_3d" ` . If junctions you expect to be connected are not, switch to ` "shell_3d" ` ,
104- increase ` mm_thickness ` , or check the mito segmentation.
159+ summaries. ` membrane_mode="slice_2d" ` (the default) caps true (non‑clipped) Z‑ends but its XY shell
160+ can still fragment across slices more readily than ` "shell_3d" ` . If junctions you expect to be
161+ connected are not, switch to ` "shell_3d" ` , increase ` mm_thickness ` , or check the mito segmentation.
105162
106163** What shifts with the membrane settings.** The membrane‑dependent outputs — ` crista_junction_count ` ,
107164the contact counts, ` avg_crista_to_membrane_nm ` , and the junction distances — change with
@@ -126,7 +183,7 @@ settings.
126183| Step | Function | File |
127184| ---| ---| ---|
128185| Membrane shell (` slice_2d ` / ` shell_3d ` ) | ` approximate_membrane ` | ` synapse_net/cristae_analysis.py ` |
129- | Crista orientation anisotropy | ` compute_crista_orientation ` (` bioimage_cpp.filters.structure_tensor_eigenvalues ` , NumPy fallback ; ` _downsampled_orientation_anisotropy ` for ` fast ` ) | ` synapse_net/cristae_analysis.py ` |
186+ | Crista orientation anisotropy | ` compute_crista_orientation ` (` bioimage_cpp.filters.structure_tensor_eigenvalues ` ; ` _downsampled_orientation_anisotropy ` for ` fast ` ) | ` synapse_net/cristae_analysis.py ` |
130187| Crista→membrane proximity | ` compute_crista_proximity ` | ` synapse_net/cristae_analysis.py ` |
131188| Junctions (crista ∩ membrane) | ` detect_contact_sites ` | ` synapse_net/cristae_analysis.py ` |
132189| Junction geodesic distances | ` compute_junction_distances ` → ` _junction_matrix_mesh ` → ` bioimage_cpp.distance.geodesic_distances_mesh ` | ` synapse_net/cristae_analysis.py ` |
@@ -171,12 +228,12 @@ labels layer.)
171228---
172229
173230## Libraries used
174- - ** NumPy** — arrays, gradients, linear algebra ( ` eigvalsh ` ) .
175- - ** SciPy** — ` scipy.ndimage ` (` binary_erosion ` , ` distance_transform_edt ` , ` gaussian_filter ` ,
176- ` label ` , ` center_of_mass ` ); ` scipy.spatial.cKDTree ` (snap junction centroids to mesh vertices).
177- - ** bioimage‑cpp** — ` bioimage_cpp.distance.geodesic_distances_mesh ` (≥ 0.6.0) for junction surface
178- geodesics, and ` bioimage_cpp.filters.structure_tensor_eigenvalues ` for the crista orientation
179- anisotropy (NumPy fallback if unavailable) .
231+ - ** NumPy** — arrays and numerics .
232+ - ** SciPy** — ` scipy.ndimage ` (` binary_erosion ` , ` distance_transform_edt ` , ` label ` ,
233+ ` center_of_mass ` ); ` scipy.spatial.cKDTree ` (snap junction centroids to mesh vertices).
234+ - ** bioimage‑cpp** (≥ 0.6.0, required) — ` bioimage_cpp.distance.geodesic_distances_mesh ` for junction
235+ surface geodesics, and ` bioimage_cpp.filters.structure_tensor_eigenvalues ` for the crista
236+ orientation anisotropy .
180237- ** scikit‑image** — ` measure.marching_cubes ` , ` measure.mesh_surface_area ` , ` measure.regionprops ` ;
181238 ` morphology.disk ` , ` morphology.local_maxima ` .
182239- ** pandas** — results table. ** tqdm** — progress. ** napari** / ** qtpy** — the widget/UI.
0 commit comments