topology-data-analysis

Topological Data Analysis

A skill for applying topological data analysis (TDA) methods to research data. Covers persistent homology, Vietoris-Rips complexes, persistence diagrams, the Mapper algorithm, and vectorization methods for integrating topological features into machine learning pipelines.

Core Concepts

Simplicial Complexes from Data

TDA extracts topological features (connected components, loops, voids) from data by building simplicial complexes at multiple scales:

Complex	Construction	Computational Cost
Vietoris-Rips	Edge if distance < epsilon	O(n^d) for d-simplices
Cech	Ball intersection (exact)	Computationally expensive
Alpha	Delaunay-based (exact in low dim)	Efficient in R^2, R^3
Cubical	Grid-based (for images)	Linear in pixels

Filtration and Persistence

Scale epsilon:  0.1    0.3    0.5    0.7    1.0
                |------|------|------|------|------|
Components:      10      6      3      2      1
  (H0 features born at 0, die at merging scale)

Loops:           0      0      1      2      0
  (H1 features born when loop forms, die when filled)

A feature that persists across many scales is a genuine topological signal; short-lived features are noise.

Persistent Homology with Ripser

Computing Persistence Diagrams

import numpy as np
from ripser import ripser
from persim import plot_diagrams

def compute_persistence(point_cloud: np.ndarray,
                         max_dim: int = 2,
                         max_edge: float = 2.0) -> dict:
    """
    Compute persistent homology of a point cloud.
    point_cloud: (n_points, n_dimensions) array
    max_dim: maximum homology dimension to compute
    max_edge: maximum edge length in Rips complex
    Returns persistence diagrams for each dimension.
    """
    result = ripser(
        point_cloud,
        maxdim=max_dim,
        thresh=max_edge,
    )

    diagrams = result["dgms"]
    summary = {}

    for dim, dgm in enumerate(diagrams):
        # Filter out infinite death times for H0
        finite = dgm[dgm[:, 1] < np.inf] if len(dgm) > 0 else dgm
        lifetimes = finite[:, 1] - finite[:, 0] if len(finite) > 0 else np.array([])

        summary[f"H{dim}"] = {
            "n_features": len(finite),
            "max_persistence": float(lifetimes.max()) if len(lifetimes) > 0 else 0,
            "mean_persistence": float(lifetimes.mean()) if len(lifetimes) > 0 else 0,
            "birth_death_pairs": finite.tolist(),
        }

    return summary

# Example: torus point cloud
def sample_torus(n=1000, R=3.0, r=1.0, noise=0.1):
    """Sample points from a torus in R^3."""
    theta = np.random.uniform(0, 2 * np.pi, n)
    phi = np.random.uniform(0, 2 * np.pi, n)
    x = (R + r * np.cos(phi)) * np.cos(theta) + np.random.normal(0, noise, n)
    y = (R + r * np.cos(phi)) * np.sin(theta) + np.random.normal(0, noise, n)
    z = r * np.sin(phi) + np.random.normal(0, noise, n)
    return np.column_stack([x, y, z])

torus = sample_torus(500)
persistence = compute_persistence(torus, max_dim=2)
# Expected: H0 has 1 long-lived component,
#           H1 has 2 prominent loops (the two fundamental cycles),
#           H2 has 1 prominent void (the cavity)