spectralbrain.core.pointclouds#

Point-cloud representation and point-cloud-specific analysis.

Provides BrainPointCloud, the mesh-free geometric container for brain structures extracted from volumetric segmentations, tractography, or surface vertex sets.

Design philosophy#

Point clouds are the atlas-free, mesh-free pathway. No marching cubes, no surface reconstruction. A volumetric segmentation label becomes a set of (x, y, z) coordinates in world space — full stop. The Laplace–Beltrami operator is approximated via graph Laplacians (kNN, Belkin–Niyogi, or Sharp–Crane robust), and spectral descriptors are computed identically to the mesh pathway.

Acquisition pathways#

┌─ FreeSurfer segmentation (aseg.mgz, thalamic_nuclei) ──► voxel → RAS
├─ FreeSurfer surface (.white, .pial) ──────────────────► vertices as points
├─ HippUnfold surface (.surf.gii, 8k/18k) ─────────────► vertices as points
├─ Raw T1w/T2w/FLAIR (via containers) ──────────────────► segment → voxel → RAS
├─ Atlas label volume (Schaefer, Julich, any .nii.gz) ──► voxel → RAS
├─ TractSeg tract masks (.nii.gz per tract) ────────────► voxel → RAS
└─ Tractography streamlines (.trk, .tck) ───────────────► streamline points

All converge to a single BrainPointCloud(points) that enters the spectral pipeline.

Classes

BrainPointCloud(points[, metadata])

Unstructured 3D point cloud for brain structures.

class spectralbrain.core.pointclouds.BrainPointCloud(points, metadata=None)[source]#

Bases: object

Unstructured 3D point cloud for brain structures.

Implements GeometricObject. No face connectivity — the Laplacian is built from a neighbourhood graph (kNN, ε-ball, or Belkin–Niyogi).

Parameters:
  • points (ndarray, shape (N, 3)) – Coordinates in world space (mm, RAS).

  • metadata (dict, optional) – Provenance (subject, hemisphere, structure, source, …).

Examples

>>> pc = BrainPointCloud.from_freesurfer_seg(
...     "aseg.mgz", label_id=17, jitter=True)
>>> decomp = pc.decompose(k=50, method="robust")
classmethod from_atlas_volume(atlas_path, label_id, *, jitter=True, seed=None, subsample=None, metadata=None)[source]#

Extract a point cloud from an arbitrary atlas NIfTI.

The atlas volume must be in the same space as the subject (template space or registered native).

Parameters:
Returns:

BrainPointCloud

Return type:

BrainPointCloud

classmethod from_freesurfer_seg(seg_path, label_id, *, jitter=True, jitter_scale=0.25, seed=None, subsample=None, metadata=None)[source]#

Load a FreeSurfer segmentation and extract a structure.

Handles aseg.mgz, aparc+aseg.mgz, ThalamicNuclei.v13.T1.mgz, amygNucVolumes.v22.mgz, hippoSfVolumes-T1.v22.mgz, or any integer label volume.

Parameters:
Returns:

BrainPointCloud

Return type:

BrainPointCloud

classmethod from_raw(raw_path, label_id, *, output_dir=None, gpu=None, jitter=True, seed=None, subsample=None, metadata=None)[source]#

End-to-end: raw T1w/T2w/FLAIR → segment → point cloud.

Chains skull-stripping and segmentation via Singularity containers, then extracts the target structure as a point cloud. No DL dependencies on host.

Parameters:
Returns:

BrainPointCloud

Return type:

BrainPointCloud

classmethod from_surface(surface_path, *, subsample=None, seed=None, metadata=None)[source]#

Load a surface file and use its vertices as a point cloud.

Works with FreeSurfer surfaces (.white, .pial, .inflated), GIfTI (.surf.gii), HippUnfold outputs, or any mesh format.

Parameters:
  • surface_path (PathLike)

  • subsample (int, optional) – FPS subsampling target.

  • seed (int, optional)

  • metadata (dict, optional)

Returns:

BrainPointCloud

Return type:

BrainPointCloud

classmethod from_tract_mask(mask_path, *, threshold=0.5, jitter=True, seed=None, subsample=None, metadata=None)[source]#

Extract a point cloud from a TractSeg binary tract mask.

TractSeg outputs one NIfTI per tract (72 tracts). Each non-zero voxel becomes a point.

Parameters:
Returns:

BrainPointCloud

Return type:

BrainPointCloud

classmethod from_tractography(trk_path, *, sampling='uniform', points_per_streamline=10, subsample=None, seed=None, metadata=None)[source]#

Extract a point cloud from tractography streamlines.

Parameters:
  • trk_path (PathLike) – .trk (TrackVis) or .tck (MRtrix) file.

  • sampling (str) – "all" — every point from every streamline. "endpoints" — first + last point per streamline. "midpoints" — middle point per streamline. "uniform" — uniformly resample each streamline to points_per_streamline points.

  • points_per_streamline (int) – For sampling="uniform": points per streamline.

  • subsample (int, optional) – FPS subsampling after extraction.

  • seed (int, optional)

  • metadata (dict, optional)

Returns:

BrainPointCloud

Return type:

BrainPointCloud

classmethod from_volume(label_volume, affine, label_id, *, jitter=True, jitter_scale=0.25, seed=None, subsample=None, metadata=None)[source]#

Extract a point cloud from a volumetric label map.

The most direct pathway: mgz np.array point cloud. No marching cubes, no mesh reconstruction. Voxel centroids in world space become points.

Parameters:
  • label_volume (ndarray, shape (X, Y, Z)) – Integer label volume (e.g. from aseg.mgz).

  • affine (ndarray, shape (4, 4)) – Voxel-to-world affine.

  • label_id (int) – Target label (e.g. 17 = left hippocampus in aseg).

  • jitter (bool) – Add sub-voxel Gaussian noise to break grid artefacts.

  • jitter_scale (float) – Jitter σ in voxel units.

  • seed (int, optional) – RNG seed for jitter.

  • subsample (int, optional) – If set, apply FPS to reduce to this many points.

  • metadata (dict, optional)

Returns:

BrainPointCloud

Return type:

BrainPointCloud

Examples

>>> data, affine = sb.io.load_nifti("aseg.mgz")
>>> pc = BrainPointCloud.from_volume(data, affine, label_id=17)
cluster(method='hdbscan', *, n_clusters=None, min_cluster_size=50, eps=None, **kwargs)[source]#

Atlas-free spatial clustering of the point cloud.

Assigns each point to a cluster without relying on any atlas or parcellation. Useful for data-driven subregion discovery.

Parameters:
  • method (str) – "hdbscan" — density-based, handles irregular shapes, auto-detects number of clusters (recommended). "spectral" — spectral clustering on kNN graph. "kmeans" — simple, fast, spherical assumption. "dbscan" — density-based, fixed epsilon.

  • n_clusters (int, optional) – For methods that require it (spectral, kmeans).

  • min_cluster_size (int) – For HDBSCAN.

  • eps (float, optional) – For DBSCAN.

  • **kwargs – Passed to the underlying clustering algorithm.

Returns:

labels (ndarray, shape (N,), int) – Cluster assignment per point. -1 = noise (HDBSCAN / DBSCAN).

Return type:

ndarray

compute_laplacian(method='robust', *, k=30, epsilon=None, sigma=None, robust_mollify=1e-05)[source]#

Construct a graph Laplacian on the point cloud.

Parameters:
  • method (str) – "knn" — Gaussian-weighted kNN graph Laplacian. "belkin_niyogi" — heat-kernel weighted Laplacian with provable convergence to the continuous LBO (Belkin & Niyogi, JCSS 2008). "robust" — Sharp & Crane tufted Laplacian (robust to noise and non-uniform density).

  • k (int) – Number of neighbours for kNN and Belkin–Niyogi.

  • epsilon (float, optional) – Bandwidth for ε-ball graph. None = auto from mean kNN distance.

  • sigma (float, optional) – Gaussian kernel bandwidth. None = auto (median distance heuristic).

  • robust_mollify (float) – Mollification for the robust method.

Returns:

  • L (SparseMatrix, shape (N, N))

  • M (MassMatrix, shape (N, N))

Return type:

tuple[spmatrix, spmatrix]

compute_normals(k=15, *, orient_to_centroid=True)[source]#

Estimate per-point normals via local PCA.

The normal at each point is the eigenvector of the local covariance matrix corresponding to the smallest eigenvalue (the direction of least variance in the neighbourhood).

Parameters:
  • k (int) – Neighbours for local PCA.

  • orient_to_centroid (bool) – Flip normals to point away from the cloud centroid (heuristic for outward orientation).

Returns:

normals (ndarray, shape (N, 3))

Return type:

ndarray[tuple[Any, …], dtype[floating]]

decompose(k=50, *, laplacian_method='robust', backend=None, **kwargs)[source]#

Compute the spectral decomposition.

Parameters:
Returns:

SpectralDecomposition

Return type:

SpectralDecomposition

denoise(k=10, threshold_sigma=2.5)[source]#

Remove density outliers (statistical outlier removal).

Parameters:
  • k (int) – Neighbours for density estimation.

  • threshold_sigma (float) – Points beyond this many σ from mean log-density are removed.

Returns:

BrainPointCloud – Cleaned copy.

Return type:

BrainPointCloud

local_curvature(k=15)[source]#

Estimate per-point curvature from local PCA eigenvalues.

The ratio λ₀ / (λ₀ + λ₁ + λ₂) of the smallest eigenvalue to the sum is a measure of local planarity — smaller values mean flatter (lower curvature), larger values mean more curved.

Parameters:

k (int) – Neighbours for local PCA.

Returns:

  • curvature (ndarray, shape (N,)) – Normalised planarity measure (0 = flat, 1 = isotropic).

  • linearity (ndarray, shape (N,)) – (λ₀ - λ₁) / λ₀ — high for tube-like structures.

Return type:

tuple[ndarray[tuple[Any, …], dtype[floating]], ndarray[tuple[Any, …], dtype[floating]]]

subsample(n, *, seed=None)[source]#

Return a subsampled copy via farthest-point sampling.

Parameters:
  • n (int) – Target number of points.

  • seed (int, optional)

Returns:

BrainPointCloud

Return type:

BrainPointCloud

surface_area()[source]#

Approximate surface area via convex hull.

For point clouds without face topology, the convex hull area is a rough proxy. Used for eigenvalue normalisation.

Return type:

float

property coordinates: ndarray[tuple[Any, ...], dtype[floating]]#

Return the point coordinates array.

property kd_tree: cKDTree#

Cached spatial index.

property n_points: int#

Return the number of points in the cloud.