Simplicial Complexes

Overview & Motivation

Take 200 points sampled from a ring (an annulus) and 200 points sampled from a Gaussian blob. Match their means and variances exactly. Now try to tell them apart: a $t$-test sees nothing, PCA sees nothing, and the covariance matrices are indistinguishable. Yet the ring has a hole through its center and the blob doesn’t — and that topological difference can be the signal that matters, whether you’re classifying dynamical regimes, detecting sensor network anomalies, or finding structure in high-dimensional embeddings.

The problem is that point clouds are discrete — they don’t have “holes” in any obvious sense. To detect shape, we need to build a continuous structure on top of the raw data. That structure is a simplicial complex: a systematic way of connecting nearby points with edges, triangles, and their higher-dimensional analogues to create a combinatorial scaffolding that does have well-defined topology.

Once you have a simplicial complex, the full machinery of algebraic topology becomes available — you can count connected components ($\beta_0$), detect loops ($\beta_1$), find enclosed voids ($\beta_2$), and quantify how these features persist across scales. If you’ve used ripser, gudhi, or scikit-tda, you’ve been building simplicial complexes without thinking about it. This topic makes the implicit explicit.

---

Formal Framework

The Building Blocks: Simplices

In concrete terms:

• A 0-simplex is a point
• A 1-simplex is an edge (line segment between two points)
• A 2-simplex is a filled triangle
• A 3-simplex is a filled tetrahedron
• A $k$-simplex for $k \geq 4$ follows the same pattern — it’s the “solid” version of the shape, not just the boundary

From Simplices to Complexes

The first condition says “if you include a triangle, you must include all its edges and vertices.” The second says “two simplices can only intersect along a shared face — no partial overlaps.” Together, these ensure the complex is combinatorially well-behaved.

Abstract vs. Geometric

The distinction matters: a geometric simplicial complex lives in $\mathbb{R}^n$ with coordinates and convex hulls. An abstract simplicial complex is purely combinatorial — it’s a list of which vertex subsets are “connected.” For TDA, we work with abstract complexes (defined by distance thresholds) and use geometric realizations only for visualization.

The Vietoris-Rips Construction

This is the construction that matters for data analysis. Given a point cloud and a scale parameter, it builds a simplicial complex by connecting points that are “close enough.”

---

Homology Preview

A simplicial complex has well-defined shape — but how do we quantify that shape? The answer is homology: an algebraic machine that assigns integer invariants to a complex, counting its connected components, loops, voids, and higher-dimensional “holes.”

We give the key definitions here at an intuitive level. The full algebraic development — including the matrix computations that make homology practical — is developed on the Persistent Homology page.

Working over $\mathbb{Z}_2$ means we don’t need to worry about orientations: coefficients are 0 or 1, and $1 + 1 = 0$. This is the standard choice in TDA software.

The boundary of an edge $[v_0, v_1]$ is $[v_0] + [v_1]$ — its two endpoints. The boundary of a triangle $[v_0, v_1, v_2]$ is $[v_0, v_1] + [v_0, v_2] + [v_1, v_2]$ — its three edges. The fundamental property is that the boundary of a boundary is always zero: $\partial_{k-1} \circ \partial_k = 0$ (proved on the Persistent Homology page, Theorem 1).

In geometric terms:

• $\beta_0$ = number of connected components
• $\beta_1$ = number of independent loops (cycles not bounding a filled region)
• $\beta_2$ = number of enclosed voids (cavities)

This is the core mechanism: the same cycle can be a hole or not, depending on whether something fills it in. Persistent homology tracks exactly when filling events happen as the scale parameter $\varepsilon$ increases.

---

The Euler Characteristic

The Euler characteristic is the oldest and simplest topological invariant — Euler discovered it for polyhedra in 1758, long before homology was invented. It turns out to be a shadow of the Betti numbers.

For a surface, this reduces to the familiar $V - E + F$ formula. For a convex polyhedron, Euler proved that $V - E + F = 2$ — regardless of whether the polyhedron is a cube, a dodecahedron, or any other convex shape. This is a topological invariant: it depends only on the shape’s topology (homeomorphic to a sphere), not on the specific triangulation.

The explorer below lets you build a simplicial complex by adding and removing simplices, watching the Euler characteristic update in real time:

Euler Characteristic Explorer

Interactive visualization coming soon — build a simplicial complex and watch χ = V − E + F update in real time.

---

Other Complex Constructions

The Vietoris-Rips complex is the workhorse of TDA, but it is not the only way to build topology from data. Two important alternatives trade off geometric fidelity against computational cost.

The alpha complex has a crucial advantage: its size is bounded by $O(n^{\lceil d/2 \rceil})$ regardless of $\varepsilon$, because the Delaunay triangulation has at most that many simplices. In low dimensions ($d = 2, 3$), this is nearly linear in $n$ — far smaller than the $\binom{n}{k+1}$ potential simplices of the VR complex. Moreover, the Nerve Theorem applies to the alpha complex (it is a subcomplex of the Čech complex), so its homology faithfully captures the topology of the union of balls.

---

Computational Notes

This combinatorial explosion is why no one actually constructs the full VR complex. The state-of-the-art approach, implemented in Ripser (Bauer, 2021), avoids explicit construction entirely:

1. Implicit representation: simplices are never stored — they are generated on-the-fly from the distance matrix during the matrix reduction algorithm.
2. Clearing optimization: when a simplex $\sigma$ creates a birth event, the simplex that kills it (its “partner”) can be identified without processing the remaining columns — roughly halving the work.
3. Cohomology: Ripser computes cohomology (the dual of homology) because the coboundary matrix is sparser than the boundary matrix for VR complexes, leading to faster column reduction.
4. Apparent pairs: many birth-death pairs can be identified from the combinatorial structure of the filtration without any matrix reduction at all.

These optimizations combine to make Ripser orders of magnitude faster than naive implementations. For 1000 points in $\mathbb{R}^{10}$, computing $H_0$ and $H_1$ persistence takes seconds, not hours.

---

Visual Intuition

The visualization below shows a Vietoris-Rips complex built from points sampled near a circle. Drag the $\varepsilon$ slider to watch the complex form:

• At small $\varepsilon$: only isolated vertices — the complex captures no structure
• At moderate $\varepsilon$: edges connect nearby points, and you can see the circular arrangement emerge
• At large $\varepsilon$: triangles fill in, eventually collapsing the loop into a solid disk

ε = 0.30

What to notice: There’s a “sweet spot” for $\varepsilon$ where the edges form a cycle around the circle but the interior hasn’t filled in yet. At that scale, the complex has a one-dimensional hole — this is the $H_1$ feature that persistent homology will detect. Finding that sweet spot automatically, across all possible scales simultaneously, is exactly what persistent homology does.

The Betti number tracker below shows $\beta_0$ and $\beta_1$ updating as the scale parameter changes — watch the loop appear and then get filled in:

Betti Number Tracker

Interactive visualization coming soon — drag ε and watch β₀, β₁ update as the complex evolves.

---

Working Code

Building the Vietoris-Rips Complex

Building the Vietoris-Rips complex from scratch in Python:

import numpy as np
from itertools import combinations
from scipy.spatial.distance import pdist, squareform

def vr_complex(points, epsilon):
    """Build Vietoris-Rips complex at scale epsilon.

    Parameters
    ----------
    points : ndarray of shape (n, d)
        Point cloud in R^d.
    epsilon : float
        Scale parameter. Points within this distance are connected.

    Returns
    -------
    vertices : list of int
    edges : list of (int, int)
    triangles : list of (int, int, int)
    """
    dist_matrix = squareform(pdist(points))
    n = len(points)
    vertices = list(range(n))
    edges = [(i, j) for i, j in combinations(range(n), 2)
             if dist_matrix[i, j] <= epsilon]
    triangles = [(i, j, k) for i, j, k in combinations(range(n), 3)
                 if all(dist_matrix[a, b] <= epsilon
                        for a, b in combinations([i, j, k], 2))]
    return vertices, edges, triangles

# Example: 12 points on a circle
theta = np.linspace(0, 2 * np.pi, 12, endpoint=False)
points = np.column_stack([np.cos(theta), np.sin(theta)])

for eps in [0.3, 0.8, 1.2, 1.8]:
    v, e, t = vr_complex(points, eps)
    print(f"ε={eps:.1f}: {len(v)} vertices, {len(e)} edges, {len(t)} triangles")

ε=0.3: 12 vertices, 0 edges, 0 triangles
ε=0.8: 12 vertices, 12 edges, 0 triangles
ε=1.2: 12 vertices, 24 edges, 12 triangles
ε=1.8: 12 vertices, 42 edges, 64 triangles

At $\varepsilon = 0.8$, we get exactly 12 edges forming a cycle — the complex “sees” the circle. By $\varepsilon = 1.8$, the complex is nearly complete (66 possible edges, 220 possible triangles) and the topological signal is lost in combinatorial noise. This is why we need persistent homology: it tracks which features are real (they persist across a range of scales) and which are artifacts.

For real-world analysis, use ripser instead of this brute-force construction:

from ripser import ripser

# ripser builds the VR filtration and computes persistence in one step
result = ripser(points, maxdim=1)
print(result['dgms'][0][:5])  # H₀ intervals (connected components)
print(result['dgms'][1][:5])  # H₁ intervals (loops)

Betti Numbers from Boundary Matrices

Computing Betti numbers directly from the boundary matrices, using the rank-nullity formula $\beta_k = \dim(\ker \partial_k) - \text{rank}(\partial_{k+1})$:

import numpy as np

def boundary_matrix_1(vertices, edges):
    """Build the boundary operator ∂₁ over Z₂.
    Rows = vertices, columns = edges. Entry (v, e) = 1 if v is an endpoint of e."""
    n_v, n_e = len(vertices), len(edges)
    D1 = np.zeros((n_v, n_e), dtype=int)
    v_idx = {v: i for i, v in enumerate(vertices)}
    for j, (a, b) in enumerate(edges):
        D1[v_idx[a], j] = 1
        D1[v_idx[b], j] = 1
    return D1 % 2  # Z₂ coefficients

def boundary_matrix_2(edges, triangles):
    """Build the boundary operator ∂₂ over Z₂.
    Rows = edges, columns = triangles. Entry (e, t) = 1 if e is a face of t."""
    e_set = {tuple(sorted(e)): i for i, e in enumerate(edges)}
    n_e, n_t = len(edges), len(triangles)
    D2 = np.zeros((n_e, n_t), dtype=int)
    for j, (a, b, c) in enumerate(triangles):
        for face in [(a,b), (a,c), (b,c)]:
            key = tuple(sorted(face))
            if key in e_set:
                D2[e_set[key], j] = 1
    return D2 % 2

def gf2_rank(matrix):
    """Compute matrix rank over Z₂ via Gaussian elimination (XOR)."""
    A = np.array(matrix, dtype=np.uint8, copy=True) % 2
    n_rows, n_cols = A.shape
    rank = 0
    pivot_row = 0
    for col in range(n_cols):
        found = None
        for row in range(pivot_row, n_rows):
            if A[row, col]:
                found = row
                break
        if found is None:
            continue
        if found != pivot_row:
            A[[pivot_row, found]] = A[[found, pivot_row]]
        for row in range(pivot_row + 1, n_rows):
            if A[row, col]:
                A[row] ^= A[pivot_row]
        rank += 1
        pivot_row += 1
        if pivot_row == n_rows:
            break
    return rank

def betti_numbers(vertices, edges, triangles):
    """Compute β₀, β₁ over Z₂ via rank-nullity."""
    D1 = boundary_matrix_1(vertices, edges)
    D2 = boundary_matrix_2(edges, triangles)
    rank_D1 = gf2_rank(D1)
    rank_D2 = gf2_rank(D2)
    beta_0 = len(vertices) - rank_D1
    beta_1 = (len(edges) - rank_D1) - rank_D2
    return beta_0, beta_1

# Wireframe triangle (no fill): should give β₀=1, β₁=1
verts = [0, 1, 2]
edgs = [(0,1), (0,2), (1,2)]
print(f"Wireframe triangle: β₀={betti_numbers(verts, edgs, [])[0]}, "
      f"β₁={betti_numbers(verts, edgs, [])[1]}")

# Filled triangle: should give β₀=1, β₁=0
print(f"Filled triangle:    β₀={betti_numbers(verts, edgs, [(0,1,2)])[0]}, "
      f"β₁={betti_numbers(verts, edgs, [(0,1,2)])[1]}")

Wireframe triangle: β₀=1, β₁=1
Filled triangle:    β₀=1, β₁=0

Alpha Complexes with GUDHI

For low-dimensional data, the alpha complex is dramatically more efficient than VR:

import gudhi
import numpy as np

# 100 points on a noisy circle in R²
theta = np.random.uniform(0, 2 * np.pi, 100)
r = 1.0 + np.random.normal(0, 0.1, 100)
circle_pts = np.column_stack([r * np.cos(theta), r * np.sin(theta)])

# Alpha complex: size bounded by O(n) in R²
alpha = gudhi.AlphaComplex(points=circle_pts)
stree = alpha.create_simplex_tree()
print(f"Alpha complex: {stree.num_simplices()} simplices")

# Compare with VR complex at a similar scale
vr = gudhi.RipsComplex(points=circle_pts, max_edge_length=0.8)
stree_vr = vr.create_simplex_tree(max_dimension=2)
print(f"VR complex:    {stree_vr.num_simplices()} simplices")

# Compute persistence — both should detect the circle's H₁ feature
stree.persistence()
dgm_alpha = stree.persistence_intervals_in_dimension(1)
print(f"Alpha H₁ features: {len(dgm_alpha)}, max persistence: "
      f"{max(d-b for b, d in dgm_alpha if d < float('inf')):.3f}")

Euler Characteristic Verification

Verify the Euler–Poincaré formula on the VR complex at various scales:

for eps in [0.3, 0.8, 1.2, 1.8]:
    v, e, t = vr_complex(points, eps)
    chi_simplex = len(v) - len(e) + len(t)
    b0, b1 = betti_numbers(list(range(len(v))), e, t)
    chi_betti = b0 - b1  # β₂ = 0 for dim ≤ 2
    print(f"ε={eps:.1f}: χ(simplex)={chi_simplex}, χ(Betti)={chi_betti}, "
          f"match={chi_simplex == chi_betti}")

---

Connections & Cross-Track Bridges

Simplicial complexes are the combinatorial substrate for nearly everything in topological data analysis and connect deeply to several other areas of mathematics on this site.

Within the Topology & TDA track:

• Persistent Homology tracks how the homology of $\text{VR}(X, \varepsilon)$ changes as $\varepsilon$ increases. The filtration — the nested sequence of VR complexes at increasing scales — is the input to the persistence algorithm.
• Čech Complexes & Nerve Theorem develops the geometrically exact alternative to VR, where the Nerve Theorem guarantees that the complex faithfully captures the topology of the union of balls.
• The Mapper Algorithm constructs a 1-dimensional simplicial complex (a graph) that captures the data’s topological skeleton via the nerve of a cover. The Mapper graph is the 1-skeleton of the nerve of a refined open cover.
• Sheaf Theory defines cellular sheaves on the face poset of a simplicial complex — assigning vector spaces to simplices and linear maps to face inclusions. The coboundary map $\delta_0$ in sheaf cohomology is a generalization of the boundary operator $\partial_k$ developed here.

Across tracks:

• Graph Laplacians & Spectrum: The 1-skeleton of a simplicial complex is a graph, and its Laplacian spectrum governs spectral clustering. The boundary operator $\partial_1$ and the graph Laplacian $L = \partial_1^T \partial_1$ (over $\mathbb{R}$) are directly related — the Laplacian measures how “boundary-like” a 0-chain is, and $\ker(L) = \ker(\partial_1^T)$ counts connected components.
• Categories & Functors: Simplicial complexes form a category $\textbf{SComp}$ with simplicial maps as morphisms. Homology $H_k(-; \mathbb{F})$ is a functor from $\textbf{SComp}$ to the category of vector spaces — a simplicial map $f: K \to L$ induces a linear map $f_*: H_k(K) \to H_k(L)$ that respects composition. This functoriality is what makes persistent homology well-defined: the inclusion maps $K_i \hookrightarrow K_j$ in a filtration induce the maps between homology groups that the persistence algorithm tracks.
• Message Passing & GNNs: The message-passing framework on graphs generalizes naturally to simplicial complexes. Simplicial neural networks pass messages along triangles and higher simplices — not just edges — enabling them to capture higher-order interactions that graph neural networks miss.

---

Exercises

Foundational

Exercise 1. A 3-simplex (tetrahedron) $[v_0, v_1, v_2, v_3]$ has faces of dimensions 0, 1, 2, and 3. Count the number of faces of each dimension and verify that the total number of faces is $2^4 - 1 = 15$. Hint: A $j$-face corresponds to a choice of $j+1$ vertices from 4.

Exercise 2. The octahedron can be triangulated as a simplicial complex with 6 vertices, 12 edges, and 8 triangular faces (it is homeomorphic to $S^2$). Verify the Euler–Poincaré formula: compute $\chi = V - E + F$ directly, then confirm that $\beta_0 = 1$, $\beta_1 = 0$, $\beta_2 = 1$ gives $\chi = 1 - 0 + 1 = 2$.

Intermediate

Exercise 3. Prove that the Vietoris-Rips complex is a flag complex — that is, it is completely determined by its 1-skeleton (the graph of edges). Specifically, show that a subset $\sigma \subseteq X$ is a simplex of $\text{VR}(X, \varepsilon)$ if and only if every pair of vertices in $\sigma$ is connected by an edge. Hint: This follows immediately from Definition 6, but spell out why the face closure property is also satisfied.

Exercise 4. Let $P$ be 6 points equally spaced on the unit circle in $\mathbb{R}^2$. For the Vietoris-Rips complex $\text{VR}(P, \varepsilon)$:

• At what value of $\varepsilon$ do the first edges appear?
• At what value of $\varepsilon$ does $\beta_1$ first become 1 (a loop appears)?
• At what value of $\varepsilon$ does $\beta_1$ return to 0 (the loop gets filled)?

Compute these thresholds analytically from the geometry, then verify with ripser.

Advanced

Exercise 5. Implement the witness complex from scratch for a set of $n = 500$ points on a noisy torus in $\mathbb{R}^3$ with $|L| = 50$ landmarks chosen via farthest-point sampling. Compare the persistence diagram (up to $H_1$) with the VR diagram on the full dataset. How many landmarks are needed for the witness complex to detect both independent loops ($\beta_1 = 2$) of the torus?

Exercise 6. The Čech-to-VR interleaving inequality states $\check{C}_\varepsilon(P) \subseteq \text{VR}_\varepsilon(P) \subseteq \check{C}_{2\varepsilon}(P)$. Construct a specific point configuration in $\mathbb{R}^2$ where the VR complex contains a 2-simplex (filled triangle) that the Čech complex does not, at the same scale parameter. Compute the Betti numbers of both complexes and show that they disagree on $\beta_1$. This is a concrete instance of the phantom cycle problem.