Improved approximate rips filtrations with shifted integer lattices and cubical complexes

Rips complexes are important structures for analyzing topological features of metric spaces. Unfortunately, generating these complexes is expensive because of a combinatorial explosion in the complex size. For n points in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\mathbb {R}^d$$\end{document}Rd, we present a scheme to construct a 2-approximation of the filtration of the Rips complex in the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$L_\infty $$\end{document}L∞-norm, which extends to a \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$2d^{0.25}$$\end{document}2d0.25-approximation in the Euclidean case. The k-skeleton of the resulting approximation has a total size of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$n2^{O(d\log k +d)}$$\end{document}n2O(dlogk+d). The scheme is based on the integer lattice and simplicial complexes based on the barycentric subdivision of the d-cube. We extend our result to use cubical complexes in place of simplicial complexes by introducing cubical maps between complexes. We get the same approximation guarantee as the simplicial case, while reducing the total size of the approximation to only \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$n2^{O(d)}$$\end{document}n2O(d) (cubical) cells. There are two novel techniques that we use in this paper. The first is the use of acyclic carriers for proving our approximation result. In our application, these are maps which relate the Rips complex and the approximation in a relatively simple manner and greatly reduce the complexity of showing the approximation guarantee. The second technique is what we refer to as scale balancing, which is a simple trick to improve the approximation ratio under certain conditions.


Introduction
Context. Persistent homology (Carlsson 2009;Edelsbrunner and Harer 2010;Edelsbrunner et al. 2002) is a technique to analyze data sets using topological invariants. The idea is to build a multi-scale representation of data sets and to track its homological changes across the scales.
A standard construction for the important case of point clouds in Euclidean space is the Vietoris-Rips complex (usually abbreviated as simply the Rips complex): for a scale parameter α ≥ 0, it is the collection of all subsets of points with diameter at most α. When α increases from 0 to ∞, the Rips complexes form a filtration, an increasing sequence of nested simplicial complexes whose homological changes can be computed and represented in terms of a barcode.
The computational drawback of Rips complexes is their sheer size: the k-skeleton of a Rips complex (that is, where only subsets of size at most k + 1 are considered) for n points consists of (n k+1 ) simplices because every (k +1)-subset joins the complex for a sufficiently large scale parameter. This size bound makes barcode computations for large point clouds infeasible even for low-dimensional homological features 1 . This difficulty motivates the question of what we can say about the barcode of the Rips filtration without explicitly constructing all of its simplices.
We address this question using approximation techniques. The space of barcodes forms a metric space: two barcodes are close if similiar homological features occur on roughly the same range of scales. More precisely, the bottleneck distance is used as a distance metric between barcodes. The first approximation scheme by Sheehy (2013) constructs a (1 + ε)-approximation of the k-skeleton of the Rips filtration using only n( 1 ε ) O(λk) simplices for arbitrary finite metric spaces, where λ is the doubling dimension of the metric. Further approximation techniques for Rips complexes (Dey et al. 2014) and the closely relatedČech complexes (Botnan and Spreemann 2015;Cavanna et al. 2015;Kerber and Sharathkumar 2013) have been derived subsequently, all with comparable size bounds. More recently, we constructed an approximation scheme (Choudhary et al. 2019) for theČech filtrations of n points in R d that had size n 1 ε O(d) 2 O(d log d+dk) for the k-skeleton, improving the size bound from previous work.
In Choudhary et al. (2017b), we constructed an approximation scheme for Rips filtration in Euclidean space that yields a worse approximation factor of only O(d), but uses only n2 O(d log k+d) simplices. In Choudhary et al. (2017b), we also show a lower bound result on the size of approximations: for any ε < 1/ log 1+c n with some constant c ∈ (0, 1), any ε-approximate filtration has size n (log log n) .
There has also been work on using cubical complexes to compute persistent homology, such as in Wagner et al. (2012). Cubical complexes are typically smaller than their simplicial counterparts, simply because they avoid triangulations. However, to our knowledge, there has been no attempt to utilize them in computing approximations of filtrations. Also, while there are efficient methods to compute persistence for simplicial complexes connected with simplicial maps (Dey et al. 2014; Schreiber 2017), we are not aware of such counterparts for cubical complexes. Our contributions. For the Rips filtration of n points in R d with distances taken in the L ∞ -norm, we present a 2-approximation whose k-skeleton has size at most n6 d−1 (2k+ denotes Stirling numbers of the second kind. This translates to a 2d 0.25 -approximation of the Rips filtration in the Euclidean metric and hence improves the asymptotic approximation quality of our previous approach (Choudhary et al. 2017b) with the same size bound. Our scheme gives the best size guarantee over all previous approaches. On a high level, our approach follows a straightforward approximation scheme: given a scaled and appropriately shifted integer grid on R d , we identify those grid points that are close to the input points and build an approximation complex using these grid points. The challenge lies in how to connect these grid points to a simplicial complex such that close-by grid points are connected, while avoiding too many connections to keep the size small. Our approach first selects a set of active faces in the cubical complex defined over the grid, and defines the approximation complex using the barycentric subdivision of this cubical complex.
We also describe an output-sensitive algorithm to compute our approximation. By randomizing the aforementioned shifts of the grids, we obtain a worst-case running time of n2 O(d) log + 2 O(d) M in expectation, where is the spread of the point set (that is, the ratio of the diameter to the closest distance of two points) and M is the size of the approximation.
Additionally, this paper makes the following technical contributions: -We follow the standard approach of defining a sequence of approximation complexes and establishing an interleaving between the Rips filtration and the approximation. We realize our interleaving using chain maps connecting a Rips complex at scale α to an approximation complex at scale cα, and vice versa, with c ≥ 1 being the approximation factor. Previous approaches (Choudhary et al. 2017b;Dey et al. 2014;Sheehy 2013) used simplicial maps for the interleaving, which induce an elementary form of chain maps and are therefore more restrictive.
The explicit construction of such maps can be a non-trivial task. The novelty of our approach is that we avoid this construction by the usage of acyclic carriers (Munkres 1984). In short, carriers are maps that assign subcomplexes to subcomplexes under some mild extra conditions. While they are more flexible, they still certify the existence of suitable chain maps, as we exemplify in Sect. 2. We believe that this technique is of general interest for the construction of approximations of cell complexes.
-We exploit a simple trick that we call scale balancing to improve the quality of approximation schemes. In short, if the aforementioned interleaving maps from and to the Rips filtration do not increase the scale parameter by the same amount, one can simply multiply the scale parameter of the approximation by a constant. Concretely, given maps interleaving the Rips complex R α and the approximation complex X α , we can define X α := X α/ √ c and obtain maps which improves the interleaving from c to √ c. While it has been observed that the same trick can be used for improving the worst-case distance between Rips anď Cech filtrations, 2 our work seems to be the first to make use of it in the context of approximations.
-We extend our approximation scheme to use cubical complexes instead of simplicial complexes, thereby achieving a marked reduction in size complexity. To connect the cubical complexes at different scales, we introduce the notion of cubical maps, which is a simple extension of simplicial maps to the cubical case. While we do not know of an algorithm that can compute persistence for the case of cubical complexes with cubical maps, we believe that this is a first step towards advocating the use of cubical complexes as approximating structures.
Our technique can be combined with dimension reduction techniques in the same way as in Choudhary et al. (2017b) (see Theorems 19, 21, and 22 therein), with improved logarithmic factors. We state the main results in the paper, while omitting the technical details. Updates from the conference version. An earlier version of this paper appeared at the 25th European Symposium on Algorithms (Choudhary et al. 2017a). In that version, we achieved a 3 √ 2-approximation of the L ∞ Rips filtration and correspondingly, a 3 √ 2d 0.25 -approximation of the L 2 case. In this version, we improve the weak interleaving of Choudhary et al. (2017a) to a strong interleaving to get improved approximation factors. We expand upon the details of scale balancing, among other proofs that were missing from the conference version. We add the case of cubical complexes in this version.
There is a subtle yet important distinction between the approximation complexes used in the conference version and the current result. In the conference version, our simplicial complex was built using only active faces, while the current version uses both active and secondary faces (please see Sect. 4 for definitions). This makes it easier to relate the simplicial and the cubical complexes in the current version. On the other hand the complexes are different, hence the associated proofs have been adapted accordingly.
Outline. We start by explaining the relevant topological concepts in Sect. 2. We give details of the integer grids that we use in Sect. 3. In Sect. 4 we present our approximation scheme that uses the barycentric subdivision, and present the computational aspects in Sect. 5. The extension to cubical complexes is presented in Sect. 6. We discuss practical aspects of our scheme and conclude in Sect. 7. Some details of the strong interleaving from Sect. 4 are deferred to Appendix A.

Preliminaries
We briefly review the essential topological concepts needed. More details are available in standard references (see Bubenik et al. 2015;Chazal et al. 2009;Edelsbrunner and Harer 2010;Hatcher 2002;Munkres 1984). Simplicial complexes. A simplicial complex K on a finite set of elements S is a collection of subsets {σ ⊆ S} called simplices such that each subset τ ⊂ σ is also in K . The dimension of a simplex σ ∈ K is k := |σ | − 1, in which case σ is called a k-simplex. A simplex τ is a sub-simplex of σ if τ ⊆ σ . We remark that, commonly a sub-simplex is called a "face" of a simplex, but we reserve the word "face" for a different structure. For the same reason, we do not introduce the common notation of of "vertices" and "edges" of simplicial complexes, but rather refer to 0-and 1-simplices throughout. The k-skeleton of K consists of all simplices of K whose dimension is at most k. For instance, the 1-skeleton of K is a graph defined by its 0-simplices and 1-simplices.
Given a point set P ⊂ R d and a real number α ≥ 0, the (Vietoris-)Rips complex on P at scale α consists of all simplices σ = ( p 0 , · · · , p k ) ⊆ P such that diam(σ ) ≤ α, where diam denotes the diameter. In this work, we write R α for the Rips complex at scale 2α with the Euclidean metric, and R ∞ α when using the metric of the L ∞ -norm. In either way, a Rips complex is an example of a flag complex, which means that whenever a set { p 0 , · · · , p k } ⊆ P has the property that every 1-simplex { p i , p j } is in the complex, then the k-simplex { p 0 , · · · , p k } is also in the complex.
A related complex is theČech complex of P at scale α, which consists of simplices of P for which the radius of the minimum enclosing ball is at most α. We do not studyČech complexes in this paper, but we mention them briefly while showing a connection with the Rips complex later in this section.
Let L be a simplicial complex. Letφ be a map which assigns a vertex of L to each vertex of K . A simplicial map is a map ϕ : K → L induced by a vertex mapφ, such that for every simplex { p 0 , · · · , p k } in K , the set {φ( p 0 ), · · · ,φ( p k )} is a simplex of L. For K a subcomplex of K , the inclusion map inc : K → K is an example of a simplicial map. A simplicial map is completely determined by its action on the 0-simplices of K . Chain complexes. A chain complex C * = (C p , ∂ p ) with p ∈ Z is a collection of abelian groups C p and homomorphisms ∂ p : C p → C p−1 such that ∂ p−1 • ∂ p = 0. A simplicial complex K gives rise to a chain complex C * (K ) for a fixed base field F: define C p for p ≥ 0 as the set of formal linear combinations of p-simplices in K over F, and C −1 := F. The boundary of a k-simplex with k ≥ 1 is the (signed) sum of its sub-simplices of co-dimension one 3 ; the boundary of a 0-simplex is simply set to 1. The homomorphisms ∂ p are then defined as the linear extensions of this boundary operator. Note that C * (K ) is sometimes called augmented chain complex of K , where the augmentation refers to the addition of the non-trivial group C −1 .
A chain map φ : C * → D * between chain complexes C * = (C p , ∂ p ) and D * = (D p , ∂ p ) is a collection of group homomorphisms φ p : For simplicial complexes K and L, we call a chain map φ : C * (K ) → C * (L) augmentation-preserving if φ −1 is the identity. A simplicial map ϕ : K → L between simplicial complexes induces an augmentation-preserving chain mapφ : C * (K ) → C * (L) between the corresponding chain complexes. This construction is functorial, meaning that for ϕ the identity function on a simplicial complex K ,φ is the identity function on C * (K ), and for composable simplicial maps ϕ, ϕ , we have that ϕ • ϕ =φ •φ . Homology. The p-th homology group H p (C * ) of a chain complex is defined as ker ∂ p /im ∂ p+1 . The p-th homology group of a simplicial complex K , H p (K ), is the p-th homology group of its induced chain complex C * (K ). Note that this definition is commonly referred to as reduced homology, but we ignore this distinction and consider reduced homology throughout. H p (C * ) is an F-vector space because we have chosen our base ring F as a field. Intuitively, when the chain complex is generated from a simplicial complex, the dimension of the p-th homology group counts the number of p-dimensional holes in the complex. We write H (C * ) for the direct sum of all H p (C * ) for p ≥ 0.
A chain map φ : C * → D * induces a linear map φ * : H (C * ) → H (D * ) between the homology groups. Again, this construction is functorial, meaning that it maps identity maps to identity maps, and it is compatible with compositions. Acyclic carriers. We call a simplicial complex K acyclic, if K is connected and all homology groups H p (K ) are trivial. For simplicial complexes K and L, an acyclic carrier is a map that assigns to each simplex σ in K , a non-empty acyclic subcomplex (σ ) ⊆ L, and whenever τ is a sub-simplex of σ , then (τ ) ⊆ (σ ). We say that a chain c ∈ C p (K ) is carried by a subcomplex K , if c takes value 0 except for p-simplices in K . A chain map φ : C * (K ) → C * (L) is carried by , if for each simplex σ ∈ K , φ(σ ) is carried by (σ ). We state the acyclic carrier theorem (Munkres 1984, Thm 13.3), adapted to our notation: Theorem 1 Let : K → L be an acyclic carrier. Then, -There exists an augmentation-preserving chain map φ : C * (K ) → C * (L) carried by . -If two augmentation-preserving chain maps φ 1 , φ 2 : C * (K ) → C * (L) are both carried by , then φ * 1 = φ * 2 . 4 We remark that "augmentation-preserving" is crucial in the statement: without it, the trivial chain map (that maps everything to 0) turns the first statement trivial and easily leads to a counter-example for the second claim.
Filtrations and towers. Let I ⊆ R be a set of real values which we refer to as scales. A filtration is a collection of simplicial complexes (K α ) α∈I such that K α ⊆ K α for all α ≤ α ∈ I . For instance, (R α ) α≥0 is a filtration which we call the Rips filtration. A (simplicial) tower is a sequence (K α ) α∈J of simplicial complexes with J being a discrete set (for instance J = {2 k | k ∈ Z}), together with simplicial maps ϕ α : K α → K α between complexes at consecutive scales. For instance, the Rips filtration can be turned into a tower by restricting to a discrete range of scales, and using the inclusion maps as ϕ. The approximation constructed in this paper will be another example of a tower.
We say that a simplex σ is included in the tower at scale α , if σ is not in the image of the map ϕ α : K α → K α , where α is the scale preceding α in the tower. The size of a tower is the number of simplices included over all scales. If a tower arises from a filtration, its size is simply the size of the largest complex in the filtration (or infinite, if no such complex exists). However, this is not true in general for simplicial towers, because simplices can collapse in the tower and the size of the complex at a given scale may not take into account the collapsed simplices which were included at earlier scales in the tower. Barcodes and Interleavings. A collection of vector spaces (V α ) α∈I connected with linear maps λ α 1 ,α 2 : V α 1 → V α 2 is called a persistence module, if λ α,α is the identity on V α and λ α 2 ,α 3 • λ α 1 ,α 2 = λ α 1 ,α 3 for all α 1 ≤ α 2 ≤ α 3 ∈ I for the index set I .
We generate persistence modules using the previous concepts. Given a simplicial tower (K α ) α∈I , we generate a sequence of chain complexes (C * (K α )) α∈I . By functoriality, the simplicial maps ϕ of the tower give rise to chain maps ϕ between these chain complexes. Using functoriality of homology, we obtain a sequence (H (K α )) α∈I of vector spaces with linear maps ϕ * , forming a persistence module. The same construction applies to filtrations as a special case.
Persistence modules admit a decomposition into a collection of intervals of the form [α, β] (with α, β ∈ I ), called the barcode, subject to certain tameness conditions. The barcode of a persistence module characterizes the module uniquely up to isomorphism. If the persistence module is generated by a simplicial complex, an interval [α, β] in the barcode corresponds to a homological feature (a "hole") that comes into existence at complex K α and persists until it disappears at K β .
Two persistence modules (V α ) α∈I and (W α ) α∈I with linear maps φ ·,· and ψ ·,· are said to be weakly (multiplicatively) c-interleaved with c ≥ 1, if there exist linear maps γ α : V α → W cα and δ α : W α → V cα , called interleaving maps, such that the diagram (1) commutes, that is, ψ = γ •δ and φ = δ •γ for all {. . . , α/c 2 , α/c, α, cα, . . . } ∈ I (we have skipped the subscripts of the maps for readability). In such a case, the barcodes of the two modules are 3c-approximations of each other in the sense of Chazal et al. (2009). We say that two towers are c-approximations of each other if their persistence modules are c-approximations.
Under the more stringent conditions of strong interleaving, the approximation ratio can be improved. Two persistence modules (V α ) α≥0 and (W α ) α≥0 with respective linear maps φ ·,· and ψ ·,· are said to be (multiplicatively) strongly c-interleaved if there exist a pair of families of linear maps γ α : V α → W cα and δ α : W α → V cα for c > 0, such that Diagram (2) commutes for all 0 ≤ α ≤ α (the subscripts of the maps are excluded for readability). In such a case, the persistence barcodes of the two modules are said to be c-approximations of each other in the sense of Chazal et al. (2009).
Finally, we mention a special case that relates equivalent persistence modules (Carlsson and Zomorodian 2005;Goodman et al. 2017). Two persistence modules V = (V α ) α∈I and W = (W α ) α∈I that are connected through linear maps φ, ψ respectively are isomorphic if there exists an isomorphism f α : V α → W α for each α ∈ I for which the following diagram commutes for all α ≤ β ∈ I : Isomorphic persistence modules have identical persistence barcodes. Scale balancing. Let V = (V α ) α∈I and W = (W α ) α∈I be two persistence modules with linear maps f v , f w , respectively. Let there be linear maps φ : V α/ε 1 → W α and ψ : W α → V αε 2 for 1 ≤ ε 1 , ε 2 such that all α, α/ε 1 , αε 2 ∈ I . Suppose that the following diagram commutes, for all α ∈ I .
We define a new vector space V cα := V α , where c = ε 1 ε 2 and cα ∈ I . This gives rise to a new persistence module, V = (V cα ) α∈I . The maps φ and ψ can then be interpreted as φ : Then, Diagram (4) can be re-interpreted as . V and V have the same barcode up to a scaling factor. This scaling trick also works when V and W are strongly interleaved. If we have the following commutative diagrams: (where we have skipped the maps for readability): then V and W are max(ε 1 , ε 2 )-approximations of each other. By defining V as before, the following diagrams commute for d = cε 2 = √ ε 1 ε 2 , so we can improve a max(ε 1 , ε 2 )-approximation to an √ ε 1 ε 2 -approximation. We end the section by discussing a basic but important relation betweenČech and Rips filtrations. It is well-known that for any α ≥ 0, C α ⊆ R α ⊆ C √ 2α (Edelsbrunner and Harer 2010). This gives a strong interleaving between the towers (C α ) α≥0 and (R α ) α≥0 with ε 1 = 1 and ε 2 = √ 2. Applying the scale balancing technique, we get that Lemma 1 The scaledČech persistence module (H (C 4 √ 2α )) α≥0 and the Rips persistence module (H (R α )) α≥0 are 4 √ 2-approximations of each other.

Shifted integer lattices
In this section, we take a look at simple modifications of the integer lattice. We denote by I := {α s := λ2 s | s ∈ Z} with λ > 0, a discrete set of scales. For each scale in I , we define grids which are scaled and translated (shifted) versions of the integer lattice.
Definition 1 (Scaled and shifted grids) For each scale α s ∈ I , we define the scaled and shifted grid G α s inductively as: -For s = 0, G α s is simply the scaled integer grid λZ d , where each basis vector has been scaled by λ.
where the signs of the components of the last vector are chosen independently and uniformly at random (and the choice is independent for each s).
where the last vector is chosen as in the case of s ≥ 0.
Equations (8) and (9) are consistent at s = 0. A simple example of the above construction is the sequence of grids with G α s := α s Z d for even s, and G α s := α s Z d + α s−1 2 (1, · · · , 1) for odd s. Next, we motivate the shifting of the grids. Let Vor G s (x) denote the Voronoi cell of any point x ∈ G s with respect to the point set G s . It is clear that the Voronoi cell is a cube of side length α s centered at x. The shifting of the grids ensures that each x ∈ G α s lies in the Voronoi region of a unique y ∈ G α s+1 . Using an elementary calculation, we show a stronger statement: Proof Without loss of generality, we can assume that α s = 2 and x is the origin, using an appropriate translation and scaling. Also, we assume for the sake of simplicity that G α s+1 = 2G α s + (1, · · · , 1); the proof is analogous for any other translation vector.
In that case, it is clear that y The claim follows. For an example look to Fig. 1.

Cubical complex of Z d
The integer grid Z d naturally defines a cubical complex, where each element is an axis-aligned, k-dimensional cube with 0 ≤ k ≤ d. To define it formally, let denote the set of all integer translates of faces of the unit cube [0, 1] d , considered as a convex polytope in R d . We call the elements of faces of Z d . Each face has a dimension k; the 0-faces, or vertices are exactly the points in Z d . The facets of a k-face E are the (k − 1)-faces contained in E. We call a pair of facets of E opposite facets, if they are disjoint. Naturally, these concepts carry over to scaled and shifted versions of Z d , so we define α s as the cubical complex defined by G α s .
We define a map g α s : α s → α s+1 as follows: for vertices of α s , we assign to }; the next lemma shows that this is a well-defined map. In this paper, we sometimes call g α s a cubical map, since it is a counterpart of simplicial maps for cubical complexes.

Lemma 3 Let f be k-face of α s with vertices
-if e 1 , e 2 are any two opposite facets of e, then there exists a pair of opposite facets f 1 , f 2 of f such that g α s ( f 1 ) = e 1 and g α s ( f 2 ) = e 2 .
Proof First claim: We prove the first claim by induction on the dimension of faces of G α s . Base case: for vertices, the claim is trivial using Lemma 2. Induction case: let the claim hold true for all (k − 1)-faces of G α s . We show that the claim holds true for all k-faces of G α s . Let f be a k-face of G α s . Let f 1 and f 2 be opposite facets of f , along the m-th coordinate. Let us denote the vertices of f 1 by ( p 1 , · · · , p 2 k−1 ) and those of f 2 by ( p 2 k−1 +1 , · · · , p 2 k ) taken in the same order, that is, p j and p 2 k−1 + j differ in only the m-th coordinate for all 1 ≤ j ≤ 2 k−1 . By definition, all vertices of f 1 share the m-th coordinate, and we denote coordinate of these vertices by z. Then, the m-th coordinate of all vertices of f 2 equals z + α s . Then g α s ( p j ) and g α s ( p 2 k−1 + j ) have the same coordinates, except possibly the m-th coordinate. By induction hypothesis, e 1 = g α s ( f 1 ) and e 2 = g α s ( f 2 ) are two faces of G s+1 . This implies that e 2 is a translate of e 1 along the m-th coordinate.
There are two cases: if e 1 and e 2 share the m-th coordinate, then e 1 = e 2 and therefore g α s ( f ) = e 1 = e 2 = e, so the claim follows. On the other hand, if e 1 and e 2 do not share the m-th coordinate, then they are two faces of α s+1 which differ in only one coordinate by α s+1 . So they are opposite facets of a co-dimension one face e of G α s+1 . Using induction, the claim follows.
Second claim: We prove the claim by induction over the dimension of e 1 . Base case: e 1 is a vertex. The vertices of f in Voronoi region of e 1 form f 1 . Since f is an axis parallel face and the Voronoi region is also axis-parallel, it is immediate that f 1 is a face of f . Assume that the claim is true up to dimension i. For e 1 a face of dimension i + 1, consider opposite facets e a and e b of e. By the induction claim, there exist faces would be common to both e a and e b , a contradiction. If e a is a translate of e b along the m-th coordinate, then f a is also a translate of f b along the same coordinate. Therefore f a and f b are opposite faces of a face f 1 and g α s ( f 1 ) = e 1 .
Third claim: Without loss of generality, assume that x 1 is the direction in which e 2 is a translate of e 1 . Using the second claim, let h denote the maximal face of f such that g α s (h) = e 1 . Clearly, h = f , since that would imply g α s ( f ) = e 1 = e, which is a contradiction.
Suppose h has dimension less than k − 1. Let h be the facet of f that contains h and has the same x 1 coordinates for all vertices. Then g α s (h ) = e 1 , which contradicts the maximality of h.
Therefore, the only possibility is that h is a facet f 1 of f such that g α s ( f 1 ) = e 1 . Let f 2 be the opposite facet of f 1 . From the proof of the first claim, it is easy to see that g α s ( f 2 ) = e 2 . The claim follows.

Barycentric subdivision
We discuss a special triangulation of The barycentric subdivision of α s , denoted by sd α s , is the (infinite) simplicial complex whose simplices are the flags of α s (Munkres 1984).
In particular, the 0-simplices of sd α s are the faces of α s . An equivalent geometric description of sd α s can be obtained by defining the 0-simplices as the barycenters of the faces in sd α s , and introducing a k-simplex between (k + 1) barycenters if the corresponding faces form a flag. For a simple example, see Figs. 2 and 3. It is easy to see that sd α s is a flag complex. Given a face f in α s , we write sd( f ) for the subcomplex of sd α s consisting of all flags that are formed only by faces contained in f .

Approximation scheme with simplicial complexes
We define our approximation complex for a finite set of points in R d . Recall from Definition 1 that we can define a collection of scaled and shifted integer grids G α s over a collection of scales I := {α s = 2 s | s ∈ Z} in R d . To make the exposition simple, we define our complex in a slightly generalized form.

Barycentric spans
Fix some s ∈ Z and let V denote any non-empty subset of G α s . Vertex span. We say that a face f ∈ α s is spanned by V , if the set of vertices -is non-empty, and -not contained in any facet of f .
Trivially, the vertices of α s which are spanned by V are precisely the points in V . Any face of α s which is not a vertex must contain at least two vertices of V in order to be spanned. We point out that the set of spanned faces of α s is not closed under taking sub-faces. For instance, if V consists of two antipodal points of a d-cube, the only faces spanned by V are the d-cube and the two vertices; all other faces of the d-cube contain at most one vertex and hence are not spanned.
It is simple to test whether any given k-face f ∈ α s is spanned by the set of Furthermore, for any non-empty subset W ⊆ V , the faces of α s that are spanned by W are also spanned by V . Consequently, the barycentric span of W is a subcomplex of the barycentric span of V .

Approximation complex
We denote by P ⊂ R d a finite set of points. We define two maps: a α s : P → G α s : for each point p ∈ P, we let a α s ( p) denote the grid point in G α s that is closest to p, that is, p ∈ Vor G αs (a α s ( p)). We assume for simplicity that this closest point is unique, which can be ensured using well-known methods (Edelsbrunner and Mücke 1990). We define the active vertices of G α s as that is, the set of grid points that have at least one point of P in their Voronoi cells. b α s : V α s → P: the map b α s takes an active vertex of G α s to its closest point in P.
By taking an arbitrary total order on P to resolve multiple assignments, we ensure that this assignment is unique.
Recall that the map g α s : α s → α s+1 takes grid points of G α s to grid points of G α s+1 . Using Lemma 2, it follows at once that: Recall that R ∞ α denotes the Rips complex at scale α for the L ∞ -norm. The next statement is a direct application of the the triangle inequality; let diam ∞ () denote the diameter in the L ∞ -norm.
Lemma 5 Let Q ⊆ P be a non-empty subset such that diam ∞ (Q) ≤ α s . Then, the set of grid points a α s (Q) is contained in a face of α s .
Proof We prove the claim by contradiction. Suppose that the set of active vertices a α s (Q) is not contained in a face of α s . Then, there exists at least one pair of points {x, y} ∈ Q such that a α s (x), a α s (y) are not in a common face of α s . By the definition of the grid G α s , the grid points a α s (x), a α s (y) therefore have L ∞ -distance at least 2α s . Moreover, x has L ∞ -distance less than α s /2 from a α s (x), and the same is true for y and a α s (y). By the triangle inequality, the L ∞ -distance of x and y is more than α s , which is a contradiction to the fact that diam ∞ (Q) ≤ α s .
We now define our approximation tower. For any scale α s , we define X α s as the barycentric span of the active vertices V α s ⊂ G α s . See Figs. 4, 5 and 6 for a simple illustration.
To simplify notation, we denote -the faces of α s spanned by V α s as active faces, and -the faces of active faces that are not spanned by V α s as secondary faces.
To complete the description of the approximation tower, we need to define simplicial maps of the formg α s : X α s → X α s+1 , which connect the simplicial complexes at consecutive scales. We show that such maps are induced by g α s .
Lemma 6 Let f be any active face of α s . Then, g α s ( f ) is an active face of α s+1 .
Proof Using Lemma 3, e := g α s ( f ) is a face of α s . If e is a vertex, then it is active, because f contains at least one active vertex v, and g α s (v) = e in this case. If e is not a vertex, we assume for a contradiction that it is not active. Then, it contains a facet e 1 that contains all active vertices in e. Let e 2 denote the opposite facet of e 1 in e. By Lemma 3, f contains opposite facets f 1 , f 2 such that g α s ( f 1 ) = e 1 and g α s ( f 2 ) = e 2 .
Since f is active, both f 1 and f 2 contain active vertices; in particular, f 2 contains an active vertex v. But then the active vertex g α s (v) must lie in e 2 , contradicting the fact that e 1 contains all active vertices of e.  The generated approximation complex, whose vertices consist of those of the cubical complex and the blue vertices (small dots), which are the barycenters of active and secondary faces As a result, g is well defined for each face e ∈ α s , since there exists some active face e ∈ α s with e ⊆ e , and g(e) ⊆ g(e ). By definition, a simplex σ ∈ X α s is a flag ( f 0 ⊆ · · · ⊆ f k ) of faces in α s . We set where (g α s ( f 0 ) ⊆ · · · ⊆ g α s ( f k )) is a flag of faces in α s+1 by Lemma 6, and hence is a simplex in X α s+1 . It follows thatg s : X α s → X α s+1 is a simplicial map. This completes the description of the simplicial tower X α s s∈Z .

Interleaving with the Rips module
First, we show that our tower is a constant-factor approximation of the the L ∞ -Rips filtration of P. We then show the relation between our approximation tower and the Euclidean Rips filtration of P.
We start by defining two acyclic carriers. First, we set λ = 1 and abbreviate α := α s = 2 s to simplify notation.
The barycentric span of any subset of U is a subcomplex of the barycentric span of U , so C α 1 is a carrier. Therefore, C α 1 is an acyclic carrier. -C α 2 : X α → R ∞ α : let σ be any flag of X α and let E be the smallest active face of α that contains σ (we break ties by making use of an arbitrary global order on P) 5 . We collect all the points of P that map to vertices of E under the map a α and set C α 2 (σ ) as the simplex on this set of points. By an application of the triangle inequality, we see that the L ∞ -diam of C α 2 (σ ) is at most 2α, so C α 2 (σ ) ∈ R ∞ α and is acyclic. It is also clear that C α 2 (τ ) ⊆ C α 2 (σ ) for each τ ⊆ σ , so C α 2 is an acyclic carrier.
Using the acyclic carrier theorem (Theorem 1), there exist augmentation-preserving chain maps between the chain complexes, which are carried by C α 1 and C α 2 respectively, for each α ∈ I . We obtain the following diagram of augmentation-preserving chain maps: where inc corresponds to the chain map for inclusion maps, andg denotes the chain map for the corresponding simplicial map g (we removed indices of the maps for readability). The chain complexes give rise to a diagram of the corresponding homology groups, connected by the induced linear maps c * 1 , c * 2 , inc * ,g * : Lemma 7 For all α ∈ I , the linear maps in the lower triangle of Diagram (11) commute, that is,g * = c * 1 • c * 2 .
Proof We look at the corresponding triangle in Diagram (10). We show that the (augmentation-preserving) chain mapsg and c 1 • c 2 are both carried by an acyclic carrier D : X α → X 2α . The claim then follows from the acyclic carrier theorem. Let σ ∈ X α be any flag and let E ∈ α denote the minimal active face containing σ . Let {q 1 , . . . , q k } be the active vertices of E. Let { p 1 , . . . , p m } be the set of points of P that map to {q 1 , . . . , q k } under the map a α . Since the L ∞ -diameter of { p 1 , . . . , p m } is at most 2α, using Lemma 5 we see that {a 2α ( p 1 ), . . . , a 2α ( p m )} is a face of 2α . We set D(σ ) as the barycentric span of {a 2α ( p 1 ), . . . , a 2α ( p m )}. It follows that D is an acyclic carrier.

Lemma 8 For all α ∈ I , the linear maps in the upper triangle of Diagram
Proof The proof technique is analogous to the proof of Lemma 7. We define an acyclic carrier D : R ∞ α → R ∞ 2α which carries inc and c 2 •c 1 , both of which are augmentationpreserving.
Let σ = ( p 0 , · · · , p k ) ∈ R ∞ α be any simplex. The set of active vertices lie in a face f of G 2α , using Lemma 5. We can assume that f is active, as otherwise, we argue about a facet of f that contains U . We set D(σ ) as the simplex on the subset of points in P, whose closest grid point in G 2α is any vertex of f . Using the triangle inequality we see that D(σ ) ∈ R ∞ 2α , so D is an acyclic carrier. The vertices of σ are a subset of D(σ ), so D carries the map inc. Showing that D carries c 2 • c 1 requires further explanation.
Using Lemmas 7 and 8, we see that the two persistence modules H (X α s ) s∈Z and H (R ∞ α ) α≥0 are weakly 2-interleaved. With elementary modifications in the definition of X andg, we can get a tower of the form (X α ) α≥0 . Furthermore, with minor changes in the interleaving arguments, we show that the corresponding persistence module is strongly 4-interleaved with the L ∞ -Rips module. Using scale balancing, this result improves to a strong 2-interleaving (see Lemma 16). Since the techniques used in the proof are very similar to the concepts used in this section, for the sake of brevity we defer all further details to Appendix A.
Using the strong stability theorem for persistence modules and taking scale balancing into account, we immediately get that: Theorem 2 The scaled persistence module H (X 2α ) α≥0 and the L ∞ -Rips persistence module H (R ∞ α ) α≥0 are 2-approximations of each other.
For any pair of points p, p ∈ R d , it holds that This in turn shows that the L 2 -and the L ∞ -Rips filtrations are strongly √ d-interleaved. Using the scale balancing technique for strongly interleaved persistence modules, we get:
Using Theorem 2, Lemma 9 and the fact that interleavings satisfy the triangle inequality (Bubenik and Scott 2014, Theorem 3.3), we see that the module (H (X 2α )) α≥0 is strongly 2d 0.25 -interleaved with the scaled Rips persistence module (H (R α/d 0.25 )) α≥0 . We can remove the scaling in the Rips filtration simply by multiplying the scales on both sides with d 0.25 and obtain our final approximation result:

Computational complexity
In this section, we discuss the computational aspects of constructing the approximation tower. In Sect. 5.1 we discuss the size complexity of the tower. An algorithm to compute the tower efficiently is presented in Sect. 5.2. Range of relevant scales. Set n := |P| and let C P(P) denote the closest pair distance of P. At scale α 0 := C P(P) 3d and lower, no two active vertices lie in the same face of the grid, so the approximation complex consists of n isolated 0-simplices. At scale α m := diam(P) and higher, points of P map to active vertices of a common face (by Lemma 5), so the generated complex is acyclic. We inspect the range of scales [α 0 , α m ] to construct the tower, since the barcode is explicitly known for scales outside this range. For this, we set λ = α 0 in the definition of the scales. The total number of scales is log 2 α m /α 0 = log 2 diam(P)3d C P(P) = log 2 + log 2 3d = O(log + log d), where = diam(P) C P(P) is the spread of the point set.

Size of the tower
The size of a tower is the number of simplices that do not have a preimage, that is, the number of simplex inclusions in the tower. We start by counting the number of active faces used in the tower.

Lemma 10
The number of active faces without pre-image in the tower is at most n3 d .
Proof At scale α 0 , there are n inclusions of 0-simplices in the tower, due to n active vertices. Using Lemma 2, g is surjective on the active vertices of (for any scale). Hence, no further active vertices are added to the tower. It remains to count the maximal active faces of dimension ≥ 1 without preimage. We will use a charging argument, charging the existence of such an active face to one of the points in P. We show that each point of P is charged at most 3 d − 1 times, which proves the claim. For that, we first fix an arbitrary total order ≺ on P. Each active vertex on any scale has a non-empty subset of P in its Voronoi region; we call the maximal such point with respect to the order ≺ the representative of the active vertex.
For each active face f of dimension at least one, we define the signature of f as the set of representatives of the active vertices of f . If for any set of active vertices u 1 , . . . , u k we have that v = g(u 1 ) = · · · = g(u k ), then the representative of v is one of the representatives of u 1 , . . . , u k , using Lemma 2. Therefore, the signatures of the active faces that are images of f under g are subsets of the signature of f . This implies that each maximal active face that is included has a unique maximal signature. We bound the number of maximal signatures to get a bound on the number of maximal active face inclusions. We charge the addition of each maximal signature to the lowest ordered point according to ≺.
Each signature contains representatives of active vertices from a face of α . Since each active vertex v has 3 d − 1 neighboring vertices in the grid that lie in a common face, the representative p of v can be charged 3 d − 1 times. There is a canonical isomorphism between the neighboring vertices of v at each scale. Then, for p to be charged more times, the image of v and some neighboring vertex u must be identical under g at some scale. But then, the representative of g(v) = g(u) is not p anymore, since p was the lowest ranked point in its neighborhood, hence the representative changes when the Voronoi regions are combined. So, p could not have been charged in such a case. Therefore, each point p ∈ P is indeed charged at most 3 d − 1 times.
There are n active faces of dimension 0 and at most n(3 d − 1) active faces of higher dimension. The upper bound is n + n(3 d − 1) = n3 d , as claimed.

Theorem 4 The k-skeleton of the tower has size at most
where a b denotes the Stirling number of the second kind.
Proof Each k-simplex that is included in the tower at any given scale α is a part of the barycentric subdivision of an active face that is also included at α. Therefore, we can account for the inclusion of this simplex by including the barycentric subdivision of its parent active face. From Lemma 10 at most n3 d active faces are included in the tower over all dimensions. We bound the number of k-simplices in the barycentric subdivision of a d-cube. Multiplying with n3 d gives the required bound.
Let c be any d-cube of α . To count the number of flags of length (m +1) contained in c that start with some vertex and end with c, we use similar ideas as in Edelsbrunner and Kerber (2012): first, we fix any vertex v of c and count the flags of the form v ⊆ · · · ⊆ c. Every -face in c incident to v corresponds to a subset of coordinate indices, in the sense that the coordinates not chosen are fixed to the coordinates of v for the face. With this correspondence, a flag from v to c of length (m +1) corresponds to an ordered m-partition of {1, · · · , d}. The number of such partitions is known as m! times the quantity d m , which is the Stirling number of second kind (Rennie and Dobson 1969), and is upper bounded by 2 O(d log m) . Since c has 2 d vertices, the total number of flags v ⊆ · · · ⊆ c of length (m + 1) with any vertex v is hence 2 d m! d m .
We now count the number of flags of length k + 1. Each such flag is (k + 1)-subset of some flag of length m = k + 3 that start with a vertex and end with c. There such flags and each of them has k+3 k+1 = (k + 3)(k + 2)/2 subsets of size (k + 1). The number of (k + 1)-flags is upper bounded by

Computing the tower
From Sect. 3, we know that G α s+1 is built from G α s by making use of an arbitrary translation vector (±1, · · · , ±1) ∈ Z d . In our algorithm, we pick the components of this translation vector uniformly at random from {+1, −1}, and independently for each scale. The choice behind choosing this vector randomly becomes more clear in the next lemma. From the definition, the cubical maps g α s : α s → α s+1 can be composed for multiple scales. For a fixed α s , we denote by g ( j) : α s → α s+ j the j-fold composition of g, that is, for j ≥ 1.

Lemma 11 For any k-face f ∈ α s with 1 ≤ k ≤ d, let Y denote the minimal integer j such that g ( j) ( f ) is a vertex, for a given choice of the randomly chosen translation vectors. Then, the expected value of Y satisfies
which implies that no face of α s survives more than 3 log d scales in expectation.
Proof Without loss of generality, assume that the grid under consideration is Z d and f is the k-face spanned by the vertices {{0, 1}, · · · , {0, 1} k , 0, · · · , 0}, so that the origin is a vertex of f . The proof for the general case is analogous.
We say that the x 1 -coordinate collapses in the first case and survives in the second. Both events occur with the same probability 1/2. Because the shift is chosen uniformly at random for each scale, the probability that x 1 did not collapse after j iterations is 1/2 j . f spans k coordinate directions, so it must collapse along each such direction to contract to a vertex. Once a coordinate collapses, it stays collapsed at all higher scales. As the random shift is independent for each coordinate direction, the probability of a collapse is the same along all coordinate directions that f spans. Using the union bound, the probability that g j ( f ) has not collapsed to a vertex is at most k/2 j . With Y as in the statement of the lemma, it follows that Hence, As a consequence of the lemma, the expected "lifetime" of k-simplices in our tower with k > 0 is rather short: given a flag e 0 ⊆ · · · ⊆ e , the face e will be mapped to a vertex after O(log d) steps, and so will be all its sub-faces, turning the flag into a vertex. It follows that summing up the total number of k-simplices with k > 0 over X α for all α ≥ 0 yields an upper bound of n2 O(d log k+d) as well.

Algorithm description
Recall that a simplicial map can be written as a composition of simplex inclusions and contractions of vertices (Dey et al. 2014;Kerber and Schreiber 2017). That means, given the complex X α s , to describe the complex at the next scale α s+1 , it suffices to specify -which pairs of vertices in X α s map to the same image underg, and -which simplices in X α s+1 are included at scale X α s+1 .
The input is a set of n points P ⊂ R d . The output is a list of events, where each event is of one of the three following types: -A scale event defines a real value α and signals that all upcoming events happen at scale α (until the next scale event). -An inclusion event introduces a new simplex, specified by the list of vertices on its boundary (we assume that every vertex is identified by a unique integer). -A contraction event is a pair of vertices (i, j) from the previous scale, and signifies that i and j are identified as the same from that scale.
In a first step, we estimate the range of scales that we are interested in. We compute a 2-approximation of diam(P) by taking any point p ∈ P and calculating max q∈P p − q . Then we compute C P(P) using a randomized algorithm in n2 O(d) expected time (Khuller and Matias 1995).
Next, we proceed scale-by-scale and construct the list of events accordingly. On the lowest scale, we simply compute the active vertices by point location for P in a cubical grid, and enlist n inclusion events (this is the only step where the input points are considered in the algorithm).
For the data structure, we use an auxiliary container S and maintain the invariant that whenever a new scale is considered, S consists of all simplices of the previous scale, sorted by dimension. In S, for each vertex, we store an id and a coordinate representation of the active face to which it corresponds. Every -simplex with > 0 is stored just as a list of integers, denoting its boundary vertices. We initialize S with the n active vertices at the lowest scale.
Let α < α be any two consecutive scales with , the respective cubical complexes and X , X the approximation complexes, withg : X → X being the simplicial map connecting them. Suppose we have already constructed all events at scale α.
-First, we enlist the scale event for α . -Then, we enlist the contraction events. For that, we iterate through the vertices of X and compute their value under g, using point location in a cubical grid. We store the results in a list S (which contains the simplices of X ). If for a vertex j, g( j) is found to be equal to g(i) for a previously considered vertex i, we choose the minimal such i and enlist a contraction event for (i, j). -We turn to the inclusion events: -We start with the case of vertices. Every vertex of X is either an active face or a secondary face of . Each active face must contain an active vertex, which is also a vertex of X . We iterate through the elements in S . For each active vertex v encountered, we go over all faces of the cubical complex that contain v as a vertex, and check whether they are active. For every active face E encountered that is not in S yet, we add it to S and enlist an inclusion event of a new 0-simplex. Additionally, we go over each face of E, add it to S and enlist a vertex inclusion event, thereby enumerating the secondary faces that are in E. At termination, all vertices of X have been detected. -Next, we iterate over the simplices of S of dimension ≥ 1, and compute their image underg using the pre-computed vertex map; we store the result in S . -To find the simplices of dimension ≥ 1 included at X , we exploit our previous insight that they contain at least one vertex that is included at the same scale (see the proof of Theorem 4). Hence, we iterate over the vertices included in X and find the included simplices inductively in dimension.
Let v be the current vertex under consideration; assume that we have found all ( p − 1)-simplices in X that contain v. Each such ( p − 1)-simplex σ is a flag of length p in . We iterate over all faces e that extend σ to a flag of length p + 1. If e is active, we have found a p-simplex in X incident to v. If this simplex is not in S yet, we add it and enlist an inclusion event for it. We also enqueue the simplex in our inductive procedure, to look for ( p + 1)-simplices in the next round. At the end of the procedure, we have detected all simplices in X without preimage, and S contains all simplices of X . We set S ← S and proceed to the next scale.
This ends the description of the algorithm.  O(d log k+d) and the space is bounded by n2 O(d log k+d) .
Proof In the analysis, we ignore the costs of point locations in grids, checking whether a face is active, and searches in data structures S, since all these steps have negligible costs when appropriate data structures are chosen.
Computing the image of a vertex of X costs O(2 d ) time. Moreover, there are at most n2 O(d) vertices altogether in the tower in expectation (using Lemma 10), so this bound in particular holds on each scale. Hence, the contraction events on a fixed scale can be computed in n2 O(d) time. Finding new active vertices requires iterating over the cofaces of a vertex in a cubical complex. There are 3 d such cofaces for each vertex. This has to be done for a subset of the vertices in X , so the running time is also n2 O(d) . Further, for each new active face, we go over its 2 O(d) faces to enlist the secondary faces, so this step also consumes n2 O(d) time. Since there are O(log + log d) scales considered, these steps require n2 O(d) log over all scales.
Computing the image ofg for a fixed scale costs at most O(2 d |X |). M is the size of the tower, that is, the simplices without preimage, and I is the set of scales considered. The expected bound for α∈I |X α | = O(log d M), because every simplex has an expected lifetime of at most 3 log d by Lemma 11. Hence, the cost of these steps is bounded by 2 O(d) M.
In the last step of the algorithm, we find the simplices of X included at α . We consider a subset of simplices of X , and for each, we iterate over a collection of faces in the cubical complex of size at most 2 O(d) . Hence, this step is also bounded by 2 O(d) |X | per scale, and hence bounded 2 O(d) M as well.
For the space complexity, the auxiliary data structure S gets as large as X , which is clearly bounded by M. For the output complexity, the number of contraction events is at most the number of inclusion events, because every contraction removes a vertex that has been included before. The number of inclusion events is the size of the tower. The number of scale events as described is O(log + log d). However, it is simple to get rid of this factor by only including scale events in the case that at least one inclusion or contraction takes place at that scale. The space complexity bound follows.

Dimension reduction
When the ambient dimension d is large, our approximation scheme can be combined with dimension reduction techniques to reduce the final complexity, very similar to the application in Choudhary et al. (2017b). For a set of n points P ⊂ R d , we apply the dimension reduction schemes of Johnson-Lindenstrauss (JL) (Johnson et al. 1986), Matoušek (MT) (Matoušek 1990), and Bourgain's embedding (BG) (Bourgain 1985). We then compute the approximation on the lower-dimensional point set. We only state the main results in Table 1, leaving out the proofs since they are very similar to those from Choudhary et al. (2017b).

Approximation scheme with cubical complexes
We extend our approximation scheme to use cubical complexes in place of simplicial complexes. We start by detailing a few aspects of cubical complexes.

Cubical complexes
We now briefly describe the concept of cubical complexes, essentially expanding upon the contents of Sect. 3.1. For a detailed overview of cubical homology, we refer to Kaczynski et al. (2004).

Definition
We define cubical complexes over the grids G α s . For any fixed α s , the grids G α s defines a natural collection of cubes. An elementary cube γ is a product of intervals γ = I 1 × I 2 × · · · × I d , where each interval is of the form I j = (x j , x j + m j ), such that the vertex (x 1 , · · · , x m ) ∈ G α s and each m j is either 0 or α s . That means, an (elementary) cube is simply a face of a d-cube of the grid. An interval I j is said to be degenerate if m j = 0. The dimension of γ is the number of non-degenerate intervals that defines it. We define the boundary of any interval as the two degenerate intervals that form its endpoints and denote this by ∂( Taking the boundary of any fixed subset of the intervals defining γ consecutively gives a sum of faces of γ . A cubical complex of G α s is a finite collection of cubes of G α s . We define chain complexes for the cubical case in the same way as in simplicial complexes. The chain complexes are connected by boundary homomorphisms, where the boundary of a cube is defined as: where (I 1 × · · · × ∂(I j ) × · · · × I d ) denotes the sum It can be quickly verified that for each cube γ , ∂ • ∂(γ ) = 0 since each term appears twice in the expression and the addition is over Z 2 .

Cubical maps and induced homology
Let T α s and T α t denote the cubical complexes defined by the grids G α s and G α t , respectively, for s ≤ t. We use the vertex map g : G α s → G α t to define a map between the cubical complexes. Note that if (a, b) are vertices of a cube of T α s that differ in one coordinate, then (g(a), g(b)) are vertices of a cube of T α t that differ in at most one coordinate. A cubical map is a map f : T α s → T α t defined using g, such that for each spans a cube of T α t . The cubical map can also be restricted to sub-complexes of T α s and T α t , provided that the image f (γ ) is well-defined.
Each cubical map also defines a corresponding continuous map between the underlying spaces of the respective complexes. Let x ∈ |γ | be a point in γ . Then, the coordinates of x can be uniquely written as The cubical map f gives rise to a chain map f # : C p (T α s ) → C p (T α t ) between the p-th chain groups of the complexes, for each p ∈ [0, · · · , d]. For each cube γ , so this gives a homomorphism between the chain groups.
Moving to the homology level, we get the respective homology groups H (T α s ) and H (T α t ) and the chain map from above induces a linear map between them. The concept of reduced homology and augmentation maps is also applicable to the cubical chain complexes. For a sequence of cubical complexes connected with cubical maps, this generates a persistence module.
Cubical filtrations and towers are defined in a similar manner to the simplicial case. A cubical filtration is a collection of cubical complexes (T α ) α∈I such that T α ⊆ T α for all α ≤ α ∈ I . A (cubical) tower is a sequence (T α ) α∈J of cubical complexes with J being an index set together with cubical maps between complexes at consecutive scales. A cubical tower can be written as a sequence of inclusions and contractions, where an inclusion refers to the addition of a cube and a contraction refers to collapsing a cube along a coordinate direction to either of the endpoints of the interval.

Description
We choose the simplest possible cubical complex to define our approximation cubical tower: for each scale α s , we define the cubical complex U α s as the set of active faces and secondary faces spanned by V α s . Hence the cubical complex is closed under taking faces and is well-defined. See Fig. 5 for a simple example.
Recall from Sect. 4 that for each s ∈ Z, U α s and U α s+1 are related by a cubical map g α s , which gives rise to the cubical tower U α s s∈Z .
We extend this to a tower (U α ) α≥0 by using techniques from Appendix A. In Sect. 4 we saw that the tower (X α ) α≥0 gives an approximation to the Rips filtration. The relation between the simplicial and cubical towers is trivial: X α s is simply a triangulation of |U α s |. Hence X α s and U α s have the same homology (Munkres 1984). Moreover, the simplicial map is derived from an application of the cubical map. In particular, the continuous versions of both maps are the same. For any 0 ≤ α ≤ β, let f 1 : H * (U α ) → H * (U β ) denote the homomorphism induced by the cubical map, f 2 : H * (X α ) → H * (X β ) denote the homomorphism induced by the simplicial map, and f 0 : H * (|X α | = |U α |) → H * (|X β | = |U β |) denote the homomorphism induced by the common continuous map.
It is well-established that f 1 = f 0 (Kaczynski et al. 2004, Chapter. 6) and f 2 = f 0 (Munkres 1984, Chapter. 2 To compute the cubical tower, we simply re-use the algorithm for the simplicial case, with small changes: -In the simplicial case, we used a container S to hold the simplices from the previous scale. We alter S to store the cubes from the previous scale. For each interval, we store an id and its coordinates. Each cube is stored as the set of ids of the intervals that define it. -At each scale, we enumerate the image of the cubical map by computing the image of each interval, and then use this pre-computed map to compute the image of (≥ 1)-dimensional cubes. -For the inclusions, we find all the active and secondary faces but do not compute the simplices. The inclusions in the cubical tower correspond exactly to the inclusions of active and secondary faces in the simplicial tower, so this enumerates all inclusions correctly.

Practicality
We now touch upon the practical aspects of our constructions. An implementation of our approximation scheme would be a tool that computes the (approximate) persistence barcode for any input data set. For any scheme to be useful in practice, it should be able to compute sufficiently close approximations using a reasonable amount of resources.
Our cubical tower consists of cubical complexes connected via cubical maps. To our knowledge, there are no algorithms to compute barcodes in this setting where the cubical maps are more than just trivial inclusions. As such, although our cubical scheme has exponentially lower theoretical guarantees compared to the simplicial tower, we can not hope to test it in practice unless the appropriate primitives are available. It could be an interesting research direction to develop this primitive and in particular investigate whether the techniques used in computing persistence barcodes for a simplicial tower allow a generalization to the cubical case.
It makes more sense to inspect the simplicial tower. We saw in Theorem 4 that the size of the tower is n6 d−1 (2k + 4)(k + 3)! d k + 2 . Unfortunately, this bound is already too large so that the storage requirement of the Algorithm (Theorem 7) explodes exponentially. Let us assume a conservative bound of 1 Byte of memory requirement per simplex. For a point set in d = 8 dimensions and k = 4, the complexity bound is already at least 4000 Terabytes, before factoring in n. For a point set in d = 10 dimensions and k = 5, this explodes to 10 20 Terabytes. While these are upper bounds, in practice the complexity will still need to be many orders of magnitude smaller to be feasile, which is unlikely. Even with conservative estimates our storage requirement is impractical. Therefore we are not very hopeful that implementing the scheme in its current state will provide any useful insight for high dimensional approximations. Making it implementation-worthy demands more optimizations and tools at the algorithmic level. This is worth another Algorithmic engineering project in its own right. We plan to pursue this line of research in the future. Since our focus in this paper was geared towards theoretical aspects of approximations, we exclude experimental results in the current work. We hope that a more careful implementation-focussed approach may prove more practical.
On the other hand, the upper bound for the cubical case is simply n6 d . Even for d = 10, the storage requirement would be less than 100 Megabytes before factoring in n. This is far more attractive than the simplicial case. As such, it may make more sense to invest time and effort in developing tools to compute barcodes in the cubical setup.

Summary
We presented an approximation scheme for the Rips filtration, with improved approximation ratio, size and computational complexity than previous approaches for the case of high-dimensional point clouds. In particular, we are able to achieve a marked reduction in the size of the approximation by using cubical complexes in place of simplicial complexes. This is in contrast to all other previous approaches that used simplicial complexes as approximating structures.
An important technique that we used in our scheme is the application of acyclic carriers to prove interleaving results. An alternative would to be explicitly construct chain maps between the Rips and the approximation towers; unfortunately, this make the interleaving analysis significantly more complex. While the proof of the interleaving in Sect. 4.3 is still technically challenging, it greatly simplifies by the usage of acyclic carriers. There is also no benefit in knowing the interleaving maps because they are only required for the analysis of the interleaving, and not for the actual computation of the approximation tower. We believe that this technique is of general interest for the construction of approximations of cell complexes.
Our simplicial tower is connected by simplicial maps; there are (implemented) algorithms to compute the barcode of such towers (Dey et al. 2014;Kerber and Schreiber 2017). It is also quite easy to adapt our tower construction to a streaming setting (Kerber and Schreiber 2017), where the output list of events is passed to an output stream instead of being stored in memory.
-for all α ∈ [α s , α s+1 ), we set X α = X α s , for any α s ∈ I . That means, the complex stays the same in the interval between any two scales of I , so we defineg as the identity within this interval.
These give rise to the tower (X α ) α≥0 , that is connected with the simplicial mapg. This modification helps in improving the interleaving with the Rips persistence module. First, we extend the acyclic carriers C 1 and C 2 from before to the new case: -C α 1 : R ∞ α → X 4α , α > 0: we define C 1 as before, simply changing the scales in the definition. It is straightforward to see that C 1 is still a well-defined acyclic carrier.
-C α 2 : X α → R ∞ α , α ≥ 0: this stays the same as before. It is simple to check that C 2 is still a well-defined acyclic carrier.
These give rise to augmentation-preserving chain maps between the chain complexes: using the acyclic carrier theorem as before (Theorem 1).

Lemma 12 The diagram
commutes on the homology level, for all 0 ≤ α ≤ α .
Proof Consider the acyclic carrier C 1 • inc : R ∞ α → X 4α . It is simple to verify that this carrier carries both c 1 • inc andg • c 1 , so the induced diagram on the homology groups commutes, from Theorem 1.

Lemma 13 The diagram
commutes on the homology level, for all 0 ≤ α ≤ α .
Proof We construct an acyclic carrier D : X α → R ∞ α which carries inc•c 2 and c 2 •g, thereby proving the claim (Theorem 1).
Consider any simplex σ ∈ X α and let E ∈ α be the minimal active face of containing σ . We set D(σ ) as the simplex on the set of input points of P, which lie in the Voronoi regions of the vertices of g(E). By the triangle inequality, D(σ ) is a simplex of R ∞ α , so that D is a well-defined acyclic carrier. It is straightforward to verify that D carries both c 2 •g and inc • c 2 .

Lemma 14 The diagram
commutes on the homology level, for all 0 ≤ α ≤ α .
Proof The diagram is essentially the same as the lower triangle of Diagram 10, with a change in the scales. As a result, the proof of Lemma 7 also applies for our claim directly.

Lemma 15
The diagram commutes on the homology level, for all 0 ≤ α ≤ α .
Proof The diagram can be re-interpreted as: The modified diagram is essentially the same as the upper triangle of Diagram 10, with a change in the scales and a replacement of c 1 withg • c 1 , that is equivalent to the chain map at the scale α . Hence, the proof of Lemma 8 also applies for our claim directly.