Paris, France Institute of Geometry, TU Graz<EMAIL_ADDRESS>Institute of Geometry, TU Graz<EMAIL_ADDRESS>Mickaël Buchet, Bianca
B. Dornelas, and Michael Kerber Theory of computation Computational geometry
Theory of computation Computational complexity and cryptography This research
has been supported by the Austrian Science Fund (FWF), grant numbers P 33765-N
and W1230.

###### Acknowledgements.

We thank Alexander Rolle for his valuable input and encouragement of this
work.

# Sparse Higher Order Čech Filtrations

Mickaël Buchet  Bianca B. Dornelas  Michael Kerber

###### Abstract

For a finite set of balls of radius $r$, the $k$-fold cover is the space
covered by at least $k$ balls. Fixing the ball centers and varying the radius,
we obtain a nested sequence of spaces that is called the $k$-fold filtration
of the centers. For $k=1$, the construction is the union-of-balls filtration
that is popular in topological data analysis. For larger $k$, it yields a
cleaner shape reconstruction in the presence of outliers. We contribute a
sparsification algorithm to approximate the topology of the $k$-fold
filtration. Our method is a combination and adaptation of several techniques
from the well-studied case $k=1$, resulting in a sparsification of linear size
that can be computed in expected near-linear time with respect to the number
of input points.

###### keywords:

Sparsification, $k$-fold cover, Higher order Čech complexes

##  1 Introduction

_Persistent homology_ [8, 20, 21] is a major branch of topological data
analysis with applications, for instance, in shape recognition [4], material
science [31] and biology [27, 32]. It studies the homological properties of
sequences of topological spaces. A standard construction is to take the
homogeneous union of balls, with increasing radius, centered at finitely many
points of $\mathbb{R}^{d}.$ We call these points _sites_ and refer to that
filtration as the _union-of-balls filtration_. For computational purposes, one
considers the homologically equivalent _Čech filtration_ , which is a sequence
of simplicial complexes that captures the intersection patterns of the balls
in the union-of-balls filtration [20, Chap.3; 28].

The drawback of the Čech filtration (as well as of the closely-related
Vietoris-Rips filtration) is that for $n$ sites, it consists of up to
$\binom{n}{m+1}$ $m$-simplices because every $(m+1)$-subset of balls
intersects at a sufficiently large radius. A technique to overcome this large
size is to _approximate_ the Čech (or Vietoris-Rips) filtration with another,
much smaller simplicial filtration with similar topological properties.
Technically, that means that the persistence modules induced by the homology
of the Čech filtration and its approximation are $\epsilon$-interleaved for an
arbitrary $\epsilon>0$ [10]. Several strategies have been devised to construct
such approximations with total size linear in $n$ for any fixed $\epsilon$
(see related work). Many of these approaches work by selecting only a subset
of the simplices of the Čech filtration, in which case we refer to the
approximation as a _sparsification_.

The union-of-balls filtration is a special case of the _$k$ -fold filtration_
built upon the _$k$ -fold cover_. For $n$ sites in $\mathbb{R}^{d}$ and $k\geq
1$ fixed, the $k$-fold cover is the subset of $\mathbb{R}^{d}$ consisting of
points contained in at least $k$ balls of radius $r$ centered at the sites.
Besides being a natural extension, $k$-fold filtrations are tightly related
with the $k$th neighbor distance that arises in the context of outlier removal
and processing of non-homogeneous data densities [11, 23, 39, 41]. For that
reason, they have received increased attention recently, both with regards to
computational [15, 22] and structural aspects [3].

For fixed $k$, the $k$-fold filtration can be equivalently expressed by its
nerve, which is a simplicial filtration called the _$k$ th order Čech
filtration_. It captures the intersection patterns of all $k$-wise
intersections of balls, which we call _lenses_. The aforementioned size issue
for Čech filtrations is even more important in the $k$th order case: the
filtration is defined over $\binom{n}{k}$ vertices (one for each $k$-subset of
sites) and consequently consists of
$\left(\genfrac{}{}{0.0pt}{}{\binom{n}{k}}{m+1}\right)$ $m$-simplices, making
it unrealistic to compute even for small values of $n$. Therefore we need to
reduce its size considerably while maintaining a good approximation quality.

### Contributions.

We propose the first sparsification of the $k$-fold filtration for a fixed
$k$. It is a simplicial filtration that, for $n$ sites in $\mathbb{R}^{d}$
(with constant $d$) and a given parameter $\epsilon>0$, is (multiplicatively)
$(1+\epsilon)$-interleaved with the $k$-fold filtration. Moreover, the number
of $m$-simplices in our sparsification is

$\mathcal{O}\left(nk^{k(m+1)}\left(\dfrac{96}{\epsilon}\right)^{\delta
k(m+1)}\right),$ (1)

where $\delta$ is the doubling dimension of $\mathbb{R}^{d}$. We point out
that for constant $k$ and $\epsilon$, the size of the filtration is linear in
the number of sites. This is remarkable because the $k$th order Čech
filtration, which captures the $k$-fold filtration exactly, already contains
$\binom{n}{k}$ vertices. Hence our construction avoids including the vast
majority of lenses into the sparsification.

We give an output-sensitive algorithm to compute our sparsification up to
dimension $m_{max}$ in

$\mathcal{O}\left(nk\log
n\log\Phi+Xk^{k+1}\left(\dfrac{96}{\epsilon}\right)^{k\delta}\cdot
m_{max}\right)$

expected time. Here $\Phi$ is the spread of the point set (i.e., the ratio of
diameter and smallest distance of two distinct points) and $X$ is the size of
the output complex, upper bounded by (1) with $m$ replaced by $m_{max}$. Again
considering everything but $n$ as constant, we get a running time of
$\mathcal{O}(n\log n)$.

### Techniques and related work.

The seminal work by Sheehy [42] was the first one to introduce a
sparsification technique for Vietoris-Rips filtrations yielding linear size
and $\mathcal{O}(n\log n)$ running time (assuming all other parameters as
constant). His technique extends to Čech complexes as well with minor
adaptations. Subsequent work [5, 7, 9, 17, 43] introduces several extensions,
variations, and simplifications of Sheehy’s original sparsification; all these
works share essentially the same size and complexity bounds.

Our results are achieved by combining several of these techniques used for
approximating in the case $k=1$, which required non-trivial adaptation for
larger values of $k$. The main idea is that for every site $p$, we define a
_removal radius_ such that, for radii larger than this removal radius, all
lenses involving $p$ are ignored. That means, for larger and larger radii, we
construct simplicial complexes with fewer and fewer sites to keep the size
small. To determine the removal radii of sites, we introduce the _$k$
-distance permutation_ which is an ordering of the sites based on the distance
to the $k$th closest neighbor. The $k$-distance permutation is a
generalization of the farthest point sampling [29] used in some sparsification
schemes [17, 43] and induces covering and packing properties analogous to
those of nets.

Although we opted to extend sparsification techniques, there is an alternative
line of research by Choudhary et al. [12, 13, 14] that defines approximations
of Čech complexes which are not sparsifications. They arrive at the same size
bound and complexity as the sparsifications for $k=1$. Approximate filtrations
are also actively researched in practice [2, 6, 18, 34, 38].

The $k$-fold cover and the higher order Čech complexes are also studied with
relation to _multiparameter persistence_ : considering the order $k$ as a
second varying parameter, we obtain the _multicover bifiltration_. Blumberg
and Lesnick [3] survey different multiparameter persistence approaches and
show a particularly strong stability result for multicovers. Sheehy [41]
introduces the barycentric bifiltrations, which is equivalent to the
multicover but whose size is prohibitively large. The question of computing
the multicover bifiltration exactly has been studied by Edelsbrunner and Osang
[23], whose results have been refined by Corbet et al. [15]. The latter
authors obtain an equivalent bifiltration to the multicover one but has total
size (over all choices of $k$) $\mathcal{O}(n^{d+1})$ for $n$ points in
$\mathbb{R}^{d}$ [15, Prop. 5]. Their construction rely on using _higher order
Voronoi diagrams_ and _Delaunay complexes_ [24]. That reduces the $n^{n}$ size
factor of Čech complexes, but cannot lead to linear size without further
improvements: the Delaunay filtration’s $d$-skeleton is of size
$\mathcal{O}(n^{\lceil d/2\rceil})$ [40], which is a substantial improvement
over the $\mathcal{O}(n^{d})$ size of the Čech $d$-skeleton, but still super-
linear for $d\geq 3$. Our approximation foregoes those constructions to reduce
the size dependency on $n$ further, with the trade-off that we get an
exponential dependency on $k$.

The size reduction in our construction is a consequence of ignoring lenses
after their removal radius. The idea of removing a lens beyond a certain
radius is justified geometrically by the fact that the remaining lenses cover
its entire area after a certain radius. This is only true, however, if we
_freeze_ a lens before removing it, that is, keep it unchanged for a short
time while the surrounding lenses keep growing. This concept was already
introduced in [9], from where we also adapt the elegant technique of lifting
the lenses to convex _cones_ in $\mathbb{R}^{d}\times\mathbb{R}.$ The
additional dimension, which is the radius $r,$ is needed because removing
simplices is not possible in filtrations.

The major geometric predicate for our computation is whether a set of balls is
intersecting, which can be dualized to computing the radius of the minimal
enclosing ball of the ball’s centers [19, 25, 26, 35]. However, the
aforementioned freezing of lenses makes this problem technically more
challenging. This question seems to be unaddressed in previous work, and we
give an efficient solution in the Euclidean setting.

### Outline.

Section 2 provides background definitions and results. Section 3 defines a
$k$-distance and uses it to construct the $k$-distance permutation of a point
set $P$. In Section 4 the permutation is used to define a sparse lens
filtration that approximates the $k$-fold cover. That results in a nerve
filtration that approximates the $k$th order Čech complex, as shown in Section
5. The size bound of that filtration is given in Section 6. Section 7 provides
an algorithm for computing the discrete sparse Čech filtration. We conclude
with Section 8.

##  2 Background

### Lenses and $k$-fold covers.

Given a point set $P\subseteq\mathbb{R}^{d}$ and a fixed $k\in\mathbb{N},$ an
element $p\in P$ is called a _site_ and a _$k$ -subset_ of $P$ is a subset
with $k$ sites. Let $\binom{P}{k}$ be the collection of all $k$-subsets of $P$
and $A\in\binom{P}{k}.$ Let also $B_{r}(a)$ denote the closed ball centered at
$a$ of radius $r$. The _lens_ corresponding to the $k$-subset $A$ at scale $r$
is

$L_{r}(A):=\underset{a\in A}{\bigcap}B_{r}(a).$

The _$k$ -fold cover_ of $P$ at scale $r$ is the union of lenses at scale $r$
over all $k$-subsets:

$\operatorname{L}(r,P):=\underset{A\in\left(\genfrac{}{}{0.0pt}{}{P}{k}\right)}{\bigcup}L_{r}(A).$

See Figure 1 for an example. Note the omission of $k$, which we consider
fixed, in the notation. When $P$ is clear from context, we drop it from the
notation as well and write $L_{r}$ instead.

Figure 1: Example of $2$\- (left) and $3$-fold (right) covers for a fixed
radius.

### Nerves.

We assume that the reader is familiar with (abstract) simplicial complexes
[20, Chap. 3]. For a finite collection $C$ of subsets of $\mathbb{R}^{d}$, we
can define a simplicial complex with vertex set $C$, called the _nerve_ of
$C$, as the set of all subsets of $C$ that have a non-empty mutual
intersection. Note that the nerve can contain simplices of larger dimension
than $d$. The nerve of the set of all lenses of $k$-subsets of $P$ at scale
$r$ is called the _$k$ th order Čech complex with radius $r$_ over $P$,
denoted by $\operatorname{\check{C}ech}_{r}\left(P,k\right)$.

### Filtrations and equivalence.

A collection of topological spaces (e.g., subsets of $\mathbb{R}^{d}$)
$\mathcal{C}=\\{C_{r}\\}_{r\geq 0}$ is called a _filtration_ if for all $r\leq
r^{\prime},$ it holds that $C_{r}\subseteq{C}_{r^{\prime}}$. The letter $r$
denotes the _scale parameter_ of the filtration. For $P$ and $k$ fixed, the
previous concepts yield two different ways of obtaining filtrations. On the
one hand, since $L_{r}\subseteq L_{r^{\prime}}$ for $r\leq r^{\prime}$, we get
the _$k$ -fold filtration_ $\mathcal{L}:=\\{L_{r}\\}_{r\geq 0}$. On the other
hand, we observe that $\operatorname{\check{C}ech}_{r}\left(P,k\right)$ is a
subcomplex of $\operatorname{\check{C}ech}_{r^{\prime}}\left(P,k\right)$ for
$r\leq r^{\prime}$ and hence we get the _$k$ th order Čech filtration_
$\\{\operatorname{\check{C}ech}_{r}\left(P,k\right)\\}_{r\geq 0}$.

Let $\mathcal{C}$ and $\mathcal{D}$ be two filtrations. We say that
$\mathcal{C}$ is _(homotopy) equivalent_ [37, Chap. 9] to $\mathcal{D}$ if
there exists a family of maps $\\{f_{r}\colon C_{r}\rightarrow D_{r}\\}_{r\geq
0}$ that are homotopy equivalences of spaces and additionally commute with the
inclusion maps of $\mathcal{C}$ and $\mathcal{D}$.

### Interleaving and approximations.

Let $\epsilon\geq 0$. Two filtrations $\mathcal{C}$ and $\mathcal{D}$ are
_(multiplicatively) $(1+\epsilon)$-interleaved_ if there exist families of
linear maps $f_{\cdot},g_{\cdot}$ such that the diagram

${C_{r}}$${{C}_{r(1+\epsilon)^{2}}}$${{C}_{r(1+\epsilon)^{4}}}$${{D}_{r(1+\epsilon)}}$${{D}_{r(1+\epsilon)^{3}}}$$\scriptstyle{f_{r}}$$\scriptstyle{f_{r(1+\epsilon)^{2}}}$$\scriptstyle{g_{r(1+\epsilon)}}$$\scriptstyle{g_{r(1+\epsilon)^{3}}}$

commutes for all $r$. Informally, interleaved filtrations with small
$\epsilon$ are good approximations of each other because every $D_{r}$ sits in
between two instances of $\mathcal{C}$ with close-by scale parameters.

If $\\{C_{r}\\}_{r\geq 0}$, $\\{D_{r}\\}_{r\geq 0}$ are
$(1+\epsilon_{1})$-interleaved and $\\{D_{r}\\}_{r\geq 0}$,
$\\{E_{r}\\}_{r\geq 0}$ are $(1+\epsilon_{2})$-interleaved, then
$\\{C_{r}\\}_{r\geq 0}$ and $\\{E_{r}\\}_{r\geq 0}$ are
$(1+\epsilon_{1})(1+\epsilon_{2})$-interleaved. Moreover, if $C_{r}\subseteq
D_{r}\subseteq C_{r(1+\epsilon)}$ for all $r\geq 0$, then $\mathcal{C}$ and
$\mathcal{D}$ are $(1+\epsilon)$-interleaved.

A filtration $\mathcal{C}$ is a _$(1+\epsilon)$ -approximation_ of another
filtration $\mathcal{D}$ if there exist filtrations $\mathcal{C}^{\prime}$ and
$\mathcal{D}^{\prime}$ such that $\mathcal{C}^{\prime}$ is equivalent to
$\mathcal{C}$, $\mathcal{D}^{\prime}$ is equivalent to $\mathcal{D}$ and
$\mathcal{C}^{\prime}$ and $\mathcal{D}^{\prime}$ are
$(1+\epsilon)$-interleaved. This is a symmetric relationship, so we can say
that $\mathcal{C}$ and $\mathcal{D}$ are $(1+\epsilon)$-approximate. If
additionally $C_{r}\subseteq D_{r}$ for all $r\geq 0$ we call $\mathcal{C}$ a
_$(1+\epsilon)$ -sparsification_ of $\mathcal{D}$. We point out that an
approximation between two filtrations implies interleaved _persistence
modules_ (see [20, Chap. 7]) in the sense of [10].

### The Persistent Nerve Theorem.

Consider a finite index set $I$ and a family of filtrations
$\\{U^{(i)}_{r}\\}_{r\geq 0}$ over $\mathbb{R}^{d}$, one for each $i\in I$.
The union filtration is $\\{U_{r}\\}_{r\geq 0}$, where $U_{r}:=\bigcup_{i\in
I}U^{(i)}_{r}$, and the nerve filtration is $\\{N_{r}\\}_{r\geq 0}$, where
$N_{r}$ is the nerve of $U_{r}$. The Persistent Nerve Theorem [1, Thm. 3.9]
states that if every $U^{(i)}_{r}$ is closed and convex, then
$\\{U_{r}\\}_{r\geq 0}$ and $\\{N_{r}\\}_{r\geq 0}$ are equivalent. As a
consequence, the Persistent Nerve Theorem implies that the $k$-fold and the
$k$th order Čech filtrations are equivalent: choose $I$ as the set of all
$k$-subsets of $P$ and $U^{(i)}_{r}$ as the lens indexed by $i$ at radius $r$,
which is a closed and convex set.

### Doubling dimension.

The _doubling constant_ $\Delta$ of $\mathbb{R}^{d}$ is such that any ball of
radius $r$ can be covered with at most $\Delta$ balls of radius $r/2$, for all
$r\geq 0$. The _doubling dimension_ of $\mathbb{R}^{d}$ is
$\delta:=\log_{2}\Delta,$ which is of order $\Theta(d)$ and hence constant for
this paper. Note that for finite point sets in $\mathbb{R}^{d}$ the doubling
dimension can be significantly smaller than $d$, for instance if the points
all lie close to a low-dimensional subspace.

To cover a ball $B$ of radius $r$ with balls of radius $r/4,$ one needs at
most $\Delta^{2}$ balls; with balls of radius $r/8$ one needs $\Delta^{3}$
balls and so on. Thus, to cover $B$ with balls of radius $r^{\prime},$ we have
to find the smallest $t$ such that $r/2^{t}\leq r^{\prime}.$ That is
$t=\lceil\log_{2}r/r^{\prime}\rceil.$ Then,
$\Delta^{t}\leq\Delta^{\log_{2}r/r^{\prime}+1}=2^{\delta}\left({r}/{r^{\prime}}\right)^{\delta}$
and $(2r/r^{\prime})^{\delta}$ balls of radius $r^{\prime}$ are sufficient to
cover $B.$

### Quadtreaps.

A _quadtreap_ [36] is a dynamic data structure for spherical range search. We
summarize its properties in a simplified form suitable for us: for a set $X$
of $n$ points in $\mathbb{R}^{d}$ (with $d$ constant), it can be built in
$\mathcal{O}(n\log n)$ expected time. It supports deletions of points in $X$
in expected $\mathcal{O}(\log n)$ time. Moreover, given a query point $q$ and
a radius $r$, it returns a list $S\subseteq X$ which is guaranteed to contain
all points in $X$ of distance $\leq r$ from $q$, and is guaranteed not to
contain any point in $X$ of distance $\geq 2r$ from $q$. The running time for
such a query is $\mathcal{O}(\log n+|S|)$.

##  3 $k$-distance permutation

Given some integer $k\geq 1$ and a finite data set $P\subseteq\mathbb{R}^{d}$
of $n\geq k$ sites, we define an order on the points in $P$ in which the sites
are denoted by $p_{1},\ldots,p_{n}$. Writing
$P_{i}:=\\{p_{1},\ldots,p_{i}\\}$, our order ensures that the $k$-fold cover
over $P_{i}$ approximates the $k$-fold cover over $P$, with increasing
approximation quality when $i$ increases.

The _$k$ -distance_ of $x\in\mathbb{R}^{d}$ to $P$, denoted by
$\operatorname{d}^{k}(x,P)$, is the distance from $x$ to its $k$th closest
neighbor in $P$. We define the _$k$ -distance permutation_ incrementally as
follows: we choose $p_{1},\ldots,p_{k}$ as arbitrary, pairwise distinct sites
from $P$. If $p_{1},\ldots,p_{i-1}$ are chosen for $k<i\leq n$, we set

$p_{i}:=\underset{q\in P\setminus
P_{i-1}}{\operatorname{argmax}}\operatorname{d}^{k}(q,P_{i-1}).$

Note that for $k=1$, we obtain the well-known farthest point sampling. We also
define

$\lambda_{i}:=\operatorname{d}^{k}(p_{i},P_{i-1})$

for $k+1\leq i\leq n$ and set $\lambda_{1},\ldots,\lambda_{k}$ to $\infty$, so
that the sequence $(\lambda_{1},\lambda_{2},\cdots,\lambda_{n})$ is non-
increasing. The next two properties of the $k$-distance permutation are
reminiscent of the packing and covering properties of $\epsilon$-nets [44,
Chap. 14].

######  Lemma 3.1 (Covering).

For all $k\leq i\leq n-1$, we have $\operatorname{L}(r,P_{i})\subseteq
L(r,P)\subseteq\operatorname{L}(r+\lambda_{i+1},P_{i})$.

######  Proof 3.2.

Recall the notation $L_{r}=L(r,P)$. $P_{i}\subseteq P$ immediately implies
$\operatorname{L}(r,P_{i})\subseteq L_{r}.$ Consider $x\in L_{r}.$ Then, $x\in
L_{r}(A)$ for some $A=\\{a_{1},a_{2},\ldots,a_{k}\\}\subseteq P$. If
$A\subseteq P_{i}$, the result follows. Otherwise, without loss of generality
let $a_{1}\notin P_{i}$. By definition of $\lambda_{i+1}$,
$\operatorname{d}^{k}(a_{1},P_{i})\leq\lambda_{i+1}$ and hence there are sites
$b_{1},b_{2},\ldots,b_{k}\in P_{i}$ with
$\operatorname{d}(a_{1},b_{j})\leq\lambda_{i+1}$ for all $1\leq j\leq k$.
Consequently,
$\operatorname{d}(x,b_{j})\leq\operatorname{d}(x,a_{1})+\operatorname{d}(a_{1},b_{j})\leq
r+\lambda_{i+1}$ and the $k$ closest sites to $x$ in $P_{i}$ are within
distance $r+\lambda_{i+1}$ of $x$, implying
$x\in\operatorname{L}(r+\lambda_{i+1},P_{i})$.

######  Lemma 3.3 (Packing).

For all $k+1\leq i\leq n$, each $p\in P_{i}$ has
$\operatorname{d}^{k}(p,P_{i}\setminus\\{p\\})\geq{\lambda_{i}}/{2}$.

######  Proof 3.4.

We do induction on $i$. For $i=k+1$, let $q$ be the $k$th closest neighbor of
$p_{k+1}$ in $P_{k}$. We have
$\operatorname{d}^{k}(p_{k+1},P_{k+1}\setminus\\{p_{k+1}\\})=\lambda_{k+1}\geq\lambda_{k+1}/2$
and, for any $p\in P_{k+1}\setminus\\{p_{k+1}\\}$,

$\operatorname{d}^{k}(p,P_{k+1}\setminus\\{p\\})=\max\limits_{p^{\prime}\in
P_{k+1}\setminus\\{p\\}}\operatorname{d}(p,p^{\prime})\geq\dfrac{\operatorname{d}(p,q)+\operatorname{d}(p,p_{k+1})}{2}\geq\dfrac{\operatorname{d}(q,p_{k+1})}{2}=\dfrac{\lambda_{k+1}}{2}.$

Hence the statement is true for $i=k+1$. Next we assume, for some $i\geq k+1$,
that for every $p\in P_{i},$
$\operatorname{d}^{k}(p,P_{i}\setminus\\{p\\})\geq\lambda_{i}/2$, and show the
statement for $i+1$.

For $p_{i+1}$, we have
$\operatorname{d}^{k}(p_{i+1},P_{i+1}\setminus\\{p_{i+1}\\})=\lambda_{i+1}\geq\lambda_{i+1}/2$
and the statement follows. Consider $p\in P_{i+1}\setminus\\{p_{i+1}\\}$. If
$p_{i+1}$ is not among the $k$ nearest neighbors of $p$ in $P_{i+1}$, then

$\operatorname{d}^{k}(p,P_{i+1}\setminus\\{p\\})=\operatorname{d}^{k}(p,P_{i}\setminus\\{p\\})\geq\dfrac{\lambda_{i}}{2}\geq\dfrac{\lambda_{i+1}}{2}$

by the induction hypothesis and because the $\lambda$-values are non-
increasing. Otherwise, $p_{i+1}$ is among the $k$ nearest neighbors of $p$ in
$P_{i+1}\setminus\\{p\\}$ and
$\operatorname{d}^{k}(p,P_{i+1}\setminus\\{p\\})\geq\operatorname{d}(p,p_{i+1}).$

If $\operatorname{d}(p,p_{i+1})\geq\lambda_{i+1}/2$, the claim follows.
Otherwise, every site at distance smaller than $\lambda_{i+1}/2$ of $p$ is at
distance smaller than $\lambda_{i+1}$ of $p_{i+1}$. Since
$\lambda_{i+1}=\operatorname{d}^{k}(p_{i+1},P_{i})$, there can be at most
$k-2$ sites of $P_{i}\setminus\\{p\\}$ at distance smaller than
$\lambda_{i+1}$ of $p_{i+1}$. Thus, counting $p_{i+1}$ as well, there can be
at most $k-1$ sites of $P_{i+1}\setminus\\{p\\}$ at distance smaller than
$\lambda_{i+1}/2$ of $p$ and it follows that
$\operatorname{d}^{k}(p,P_{i+1}\setminus\\{p\\})\geq\lambda_{i+1}/2$.

### Computation.

We give a simple algorithm for computing the $k$-distance permutation that has
quadratic running time in the number of input points and discuss an approach
for improving it. We call a site _ordered_ if it has already been assigned its
index in the $k$-distance permutation and _unordered_ otherwise.

The simple approach is the following. Pick $k$ sites $p_{1},\ldots,p_{k}$ and
compute, for each $y\in P\setminus P_{k}$, the distances from $y$ to $p_{i}$,
$1\leq i\leq k$. Store them in a max-heap $T_{y}$ that also has a fixed entry
identifying $y$. Up until this point we need $\mathcal{O}(nk)$ time. The next
steps are repeated iteratively. For all unordered $y$, group the $T_{y}$ in a
list $L$. When $p_{1},\ldots,p_{i-1}$ are chosen, the algorithm picks $p_{i}$
by scanning over $L$ and choosing the point with largest distance to its $k$th
nearest ordered neighbor, which takes $\mathcal{O}(n)$ time. When $p_{i}$ is
picked and becomes ordered, remove its entry from $L$. Then, by traversing all
remaining elements in $L$, identify each unordered $y$ with $p_{i}$ among
$y$’s $k$ nearest neighbors in $P_{i}$ and insert $\operatorname{d}(p_{i},y)$
to $T_{y}$. The $(k+1)$-distance from $y$ to the ordered sites, which was a
previous entry in $T_{y}$, is removed. This takes $\mathcal{O}(\log k)$ time
per element of $L$ and hence $\mathcal{O}(n\log k)$ per iteration. Since there
are $n$ iterations, the total running time is of $\mathcal{O}(n^{2}\log k)$.

This simple algorithm can be improved with the main insight that when $p_{i}$
is determined, the $k$th nearest ordered neighbor of all remaining unordered
sites is at most $\lambda_{i}$ away. Hence, unordered sites further than
$\lambda_{i}$ away from $p_{i}$ do not have to be updated. Whenever a site
$p_{i}$ is ordered, we can employ a quadtreap (Section 2) to only update the
unordered sites within distance $\lambda_{i}$. This last step can also be done
in general metric spaces with elementary but rather tedious techniques; see
[30, Sec. 3.1]. Using the packing property, the total number of updates
reduces to $\mathcal{O}(nk\log\Phi)$ (with a constant that depends
exponentially on the doubling dimension of the point set). Finally, we replace
the list $L$ by a max-heap to avoid the linear scan to search for the next
ordered points. Appendix A provides further details on how to achieve this
improvement, which results in the next theorem.

######  Theorem 3.5.

The $k$-distance permutation can be computed in expected time
$\mathcal{O}(nk\log n\log\Phi)$ with $\Phi$ the spread of the point set.

##  4 A sparse union of lenses

We define several spaces in this section and the following. Figure 2 has an
overview.

Figure 2: Schematical view of the different filtrations introduced in Sections
4 and 5. We consider $r>(1+\epsilon)$ and $r\in($.

Recall that the $k$-fold cover is defined as the union of all lenses at radius
$r$, where every lens is given by $k$ sites. For large values of $r$, most of
these lenses intersect, yielding a size explosion in its nerve, the $k$th
order Čech complex. At the same time, many lenses are eventually covered by
the union of other lenses and so may be removed from consideration.

To define the precise threshold for removal of a lens, recall that in Section
3 we ordered the sites as $p_{1},\ldots,p_{n}$ and obtained values
$\lambda_{1},\lambda_{2},\ldots,\lambda_{k}=\infty,\,\lambda_{k+1}\geq\ldots\geq\lambda_{n}$.
Fix $\epsilon\in(0,1]$. Since it is fixed, we drop $\epsilon$ from the
upcoming notation. The _freezing radius_ of a site $p_{i}$ is

We extend the definition to lenses by setting Then at radius $r$ we only
consider lenses whose freezing radius is at least $r$ and set

$U_{r}:=\bigcup L_{r}(A).$

Notice that $\mathcal{U}:=\\{U_{r}\\}_{r\geq 0}$ is _not_ a filtration: Figure
3 illustrates that $U_{r}$ might not be a subset of $U_{r^{\prime}}$ for
$r<r^{\prime}$. Even so, some useful inclusions hold as we see on the next
lemma.

$c$$b$$a$

$c$$a$$b$

Figure 3: Example in $\mathbb{R}^{2}$ with $k=2$. Left: $U_{r}$ at radius $r=$
Right: $U_{r^{\prime}}$ at radius immediately after . Even though
$r<r^{\prime},$ $U_{r}\nsubseteq U_{r^{\prime}}.$

######  Lemma 4.1.

$U_{r}\subseteq L_{r}\subseteq U_{(1+\epsilon)r}.$

######  Proof 4.2.

The first inclusion is clear. For the second inclusion, consider $x\in L_{r}$
and let $i$ be the maximal index such that

$r\leq\frac{\lambda_{i}}{\epsilon}.$

If $i=n$, then there is $A\subseteq P_{n}$ with $x\in L_{r}(A)$ because $x\in
L_{r}$ and $P=P_{n}$. By definition of the freezing radius and inequality
$(\ast)$, and thus $L_{r(1+\epsilon)}(A)\subseteq U_{(1+\epsilon)r}$. Since
$L_{r}(A)\subseteq L_{(1+\epsilon)r}(A)$, the result follows.

For $i<n$, notice that the Covering Property (Lemma 3.1) guarantees that $x$
is contained in a lens $L_{r+\lambda_{i+1}}(A)$ for some $A\subseteq P_{i}$.
Since $i$ is maximal, $\lambda_{i+1}/\epsilon<r$ and so
$L_{r+\lambda_{i+1}}(A)\subseteq L_{(1+\epsilon)r}(A)$. Moreover, $A\subseteq
P_{i}$ and inequality $(\ast)$ imply . Hence the lens of $A$ contributes to
$U_{(1+\epsilon)r}$ and as it contains $x$, the statement follows.

This lemma suggests that a $(1+\epsilon)$-interleaving with the filtration
$\mathcal{L}=\\{L_{r}\\}_{r\geq 0}$ consisting of all lenses should be
possible if we adjust $\mathcal{U}$ to obtain an actual filtration. To do that
we slightly delay the removal of lenses. More precisely, we define

$\tilde{L}_{r}(A):=\begin{cases}L_{r}(A)&r<\crit[A]\\\\[6.45831pt]
L(A)&\crit[A]\leq r\leq(1+\epsilon)\crit[A]\\\\[6.45831pt]
\emptyset&(1+\epsilon)\crit[A]<r.\end{cases}$

One can visualize the evolution of a lens as a continuous process for
increasing $r$: the lens $\tilde{L}_{r}$ grows until it reaches its freezing
radius and remains unchanged (it is “frozen”) for the interval $[$. Afterwards
it completely disappears. We call $(1+\epsilon)$ the _removal radius_ of $A$.
The construction is an adaptation of a similar one by Sheehy [42].

We write $\tilde{L}_{r}$ for the union of $\tilde{L}_{r}(A)$ over all
$A\in\binom{P}{k}$ and $\tilde{\mathcal{L}}:=\\{\tilde{L}_{r}\\}_{r\geq 0}$.
We show next that $\tilde{\mathcal{L}}$ is a filtration.

######  Lemma 4.3.

$\tilde{\mathcal{L}}$ is a filtration, i.e., for any $r\leq r^{\prime},$
$\tilde{L}_{r}\subseteq\tilde{L}_{r^{\prime}}$.

######  Proof 4.4.

If the interval $(r,r^{\prime}]$ does not contain any removal radius,
$\tilde{L}_{r}\subseteq\tilde{L}_{r^{\prime}}$ because the inclusions hold
lens-wise. Since the number of different removal radii is bounded by the
number of sites and hence finite, it suffices to show that at a removal radius
$s$, any lens that is removed is already covered by lenses that are not being
removed at $s$. In fact, we show that such a lens is covered by lenses that
are not yet frozen at $s$.

Let $A$ be the $k$-subset associated with a lens being removed at $s$ and
$x\in\tilde{L}_{s}(A)$. By definition, $s=(1+\epsilon)t$ with $t=$. This
implies that $x\in\tilde{L}_{t}(A)=L_{t}(A)\subseteq L_{t}$ because the lens
is frozen from radius $t$ on. By Lemma 4.1, it follows that $x\in
U_{(1+\epsilon)t}=U_{s}$, and therefore $x$ is contained in a lens $L_{s}(B)$
with . Thus $\tilde{L}_{s}(A)$ is covered by lenses $L_{s}(B)$ with and
$\tilde{\mathcal{L}}$ is a filtration.

######  Lemma 4.5.

$\tilde{\mathcal{L}}$ and $\mathcal{L}$ are $(1+\epsilon)$-interleaved.

######  Proof 4.6.

We show that $\tilde{L}_{r}\subseteq L_{r}\subseteq\tilde{L}_{(1+\epsilon)r}$.
Note that by definition, $\tilde{L}_{r}\subseteq L_{r}$. For the second
inclusion, observe that $U_{r}\subseteq\tilde{L}_{r}$ follows directly from
their definition. Then, Lemma 4.1 yields $L_{r}\subseteq
U_{(1+\epsilon)r}\subseteq\tilde{L}_{(1+\epsilon)r}$.

##  5 A sparse simplicial filtration

Since $\tilde{L}_{r}(A)$ is closed and convex for every $r$, the nerve of all
(non-empty) $\tilde{L}_{r}(A)$ yields a simplicial complex with the same
homotopy type as the $k$-fold cover at radius $r$. However, the collection of
simplicial complexes obtained when varying $r$ does not form a filtration
because simplices disappear from the nerve when passing a removal radius. To
overcome this problem, we adapt a construction of Cavanna et al. [9] that is
similar to a function’s graph.

### Cones.

The idea is to “stack-up” the lenses $\tilde{L}_{r}(A)$ for all radii: the
_cone_ of $A$ at radius $r$ is

$C_{r}(A):=\underset{\alpha\in[0,r]}{\bigcup}\left(\tilde{L}_{\alpha}(A)\times\\{\alpha\\}\right)\subseteq\mathbb{R}^{d}\times\mathbb{R}.$

We write $C_{r}$ for the union of $C_{r}(A)$ over all $A\in\binom{P}{k}$.
Figure 4 shows one cone.

Figure 4: A cone representing the evolution of a lens for $k=2.$ At first the
lens grows, until it is frozen. Then, the lens has static size and afterwards
it disappears. At the radius where the lens disappears, it is completely
covered by other lenses (which are not displayed in the figure).

######  Lemma 5.1.

The filtrations $\mathcal{C}=\\{C_{r}\\}_{r\geq 0}$ and $\tilde{\mathcal{L}}$
are equivalent.

######  Proof 5.2.

In the same sense as above $C_{r}$ is a stacked-up version of
$\tilde{L}_{\alpha}$ for all $\alpha\leq r$, and we can consider
$\tilde{L}_{r}$ as a subspace of $C_{r}$ via the map $x\mapsto(x,r)$ for
$x\in\tilde{L}_{r}$. Since $\tilde{\mathcal{L}}$ is a filtration, there is a
strong deformation retraction $R$ from $C_{r}$ to $\tilde{L}_{r}$, given by
$R((x,\alpha),t)=(x,(1-t)\alpha+tr)$, which naturally commutes with the
canonical inclusions. The result follows.

We define the nerve of the cones as the _sparse $k$th order Čech complex_,

$S_{r}:=\operatorname{Nrv}\left\\{C_{r}(A)\mid
A\in\left(\genfrac{}{}{0.0pt}{}{P}{k}\right)\right\\}.$

$S_{r}$ is a subcomplex of the $k$th order Čech complex for every $r$ because
$\tilde{L}_{r}\subseteq L_{r}$. Moreover, if $r\leq r^{\prime}$,
$C_{r}(A)\subseteq C_{r^{\prime}}(A)$ for all $A$. Hence
$\mathcal{S}=\\{S_{r}\\}_{r\geq 0}$ is a filtration. By the Persistent Nerve
Theorem [1, Thm. 3.9] and Lemmas 4.5 and 5.1 we obtain:

######  Lemma 5.3.

The filtrations $\mathcal{S}$ and $\mathcal{C}$ are equivalent. As a
consequence, $\mathcal{S}$ is a $(1+\epsilon)$-approximation of the $k$-fold
filtration.

### Discretization of the radius.

The filtration $\mathcal{S}$ is challenging to compute, due to the freezing of
lenses. We elaborate on these issues in Section 7. We now define a variant of
$\mathcal{S}$ which is easier to compute and also is
$(1+\epsilon)$-interleaved with the $k$-fold filtration.

Recall that the filtrations $\tilde{\mathcal{L}}$, $\mathcal{C}$ and
$\mathcal{S}$ are defined based on the freezing radii of sites, which depend
on a parameter $\epsilon>0$. To obtain a $(1+\epsilon)$-approximation for
$\epsilon\in(0,1]$ in the end, we consider the above construction of
$\mathcal{S}$ with parameter $\epsilon^{\prime}=\frac{\epsilon}{3}$, obtaining
a $(1+\frac{\epsilon}{3})$-approximation of the $k$-fold filtration.

Next, for every $r\geq 0$, let $z\in\mathbb{Z}$ be such that
$(1+\frac{\epsilon}{3})^{z}\leq r<(1+\frac{\epsilon}{3})^{z+1}$ and define

$D_{r}:=S_{(1+\epsilon/3)^{z}}.$

We call $\mathcal{D}:=\\{D_{r}\\}_{r\geq 0}$ the _discrete sparse $k$th order
Čech filtration_. It is formed by a discrete set of snapshots of $\mathcal{S}$
and kept unchanged except when passing over a snapshot radius (this is also
referred to as the _Left Kan extension_ of a discrete filtration [33, Chap.
10]).

######  Theorem 5.4.

$\mathcal{D}$ is a $(1+\epsilon)$-approximation of the $k$-fold filtration.

######  Proof 5.5.

From the definition, $D_{r}\subseteq S_{r}\subseteq D_{(1+\epsilon/3)r}$. This
interleaving implies that $\mathcal{D}$ is a
$(1+\frac{\epsilon}{3})$-approximation of $\mathcal{S}$. Since $\mathcal{S}$
is a $(1+\frac{\epsilon}{3})$-approximation of the $k$-fold filtration, by
transitivity, we get that $D_{r}$ is a
$(1+\frac{\epsilon}{3})^{2}$-approximation of the $k$-fold filtration. The
result follows by noting that $(1+\frac{\epsilon}{3})^{2}\leq 1+\epsilon$ for
all $\epsilon\in(0,1]$.

##  6 Size analysis

We bound the size of $\mathcal{D}$, i.e., the number of simplices it contains.
Since $D_{r}\subseteq S_{r}$ for all $r\geq 0,$ it is enough to bound the size
of $\mathcal{S}$ with parameter $\epsilon^{\prime}=\epsilon/3.$ Let

$C_{\infty}(A):=\bigcup_{r\geq 0}C_{r}(A)$

be the _cone of $A$_ (without dependence on a radius). Then, the size of
$\mathcal{S}$ equals the size of the nerve of the cones $C_{\infty}(A)$, where
$A$ ranges over all $k$-subsets of sites. However, the number of vertices is
not necessarily $\binom{n}{k}$ because many cones are empty: this happens in
particular when the smallest radius for which the balls around the sites of
$A$ intersect is larger than the removal radius of the lens. In fact, our
argument shows that this is the case for the vast majority of cones. The proof
for vertices extends readily to the case of $m$-wise intersections of cones,
i.e. for $(m-1)$-simplices, without change and thus we treat the general case
directly.

For fixed $m\geq 1,$ we derive an upper bound for the number of sets
$\\{A_{1},\ldots,A_{m}\\}$ such that the cones
$C_{\infty}(A_{1}),\ldots,C_{\infty}(A_{m})$ intersect. Such sets are in one-
to-one correspondence to the $(m-1)$-simplices of the sparse $k$th order Čech
filtration, hence we refer to these sets as $(m-1)$-simplices. Let
$\sigma=\\{C_{\infty}(A_{1}),\ldots,C_{\infty}(A_{m})\\}$ be a $(m-1)$-simplex
and the set of sites $P$ be ordered according to the $k$-distance permutation.
We call a site $p_{i}$ _involved_ in $\sigma$ if $p_{i}$ belongs to one of the
sets $A_{1},\ldots,A_{m}$. Note that there are at most $km$ sites involved in
$\sigma$. We say that $\sigma$ is _associated_ to a site $p_{i}$ if $p_{i}$ is
involved in $\sigma$ and all other involved sites have index smaller than $i$.
Our strategy is to upper bound the number of $(m-1)$-simplices associated to
an arbitrary $p_{i}$. We only need to consider simplices associated to $p_{i}$
that appear in the filtration, i.e., simplices whose defining cones intersect.

Fix $p_{i}$ and $\omega_{i}:=(1+\epsilon^{\prime})$. Let $B$ denote the ball
of radius $2\omega_{i}$ centered at $p_{i}$.

######  Lemma 6.1.

If $\sigma:=\\{C_{\infty}(A_{1}),\ldots,C_{\infty}(A_{m})\\}$ is an
$(m-1)$-simplex associated to $p_{i}$ whose cones intersect, then all sites
involved in $\sigma$ are contained in $B$.

######  Proof 6.2.

Let $\alpha$ denote the minimal radius such that all the balls around sites
involved in $\sigma$ intersect. This is the radius of the minimum enclosing
ball of the involved sites. Any common intersection of the cones must happen
at scale $r\geq\alpha.$

On the other hand, assume wlog that $p_{i}\in A_{1}$. Since $p_{i}$ has
maximal index in $A_{1}$, we have that . Hence the removal radius of $A_{1}$
is equal to $\omega_{i}$ and it follows that the cone of $A_{1}$ is empty for
all radii greater than $\omega_{i}$. Therefore any common intersection of the
cones of $\sigma$ must happen at scale $r\leq\omega_{i}$. Hence, as we assume
that the cones do intersect, we must have that $\alpha\leq\omega_{i}$.

Then, since $\operatorname{d}(q,p_{i})\leq 2\alpha$, any involved site $q$
lies within distance $2\omega_{i}$ from $p_{i}$.

Hence the involved sites of $\sigma$ are close to $p_{i}$ in the sense of the
lemma. We can furthermore guarantee that the points of $P_{i}$ are not too
densely packed in $B$ using the packing property of the $k$-distance
permutation.

######  Lemma 6.3.

The ball $B$ contains at most $\hypertarget{link_sites-
bound}{\Gamma}:=k\left(\dfrac{96}{\epsilon}\right)^{\delta}$ sites of
$P_{i}=\\{p_{1},\ldots,p_{i}\\}$, where $\delta$ is the doubling dimension of
$\mathbb{R}^{d}$.

######  Proof 6.4.

We cover $B$ by balls of radius ${\lambda_{i}}/{4}$. That can be done with at
most
$\zeta=\left(\dfrac{16(1+\epsilon^{\prime})^{2}}{\epsilon^{\prime}}\right)^{\delta}$
balls (see _Doubling Dimension_ in Section 2). By the Packing Lemma 3.3, each
open ball of radius $\lambda_{i}/4$ contains at most $k$ sites of $P_{i}$,
thus the total number of sites in $B$ is at most $k\zeta$. The bound follows
because $\epsilon^{\prime}=\frac{\epsilon}{3}$ and $\epsilon\leq 1$, hence
$(1+\epsilon^{\prime})^{2}\leq\frac{16}{9}<2$.

Bounding the number of non-empty $m$-intersections of cones is now a matter of
simple combinatorics. Recall the $\Gamma$ notation from Lemma 6.3.

######  Theorem 6.5.

The number of $(m-1)$-simplices of the sparse $k$th order Čech complex with
non-empty cone intersection is at most

$n\cdot\Gamma^{km}=nk^{km}\left(\dfrac{96}{\epsilon}\right)^{\delta km}.$

######  Proof 6.6.

Fix $p_{i}$. Every $(m-1)$-simplex associated to $p_{i}$ with non-empty cone
intersection has up to $km$ involved sites, which all lie in $B$ by Lemma 6.1.
Moreover, all involved sites are in $P_{i}$ and, by Lemma 6.3, there are at
most $\Gamma$ of those sites in $B$ to choose from. It follows that there are
at most $\Gamma^{km}$ different choices possible. This upper bound holds for
every $p_{i}$, so multiplying by the number of sites $n$ yields the result.

We remark that the bounds on this section are not tight and slightly better
ones could be easily achieved, by keeping binomials in place or avoiding some
approximations. However the improvements would be minor.

##  7 Computation

We now present an algorithm to construct the discrete sparse Čech filtration.
As in the previous sections, let us fix a finite set
$P\subseteq\mathbb{R}^{d}$, an integer $k>0$ and $\epsilon\in(0,1]$. Assume
that $P=\\{p_{1},\ldots,p_{n}\\}$ has the indices ordered with respect to the
$k$-distance permutation, and that we have computed the corresponding values
$\lambda_{1},\ldots,\lambda_{n}$ as discussed in Section 3. The algorithm
outputs the discrete sparse $k$th order Čech filtration $\\{D_{r}\\}_{r\geq
0}$ as a list of simplices with their corresponding _critical value_ , i.e.,
the smallest parameter value $r$ for which the simplex is part of the
filtration. Note that by definition of the discrete sparse Čech filtration,
every critical value is of the form $(1+\epsilon/3)^{z}$ for some integer $z$.

### Friends.

Our algorithm follows the approach and notation of Section 6. For every
$p_{i}$, we compute all simplices associated to $p_{i}$ in the filtration
together with their critical value. To do so, we first find, among
$p_{1},\ldots,p_{i-1}$, all sites of distance at most $2\omega_{i}$ from
$p_{i}$, where $\omega_{i}=(1+\epsilon^{\prime})$ (compare Lemma 6.1). We call
these points _friends_ of $p_{i}$. We compute friends using a quadtreap data
structure, as introduced in Section 3, which we query for every $p_{i}$ at
$2\omega_{i}$. $p_{i}$ is added to the quadtreap after the $i$th iteration
(adding an element costs $\mathcal{O}(\log n)$ in expectation as well). Hence
the expected running time for this loop is $\mathcal{O}(n\log n+\Sigma)$,
where $\Sigma$ is the number of reported points. These reported points have
distance at most $4\omega_{i}$ from the respective $p_{i}$ (since the queries
are approximate), and by the same argument as in Lemma 6.3, the number of
sites reported for $p_{i}$ is at most
$k\left(\dfrac{192}{\epsilon}\right)^{\delta}=\mathcal{O}(k(1/\epsilon)^{\delta})$.
We traverse the list and remove all “false friends” of distance more than
$2\omega_{i}$. Thus we get the friends of $p_{i}$ for all sites $p_{i}$ in
expected time

$\mathcal{O}(n\log n+nk(1/\epsilon)^{\delta}).$ (1)

Note that the number of friends is bounded by Γ as defined in Lemma 6.3.

### Vertices.

Next, we compute the vertices of the filtration associated to $p_{i}$, for
each $i\geq k$. We proceed by brute-force, just enumerating all $k$-tuples
formed by $p_{i}$ and $k-1$ of its friends and checking for every $k$-tuple
whether their cone is non-empty. The last condition is simple to check, as the
cone is non-empty if and only if the radius $\alpha$ of the minimum enclosing
ball of the $k$ sites is at most . In this case, the critical radius of the
vertex is set to $(1+\epsilon/3)^{z}$, where $z$ is the smallest integer such
that $(1+\epsilon/3)^{z}\geq\alpha$. Computing $\alpha$ and $z$ per vertex
requires $\mathcal{O}(k)$ in expectation [16]. Hence we calculate all vertices
in $\mathcal{O}((n-k)k\Gamma^{k-1})$ time.

### Simplices.

For higher-dimensional simplices associated to $p_{i}$, we proceed inductively
by dimension, up to a maximal dimension $m_{max}$. Fix a $(m-1)$-simplex
$\sigma=\\{A_{1},\ldots,A_{m}\\}$ associated to $p_{i}$. We compute all
cofacets of $\sigma$ in the filtration, that is, all $m$-simplices that
contain $\sigma$ and one further vertex $A_{m+1}$. Notice that one could order
the $k$-subsets and avoid computing all cofacets, computing instead only
cofacets with larger index in the ordering. This would remove a $k$ factor
from the computation expected time locally, but does not change the final
bound in Theorem 7.1, which has a $k^{k}$ factor. We compute all cofacets for
simplicity. Since each element of $A_{m+1}$ is either a friend of $p_{i}$ or
$p_{i}$ itself, we enumerate all $k$-tuples consisting of these sites and
check whether they form a vertex of the filtration. This takes
$\mathcal{O}(k)$ time per vertex, just by re-doing the check from the previous
step, except that $k$-tuples associated to $p_{j}$ must be checked at radius
$\min\\{$. That is because if the $k$-tuple cone becomes non empty only after
$\omega_{i}$, then it cannot contribute to a coface of $\sigma$.

For a vertex $A_{m+1}$ of the filtration, check whether the cones of
$A_{1},\ldots,A_{m+1}$ intersect is technically challenging because the cones
might intersect for a radius where one or several cones are frozen. One cannot
resolve this question by a simple minimal enclosing ball computation. In fact,
we are not aware of an efficient way to compute the smallest intersection
radius of such cones in general. However, as demonstrated in Appendix B, given
a collection of cones $A_{1},\ldots,A_{m+1}$ and a fixed radius $r$, we can
decide whether the cones (or rather, the corresponding lenses $\tilde{L}_{r}$)
intersect at radius $r$ by a reduction to a minimum enclosing ball of balls
instance [26] in expected $\mathcal{O}(k(m+1))$ time. We use this predicate
and query whether the cones intersect at the smallest removal radius of
$A_{1},\ldots,A_{m+1}$. That decides whether the $m$-simplex is in the
filtration, and its critical radius can be computed by two more minimal
enclosing ball computations. This is only possible because we have discretized
the filtration, as discussed in Appendix B. It follows that the expected
running time spent per $(m-1)$-simplex of the filtration is
$\mathcal{O}(mk\Gamma^{k})$, and doing this over all simplices of the
filtration up to dimension $m_{max}$ yields a total expected complexity of
$\mathcal{O}(X\cdot m_{max}\cdot k\cdot\Gamma^{k})$ for this step, where $X$
is the total number of simplices in the filtration.

This concludes the description of the algorithm. Recalling (1) we obtain a
complexity of

$\mathcal{O}\left(n\log
n+\frac{nk}{\epsilon^{\delta}}+(n-k)k\Gamma^{k-1}+X\cdot m_{max}\cdot
k\cdot\Gamma^{k}\right),$ (2)

where the second and third terms are dominated by the last one because $X\geq
n-k$ (see Appendix C). Together with the algorithm from Section 3, we arrive
at the result

######  Theorem 7.1.

Given a set $P$ of $n$ points, $k\geq 0$ and $\epsilon\in(0,1]$, the discrete
sparse $k$th order Čech filtration up to dimension $m_{max}$ can be computed
in time

$\mathcal{O}\left(nk\log
n\log\Phi+Xk^{k+1}m_{max}\left(\dfrac{96}{\epsilon}\right)^{k\delta}\right),$

where $\Phi$ is the spread of $P$ and $X$ is the total number of simplices in
the filtration.

##  8 Conclusions

We introduced the first $(1+\epsilon)$-approximate filtration of the higher
order $k$-fold filtration and provided an algorithm for computing it. If $k$
and $\epsilon$ are considered as constants and the input point set has
constant spread, the algorithm runs in time $\mathcal{O}(n\log n)$ and yields
a filtration of size $\mathcal{O}(n)$, which are the same favorable properties
of the well-studied case $k=1$.

There are various avenues to strengthen and generalize our results. First of
all, our method has concentrated on the Euclidean case, but our approach
mostly generalizes to point sets in arbitrary metric spaces – the algorithm
cannot use the quadtreap data structure anymore, but there is no need for it,
since the algorithm by Har-Peled and Mendel [30, Sec. 3.1] can be adapted to
the $k$-distance case with little effort. Also, the friends of $p_{i}$
(Section 7) can be computed with a slight adaptation of their techniques; we
used quadtreaps mostly for the ease of presentation. However, the computation
of critical values of simplices described in Appendix B is for the Euclidean
case only, and the complexity of this step remains unspecified for a general
metric space. This is common in related work; see, for instance [9, Sec. 5].

Another natural goal is to remove the dependence on the spread. This
dependence is caused by the computation of the $k$-distance permutation which
is inspired by the algorithm of [30, Sec. 3.1]. In the same paper [30, Sec.
3.2–3.3], they describe an approach to remove the spread from the bound (for
$k=1$) using an approximate version of the greedy permutation. While our
construction of the sparsified filtration can be easily adapted to work with
an approximate version of the $k$-distance permutation, it seems less
straight-forward to generalize the computation to the $k$-distance, even in
the Euclidean case. We leave this for future work.

While we concentrate on the case of a single value of $k$, we pose the
question whether our methods can be used to approximate the multi-cover
bifiltration, as studied in [15]. The extension is not straight-forward
because there is no direct relation between the approximate filtrations on
level $k$ and $k+1$. We speculate that the technique of double-nerve
constructions of [15] could be useful in this context. The presence of an
exponential factor on $k$ in our size bounds suggests a restriction of our
approach’s usability to small portions of the bifiltration, for small $k$. The
exponential factor on $k$ also carries over to the expected computation time.
Reducing that dependency on $k$ is another possible line of future work. Note
that [15, Prop. 5] gives a size bound of $\mathcal{O}\left(n^{d+1}\right)$ for
the exact version, but we ask whether a polynomial bound on $k$ could be
achieved without such a blow-up in the dependency on $n$.

Finally, a natural question is the practicality of our algorithm. We remark
that even for $k=1$, while some work has been devoted to practical aspects of
computing sparsifications [2, 6, 18, 34, 38], the actual practical computation
is still an unresolved problem. We think that the natural order for a
practically efficient solution would be to first identify best practices in
the simpler $k=1$ case and subsequently try to adapt them to larger values of
$k$. So, while we would be curious about the performance of our algorithm,
such an evaluation seems to be premature at the moment.

## References

  * [1] Ulrich Bauer, Michael Kerber, Fabian Roll, and Alexander Rolle.  A unified view on the functorial nerve theorem and its variations.  arXiv preprint, 2022.  arXiv:2203.03571.
  * [2] Nello Blaser and Morten Brun.  Sparse nerves in practice.  In Andreas Holzinger, Peter Kieseberg, A M. Tjoa, and Edgar R. Weippl, editors, Machine Learning and Knowledge Extraction - CD-MAKE 2019, Canterbury, UK, August 26-29, 2019, Proceedings, volume 11713 of Lecture Notes in Computer Science, pages 272–284. Springer, 2019.  doi:10.1007/978-3-030-29726-8\\_17.
  * [3] Andrew J. Blumberg and Michael Lesnick.  Stability of 2-parameter persistent homology.  CoRR, 2020.  arXiv:2010.09628.
  * [4] Thomas Bonis, Maks Ovsjanikov, Steve Oudot, and Frédéric Chazal.  Persistence-based pooling for shape pose recognition.  In Alexandra Bac and Jean-Luc Mari, editors, Computational Topology in Image Context - 6th International Workshop, CTIC 2016, Marseille, France, June 15-17, 2016, Proceedings, volume 9667 of Lecture Notes in Computer Science, pages 19–29. Springer, 2016.  doi:10.1007/978-3-319-39441-1\\_3.
  * [5] Magnus B. Botnan and Gard Spreemann.  Approximating persistent homology in euclidean space through collapses.  Applicable Algebra in Engineering, Communication and Computing, 26(1-2):73–101, Jan 2015.  doi:10.1007/s00200-014-0247-y.
  * [6] Bernhard Brehm and Hanne Hardering.  Sparips.  CoRR, 2018.  arXiv:1807.09982.
  * [7] Mickaël Buchet, Frédéric Chazal, Steve Y. Oudot, and Donald R. Sheehy.  Efficient and robust persistent homology for measures.  Comput. Geom., 58:70–96, 2016.  doi:10.1016/j.comgeo.2016.07.001.
  * [8] Gunnar Carlsson.  Topology and data.  Bulletin of the American Mathematical Society, 46(2):255–308, 2009\.  doi:10.1090/S0273-0979-09-01249-X.
  * [9] Nicholas J. Cavanna, Mahmoodreza Jahanseir, and Donald R. Sheehy.  A geometric perspective on sparse filtrations.  In Proceedings of the 27th Canadian Conference on Computational Geometry, CCCG 2015, Kingston, Ontario, Canada, August 10-12, 2015. Queen’s University, Ontario, Canada, 2015.  URL: http://research.cs.queensu.ca/cccg2015/CCCG15-papers/01.pdf.
  * [10] Frédéric Chazal, David Cohen-Steiner, Marc Glisse, Leonidas J. Guibas, and Steve Oudot.  Proximity of persistence modules and their diagrams.  In John Hershberger and Efi Fogel, editors, Proceedings of the 25th ACM Symposium on Computational Geometry, Aarhus, Denmark, June 8-10, 2009, pages 237–246. ACM, 2009.  doi:10.1145/1542362.1542407.
  * [11] Frédéric Chazal, David Cohen-Steiner, and Quentin Mérigot.  Geometric inference for probability measures.  Found. Comput. Math., 11(6):733–751, 2011.  doi:10.1007/s10208-011-9098-0.
  * [12] Aruni Choudhary, Michael Kerber, and Sharath Raghvendra.  Improved topological approximations by digitization.  In Timothy M. Chan, editor, Proceedings of the Thirtieth Annual ACM-SIAM Symposium on Discrete Algorithms, SODA 2019, San Diego, California, USA, January 6-9, 2019, pages 2675–2688. SIAM, 2019.  doi:10.1137/1.9781611975482.166.
  * [13] Aruni Choudhary, Michael Kerber, and Sharath Raghvendra.  Polynomial-sized topological approximations using the permutahedron.  Discret. Comput. Geom., 61(1):42–80, 2019.  doi:10.1007/s00454-017-9951-2.
  * [14] Aruni Choudhary, Michael Kerber, and Sharath Raghvendra.  Improved approximate rips filtrations with shifted integer lattices and cubical complexes.  J. Appl. Comput. Topol., 5(3):425–458, 2021.  doi:10.1007/s41468-021-00072-4.
  * [15] René Corbet, Michael Kerber, Michael Lesnick, and Georg Osang.  Computing the multicover bifiltration.  Discret. Comput. Geom., 2023.  doi:10.1007/s00454-022-00476-8.
  * [16] Mark de Berg, Otfried Cheong, Marc J. van Kreveld, and Mark H. Overmars.  Computational geometry: algorithms and applications, 3rd Edition.  Springer, 2008.  URL: https://www.worldcat.org/oclc/227584184.
  * [17] Tamal K. Dey, Fengtao Fan, and Yusu Wang.  Computing topological persistence for simplicial maps.  In Siu-Wing Cheng and Olivier Devillers, editors, 30th International Symposium on Computational Geometry, SoCG 2014, Kyoto, Japan, June 08 - 11, 2014, page 345. ACM, 2014.  doi:10.1145/2582112.2582165.
  * [18] Tamal K. Dey, Dayu Shi, and Yusu Wang.  Simba: An efficient tool for approximating rips-filtration persistence via simplicial batch collapse.  ACM J. Exp. Algorithmics, 24(1):1.5:1–1.5:16, 2019.  doi:10.1145/3284360.
  * [19] Martin E. Dyer.  A class of convex programs with applications to computational geometry.  In David Avis, editor, Proceedings of the Eighth Annual Symposium on Computational Geometry, Berlin, Germany, June 10-12, 1992, pages 9–15. ACM, 1992.  doi:10.1145/142675.142681.
  * [20] Herbert Edelsbrunner and John L. Harer.  Computational topology.  American Mathematical Society, Providence, RI, 2010.  An introduction.  doi:10.1090/mbk/069.
  * [21] Herbert Edelsbrunner, David Letscher, and Afra Zomorodian.  Topological persistence and simplification.  Discret. Comput. Geom., 28(4):511–533, 2002.  doi:10.1007/s00454-002-2885-2.
  * [22] Herbert Edelsbrunner and Georg Osang.  A simple algorithm for higher-order delaunay mosaics and alpha shapes.  CoRR, 2020.  arXiv:2011.03617.
  * [23] Herbert Edelsbrunner and Georg Osang.  The multi-cover persistence of euclidean balls.  Discret. Comput. Geom., 65(4):1296–1313, 2021.  doi:10.1007/s00454-021-00281-9.
  * [24] Herbert Edelsbrunner and Raimund Seidel.  Voronoi diagrams and arrangements.  Discret. Comput. Geom., 1:25–44, 1986.  doi:10.1007/BF02187681.
  * [25] Kaspar Fischer.  Smallest enclosing balls of balls: Combinatorial structure & algorithms.  PhD thesis, Swiss Federal Institute of Technology, ETH Zürich, 2005.  URL: https://people.inf.ethz.ch/emo/DoctThesisFiles/fischer05.pdf.
  * [26] Kaspar Fischer and Bernd Gärtner.  The smallest enclosing ball of balls: combinatorial structure and algorithms.  Int. J. Comput. Geom. Appl., 14(4-5):341–378, 2004.  doi:10.1142/S0218195904001500.
  * [27] Marcio Gameiro, Yasuaki Hiraoka, Shunsuke Izumi, Miroslav Kramar, Konstantin Mischaikow, and Vidit Nanda.  A topological measurement of protein compressibility.  Japan Journal of Industrial and Applied Mathematics, 32(1):1–17, 2015.
  * [28] Robert W. Ghrist.  Elementary applied topology, volume 1.  Createspace Seattle, WA, 2014.
  * [29] Teofilo F. Gonzalez.  Clustering to minimize the maximum intercluster distance.  Theoretical Computer Science, 38:293–306, 1985.  URL: https://www.sciencedirect.com/science/article/pii/0304397585902245, doi:10.1016/0304-3975(85)90224-5.
  * [30] Sariel Har-Peled and Manor Mendel.  Fast construction of nets in low-dimensional metrics and their applications.  SIAM J. Comput., 35(5):1148–1184, 2006.  doi:10.1137/S0097539704446281.
  * [31] Yasuaki Hiraoka, Takenobu Nakamura, Akihiko Hirata, Emerson G. Escolar, Kaname Matsue, and Yasumasa Nishiura.  Hierarchical structures of amorphous solids characterized by persistent homology.  Proceedings of the National Academy of Sciences, 113(26):7035–7040, 2016.
  * [32] Lida Kanari, Pawel Dlotko, Martina Scolamiero, Ran Levi, Julian C. Shillcock, Kathryn Hess, and Henry Markram.  A topological representation of branching neuronal morphologies.  Neuroinformatics, 16(1):3–13, 2018.  doi:10.1007/s12021-017-9341-1.
  * [33] Saunders Mac Lane.  Categories for the working mathematician., volume 5.  New York, NY: Springer, 1998.
  * [34] Clément Maria, Pawel Dlotko, Vincent Rouvreau, and Marc Glisse.  Rips complex.  In GUDHI User and Reference Manual. GUDHI Editorial Board, 3.5.0 edition, 2022.  URL: https://gudhi.inria.fr/doc/3.5.0/group__rips__complex.html.
  * [35] Nimrod Megiddo.  On the ball spanned by balls.  Discret. Comput. Geom., 4:605–610, 1989.  doi:10.1007/BF02187750.
  * [36] David M. Mount and Eunhui Park.  A dynamic data structure for approximate range searching.  In David G. Kirkpatrick and Joseph S. B. Mitchell, editors, Proceedings of the 26th ACM Symposium on Computational Geometry, Snowbird, Utah, USA, June 13-16, 2010, pages 247–256. ACM, 2010.  doi:10.1145/1810959.1811002.
  * [37] James Munkres.  Topology.  Prentice Hall, 2ed edition, 2000.
  * [38] Julian B. Pérez, Sydney Hauke, Umberto Lupo, Matteo Caorsi, and Alberto Dassatti.  Giotto-ph: A python library for high-performance computation of persistent homology of vietoris-rips filtrations.  CoRR, 2021.  arXiv:2107.05412.
  * [39] Jeff M. Phillips, Bei Wang, and Yan Zheng.  Geometric inference on kernel density estimates.  In Lars Arge and János Pach, editors, 31st International Symposium on Computational Geometry, SoCG 2015, June 22-25, 2015, Eindhoven, The Netherlands, volume 34 of LIPIcs, pages 857–871. Schloss Dagstuhl \- Leibniz-Zentrum für Informatik, 2015.  doi:10.4230/LIPIcs.SOCG.2015.857.
  * [40] Raimund Seidel.  On the number of faces in higher-dimensional voronoi diagrams.  In D. Soule, editor, Proceedings of the Third Annual Symposium on Computational Geometry, Waterloo, Ontario, Canada, June 8-10, 1987, pages 181–185. ACM, 1987.  doi:10.1145/41958.41977.
  * [41] Donald R. Sheehy.  A multicover nerve for geometric inference.  In Proceedings of the 24th Canadian Conference on Computational Geometry, CCCG 2012, Charlottetown, Prince Edward Island, Canada, August 8-10, 2012, pages 309–314, 2012.  URL: http://2012.cccg.ca/papers/paper52.pdf.
  * [42] Donald R. Sheehy.  Linear-size approximations to the vietoris-rips filtration.  Discret. Comput. Geom., 49(4):778–796, 2013.  doi:10.1007/s00454-013-9513-1.
  * [43] Donald R. Sheehy.  A sparse delaunay filtration.  In Kevin Buchin and Éric Colin de Verdière, editors, 37th International Symposium on Computational Geometry, SoCG 2021, June 7-11, 2021, Buffalo, NY, USA (Virtual Conference), volume 189 of LIPIcs, pages 58:1–58:16. Schloss Dagstuhl - Leibniz-Zentrum für Informatik, 2021\.  doi:10.4230/LIPIcs.SoCG.2021.58.
  * [44] Wilson A. Sutherland.  Introduction to metric and topological spaces.  Oxford Mathematics Series. Oxford University Press, 2 edition, 2009.

##  Appendix A Details on $k$-greedy permutation computation

We change the $\mathcal{O}(n^{2}\log k)$ algorithm from Section 3. Recall that
$T_{y}$ is a max-heap containing the distances from $y$ to its $k$ closest
ordered sites and that the list $L$ groups the $T_{y}$ of all $y$ unordered.
The modifications are as follows: at initialization, set up a quadtreap with
all $n-k$ unordered sites. After $p_{i}$ is picked, query the quadtreap at
$p_{i}$ with radius $\lambda_{i}$, obtaining a set $S_{i}$. Then delete
$p_{i}$ from the quadtreap and iterate through $S_{i}$ updating the $k$
nearest neighbors with respect to $P_{i}$.

Moreover, to avoid the linear scan of $L$ to identify $p_{i}$, we replace $L$
by a max heap whose priority value is the $k$-distance to ordered sites, such
that $p_{i}$ can be determined in constant time. This requires to maintain the
heap property when the $k$th nearest neighbor of an unordered site changes, so
the running time of updating an element in $L$ worsens from $\mathcal{O}(\log
k)$ to $\mathcal{O}(\log n)$ (this worsening can be avoided by altering the
organization of the heap slightly, but we omit details for brevity). It is
then not difficult to see that the expected running time is bounded by
$\mathcal{O}(nk\log k+n\log n+\Sigma\log n)$, where $\Sigma$ is the
accumulated size of all reported sets $S_{i}$ of the quadtreap. It remains to
bound $\Sigma$.

The _spread_ $\Phi$ of a point set is the ratio of the diameter of the point
set to the shortest distance between two distinct points in the point set. We
claim that every site appears in at most $\mathcal{O}(k\log\Phi)$ sets
$S_{i}$, which results in the bound $\Sigma=\mathcal{O}(nk\log\Phi)$. To see
that consider the interval $(2^{s},2^{s+1}]$ for some integer $s$ and the
sites $p_{i},\ldots,p_{j}$ with $i\leq j$ such that
$\lambda_{i},\ldots,\lambda_{j}$ lies in that interval. Note that there are at
most $\log\Phi$ such intervals with at least one $\lambda$-value inside, so it
is enough to show that every site appears in at most $\mathcal{O}(k)$ sets
among $S_{i},\ldots,S_{j}$. Fix a site $p$. To be in the list $S_{\ell}$, it
must be at distance at most $2\cdot 2^{s+1}$ from $p_{\ell}$. In other words,
all sites among $p_{i},\ldots,p_{j}$ for which $p$ is reported lie in a ball
of radius $2^{s+2}$. This ball can be covered by $32^{\delta}=\mathcal{O}(1)$
balls of radius $2^{s-2}$ (see _Doubling Dimension_ in Section 2). None of
these balls can contain more than $k$ sites of $p_{i},\ldots,p_{j}$: if
$p_{\ell}$ is one site in such a ball with $k$ closest ordered sites in the
ball, we have that $d_{k}(p_{\ell},P_{j}\setminus\\{p_{\ell}\\})\leq
2^{s-1}=\dfrac{2^{s}}{2}<\dfrac{\lambda_{j}}{2}$, which is a contradiction to
the Packing property. The bound of $\Sigma$ follows. Hence we have proved the
run time claimed in Theorem 3.5.

##  Appendix B Critical value of simplices

We elaborate on how to determine if and when a $m$-simplex
$\sigma=(A_{0},\ldots,A_{m})$, with $m\geq 1$, appears in the filtration
$\mathcal{D}$. Recall that $\sigma$ is in the filtration if there is some
$r\geq 0$ such that the restricted lenses
$\tilde{L}_{r}(A_{0}),\ldots,\tilde{L}_{r}(A_{m})$ have a common intersection.

Let $r_{MEB}$ denote the radius of the minimum enclosing ball of the sites in
$A_{0}\cup\ldots\cup A_{m}$. We can compute this radius in expected linear
time in the number of sites (see [16]), which is at most $k(m+1)$. Let
$\Lambda$ denote the smallest freezing time among . The unrestricted lenses
$L_{r}(A_{0}),\ldots,L_{r}(A_{m})$ intersect if and only if $r\geq r_{MEB}$.
There are three cases to consider: first, if $r_{MEB}\leq\Lambda$, then
$\tilde{L}_{r_{MEB}}(A_{i})=L_{r_{MEB}}(A_{i})$ for all $i$ and they intersect
at radius $r_{MEB}$. In that case, $\sigma$ belongs to the filtration and its
critical value is the minimal value of the form $(1+\epsilon^{\prime})^{z}$
that is greater or equal to $r_{MEB}$. Second, if
$r_{MEB}>(1+\epsilon^{\prime})\Lambda$, the unrestricted lenses do not
intersect until the removal radius of one lens. Since the restricted lenses
are subsets of their unrestricted counterparts, it follows that
$\tilde{L}_{r}(A_{0}),\ldots,\tilde{L}_{r}(A_{m})$ do not intersect before one
of them becomes empty. Hence $\sigma$ does not belong to the filtration in
this case.

The final and most challenging case is if
$\Lambda<r_{MEB}\leq(1+\epsilon^{\prime})\Lambda$. To manage this case, we
first consider the following general geometric problem: given balls
$\beta_{1},\ldots,\beta_{\ell}$ with centers $x_{1},\ldots,x_{\ell}$ and radii
$r_{1},\ldots,r_{\ell}$, decide whether
$\beta_{1}\cap\ldots\cap\beta_{\ell}\neq\emptyset$ $(\star)$. Assume wlog that
all radii $r_{1},\ldots,r_{\ell}$ lie between $0$ and $1$ (otherwise, scale
the instance). Now, define $\beta_{j}^{\prime}$ as the ball with center
$x_{j}$ and radius $1-r_{j}$, for $1\leq j\leq\ell$. We use the following
geometric predicate to deal with this problem:

######  Lemma B.1.

The balls $\beta_{1},\ldots,\beta_{\ell}$ have a common intersection point if
and only if there exists a ball of radius $1$ that contains
$\beta_{1}^{\prime},\ldots,\beta_{\ell}^{\prime}$.

######  Proof B.2.

Assume there is a ball $B_{1}(c)$ that contains each $\beta_{j}^{\prime}$ .
Then $d(c,x_{j})+(1-r_{j})\leq 1$ and thus $d(c,x_{j})\leq r_{j}$ for each
$1\leq j\leq\ell$. It follows that $c\in\beta_{j}$ for every $1\leq
j\leq\ell$, implying “$\Leftarrow$”. The argument can be reversed to show the
equivalence.

The question whether a set of balls is contained in a ball of fixed radius can
be answered by computing the _minimum enclosing ball of balls_ [26]. Hence
this geometric predicate can be used to solve $(\star)$. For constant
dimension, the expected running time for computing the minimum enclosing ball
of $\ell$ balls is $\mathcal{O}\left(\ell\right)$, the same as for standard
minimum enclosing ball computations.

We now use this predicate to decide whether $\sigma=(A_{0},\ldots,A_{m})$ is
in the filtration. Write $\omega:=(1+\epsilon^{\prime})\Lambda$ for brevity.
For $0\leq j\leq m$, set $r_{j}$ as follows: if the restricted lens of $A_{j}$
is not frozen when the first lens disappears, i.e., if $\omega\leq$, set
$r_{j}:=\omega$. Otherwise, $r_{j}:=$. Next, for each $A_{j}$ we define $k$
balls with radius $r_{j}$, each centered at a different site of $A_{j}$. These
$k(m+1)$ balls intersect if and only if the restricted lenses
$\tilde{L}_{\omega}(A_{0}),\ldots,\tilde{L}_{\omega}(A_{m})$ intersect.
Furthermore, the latter is equivalent to $\sigma$ being in the filtration,
because $\omega$ is the largest radius at which the restricted lenses could
intersect (afterwards one of them is empty and so is their intersection). Thus
we decide whether $\sigma$ belongs to the filtration by creating an instance
of the ball intersection problem with $\ell=k(m+1)$, which takes time
$\mathcal{O}(km)$.

If $\sigma$ is in the filtration, it remains to find its critical value. Note
that by our discretization, the critical value is necessarily an integer power
of $(1+\epsilon^{\prime})$. Let $z$ be the smallest integer such that
$(1+\epsilon^{\prime})^{z}\geq\omega$. By construction, the cones at this
radius intersect (because they intersect at radius $\omega$), hence the
critical value is at most $(1+\epsilon^{\prime})^{z}$. On the other hand, the
critical value must be at least $(1+\epsilon^{\prime})^{z-1}$ because
$(1+\epsilon^{\prime})^{z-1}<\omega$ implies
$(1+\epsilon^{\prime})^{z-2}<\omega/(1+\epsilon^{\prime})=\Lambda$: if the
restricted lenses would intersect at radius $(1+\epsilon^{\prime})^{z-2}$ or
smaller, we would have that $r_{MEB}\leq\Lambda$, contradicting the assumption
in the case distinction above.

This suggests a simple algorithm to compute the critical value: set up another
instance of the balls intersection predicate, just as before, but using the
radius $(1+\epsilon^{\prime})^{z-1}$ instead of $\omega$. If the corresponding
balls intersect, we set the critical value of $\sigma$ to be
$(1+\epsilon^{\prime})^{z-1}$. Otherwise, we set it to
$(1+\epsilon^{\prime})^{z}$.

Considering all discussed above, the running time of deciding whether one
simplex is in $\mathcal{D}$ and computing its critical value is linear in
$km$, because it requires the computation of at most three minimum enclosing
ball instances (one for computing $r_{MEB}$, one to check for intersections at
$\omega$ and one to check for intersections at $(1+\epsilon^{\prime})^{z-1}$).
See Algorithm 1 for the pseudo-code of the entire computation of the
filtration.

We point out that we were not able to find an algorithm to compute exactly the
smallest $r$ for which $\tilde{L}_{r}(A_{0}),\ldots,\tilde{L}_{r}(A_{m})$ have
a common intersection, which is the reason why we need to introduce the
discretization of the sparse filtration. Note that the smallest $r$ can be
approximated arbitrarily closely by combining the minimum enclosing ball
method described above and a binary search, however the cost of this
approximation seems prohibitive as it does not improve significantly the
overall approximation quality of the filtration and its persistence module.

Data: $k$, $\epsilon^{\prime},$ $P$, $n=|P|$.

Result: $\\{D_{r}\\}_{r\geq 0}$ as a list of pairs (simplex, critical value).

Initialize _$D_{r}=\emptyset$._order(P); /* order the set according to
$k$-distance permutation */

1 for  _$i=n,n-1,\ldots,k+1$_ do

2
$z=\left\lfloor\log_{(1+\epsilon^{\prime})}\dfrac{\lambda_{i}(1+\epsilon^{\prime})^{2}}{\epsilon^{\prime}}\right\rfloor$;

3 $F=\emptyset$;

4 $\omega=\operatorname{frz}\left(p_{i}\right)(1+\epsilon^{\prime})$;

range search for $J=P_{i}\cap B_{2\omega}(p_{i})$; /* find friends of $p_{i}$
*/

/* add vertices with $p_{i}$ */

5 for  _all $k$-subsets $V\subseteq J$ such that $p_{i}\in V$_ do

6 if  _$r_{MEB}(V)\leq$_ then

7
$F=F\cup\left\\{\left(V,(1+\epsilon^{\prime})^{\left\lceil\log_{(1+\epsilon^{\prime})}r_{MEB}(V)\right\rceil}\right)\right\\}$;

8

9 end if

10

11 end for

12

/* add higher dimensional simplices with $p_{i}$ */

13 for  _$m=2,3,\cdots$_ do

14 $\mathcal{C}=$ all $(m-1)$-dimensional simplices in $F$;

15 if  _$\mathcal{C}=\emptyset$_ then break;

16 for  _$\sigma$ a $m$-dimensional simplex over the sites in $J$ s.t
$\exists\,\tau\in\mathcal{C}$ with $\tau\subseteq\sigma$_ do

17 if  _lenses in $\sigma$ intersect at radius $(1+\epsilon^{\prime})^{z}$_
then

18 $F=F\cup\\{(\sigma,(1+\epsilon^{\prime})^{z})\\}$;

19

20 end if

21 else if  _lenses in $\sigma$ intersect at radius
$\operatorname{frz}\left(p_{i}\right)(1+\epsilon^{\prime})$_ then

22 $F=F\cup\\{(\sigma,(1+\epsilon^{\prime})^{z+1})\\}$;

23

24 end if

25

26 end for

27

28 end for

29

/* remove $p_{i},$ update the filtration and continue */

30 $P=P\setminus\\{p_{i}\\}$;

31 $D_{r}=D_{r}\cup F$;

32

33 end for

34

Output: $D_{r}$.

Algorithm 1 Pseudo-code for computing the filtration

##  Appendix C Lower bound for output size

We show a lower bound for the number of output simplices in the sparse
filtration and elaborate on why does this lower bound guarantee that the term
$({nk}/{\epsilon^{\delta}}+(n-k)k\Gamma^{k-1})$ is
$\mathcal{O}(Xk\Gamma^{k}\cdot m_{max})$ on Equation (2).

######  Lemma C.1.

The number of output simplices is $X\geq n-k$.

######  Proof C.2.

We show that for every site $p_{i}$ with $\lambda_{i}<\infty$, there is one
vertex associated to $p_{i}$ that is inserted in the sparse filtration. Then,
the number of output simplices is at least the number of $p_{i}$s satisfying
that condition, which is $n-k$.

Recall that $\lambda_{i}=\operatorname{d}^{k}(p_{i},P_{i-1})$. Since by
definition
$\operatorname{d}^{k}(p_{i},P_{i-1})\geq\operatorname{d}_{k-1}(p_{i},P_{i-1})=:r$,
it holds that the lens $A$ defined by $p_{i}$ and its $k-1$ nearest neighbors
in $P_{i-1}$ becomes non-empty at most at radius $r\leq\lambda_{i}$.
Furthermore, since $A$ depends only on $p_{i}$ and sites with index smaller
than $i$, it is associated to $p_{i}$. Hence the vertex corresponding to $A$
is added to the filtration before $p_{i}$ is frozen.

Since $k,\Gamma\geq 1$, $m_{max}\geq 0$, the lemma immediately implies that
$(n-k)k\Gamma^{k-1}$ is $\mathcal{O}(Xk\Gamma^{k}\cdot(m_{max}+1))$. To see
that ${nk}/{\epsilon^{\delta}}$ is dominated by
$(Xk\Gamma^{k}\cdot(m_{max}+1))$ as well, we assume $n>2$ and consider two
cases. If $k\geq n/2>1$, then $\Gamma\geq 2$ and thus $k\geq n/\Gamma^{k-1}$.
In that case, $k\Gamma^{k}\geq n\Gamma$ and the result follows. If $k<n/2$,
then $n$ is $\mathcal{O}(n-k)$ and, by the previous lemma, it is
$\mathcal{O}(X)$. The result follows because $\Gamma>k/\epsilon^{\delta}$. To
finish we use that $m_{max}+1=\mathcal{O}(m_{max})$. Note that $n>2$ is not a
restrictive assumption, since sparsification becomes relevant only for large
values of $n$.