Estimating the Persistent Homology of ${ℝ}^n$-Valued Functions Using Function-Geometric Multifiltrations

In this paper we address both Q1 and Q2 using a unified persistence-based approach.

For Q2 we study the extension of the estimator $H_{\ast }({ℛ}^{\delta \to 2\delta }({𝒻|}_P))$ from [10] to vector-valued functions $𝒻:X\to {ℝ}^n$, where the codomain ${ℝ}^n$ is equipped with the product order, and we generalize its approximation guarantee to this setting. For the sake of generality we assume $𝒻$ to be $\omega$-continuous for some modulus of continuity $\omega$, that is: ${\parallel 𝒻(x)-𝒻(y)\parallel }_{\mathrm{\infty }}\le \omega (d_X(x,y))$ for any $x,y\in X$, where $\omega :{ℝ}_{\ge 0}\to {ℝ}_{\ge 0}$ is a non-decreasing subadditive function such that $\omega (\delta )$ goes to $0$ as $\delta$ does. This includes Lipschitz or Hölder continuous functions as special cases. Under this assumption, our first main result (Theorem 3) states that the estimator $H_{\ast }({ℛ}^{\delta \to 2\delta }({𝒻|}_P))$ is $\omega (2\delta )$-interleaved with its target $H_{\ast }(𝒻)$ for any choice of $\delta$ within the range $[2\epsilon ,{\varrho }_X/2)$, under the same regularity condition on $X$ and $\epsilon$-sampling condition on $P$ as before. Our proof takes inspiration from the one in [10] but uses a new, purely diagrammatic formulation that clarifies it and makes it hold for any $n\ge 1$.

For Q1 we propose two complementary approaches. The first one infers a relevant choice of $\delta$ under an appropriate statistical model. This model assumes the points $P$ are sampled i.i.d. according to some unkown probability distribution $\mu$ that is supposed to be $(a,b)$-standard. When parameters $a,b$ are known, we show that there is an explicit optimal choice ${\delta }_k$ for parameter $\delta$, as a function of the sample size $k$, that makes the sequence of approximations ${\left(H_{\ast }({ℛ}^{{\delta }_k\to 2{\delta }_k}({𝒻|}_P))\right)}_{k\in \mathbb{N}}$ a quasi-minimax estimator of $H_{\ast }(𝒻)$. When $a,b$ are unknown, as would generally be the case in practice, we show that ${\delta }_k$ can be effectively estimated by some ${\widehat{\delta }}_k$ under our statistical model, so that the sequence ${\left(H_{\ast }({ℛ}^{{\widehat{\delta }}_k\to 2{\widehat{\delta }}_k}({𝒻|}_P))\right)}_{k\in \mathbb{N}}$ is again a quasi-minimax estimator of $H_{\ast }(𝒻)$. These statements are summarized in Theorem 5.

Our second approach to addressing Q1 is meant for cases where our statistical model does not hold and $\epsilon$ cannot be estimated reliably. We follow the general philosophy underpinning persistent homology by letting $\delta$ vary over ${ℝ}_{\ge 0}$ and treating it as an extra parameter for persistence. This means that our estimator now becomes the $(n+1)$-parameter persistence module $H_{\ast }({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$ indexed over ${ℝ}_{\ge 0}\times {ℝ}^n$, where the notation $\bullet$ indicates that the Rips parameter is also a filtration parameter. In order to enable comparisons, the target $H_{\ast }(𝒻)$ must also be extended to ${ℝ}_{\ge 0}\times {ℝ}^n$, which we do by taking its left Kan extension ${\mathrm{Lan}}_{\iota }{H}_{\ast }(𝒻)$ along the poset embedding $\iota :{ℝ}^n\stackrel{\cong }{⟶}\{0\}\times {ℝ}^n\hookrightarrow {ℝ}_{\ge 0}\times {ℝ}^n$. The extension ${\mathrm{Lan}}_{\iota }{H}_{\ast }(𝒻)$ merely contains copies of $H_{\ast }(𝒻)$ in each vertical hyperplane $\{\delta \}\times {ℝ}^n$, connected by horizontal identities, as illustrated in Figure 1. Our third main result (Theorem 6) states that the $n$-dimensional interleavings happening between $H_{\ast }(𝒻)$ and $H_{\ast }({ℛ}^{\delta \to 2\delta }({𝒻|}_P))$ inside each vertical hyperplane $\{\delta \}\times {ℝ}^n$ for $\delta \in [2\epsilon ,{\varrho }_X/2)$ commute with the horizontal morphisms in $H_{\ast }({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$ and in ${\mathrm{Lan}}_{\iota }{H}_{\ast }(𝒻)$, so that, within any vertical slab $[2\epsilon ,{\delta }_0]\times {ℝ}^n$ with , they all together form a vertical $\omega (2{\delta }_0)$-interleaving (i.e., an $\omega (2{\delta }_0)$-interleaving along the direction ${(0,1,\mathrm{\cdots },1)}^{\top }$ in ${ℝ}_{\ge 0}\times {ℝ}^n$) between $H_{\ast }({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$ and ${\mathrm{Lan}}_{\iota }{H}_{\ast }(𝒻)$.

Since vertical interleavings are stronger than ordinary interleavings (Remark 1), Theorem 6 not only provides approximation guarantees for any interleaving-stable invariant of ${\mathrm{Lan}}_{\iota }{H}_{\ast }(𝒻)$ by the corresponding invariant of $H_{\ast }({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$, but it also offers refined statements for invariants that behave in a particular way under vertical interleavings, such as for instance multigraded Betti numbers as proven in the full version [1] and illustrated in Figure 2.

In practice, function values may be corrupted with measurement noise, meanwhile geodesic distances may have to be approximated from the input data. Our framework can account for these imprecisions by suitably adapting the parameters of the estimators and their approximation bounds (see the full version [1, Section 3.3]).

Still for practical purposes, we discuss the computation of our estimators, in particular we provide an algorithm for computing a free presentation of $H_{\ast }({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$, and thus also of $H_{\ast }({ℛ}^{\delta \to 2\delta }({𝒻|}_P))$ for any fixed $\delta \ge 0$, in time that is comparable to that of computing a free presentation of the persistent homology of a single multifiltration. From a free presentation, a variety of invariants, including multigraded Betti numbers, can be derived using existing software. We have implemented our algorithm and used it in our experiments.

We provide background material in Section 2. Our problem statement, estimators, and main results are detailed in Section 3. In Section 4 we discuss algorithmic aspects. Finally, in Section 5 we present experimental results that illustrate our theoretical guarantees and we investigate additional properties of our estimators. Further algorithmic details, as well as the proofs of our main results, can be found in the full version [1].

We see ${ℝ}^n$ as the product of $n$ copies of the real line, equipped with the product order noted $\le$. Thus, two points $𝒙,𝒚\in {ℝ}^n$ satisfy $𝒙\le 𝒚$ whenever ${𝒙}_i\le {𝒚}_i$ for all $1\le i\le n$. We denote by $(-\mathrm{\infty },𝒙]$ the downset $\{𝒚\in {ℝ}^n\mid 𝒚\le 𝒙\}$, and by $[𝒙,+\mathrm{\infty })$ the upset $\{𝒚\in {ℝ}^n\mid 𝒚\ge 𝒙\}$. An $n$-parameter filtration is a functor $ℱ:{ℝ}^n\to \mathrm{𝐓𝐨𝐩}$ whose constituent maps are inclusions, which means that we have a topological space ${ℱ}_{𝒙}$ for all $𝒙\in {ℝ}^n$ and an inclusion ${ℱ}_{𝒙,𝒚}:{ℱ}_{𝒙}\hookrightarrow {ℱ}_{𝒚}$ for all $𝒙\le 𝒚\in {ℝ}^n$. An $n$-parameter persistence module is a functor $M:{ℝ}^n\to \mathrm{𝐕𝐞𝐜}$, which means that we have a $𝕜$-vector space $M_{𝒙}$ for all $𝒙\in {ℝ}^n$ and a $𝕜$-linear map $M_{𝒙,𝒚}:M_{𝒙}\to M_{𝒚}$ for all $𝒙\le 𝒚\in {ℝ}^n$, such that $M_{𝒙,𝒙}={\mathrm{Id}}_{M_{𝒙}}$ and $M_{𝒙,𝒛}=M_{𝒚,𝒛}\circ M_{𝒙,𝒚}$ for all $𝒙\le 𝒚\le 𝒛\in {ℝ}^n$. A morphism of persistence modules $M\to N$ is a natural transformation between functors, i.e., a family of linear maps ${\phi }_{𝒙}:M_{𝒙}\to N_{𝒙}$ such that $N_{𝒙,𝒚}\circ {\phi }_{𝒙}={\phi }_{𝒚}\circ M_{𝒙,𝒚}$ for all $𝒙\le 𝒚\in {ℝ}^n$.

Throughout the paper, unless otherwise specified, $(X,d_X)$ is a compact geodesic metric space. This means that, for all $x,y\in X$, there is a shortest continuous path from $x$ to $y$ in $X$, of length equal to . In particular, the space is path-connected. Let $P$ be a finite set of points of $X$. Given $\epsilon >0$, we say $P$ is a geodesic $\epsilon$-sample of $X$ if , where $d_{\text{H}}$ denotes the Hausdorff distance in $(X,d_X)$. Define $B_X(x,r)$ as the open geodesic ball centered at $x\in X$ with radius $r$, given by A ball $B_X(x,r)$ is said to be convex if, for any pair of points $y$ and $y^{\prime }$ in $B_X(x,r)$, the shortest path in $X$ that connects $y$ to $y^{\prime }$ is unique and included in $B_X(x,r)$. Let $\varrho (x):=inf\{r>0\mid B(x,r)\text{ is not convex}\}$, and let ${\varrho }_X:=inf\{\varrho (x)\mid x\in X\}\ge 0$. This quantity, called the convexity radius of $X$, plays an important role because intersections of convex sets are convex and therefore contractible, so the Nerve Lemma [13, Corollary 4G.3] applies to unions of balls of radii less than ${\varrho }_X$.

Given a geodesic metric space $(X,d_X)$ and a function $𝒻:X\to {ℝ}^n$, the sublevel filtration $ℱ:{ℝ}^n\to \mathrm{𝐓𝐨𝐩}$ of $𝒻$ is defined by ${ℱ}_{𝒙}:={𝒻}^{-1}((-\mathrm{\infty },𝒙])$ for all $𝒙\in {ℝ}^n$. The persistent homology of $ℱ$ (and, by extension, of $𝒻$) is the persistence module $H_{\ast }(ℱ)$ (also written $H_{\ast }(𝒻)$) defined by $H_{\ast }{(ℱ)}_{𝒙}:=H_{\ast }({ℱ}_{𝒙})$ and $H_{\ast }{(ℱ)}_{𝒙,𝒚}:=H_{\ast }({ℱ}_{𝒙,𝒚})$ for all $𝒙\le 𝒚\in {ℝ}^n$.

Given a finite set $P$ of points in X, and a real parameter $\delta \ge 0$, we call $\delta$-offset of $P$ the union of open geodesic balls ${𝒪}^{\delta }(P):={\bigcup }_{p\in P}{B}_X(p,\delta )$. We call $\delta$-Čech Complex of $P$ the nerve of the collection of balls $\{B_X(p,\delta )\mid p\in P\}$, defined as the abstract simplicial complex ${𝒞}^{\delta }(P):=\{\mathrm{\varnothing }\ne \sigma \subseteq P\mid {\bigcap }_{p\in \sigma }{B}_X(p,\delta )\ne \mathrm{\varnothing }\}$. For , the balls $B_X(p,\delta )$ and their intersections are either empty or convex, so the Nerve Lemma [13, Corollary 4G.3] ensures that the homology groups $H_{\ast }({𝒪}^{\delta }(P))$ and $H_{\ast }({𝒞}^{\delta }(P))$ are isomorphic. We call $\delta$-Rips complex of $P$ the abstract simplicial complex . Čech and Rips complexes of $P$ are related as follows: ${𝒞}^{\delta }(P)\subseteq {ℛ}^{2\delta }(P)\subseteq {𝒞}^{2\delta }(P)\text{ for all }\delta \ge 0.$

Let $𝒗\in {ℝ}_{\ge 0}^n$ be a fixed vector, and let $\epsilon \ge 0$. Given any functors $F,G:{ℝ}^n\to 𝒞$ to some category $𝒞$, we say $F$ and $G$ are $\epsilon 𝐯$-interleaved if there are two natural transformations $f:F\to G[\epsilon 𝒗]$ and $g:G\to F[\epsilon 𝒗]$ such that $g[\epsilon 𝒗]\circ f={\phi }_F^{2\epsilon 𝒗}$ and $f[\epsilon 𝒗]\circ g={\phi }_G^{2\epsilon 𝒗}$. Here, $-[\epsilon 𝒗]$ denotes the $\epsilon 𝒗$-shift functor, defined on objects by $F{[\epsilon 𝒗]}_{𝒙}:=F_{𝒙+\epsilon 𝒗}$ and on morphisms by $F{[\epsilon 𝒗]}_{𝒙,𝒚}:=F_{𝒙+\epsilon 𝒗,𝒚+\epsilon 𝒗}$ for all $𝒙\le 𝒚\in {ℝ}^n$. Meanwhile, ${\phi }_F^{2\epsilon 𝒗}:F\to F[2\epsilon 𝒗]$ and ${\phi }_G^{2\epsilon 𝒗}:G\to G[2\epsilon 𝒗]$ denote the natural morphisms from $F$ and $G$ to their respective $2\epsilon 𝒗$-shifts, defined by ${({\phi }_F^{2\epsilon 𝒗})}_{𝒙}:=F_{𝒙,𝒙+2\epsilon 𝒗}$ and ${({\phi }_G^{2\epsilon 𝒗})}_{𝒙}:=G_{𝒙,𝒙+2\epsilon 𝒗}$ for all $𝒙\in {ℝ}^n$. The $𝐯$-interleaving distance is defined as follows:

In the following, we consider two specific types of interleavings:

1. 1.

   Ordinary interleavings: Let $𝒗=\mathrm{𝟏}:={(1,\mathrm{\cdots },1)}^T$. Then, $\epsilon \mathrm{𝟏}$-interleavings are ordinary $\epsilon$-interleavings from the TDA literature, and $d_{\mathrm{i}}^{\mathrm{𝟏}}$ is called the ordinary interleaving distance.

2. 2.

   Vertical interleavings: Let $𝒗={\mathrm{𝟏}}_0:={(0,1,\mathrm{\cdots },1)}^T.$ Then, $\epsilon {\mathrm{𝟏}}_0$-interleavings are called vertical $\epsilon$-interleavings, and $d_{\mathrm{i}}^{{\mathrm{𝟏}}_0}$ is called the vertical interleaving distance.

Suppose we know $\epsilon$ or some reasonable estimate. Then we can approximate our target $H_{\ast }(𝒻)$ by restricting our estimators to a vertical hyperplane $\{\delta \}\times {ℝ}^n\subseteq {ℝ}_{\ge 0}\times {ℝ}^n$ for a suitable choice of parameter $\delta$. We denote these restrictions respectively by $H_{\ast }({𝒪}^{\delta }({𝒻|}_P))$, $H_{\ast }({𝒞}^{\delta }({𝒻|}_P))$ and $H_{\ast }({ℛ}^{\delta \to 2\delta }({𝒻|}_P))$.

The proof of this result is given in the full version [1, Section 4]. The conclusion still holds under noisy input function values or geodesic distances, as detailed in [1, Section 3.3].

Theorem 3 can be read in two different ways. First, when the modulus of continuity $\omega$ is known, we obtain explicit bounds on the interleaving distance between our estimators and the target $H_{\ast }(𝒻)$. Second, when we only know that $𝒻$ is continuous, the Heine–Cantor Theorem implies that $𝒻$ admits some modulus of continuity $\omega$, and even without knowning that $\omega$, we are guaranteed by Theorem 3 that our estimators converge to $H_{\ast }(𝒻)$ as $\epsilon$ goes to zero under suitable choices of parameter $\delta$.

When $\epsilon$ is known exactly, the best choice for $\delta$ is the lower bound of its admissible interval, i.e., $\delta =\epsilon$ in Items i and ii and $\delta =2\epsilon$ in Item iii. When $\epsilon$ is only known approximately, the admissible interval provides some leeway for the choice of $\delta$.

An immediate consequence of Theorem 3 is that $d_{\mathrm{i}}^{\mathrm{𝟏}}$-stable invariants computed on our estimators approximate the corresponding invariants defined on the target $H_{\ast }(𝒻)$. For instance, the approximation guarantee obtained for the multigraded Betti numbers of $H_{\ast }({ℛ}^{\delta \to 2\delta }({𝒻|}_P))$ is given in the full version [1, Corollary 3.3].

Another consequence of Theorem 3 (i) is that, under our assumptions, the persistence module $H_{\ast }(𝒻)$ satisfies the following form of tameness, as proven in the full version [1].

Suppose now the points of $P$ are i.i.d. samples drawn from some unknown probability measure $\mu$ supported on $X$. In this setting, standard statistical techniques can be used to estimate $\epsilon$ with high probability, then Theorem 3 can be applied to get statistical estimates of $H_{\ast }(𝒻)$ or of any $d_{\mathrm{i}}^{\mathrm{𝟏}}$-stable invariant thereof. We summarize our results in Theorem 5. The relevant statistical background and definitions are provided in the full version [1, Section 2], together with the analysis of the statistical estimators [1, Section 6].

When $\epsilon$ is unknown and cannot be effectively estimated – for instance when the assumptions in Theorem 5 are not satisfied, we adopt the traditional approach in persistence theory which is to not fix the scale parameter $\delta$ but rather to consider it as an extra filtration parameter. This means using the full $(n+1)$-parameter modules $H_{\ast }({𝒪}^{\bullet }({𝒻|}_P))$, $H_{\ast }({𝒞}^{\bullet }({𝒻|}_P))$ and $H_{\ast }({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$ as estimators. In order to state approximation results, we extend our target $H_{\ast }(𝒻)$ to a persistence module over ${ℝ}_{\ge 0}\times {ℝ}^n$ by taking its left Kan extension ${\mathrm{Lan}}_{\iota }{H}_{\ast }(𝒻)$ along the poset embedding $\iota :{ℝ}^n\hookrightarrow {ℝ}_{\ge 0}\times {ℝ}^n$ given by $𝒙\mapsto (0,𝒙)$ for all $𝒙\in {ℝ}^n$. Since the category ${ℝ}^n$ is small and the category $\mathrm{𝐕𝐞𝐜}$ is cocomplete, ${\mathrm{Lan}}_{\iota }{H}_{\ast }(𝒻)$ is well-defined as a functor ${ℝ}_{\ge 0}\times {ℝ}^n\to \mathrm{𝐕𝐞𝐜}$, and given pointwise by the colimit formula [19]:

Moreover, ${\mathrm{Lan}}_{\iota }{H}_{\ast }(𝒻)$ has structural morphisms that are identity maps horizontally (i.e., along the first coordinate axis) and the structural morphisms of $H_{\ast }(𝒻)$ vertically (i.e., along any direction orthogonal to the first coordinate axis), as illustrated in Figure 1.

The proof is given in the full version [1, Section 5]. Similar to Theorem 3, the conclusion still holds under noisy input function values or geodesic distances – see [1, Section 3.3].

This approximation result encapsulates the previous one in the sense that restricting the estimators to the vertical hyperplane $\{\delta \}\times {ℝ}^n$ for $\delta ={\delta }_0$ recovers Theorem 3. However, the result says something deeper, namely: that the interleavings in the vertical hyperplanes $\{\delta \}\times {ℝ}^n$, for $\delta$ ranging over $[\epsilon ,{\delta }_0]$ (resp. $[2\epsilon ,{\delta }_0]$), commute with the horizontal structural morphisms of the modules, so that they all together form a vertical interleaving in the slab $[\epsilon ,{\delta }_0]\times {ℝ}^n$ (resp. $[2\epsilon ,{\delta }_0]\times {ℝ}^n$). Intuitively, this means that prominent features in $H_{\ast }(𝒻)$ are not just present in $H_{\ast }({ℛ}^{\delta \to 2\delta }({𝒻|}_P))$ at individual scales $\delta$ in the range $[2\epsilon ,{\delta }_0]$, but that they persist across large ranges of scales in $H_{\ast }({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$. In particular, they are not ephemeral in ${ℝ}_{\ge 0}\times {ℝ}^n$ and can be captured by stable invariants.

In terms of approximation accuracy, the result says that, the more one looks to the left (i.e., toward small positive values of ${\delta }_0$) in the parameter space ${ℝ}_{\ge 0}\times {ℝ}^n$, the more precisely the estimators approximate the target, until some point where ${\delta }_0$ becomes too small and the approximation breaks down. On the opposite side, the more one looks to the right (i.e., toward large values of ${\delta }_0$), the more the estimators drift away from the target, until again the approximation eventually breaks down if the estimator is $H_{\ast }({𝒞}^{\bullet }({𝒻|}_P))$ or $H_{\ast }({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$. See Figure 2 for an illustration, and Figure 3 for a comparison of Theorem 6 with Theorem 3.

Theorem 6, like its counterpart Theorem 3, implies that provably $d_{\mathrm{i}}^{\mathrm{𝟏}}$-stable invariants computed from our estimators approximate the corresponding invariants defined on their target, which this time is ${\mathrm{Lan}}_{\iota }{H}_{\ast }(𝒻)$. And since the interleaving is now vertical – hence stronger than an ordinary interleaving according to Remark 1, we can get more precise statements for some invariants. For instance, the induced bottleneck matching on multigraded Betti numbers is vertical, as proven in the full version [1] and illustrated in Figure 2.

We validate the claims in Theorems 3 and 6 through a simple example consisting of two circles with distinct radii and concentration levels, as shown in Figure 4 (a). In Figure 4 (b) we analyze the error between the estimator $H_1({ℛ}^{\delta \to 2\delta }({𝒻|}_P))$ and the target $H_1(𝒻)$ as a function of the scale parameter $\delta$. As guaranteed by Theorem 3, for $\delta$ within the range $[2\epsilon ,{\varrho }_X/2)$ the error between the estimator and the target is small, controlled by $\delta$. The approximation remains good even for $\delta$ going beyond that range, which is not explained by the theory but occurs here because the space and function are very regular in this example.

In Figure 4 (c) we show a visualization of our estimator $H_1({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$ using MMA [17]. As expected from Theorem 6, we observe two topological features in our estimator, corresponding to the two circles in the dataset. The values of $\delta$ at which these features appear then disappear are aligned with the endpoints of the interval of values of $\delta$ for which $H_1({ℛ}^{\delta \to 2\delta }({𝒻|}_P))$ is a good approximation of $H_1(𝒻)$ in Figure 4 (b).

For comparison, in Figure 4 (d) we show a visualization of the estimator $H_1({𝒞}^{\bullet }({𝒻|}_P))$ using MMA, for which we use the Euclidean distance instead of the geodesic distance in order to enable the construction of the function–Čech filtration. Again, we see two topological features corresponding to the two circles, and $H_1({𝒞}^{\delta }({𝒻|}_P))$ is a good approximation of $H_1(𝒻)$ for $\delta$ within a fairly large range.

We now add noise to the dataset sampled from the two circles, as shown in Figure 5 (a). As proven in the full version [1, Section 3.3], the trend is the same. While the theory requires a factor $\alpha >2$ in $H_1({ℛ}^{\bullet \to \alpha \bullet }({𝒻|}_P))$ to handle metric noise, Figure 5 shows that, in this example, the estimator $H_1({ℛ}^{\bullet \to 2\bullet }({𝒻|}_P))$ performs well even with $\alpha =2$.

In this section, we consider the digitized immunohistochemistry dataset from [24]. This dataset is comprised of manually annotated cell locations from three cell types, namely cytotoxic T lymphocytes, regulatory T lymphocytes, and macrophages, each characterized by the expression of one of the proteins CD8, FoxP3, and CD68, respectively – see Figure 6.

The index $i\in \{1,2,3\}$ refers to the protein in the ordered set $\{\mathrm{CD8},\mathrm{FoxP3},\mathrm{CD68}\}$, and ${𝒻}_i:{[0,1]}^2\to ℝ$ denotes the negative density estimate associated with that protein. Since the functions ${𝒻}_i$ are unknown, we estimate the density of each point cloud using kernel density estimation, which yields approximations – also called ${𝒻}_1,{𝒻}_2,{𝒻}_3$ for simplicity – that are Lipschitz continuous by construction. We then discretize the domain ${[0,1]}^2$ on a grid (using a resolution of $\mathrm{1\; 000}\times \mathrm{1\; 000}$ for Section 5.2.1 and of $500\times 500$ for Section 5.2.2) and compute the persistence module of the cubical complex induced by the estimated function on this grid, which serves as a proxy for the ground-truth target.

In the following experiments, in order to reach smaller Hausdorff distances between the sample $X_k$ and the domain $X={[0,1]}^2$, we first draw $10,000$ points uniformly from $X$, and then select $X_k$ using $k$-means farthest point sampling applied to this set. In practice, the choice of the rate ${\delta }_k$ prescribed in Theorem 5 can be improved significantly when the sample $X_k$ is well-behaved. We propose heuristics to estimate good rates in non-asymptotic regimes.

The Brownian motion is not Lipschitz continuous but almost surely almost everywhere $\omega$-continuous for any modulus of continuity $\omega$ of the form $\omega :x\mapsto c{x}^{\frac{1}{2}-\epsilon }$, with $\epsilon >0$ and $c>0$ large enough [20, Corollaries 1.20 and 5.3]. Let $𝒻:[0,1]\to ℝ$ be a realization of a standard Brownian motion, which can be approximated by a Rademacher random walk.

For each $k\in \mathbb{N}$, we build a $k$-sample $X_k$ in $[0,1]$ according to the uniform distribution, as shown in Figure 11 for several values of $k$. Then, since our goal here is to illustrate Theorem 3, we assume we have access to the sampling error ${\epsilon }_k:=d_{\text{H}}(X_k,[0,1])$ and so we use the estimator $H_0({𝒞}^{{\epsilon }_k}({𝒻|}_{X_k}))$ to approximate the target $H_0(𝒻)$ for various sample sizes $k$. The corresponding results are shown in Figure 12. A regression analysis of the experimental data indicates that $d_{\mathrm{i}}^{\mathrm{𝟏}}(H_0({𝒞}^{{\epsilon }_k}({𝒻|}_{X_k})),H_0(𝒻))\approx {10}^{0.16}\epsilon _k^{}{}_{}{}^{0.47},$ which aligns closely with the theoretical prediction from Theorem 3, stating that the interleaving distance between $H_0({𝒞}^{{\epsilon }_k}({𝒻|}_{X_k}))$ and $H_0(𝒻)$ should be at most $\omega ({\epsilon }_k)\lesssim \epsilon _k^{}{}_{}{}^{1/2}$.