Morphometry Pipeline

​

Overview

• Volume - Estimated using winding number and Monte Carlo sampling
• Surface area - Computed from mesh vertex areas
• Sphericity - Shape measure comparing volume to surface area
• Principal components - Orientation and elongation of structures
• Distance transform - Distance from interior points to surface
• Spatial location - Centroid position of each component

​

Functions

​

component_morphometry_pipeline

component_morphometry_pipeline(
    mesh,
    labels,
    select_label,
    post_synapse_mappings=None,
    vertex_features=None,
    split_laplacian="graph",
    split_threshold=None,
    split_min_size=10,
    bound_volume_threshold=None,
    verbose=False,
)

• results_df: DataFrame with one row per component and columns for all measurements
• corrected_components: Array of component labels after splitting, same length as vertices

import numpy as np
import pandas as pd
from meshmash.morphometry_pipeline import component_morphometry_pipeline
from meshmash.cave import get_synapse_mapping

# Assume you have a mesh and domain labels from the HKS pipeline
mesh = (vertices, faces)  # Your neuron mesh
labels = np.load('domain_labels.npy')  # Labels from HKS pipeline

# Mark dendritic spines as label 1 (from your classification)
spine_label = 1

# Get synapse mappings
synapse_mappings = get_synapse_mapping(
    root_id, mesh, client, side="post"
)[:, 1]  # Extract vertex indices

# Run morphometry analysis
results, corrected_labels = component_morphometry_pipeline(
    mesh,
    labels,
    select_label=spine_label,
    post_synapse_mappings=synapse_mappings,
    split_threshold=0.1,  # Split elongated components
    bound_volume_threshold=50,  # Skip components > 50 μm³ bounding box
    verbose=True
)

print(f"Analyzed {len(results)} spine components")
print(results.head())

# Examine measurements
print(f"\nVolume range: {results['size_nm3'].min():.0f} - {results['size_nm3'].max():.0f} nm³")
print(f"Mean sphericity: {results['sphericity'].mean():.3f}")
print(f"Spines with synapses: {(results['n_post_synapses'] > 0).sum()}")

​

Output Columns

​

Geometric Measurements

​

Spatial Features

​

Principal Component Analysis

​

Distance Transform

​

Topology

​

Connectivity

​

Diagnostics

​

Component Splitting

1. Computes the graph Laplacian for each component
2. Calculates the Fiedler vector (second eigenvector)
3. If the Fiedler eigenvalue is below split_threshold, bisects the component
4. Each half becomes a new component if it meets split_min_size
5. Process repeats recursively until no more splits occur

# Split components aggressively
results, corrected = component_morphometry_pipeline(
    mesh,
    labels,
    select_label=1,
    split_threshold=0.15,  # Lower = more splitting
    split_min_size=15,     # Minimum 15 vertices per component
    verbose=True
)

print(f"Split {len(np.unique(labels[labels==1]))} components into {len(results)}")

​

Volume Estimation

1. Compute axis-aligned bounding box around the component
2. Generate uniform random sample points within the box
3. Compute winding number for each point (inside if > 0.5)
4. Estimate volume as: (fraction inside) × (bounding box volume)

​

Best Practices

• Component size: Works best for structures between 100-10,000 vertices
• Bounding threshold: Set bound_volume_threshold to exclude large artifacts (typically 50-100 μm³)
• Splitting: Use split_threshold=0.1-0.2 for spines, higher values for more conservative splitting
• Synapse mapping: Provides n_post_synapses count when analyzing postsynaptic structures
• Missing values: Components that fail measurements are excluded from results

​

Performance

• Bounding box volume (determines sample points)
• Number of vertices in component
• Whether splitting occurs

• Small spines (< 1 μm³): 0.01-0.1 seconds
• Medium boutons (1-10 μm³): 0.1-0.5 seconds
• Large components (> 10 μm³): 0.5-2 seconds

​

Notes

• All spatial measurements are in nanometers (nm)
• Volume/area measurements assume vertex coordinates are in nm
• Components with fewer than 10 vertices are automatically skipped
• The algorithm filters autapses when counting synapses
• Rows with any NaN values are dropped from results
• Component IDs are preserved in the DataFrame index as component_id