Filtrations

Introduction

A filtration is a nested sequence of topological spaces built from data. For images, filtrations are constructed by thresholding intensity values, creating a series of binary images at different levels.

Filtrations are the foundation of persistent homology, allowing us to track how topological features evolve across scales.

What is a Filtration?

Formal Definition (Simplified)

A filtration is a sequence of topological spaces:

\[\emptyset = X_0 \subseteq X_1 \subseteq X_2 \subseteq \cdots \subseteq X_n = X\]

where each \(X_i\) is contained in the next.

For Images:

Each \(X_i\) corresponds to regions where intensity is below (or above) a threshold \(\tau_i\).

Intuitive Understanding

Think of a filtration as a movie showing structure appearing/disappearing:

  1. Start with nothing (empty set)

  2. Gradually add pixels/voxels based on intensity

  3. End with the complete image

Analogy: Watching a photograph develop in a darkroom - structure emerges gradually.

Building Filtrations from Images

Sublevel Filtration

Definition: Include pixels with intensity \(\leq\) threshold.

Process:

  1. Start with lowest intensity value

  2. Gradually increase threshold

  3. Add pixels as threshold increases

  4. End with all pixels

Use Cases:

  • Dark features on bright background

  • Vessels in angiography (low intensity)

  • Air in CT scans (dark regions)

In MedTDA:

from medtda import BarcodeExtractor

extractor = BarcodeExtractor(filtration_type='sublevel')
barcodes = extractor.execute(image)

Superlevel Filtration

Definition: Include pixels with intensity \(\geq\) threshold.

Process:

  1. Start with highest intensity value

  2. Gradually decrease threshold

  3. Add pixels as threshold decreases

  4. End with all pixels

Use Cases:

  • Bright features on dark background

  • Contrast-enhanced regions (high intensity)

  • Bright lesions in MRI

  • Enhanced tumors in CT

In MedTDA:

extractor = BarcodeExtractor(filtration_type='superlevel')
barcodes = extractor.execute(image)

Choosing Filtration Direction

Feature Type

Filtration

Example

Dark on bright background

Sublevel

Blood vessels

Bright on dark background

Superlevel

Contrast-enhanced lesions

Both types present

Try both

Mixed tissue types

Tip: Visualize your image - which features stand out?

Cubical Complexes

What are Cubical Complexes?

Cubical complexes are topological structures built from cubes (pixels in 2D, voxels in 3D).

Why Cubical?

  • Images are naturally defined on grids

  • Direct construction from pixels/voxels

  • No need for triangulation

  • Computationally efficient

Alternative: Simplicial complexes (triangulations) - more common in point cloud TDA but less natural for images.

Construction Methods

MedTDA supports two cubical complex construction methods:

T-Construction (Top-Cell)

Each pixel/voxel becomes a top-dimensional cell - a 2-cube (square) in 2D images or a 3-cube (cube) in 3D images. The filtration value equals the pixel/voxel intensity.

Properties:

  • More features detected

  • Finer topological detail

  • Default and recommended method

In MedTDA:

extractor = BarcodeExtractor(construction='T')

V-Construction (Vertex-Based)

Each pixel/voxel becomes a vertex (0-cube). Higher-dimensional cells are built from the vertex grid, inheriting intensities from the vertices.

Properties:

  • Fewer features detected

  • Coarser topological detail

  • Slightly faster - use when speed is critical

In MedTDA:

extractor = BarcodeExtractor(construction='V')

Comparison

Property

T-Construction

V-Construction

Cells per pixel/voxel

One top-cell

One vertex

Feature count

More features

Fewer features

Detail level

Finer

Coarser

Computation speed

Slightly slower

Slightly faster

Recommendation

Default choice

Speed-critical applications

Practical Difference:

# T-construction: More detailed
extractor_t = BarcodeExtractor(construction='T')
barcodes_t = extractor_t.execute(image)
print(f"H1 features: {len(barcodes_t['H1'])}")  # e.g., 47 features

# V-construction: Faster, less detail
extractor_v = BarcodeExtractor(construction='V')
barcodes_v = extractor_v.execute(image)
print(f"H1 features: {len(barcodes_v['H1'])}")  # e.g., 32 features

Filtration Process

As we sweep through thresholds during a sublevel filtration:

H0 (Connected Components):

  • Birth: New component appears

  • Death: Two components merge

  • One component persists forever (infinite persistence)

H1 (Loops):

  • Birth: Cycle forms

  • Death: Cycle gets filled

  • All H1 features eventually die

H2 (Voids - 3D only):

  • Birth: Enclosed void forms

  • Death: Void gets filled

  • All H2 features eventually die

Multi-Dimensional Images

2D Images (H0, H1): H0 captures separate regions/blobs; H1 captures loops and enclosed areas.

extractor = BarcodeExtractor(max_dimension=1)
barcodes = extractor.execute(image_2d)
h0 = barcodes['H0']  # Connected components
h1 = barcodes['H1']  # Loops

3D Images (H0, H1, H2): H0 captures separate regions; H1 captures tunnels/loops; H2 captures enclosed voids and cavities.

extractor = BarcodeExtractor(max_dimension=2)
barcodes = extractor.execute(image_3d)
h0 = barcodes['H0']  # Connected components
h1 = barcodes['H1']  # Loops/tunnels
h2 = barcodes['H2']  # Voids/cavities

Practical Considerations

Image Preprocessing

Filtrations work best on preprocessed images:

  1. Normalization: Consistent intensity ranges

from medtda import Preprocessor

preprocessor = Preprocessor(
    normalize=True,
    normalize_method='robust'
)
processed = preprocessor.preprocess(image)

  1. Resampling: Isotropic voxels

preprocessor = Preprocessor(spacing=(1.0, 1.0, 1.0))

  1. Masking: Focus on region of interest

preprocessor = Preprocessor(crop_to_roi=True)

Computational Complexity

Time complexity: Roughly O(n³) for n³ voxels

Memory: Stores cubical complex

Tips for large images:

  • Downsample if possible

  • Crop to ROI first

  • Use V-construction for speed

  • Process 2D slices instead of full 3D

Mathematical Foundations

Cubical Complexes

A cubical complex is a collection of elementary cubes (vertices, edges, squares, cubes, etc.) in Euclidean space. Medical images naturally form cubical complexes where pixels/voxels are the top-dimensional cubes.

Elementary Cubes: In dimension n, an elementary k-cube is a product:

\[[a_1, b_1] \times [a_2, b_2] \times \cdots \times [a_n, b_n]\]

where exactly k intervals have \(a_i < b_i\) (non-degenerate) and the rest have \(a_i = b_i\) (degenerate).

Examples:

  • 0-cube (vertex): all intervals degenerate, a single point

  • 1-cube (edge): one non-degenerate interval

  • 2-cube (square): two non-degenerate intervals

  • 3-cube (cube): three non-degenerate intervals

Boundary Operator: The boundary of a cube consists of all faces (cubes of one dimension lower) on its surface. For a 2D pixel at (i, j):

\[\partial([i,i+1] \times [j,j+1]) = \text{(four edges around the pixel)}\]

Filtered Complexes

A filtration is an increasing sequence of complexes where each complex is a subcomplex of the next:

\[X_0 \subseteq X_1 \subseteq X_2 \subseteq \cdots \subseteq X_n\]

From Images: Given a grayscale image I with intensity function f(x), filtrations are built by thresholding:

Sublevel sets: \(X_t = \{x : f(x) \leq t\}\)

Superlevel sets: \(X_t = \{x : f(x) \geq t\}\)

Lower-star filtration: For cubical complexes, each lower-dimensional cube (edge, vertex) enters the filtration at the minimum intensity of its cofaces (adjacent pixels):

\[t(\text{cube}) = \min\{f(\text{pixel}) : \text{cube} \in \partial(\text{pixel})\}\]

This ensures the filtration respects the cubical complex structure.

Cubical Homology

Chain Groups: The k-th chain group \(C_k\) consists of formal sums of k-dimensional cubes with coefficients in a field (typically \(\mathbb{Z}/2\mathbb{Z}\)).

Boundary Map: The boundary operator \(\partial_k : C_k \to C_{k-1}\) maps each k-cube to its boundary (collection of (k-1)-faces).

Homology Groups: The k-th homology group is:

\[H_k = \ker(\partial_k) / \text{im}(\partial_{k+1})\]

This quotient captures:

  • Cycles: k-dimensional structures with no boundary (\(\partial_k = 0\))

  • Boundaries: k-dimensional structures that are boundaries of (k+1)-structures

  • Homology classes: Equivalence classes of cycles that differ by boundaries

Interpretation:

  • H_0: Connected components (0-cycles not bounding anything)

  • H_1: Loops (1-cycles that don’t bound 2-dimensional regions)

  • H_2: Voids (2-cycles that don’t bound 3-dimensional volumes)

See Also