In topological data analysis (TDA), persistence diagrams have been a succesful tool.
To compare them, Wasserstein and Bottleneck distances are commonly used.
We address the shortcomings of these metrics and show a way to investigate them in a systematic way by introducing bottleneck profiles.
This leads to a notion of discrete Prokhorov metrics for persistence diagrams as a generalization of the Bottleneck distance.
They satisfy a stability result and bounds with respect to Wasserstein metrics.
We provide algorithms to compute the newly introduced quantities and end with an discussion about experiments.

§ INTRODUCTION

The field of topological data analysis (TDA) is becoming a popular tool to study the structure of complex data.
One of the major tools of TDA is persistent homology (PH) <cit.>.
Its pipeline takes a (often highly complex) point cloud in Euclidean space as input and produces a point cloud in the plane, the persistence diagram (PD), as output.
Intuitively, persistence diagrams serve as a summary of the shape of the input data.
As a consequence, one can compare different shapes indirectly, by comparing their PDs.
The need for a robust and computationally efficient notion of distance for PDs arises.
Classically, one uses the Bottleneck and Wasserstein distances to this end <cit.>.
However, the Bottleneck distance only picks up the single biggest difference and the Wasserstein distance is prone to noise, as it picks up every difference no matter how small.

This fact motivates our work to search for new metrics.
Starting from the investigation of bottlenecks, we introduce the notion of the bottleneck profile of two PDs, which is a map $\R_{\ge 0} \to \N \cup\{\infty\}$.
This tool summarizes metric information at varying scales and generalizes the Bottleneck distance.
Also the Wasserstein distance can be, in special cases, computed from the bottleneck profile; in general, it can be bounded given a bottleneck profile.

The bottleneck profiles arises naturally as one discrete version of the Prokhorov distance, which is a classical tool in probability theory.
It turns out that the Bottleneck and the Prokhorov distance are just two instances of a whole family of Prokhorov-style metrics discussed in this paper.
This family is parametrised by subclass of functions $f\colon [0,\infty) \to [0,\infty)$.
Not every function $f$ gives in fact rise to a genuine metric; we examine the conditions on $f$ in which cases it does.

In addition to theoretical developement, in this paper we discuss algorithms to compute the bottleneck profile and various Prokhorov-type distances.
We provide a runtime analysis and experiments on a number of data sets. The algorithms are provided as an open source implementation.

§ BACKGROUND

§.§ Measure Theory

Let $(X,d)$ be a metric space.
It is complete if every Cauchy sequence has a limit in $X$.
It is separable if it has a countable dense subset.
A complete separable metrizable topological space is called a polish space.
For example, all Euclidean spaces $\R^n$ are polish.
Polish spaces are a convenient setting for measure or probability theory.

In general, we endow $X$ with the Borel $\Sigma$-algebra $\mathfrak{B}(M)$ and denote the set of probability measures by $\mathcal{P}(X)$.

Let us recall an important inequality:

Let $(X,\Sigma, \mu)$ be a measure space and let $f\colon X \to \R$ be a measurable function.
Then for any $p>0$,

\[
\mu\{|f(x)|>t\} \le \frac{1}{t^p}\int\limits_{|f|\ge t}|f|^p.
\]

§.§.§ Metrics for Probablity Measures

There are various ways to compare different probability measures.

For $p\ge 1$, the $p$-Wasserstein metric is

\begin{align*}
W_p\colon \mathcal{P}(X) \times \mathcal{P}(X) & \to [0,\infty]\\
W_p(\mu,\nu) & = \inf\limits_\gamma\left( \int d(x,y)^p d\gamma(x,y)\right)^\frac{1}{p},
\end{align*}

where the infimum is taken over all couplings $\gamma$ with marginals $\mu$ and $\nu$.

The $1$-Wasserstein metric is also known as Kantorovich metric or earth mover's distance.
The latter name is motivated by the idea of thinking about $\gamma$ as a transport plan for moving a pile of earth $\mu$ into the pile $\nu$.
The cost of transportation equals the distance by which the earth is moved.

Intuitively, there are two different ways to "slightly change" a measure.
The first one is to move all the mass by a tiny distance.
The second one is to move a tiny part of the mass arbitrarily, possibly very far away.
While the Wasserstein metric is stable under perturbations of the first kind, small chaanges of the second kind can result in large differences in the metric.
The Prokhorov metric seeks to resolve this problem.
It is constructed in such a way that an $\eps$-neighborhood of a measure is characterized as follows:
One may move $\eps$ of the mass arbitrarily and the rest by at most $\eps$.

We now formalize this idea.

We denote measures by $\mu$ and $\nu$.
Couplings, that are measures on the product space with prescribed marginalization, are denoted $\gamma$.
For a Borel set $A\subset X$, the (open) $\eps$-ball around $A$ is

\[
A_\eps = \{ x\in X \colon d(x,A)<\eps\}
\]

The Prokhorov metric on $\mathcal{P}(X)$ is given by

\begin{align*}
\pi\colon\mathcal{P}(X) \times \mathcal{P}(X) & \to [0,\infty]\\
\pi(\mu, \nu) &= \inf\{\eps>0 \colon \forall A \in \mathfrak{B}(X): \mu(A) \le \nu(A_{\eps})+\eps \text{ and}\\
&\phantom{=inf\{} \nu(A) \le \mu(A_{\eps})+\eps\}.
\end{align*}

By Strassen's Theorem (cf. Remark 1.29 in <cit.>, Appendix 1.4), an alternative  characterization of the Prokhorov metric is given in terms of couplings $\gamma$ which marginalize to $\mu$ and $\nu$,

\begin{align}
\pi(\mu, \nu) &= \inf \{ \inf\{\eps>0 \colon\gamma(\{(x,y)\colon d(x,y)>\eps\})<\eps\} \colon\notag \\ & \phantom{= inf \{} \gamma \text{ has marginals } \mu \text{ and  } \nu\}.\label{eqn:MeasureProkhorov}
\end{align}

This allows for a discretization suitable for persistence diagrams, see Section <ref>.

Let $x_1, x_2 \in X$ and consider the Dirac measures $\delta_{x_1}, \delta_{x_2}$.
We claim that $\pi(\delta_{x_1}, \delta_{x_2}) = d(x_1,x_2)$.
The only coupling with correct marginals is $\delta_{(x_1,x_2)}$.
Then we have

\[
\delta_{(x_1,x_2)}(\{(x,y)\colon d(x,y)>\eps\}) = \mathbbm{1}_{\eps<d(x_1,x_2)}
\]

and consequently

\[
\inf\{\eps >0 \colon \mathbbm{1}_{\eps<d(x_1,x_2)}<\eps\} = d(x_1, x_2).
\]

In their survey <cit.>, Gibbs and Su show that

\[
\pi^2 \le W_1 \le (1+ \text{diam}(X)) \; \pi,
\]

where $\text{diam}(X)$ is the diameter of the underlying space.

For more on metrics of probability measures, see the book <cit.>; references for optimal transport include <cit.> which takes on a computational perspective.

§.§ Persistent Homology

The category $\pmod$ is the functor category

\[
\pmod = \text{Fun}(\R, k\textbf{-mod}^{fd})
\]

from the reals as a poset category to finite dimensional vector spaces. Its elements are called pointwise finite dimensional (p.f.d.) persistence modules.
A p.f.d. persistence module $\mathbb{A} = (A_t)_{t\in \R}$ comes with transition maps

\[
a_s^t \colon A_s \to A_t
\]

for $s\le t$.
If in addition, there are finitely many $t_1, \ldots t_n$ such that $a_s^t$ fails to be an isomorphism only if $s<t_i< t$ for some $i$, the p.f.d. persistence module is said to be tame.

An interval module for an interval $J\subset \R$ is a p.f.d persistence module with

\[
\mathbb{I}_{J, t} = \begin{cases}k &\textnormal{if }t\in J,\\
\end{cases}
\]

\[
\mathbb{I}_{J, s}\to \mathbb{I}_{J, t} = \begin{cases}\id &\textnormal{if } s,t\in J,\\
\end{cases}
\]

for $s<t$.
The start and endpoint of $J$ are referred to as birth time $b(J)$ and death time $d(J)$, respectively.

Note that we do not specify whether the endpoints are contained in the interval; they may be $\pm \infty$.
Interval modules are of special interest because p.f.d. persistence modules admit an interval decomposition.

Let $\mathbb{A}\in \pmod$. Then there exists a collection of intervals $\mathcal{J}$ such that

\[
\mathbb{A}\cong \bigoplus\limits_{J\in \mathcal{J}}\mathbb{I}_J.
\]

Any interval is possible in this decomposition.
The collection of birth and death times of the intervals from the decomposition form the
persistence diagram, a multiset containing

\[
\Dgm(\mathbb{A}) = \{(b(J),d(J)) \in \R^2 \colon J \in \mathcal{J}\}.
\]

and additionally all the points on the diagonal with countable multiplicities.

Since the persistence module is pointwise finite dimensional, off-diagonal points have finite multiplicities.
For a persistence diagram $X$ denote the multiplicity by $\mu_X\colon \R^2 \to \N$.
We follow the notation of <cit.>:
For a persistence diagram $X$ let $X_0$ denote its off-diagonal points and $X_0'$ their projection onto the diagonal.
We compare persistence diagrams $X,Y$ using the $p$-Wasserstein metric:

\begin{align*}
W_p(X,Y) &= \inf\limits_{\eta\colon X \to Y} \left(\sum \limits_{x \in X} d(x, \eta(x))^p\right)^\frac{1}{p}\\
&=\inf\limits_{\eta\colon X \to Y} \left(\sum \limits_{x \in X_0, \eta(x)\in Y_0} d(x, \eta(x))^p\right.\\
& \phantom{=\inf\limits_{\eta\colon X \to Y} \left(\right.}
\left.+ \sum \limits_{x \in X_0', \eta(x)\in Y_0} d(x, \eta(x))^p
+ \sum \limits_{x \in X_0, \eta(x)\in Y_0'} d(x, \eta(x))^p
\right)^\frac{1}{p},%\\
%&= \inf\limits_{\eta\colon X \to Y} \left(\sum \limits_{x \in X_0, \eta(x)\in Y_0} d(x, \eta(x))^p
% + \sum \limits_{x \in X_0', \eta(x)\in Y_0} \left(\frac{\eta(x)_2-\eta(x)_1}{2}\right)^p
% + \sum \limits_{x \in X_0, \eta(x)\in Y_0'} \left(\frac{x_2-x_1}{2}\right)^p
\end{align*}

where $\eta \colon X \to Y$ is a bijection between the persistence diagrams including diagonal points.
Notice the analogy with Definition <ref>.
We replace probability measures by counting measures and hence turn the integral into a sum.
The infimum is taken over all matchings instead of all couplings.
This observation will serve as a blueprint for the construction of the discrete Prokhorov metric for persistence diagrams in the next chapter.

The motivation to compare persistence diagrams comes from topological data analysis, where they serve as a summary statistic of topological information.

Given a finite subset $X$ of some metric space, we can consider the Vietoris-Rips complex, cf. <cit.> III.2.
This is a filtered simplicial complex; for filtration value $r>0$ it is given by

\[
VR(X)[r] = \{S \subset X \colon \diam(S)<2r \}.
\]

Applying the homology functor gives rise to a persistence module:
For $r\le s$, the inclusion

\[
VR(X)[r] \hookrightarrow VR(X)[s]
\]

induces a module homomorphism

\[
H_\ast(VR(X)[r]) \to H_\ast(VR(X)[s]).
\]

This is called the (Vietoris-Rips) persistent homology of $X$, denoted $PH_\ast(X)_t$.
Summands in its interval decomposition are interpreted as topological features which are "born" at a certain point in the filtration and "persist" for some time.
They are regarded to be more significant the longer the intervals are.
The following theorem ascertains that this is a useful tool.

Let $X,Y$ be finite metric spaces, fix some $k\ge 0$.

Then we have

\[W_\infty(\Dgm(PH_k(X)), \Dgm(PH_k(Y))) \le d_{GH}(X,Y),\]

where $d_{GH}$ is the Gromov-Hausdorff distance <cit.>.

In other words, if we change the input point cloud by $\eps$ in the Gromov-Hausdorff metric, the resulting PDs differ by at most $\eps$ in the Bottleneck distance.

§ BOTTLENECK PROFILES

The bottleneck distance $W_\infty$ has a major drawback:
It only captures the single most extreme difference between two persistence diagrams.
In other words: One and the same bottleneck distance can be attained in different ways.
For example, the bottleneck distance can be attained once or many times.
Moreover, the other distances in a matching are not taken into account.
We introduce the notion of the bottleneck profile to address the topic of secondary, tertiary,... bottlenecks and their multiplicities.

Given two persistence diagrams $X, Y$, define their bottleneck profile to be

\[
D_{X,Y}\colon \R_{\ge 0} \to [0,\infty],\qquad t \mapsto \inf\limits_{\eta\colon X \to Y} |\{x \colon d(x,\eta(x))>t\}|.
\]

For $d: \R^2 \times \R^2 \to \R_{\ge 0}$ we take a $p$-metric $d(x,y) = \|x-y\|_p$, where the choice of $p$ might depend on the setting.
For example, when comparing with the $p$-Wasserstein distance, one might like to choose this same $p$.

The bottleneck profile $D_{X,Y}$ is monotonically decrseaing.

Let $\eta\colon X \to Y$ be the matching realizing $D_{X,Y}(s)$ for some $s$.
Let now $t>s$, then every distance longer than $t$ is in particular longer than $s$ and consequently

\[
|\{x \colon d(x,\eta(x))>t\}| \le |\{x \colon d(x,\eta(x))>s\}| = D_{X,Y}(s).
\]

Taking the infimum over all matchings decreases the left hand side and yields $D_{X,Y}(t)$.

Knowing this, it is interesting when the bottleneck profile becomes zero.

The bottleneck profile $D_{X,Y}(t)$ vanishes for $t>W_\infty(X,Y)$.

By the definition, the bottleneck distance is the smallest $t>0$ such that there is a matching mapping all points within distance $t$.
In formulas,

\begin{align*}
W_\infty(X,Y) &= \inf\{t>0 \colon \inf\limits_{\eta\colon X \to Y} |\{x \colon d(x,\eta(x))>t\}| = 0\}\\
&=\inf\{t>0 \colon D_{X,Y}(t) = 0\}.\qedhere
\end{align*}

Similarly, the secondary bottleneck is $\inf\{t>0 \colon D_{X,Y}(t) = 1\}$ and so on.
This motivates the name bottleneck profile.

If we take one of the persistence diagrams to be the empty one, there is only one choice of matching: everything is paired with the diagonal.
As a consequence,

\[
D_{X,\emptyset}(t) = |\{x\colon \frac{x_2-x_1}{2}>t\}| = |\{x \colon x_1 + 2t < x_2\}|.
\]

This is also known as the stable rank function corresponding to the contour $C(a, \eps) = a + 2 \eps$, introduced in <cit.>.

The PD $X$ has bottleneck distance $3$ to each of the PDs $X,Y, Z$ (first three images).
However, it is attained with different multiplicities, which one can read off from the bottleneck profile (right-most image)

Consider some particular simple persistence diagrams.
The first three parts of Figure <ref> each show a base diagram ("Diagram $X$", in blue) with four points and perturbations of it:
The orange diagram ("Diagram $Y$") in the first image is obtained by shifting the blue one by three.
The green diagram ("Diagram $Z$") shifts the top point of $X$ by three, the next point by two, the third by one and leaves the lowest point unchanged.
For the yellow diagram ("Diagram $W$") in the third image, we only shift two points from $X$ by three and leave the other two untouched.
Clearly, the bottleneck distance between the base diagram and each of the shifted versions is three.
But the amount of shifted points is reflected in the bottleneck profile:
While $D_{X,Y}(t)$ is four, $D_{X,Z}(t)$ is two (i. e. the multiplicity of the bottleneck) for $0<t<3$.
And $D_{X,W}$ displays more steps, reflecting the fact that there are secondary and tertiary bottlenecks.

Note that the function $D$ enjoys some properties reminiscent of a metric (hence the notation $D$):
It is obviously symmetric.
The triangle inequality does not hold pointwise but in a scaled version, that is:

$D_{X,Z}(s+t)\le D_{X,Y}(s)+D_{Y,Z}(t)$.

This follows from the triangle inequality on $\R^2$.
For some $t$, let $\eta_{XY} \colon X \to Y$ and $\eta_{YZ}\colon Y \to Z$ denote optimal matchings realizing $D_{X,Y}(t)$ and $D_{Y,Z}(t)$, respectively.
Let $\eta = \eta_{Y,Z}\circ \eta_{X,Y}\colon X \to Z$ be the matching obtained by composition.
It suffices to show that

\[
|\{x\colon d(x,\eta(x))>s+t\}|\le |\{x\colon d(x,\eta_{X,Y}(x))>s\}| + |\{y\colon d(y,\eta_{Y,Z}(y))>t\}|,
\]

because the left hand side only decreases if we take the infimum over all matchings.
Hence we have to investigate what happens when a point $x$ is matched to $\eta(x)$ which is farther apart than $s+t$.
Note that $\eta(x) = \eta_{Y,Z}( \eta_{X,Y}(x))$, so we compare the distances of the matched points using the triangle inequality,

\[
s+t < d(x,\eta(x)) \le d(x,\eta_{X,Y}(x)) + d(\eta_{X,Y}(x), \eta(x)).
\]

[vertex] (v1) ;
[vertex] (v2) [right=of v1] ;
[vertex] (v3) [below=of v1] ;

[thick] (v1)–(v2)–(v3)–(v1);

[left] at (v1.west)$\eta_{X,Y}(x)$;
[right] at (v2.east)$\eta(x) = \eta_{Y,Z}( \eta_{X,Y}(x))$;
[below] at (v3.south)$x$;

The situation in the proof of Lemma <ref>

Therefore, it cannot be that both $d(x,\eta_{X,Y}(x))<s$ and $d(\eta_{Y,Z}(x), \eta(x))<t$ (compare Figure <ref>).
Consequently, we obtain an injective map

\[
\{x\colon d(x,\eta(x))>s+t\} \hookrightarrow \{x\colon d(x,\eta_{X,Y}(x))>s\} \sqcup \{y\colon d(y,\eta_{Y,Z}(y))>t\}
\]

and conclude

\[
|\{x\colon d(x,\eta(x))>s+t\}|\le |\{x\colon d(x,\eta_{X,Y}(x))>s\}| + |\{y\colon d(y,\eta_{Y,Z}(y))>t\}|.\qedhere
\]

However, $D_{X,Y}(t)=0$ for all $t>0$ implies $X=Y$ only under some finiteness assumtions.
For example, consider a converging sequence $(a_n)_{n\in \N} \subset \R^2$ above the diagonal with limit $a\not\in (a_n)$.
Set $X$ to consist of all elements of the sequence $\{a_n \colon n\in \N\}$.
Set $Y$ to be $X \cup \{a\}$.
Then for all $\eps > 0$ there exists $\eta \colon X\to Y$ such that $d(x, \eta(x))<\eps$ for all $x\in X$.
Therefore, $D_{X,Y}(t) = 0$ for every $t>0$, but $X\neq 0$.

Following <cit.>, we denote by $\Bar{\mathcal{B}}$ the set of persistence diagrams such that for each $\eps>0$ there are finitely many points of persistence $>\eps$.
The next lemma is an immediate consequence of <cit.>.

The bottleneck profile satisfies $D_{X, X}(t) = 0$ for all $t>0$.
Moreover, $D_{X, Y}(t) = 0$ for all $t>0$ implies $X=Y$ for $X, Y\in \bar{\mathcal{B}}$.

If $D_{X, Y}(t) = 0$ for all $t>0$, then $W_\infty(X,Y)=0$ by Lemma <ref>.
Now for $X,Y \in \bar{\mathcal{B}}$, this only happens if $X=Y$ by <cit.>.

§.§ Relation to  Wasserstein distances

We have already seen how the bottleneck profile is related to the bottleneck distance.
This is actually part of a more general result comparing it to $p$-Wasserstein metrics.

Let $X, Y$ be two persistence diagrams, let $p>0$. Then

\[
D_{X,Y}(t)  \leq \frac{1}{t^p}W_p(X,Y)^p.
\]

This follows from the Chebychev inequality (Lemma <ref>) for counting measures.
To spell out the details, estimate that for every bijection $\eta$

\begin{align*}
|\{x \colon d(x,\eta(x)) >t\}|
&= \sum \limits_{\{x \colon d(x,\eta(x)) >t\}\}} 1\\
&\leq \sum \limits_{\{x \colon d(x,\eta(x)) >t\}\}} \frac{d(x,\eta(x))^p}{t^p}\\
&\leq \sum \limits_{x} \frac{d(x,\eta(x))^p}{t^p}\\
&=\frac{1}{t^p}\sum \limits_{x} d(x,\eta(x))^p.
\end{align*}

Now choosing $\eta$ to minimize the right hand size, we have by definition of the Wasserstein distance an estimate for $D_{X,Y}$:

\[
D_{X,Y}(t) = |\{x \colon d(x,\eta(x)) >t\}| \leq \frac{1}{t^p}W_p(X,Y)^p.\qedhere
\]

This is illustrated by Figure <ref>.

An example for the relation between $D_{X,Y}$ and the Wasserstein distance. See Section <ref> for details about the underlying persistence diagrams.

Note that, consistently with Lemma <ref>, in the limit we have

\[
\left(\frac{W_p}{t}\right)^p \longrightarrow \begin{cases}
\infty & \text{if }t<W_\infty\\
1 & \text{if }t=W_\infty\\
0 & \text{if }t>W_\infty.
\end{cases}
\]

§.§ Algorithms

Recall the definition

\begin{align*}
D_{X,Y}(t) &= \inf\limits_{\eta} |\{x \colon d(x,\eta(x))>t\}|\\
&=\sup\limits_{\eta} |\{x \colon d(x,\eta(x))\le t\}|.
%\pi_f(X,Y) &= \inf\{t>0 \colon D_{X,Y}(t) < f(t)\}.
\end{align*}

The computation of $\sup\limits_{\eta} |\{x \colon d(x,\eta(x))\le t\}|$ is a version of the unweighted maximum cardinality bipartite matching problem.
First, set up the following notation (following <cit.>).
Denote by$X_0$ the off-diagonal points of $X$ and by $X_0'$ their projections to the diagonal (and analogously for $Y$).
Set $U = X_0 \cup Y_0'$ and $V = Y_0 \cup X_0'$ and consider the bipartite graph $G=(U \cup V, E)$ with $e=\{u,v\} \in E$ if either of the following holds:

* $u\in X_0, v\in Y_0$ and $d(u,v) \le t$,

* $u\in X_0, v \in X_0'$ is its projection to the diagonal and $d(u,v) \le t$,

* $v\in Y_0, u \in Y_0'$ is its projection to the diagonal and $d(u,v) \le t$,

* $u \in Y_0'$ and $v \in X_0'$.

Let $M\subset E$ be a matching of maximal cardinality.
Observe that such a matching corresponds to a bijection $\eta \colon X \to Y$ maximizing $|\{x \colon d(x,\eta(x))\le t\}|$.

To estimate the runtime of this algorithm, let $n = |X|+|Y|$.
We solve the unweighted maximum cardinality bipartite matching problem using the Hopcroft-Karp algorithm <cit.>.
While this takes $O(n^{2.5})$ in the worst case, we perform a variant which exploits the geometric nature of the setting, as suggested in <cit.>.
To build the layer graph, we do not perform DFS but a neighbor search based on a kD tree.
Deletion and insertion in this data structure take $O(\sqrt{n})$ (<cit.>).
Hence, finding the alternating paths in one layered graph takes $O(n^{1.5})$.
This needs to be done $O(\sqrt{n})$ times in order to find a matching of maximum cardinality <cit.>.
In conclusion, we find the following:

Let $X, Y$ be finite persistence diagrams and denote $n = |X|+|Y|$.
The value of the bottleneck profile $D_{X,Y}(t)$ can be computed in $O(n^2)$.

Instead of using Hopcroft-Karp, one can regard the matching problem as an LP.
For each $x \in X$ and $y\in Y$, we have a binary variable $f_{xy}$ indicating whether the edge from $x$ to $y$ is in the matching.
The coefficients (the cost of the edge) are given by

\[
c_{xy} = \begin{cases}
1 &\text{if } d(x,y) > t,\\
0 & \text{otherwise.}
\end{cases}
\]

The objective is

\begin{align*}
\text{minimize } &\sum \limits_{x,y} c_{xy}f_{xy}\\
\text{subject to }    &\forall x \in X\colon \sum\limits_y f_{xy} = 1, \;
\forall y \in Y\colon \sum\limits_x f_{xy} = 1.
\end{align*}

§ DISCRETE PROKHOROV METRICS FOR PERSISTENCE DIAGRAMS

A straight-forward discretization of the coupling characterization of the probabilistic Prokhorov metric (<ref>) gives the main notion of this section.

Given two persistence diagrams $X, Y$, consider matchings $\eta\colon X \to Y$ to define their Prokhorov distance as

\begin{align*}
\pi(X,Y) &= \inf\{t>0\colon D_{X,Y}(t)<t\}\\
&= \inf\{t>0\colon \inf\limits_{\eta\colon X\to Y} |\{x \colon d(x,\eta(x))>t\}|<t\}.
\end{align*}

In other words, we look at the intersection of the bottleneck profile with the diagonal.
Similarly, we have already seen that the bottleneck distance arises as the intersection of $D_{X,Y}$ with the horizontal axis.
This motivates the following definition, which generalizes both the bottleneck and the Prokhorov distance.

Given a fixed monotonically increasing function $f\colon \R_{\ge 0} \to \R_{\ge 0}$, define for any two PDs $X,Y$ their $f$-Prokhorov distance to be

\begin{align*}
\pi_f(X,Y) &= \inf\{t>0\colon D_{X,Y}(t)<f(t)\}\\
&= \inf\{t>0\colon \inf\limits_{\eta} |\{x\colon d(x,\eta(x))>t\}|<f(t)\} .
\end{align*}

Plugging in $f = id$ gives the Prokhorov distance, plugging in $f=0$ gives the Bottleneck metric.

Let $X=\{x\}$ and $Y=\{y\}$ both consist of one point each and assume that $d(x,y)<d(x,x')+d(y,y')$, where the prime denotes the projection to the diagonal.
That means that $x\mapsto y$ is an optimal matching.
Consequently, the bootleneck profile looks as follows:

\[
D_{X,Y}(t) = \begin{cases}
2 & \text{if } t=0,\\
1 & \text{if } 0<t\le d(x,y),\\
0 & \text{if } t> d(x,y).
\end{cases}
\]

It follows that

\[
\pi_f(X,Y) = \min(f^{-1}(1), d(x,y)).
\]

For two Prokhorov-close PDs, we require the measure of unmatched points to be small.
Points with small persistence get matched to the diagonal and thus do not blow up the Prokhorov distance.
Hence it is robust with respect to noise.

The triangle inequality follows from Lemma <ref>.

Fix a superadditive function $f\colon \R_{\ge 0} \to \R_{\ge 0}$. For persistence diagrams $X,Y,Z$, we have

\[
\pi_f(X,Z) \le \pi_f(X,Y)+\pi_f(Y,Z).
\]

In order to show that

\begin{align*}
\inf \{t>0 \colon D_{X,Z}(t)<f(t)\} \le \pi_f(X,Y) + \pi_f(Y,Z),
\end{align*}

it suffices to prove

\begin{align*}
D_{X,Z}(\pi(X,Y) + \pi(Y,Z))< f(\pi(X,Y)) + f(\pi(Y,Z));
\end{align*}

the assertion then follows by superadditivity.
Now by Lemma <ref>, we have

\[
D_{X,Z}(\pi_f(X,Y) + \pi_f(Y,Z)) \le D_{X,Y}(\pi_f(X,Y))+D_{Y,Z}(\pi_f(Y,Z)).
\]

Using the definition of the metric, we obtain

\begin{align*}
\end{align*}

as desired.

Non-negative superadditive functions are in particular monotonically increasing.
For example, convex functions are supperadditive.
A particularly simple family of such functions are polynomials with non-negative coefficients and absolute term $0$.
Piecewise linear functions constitute another set of examples.
As the symmetry is clear, we have shown:

Fix a superadditive function $f\colon \R_{\ge 0} \to \R_{\ge 0}$. The discrete $f$-Prokhorov metric is an extended pseudometric.

Just like for the bottleneck distance, we need some finiteness property for the $\pi_f$ to be a genuine metric.
Let $\Bar{\mathcal{B}}$ denote the persistence diagrams which for every $\eps>0$ have only finitely many points of persistence $>\eps$.
Then Lemma <ref> implies:

Let $f\colon \R_{\ge 0} \to \R_{\ge0}$ with $\lim \limits _{t\searrow 0} f(t) = 0$.
For $X,Y\in \Bar{\mathcal{B}}$, we have $\pi_f(X,Y) = 0$ only if $X=Y$.

If $\pi_f(X,Y) = 0$, then $D_{X,Y}(t)<f(t)$ for all $t>0$.
As the bottleneck profile is monotonically decreasing and $\lim \limits _{t\searrow 0} f(t) = 0$, this implies $D_{X,Y}(t)=0$ for all $t>0$.
By Lemma <ref>, this happens only if $X=Y$.

Our next task is to investigate how $\pi_f$ depends on the function $f$.
While from a metric point of view, we need to fix $f$, the context of data science suggests a different perspective:
For given training data (a fixed set of persistence diagrams) adjust $f$ to obtain a metric that performs well on it (e.g. in a classification problem, cf. section <ref>).

Let $f, g \colon \R_{\ge0} \to \R_{\ge 0}$ such that $f(t)\le g(t)$ for all $t\ge 0$.
Then for any two persistence diagrams $X,Y$, we have $\pi_g(X,Y) \le \pi_f(X,Y)$.

If $t>0$ satisfies $D_{X,Y}(t)<f(t)$, then also $D_{X,Y}(t)<g(t)$.

\[
\inf\{t>0\colon D_{X,Y}(t)<g(t)\} \le \inf\{t>0\colon D_{X,Y}(t)<f(t)\}
\]

and by definition $\pi_g(X,Y) \le \pi_f(X,Y)$.

For fixed diagrams, the Prokhorov metric is continuous with respect to the functions in supremum metric.

Fix two persistence diagrams $X,Y$. Let $f\colon \R_{\ge 0} \to \R_{\ge 0}$ be superadditive with $\lim \limits_{t\to 0} f(t) = 0$ but $f(t)>0$ for all $t>0$.

\[
\forall \eps > 0 \; \exists \delta > 0 \colon \forall g \colon \|f-g\|_\infty < \delta \Rightarrow |\pi_f(X,Y) - \pi_g(X,Y)| < \eps.
\]

Let $\eps > 0$ and choose $0<\delta < \min(f(\eps), 1)$.
Let $g\colon \R_{\ge 0} \to \R_{\ge0}$ with $\|f-g\|_\infty < \delta$.
If $D_{X,Y}(\pi_f(X,Y)) \neq D_{X,Y}(\pi_g(X,Y))$ then $f$ and $g$ intersect the bottleneck profile $D_{X,Y}$ at different integral values.
This means $\|f-g\|_\infty\geq 1$, which contradicts our choice of $\delta$.
Hence we have $D_{X,Y}(\pi_f(X,Y)) = k =D_{X,Y}(\pi_g(X,Y))$ for some $k$.
Note that we can write $\pi_f(X,Y) = f^{-1}(k)$ and $\pi_g(X,Y) = g^{-1}(k)$.
If $\pi_f(X,Y) \ge \pi_g(X,Y)$, then

\[
f(\eps) > \delta > \|f-g\|_\infty  \ge g(g^{-1}(k)) - f(g^{-1}(k)) = k - f(g^{-1}(k)).
\]

Equivalently, $f(\eps)+ f(g^{-1}(k)) > k$ and we can apply superadditivity to obtain

\[
f(\eps + g^{-1}(k)) > k.
\]

Using monotonicity of $f$, we have the following equivalences:

f(+ g^-1(k))    > k

+ g^-1(k)    > f^-1(k)

> f^-1(k) - g^-1(k).

Analogously, if $\pi_g(X,Y) \ge \pi_f(X,Y)$, then

\[
f(\eps) > \delta > \|f-g\|_\infty  = f(g^{-1}(k)) - g(g^{-1}(k)) \ge f(g^{-1}(k)) - k.
\]

Equivalently, $k + f(\eps) > f(g^{-1}(k))$ and we can again apply superadditivity to now obtain

\[
f(\eps + f^{-1}(k)) > f(g^{-1}(k)).
\]

Using monotonicity of $f$, we have the following equivalences:

f(+ f^-1(k))    > f(g^-1(k))

+ f^-1(k)    > g^-1(k)

> g^-1(k) - f^-1(k).

Together, we have shown that $|\pi_f(X,Y) - \pi_g(X,Y)| < \eps$, as desired.

§.§ Comparison with Wasserstein

Fix a persistence diagram $X$ and consider Wasserstein metrics and Prokhorov distances to some other diagram $Y$.
We can perturb $Y$ by adding more "noise".
More precisely, we add $k$ points whose distance to the diagonal is less than $\pi_f(X,Y)$ and denote this diagram by $Y_k$.
This does not affect the Prokhorov metric at all, while $W_p(X,Y_k)$ goes to infinity when $k$ does.
This is what we mean when we say that the Prokhorov metric is more robust with respect to noise compared to the Wasserstein metric.
In other (more mathematical) words, the identity map $\id \colon (Dgm, \pi_f)\to (Dgm, W_p)$, where $Dgm$ is the set of all persistence diagrams, is nowhere continuous[To avoid such problems, one usually restricts to a subset of $Dgm$ of diagrams with "finite $p$th moment" <cit.>  when using $p$-Wasserstein distances.].
This motivates the further investigation in this section.

Similarly to the proofs in <cit.> for the measure-theoretic variants, we can bound our metric in terms of the Wasserstein distance.
It turns out that the metrics $\pi_q = \pi_{t\mapsto t^q}$ are of special interest.

Let $p\ge 1, q>0, c>0$ and $f(t) = c \cdot t^q$. For two persistence diagrams $X,Y$ we have

\[
\pi_f(X,Y) \le W_p(X,Y)^{\frac{p}{p+q}}\cdot c^{\frac{-1}{p+q}}.
\]

Recall from Lemma <ref> that

\[
D_{X,Y}(t) = |\{x \colon d(x,\eta(x)) >t\}| \leq \frac{1}{t^p}W_p(X,Y)^p.
\]

We now want to find a suitable value of $t$ such that $D_{X,Y}(t)<c \cdot t^q$ to infer that $\pi_f(X,Y)\leq t$.
Plugging in $t = W_p(X,Y)^{\frac{p}{p+q}} \cdot c^{\frac{-1}{p+q}}$, one obtains

\begin{align*}
|\{x \colon d(x,\eta(x)) >W_p(X,Y)^{\frac{p}{p+q}}\cdot c^{\frac{-1}{p+q}}\}|
&\leq \frac{W_p(X,Y)^p}{W_p(X,Y)^{\frac{p^2}{p+q}}\cdot c^{\frac{-p}{p+q}}} \\
&= \left(W_p(X,Y)^{\frac{p}{p+q}}\cdot c^{\frac{-1}{p+q}}\right)^q.
\end{align*}

\[
\inf_\eta |\{x \colon d(x,\eta(x)) >W_p(X,Y)^{\frac{p}{p+q}}\}| \leq \left(W_p(X,Y)^{\frac{p}{p+q}}\right)^q
\]

and we conclude $\pi_f(X,Y) \le W_p(X,Y)^{\frac{p}{p+q}} \cdot c^{\frac{-1}{p+q}}$ as desired.

The map $\id \colon (Dgm, W_p) \to (Dgm, \pi_{t^q})$ is continuous.

Specializing to $c=1$ and $p\in\{1,\infty\}$ or $q\in\{0,1\}$, we obtain:

The following inequalities hold:

[3em]$0$   [3em]$1$   [3em]$q$
$1$     $d_B \le W_1$     $D_P \le \sqrt{W_1}$     $\pi_q \le W_1^{\frac{1}{1+q}}$
$\infty$    $d_B\le d_B$    $\pi\le d_B$    $\pi_q \le d_B$
$p$     $d_B \le W_p$     $\pi \le W_p^{\frac{p}{p+1}}$     $\pi_q \le W_p^{\frac{p}{p+q}}$

In particular, the Bottleneck Stability Theorem <cit.> implies stability for the new metrics:

Let $X,Y$ be finite metric spaces, fix some $k \in \N$.
Then we have

\[
\pi_f(\Dgm(PH_k(X)), \Dgm(PH_k(Y))) \le d_{GH}(X,Y),
\]

where $d_{GH}$ is the Gromov-Hausdorff distance.

$W_q(X,Y)^q \leq \pi_q(X,Y)^q (\max(d(X,Y))^q+|X|+|Y|-1) $

For an arbitrary bijection $\eta \colon X \to Y$ and let $t>0$ be such that $|\{d(x,\eta(x))>t\}| \le t^q$.
We estimate:

\begin{align*}
W_q(X,Y)^q &\leq \sum_{x} d(x,\eta(x))^q\\
&=\sum_{d(x,\eta(x))>t}d(x,\eta(x))^q+\sum_{d(x,\eta(x))\le t}d(x,\eta(x))^q\\
&\le |\{d(x,\eta(x))>t\}| \max(d(X,Y))^q+t^q |\{d(x,\eta(x))\le t\}|\\
&= |\{d(x,\eta(x))>t\}| \max(d(X,Y))^q+t^q (|X|+|Y|-|\{d(x,\eta(x))>t\}|)\\
&= |\{d(x,\eta(x))>t\}| (\max(d(X,Y))^q-t^q) + t^q (|X|+|Y|)\\
&\le t^q \max(d(X,Y))^q -t^{2q} +t^q (|X|+|Y|)
\end{align*}

Taking the infimum over all such $t$ we obtain the desired inequality

\[
W_q(X,Y)^q \le \pi_q(X,Y)^q(\max(d(X,Y))^q+|X|+|Y|-1).\qedhere
\]

Combining the two inequalities from Propositions <ref> and <ref>, we obtain a comparison for different Wasserstein metrics.

$W_q(X,Y)^q \le W_p(X,Y)^{\frac{pq}{p+q}}(\max(d(X,Y))^q+|X|+|Y|-1)$.
In particular, $W_1(X,Y) \le W_2(X,Y)^{\frac{2}{3}}(\max(d(X,Y))+|X|+|Y|-1).$

§.§ Metric and Topological Properties

Using the comparison with Wasserstein (Section <ref>) and the results from <cit.>, we address questions of convergence, separability and compactness. We run into similar issues as <cit.>, Theorems 4.20, 4.24, 4.25 and <cit.>, section 3.

The space of persistence diagrams with finite $p$th moment (for some $p$) endowed with the $c \cdot t^q$-Prokhorov metric is separable.

Let $\eps>0$, $X$ a persistence diagram and $p\ge 1$.
Let $S$ be a countable dense subset for the $p$-Wasserstein metric; this exists by <cit.>, Theorem 12.
Let $X_S\in S$ be a persistence diagram such that $W_p(X, X_S) < \eps^{\frac{p+q}{p}}\cdot c^{\frac{1}{p}}$.
Then by Proposition <ref>, we have

\[
\pi_{t^q}(X,X_S)\le W_p(X, X_S)^{\frac{p}{p+q}} \cdot c^{\frac{-1}{p+q}} < \eps^{\frac{p}{p+q}\frac{p+q}{p}} \cdot c^{\frac{-1}{p+q}} \cdot  c^{\frac{1}{p+q}} = \eps.
\]

The example of Theorem 4.20 in <cit.> shows that we indeed need some restriction on the persistence diagrams.

Recall that $\Bar{\mathcal{B}}$ denotes the persistence diagrams which for all $\eps>0$ have finitely many points of persistence $>\eps$.
The next Theorem is a consequence of <cit.>, which asserts that the bottleneck distance makes $\Bar{\mathcal{B}}$ into a polish space.

The space $\Bar{\mathcal{B}}$ endowed with the Prokhorov metric $\pi_f$ is polish.

Let $(X_n) \subset \Bar{\mathcal{B}}$ be a Cauchy sequence with respect to the Prokhorov metric $\pi_f$.
Let $\eps>0$.
For $\eps<f^{-1}(1)$, the inequality $\pi_f(X_m,X_n)<\eps$ implies

\[
D_{X_m, X_n}(\eps)<f(\eps)\leq f(f^{-1}(1)) = 1
\]

by monotonicity of $f$.
Hence, bottleneck and Prokhorov metric coincide for diagrams which are close to each other.
In particular, $X_n$ is a Cauchy sequence with respect to the Bottleneck distance.
By completeness of $\Bar{\mathcal{B}}$ with this metric, there is a limit diagram $X\in \Bar{\mathcal{B}}$ to which the sequence converges.
As $X_n$ approaches $X$ in the Bottleneck distance, it also does in the Prokhorov metric.

Now for separability, consider a subset $A\subset\Bar{\mathcal{B}}$ which is dense with respect to the bottleneck distance.
Let $X\in \Bar{\mathcal{B}}$ and $0<\eps \le f^{-1}(1)$.
Then there exists $Y \in A$ with $W_\infty(X,Y) < \eps$.
But by the previous considerations, $\pi_f(X,Y) = W_\infty(X,Y)$ by the choice of $\eps$.
Therefore, $A$ is dense in $\Bar{\mathcal{B}}$ with respect to $\pi_f$ as well.

§.§ Algorithms

Now we will provide an algorithm to compute $\pi_f(X,Y)$ for monotonically increasing functions $f$.
In this case, there is always just one intersection of the graphs of $D$ and $f$ so that we can find its location by bisection.

Let $f\colon [0,\infty)\to [0,\infty)$ be monotonically increasing. Assume that the values and preimages of $f$ can be computed in $O(1)$.
Then $\pi_f(X,Y)$ can be computed in $O(n^2 \log(n))$.

First, observe that $\pi_f(X,Y) \in \{d(x,y)\colon x\in X, y\in Y\} \cup f^{-1}(\N_{\le |X|+|Y|}) = T_1$.
To perform a binary search, we sort the possible solutions $T_1$ as a preprocess, which has runtime complexity $O(n^2 \log(n))$.
In each iteration of the binary search we pick the median $t\in T_i$.
Next we compute the value of the bottleneck profile $D_{X,Y}(t)$ using Proposition <ref>, taking $O(n^2)$.
Then we compute $f(t)$, which by assumption takes $O(1)$.
Now if $D_{X,Y}(t)>f(t)$ set $T_{i+1}$ to be the right half, if $D_{X,Y}(t)\le f(t)$ set $T_{i+1}$ to be the left half of $T_i$.
InputPersistence diagrams $X,Y$; function $f$Output$\pi_f(X,Y)$
$T = \{d(x,y)\colon x\in X, y\in Y\} \cup f^{-1}(\N_{\le |X|+|Y|})$
sort $T$
$L = 0; R = \text{length}(T)$
$m = \lfloor \frac{R+L}{2}\rfloor$
$t = T[m]$

The binary search to compute $\pi_f(X,Y)$

Hence we obtain a runtime of $O(n^2 \log{n})$ for the binary search as well.

In particular, if one uses a more efficient geometric data structure to improve the runtime of the matching algorithm, the sorting preprocessing dominates the runtime.
Compare <cit.>, Theorem 3.2 and the preceding discussion therein for more details and possible improvements of the runtime complexity.
Please refer to Appendix <ref> for details about our implementation and its availability.

§ EXPERIMENTS

A simple application of the bottleneck profile, based on simple synthetic persistence diagrams, was already presented in Example <ref>. We now turn to a little more sophisticated data sets to illustrate the usage and advantages of the Prokhorov distance.
In particular, we consider persistence diagrams that actually arise in applications of TDA.
We use the library <cit.> for standard machine learning algorithms (in particular $K$Neighbors) and <cit.> for the images.
For the Bottleneck and Wasserstein metrics we use <cit.> and <cit.>, respectively.
To score the different metrics, we use KMedoids clustering with the "build" initialization <cit.>, <cit.> as well as the silhouette score <cit.>.
All computations were carried out on a laptop with an Intel i5-8265U CPU with $1.60$ GHz and $8$ GB memory.

§.§ Three Circles

Consider three different shapes in $\mathbb{R}^2$: a) a big circle ($r=6$), b) a big ($r=6$) and a medium circle ($r=4$), c) a big ($r=6$), a medium ($r=4$) and small circle ($r=2$).
We take five samples with noise from each shape according to Table <ref>.

shape      number of circles     radii     samples     noise                                 colour in the figures
a          1     6           120        uniform from $[-0.2,0.2]^2$     blue

b          2     6, 4        300        uniform from $[-0.23,0.23]^2$     red

c          3     6, 4, 2     120        uniform from $[-0.2,0.2]^2$     green

The three shapes: one two and three circles.

One two and three noisy circles and their PDs for the first persistent homology.

For each point cloud we compute the first persistent homology modules of it alpha complex filtration and represent them as PDs (see Figure <ref>).
We can look at the averaged $D$-function for each pair of shapes (Figure <ref>).
After careful inspection of this figure, we come up with the choice of $f(t) = t^3\cdot 20^t$ to separate three bottleneck profiles in a most efficient way:
Between around $0.55$ and $0.65$, the averaged bottleneck profiles involving shape c) with the small circle decrease, while the one comparing a) and b) stays constant.
Intersecting with a function in this interval will provide a good choice for the Prokhorov distance:
It puts the two and three circles closest to each other and one and three circles the farthest apart.

The averaged bottleneck profile for the three circles.

Now we want to compare the Bottleneck, Prokhorov and Wasserstein distances.

MDS plots of the dataset in Section <ref>.

Distance matrices of the dataset in Section <ref>.

The bottleneck distance between shapes a) and both b) and c) is roughly the same.
This distance does not take the presence of the additional small circle in shape c).
By blowing up the sample size and the noise in shape b), the Wasserstein distance from a) and c) to it are artificially blown up (Figures <ref> and <ref>).
The Prokhorov distance is built to avoid these pitfalls and nicely captures the geometry of the setting.
The MDS plot for Prokhorov agrees with our intuition and places b) between a) and c) (Figures <ref>).

§.§ 3D Segmemtation

We adapt an example from <cit.> and <cit.>, which is based on the dataset <cit.>.
The task is to classify 3D-meshes based on the persistence diagrams of certain functions defined on them.
We divide the diagrams into 171 training points and 5529 test points.
We evaluate the classification accuracy using i) KNeighbors, ii) KMedoids, iii) silhouette score.
The results of classification are presented in the Tables <ref>, <ref> and <ref>, respectively.
On the training data, we select the function $f$ to yield a high score.
We restrict ourselves to quadratic polynomials to avoid overfitting.

The $K$-neighbors classifier is a number between zero and one indicating the percentage of points classified correctly.
While the Prokhorov distance is quickly computed, it yields the worst score.

The $K$-medoids score is a number between zero and one indicating the percentage of points that have the same label as the medoid they are assigned to.
Here, Prokhorov scores best out of all metrics.
The 1-Wasserstein distance is faster to compute, but yields a worse score.

The silhouette score is a number between minus one and one.
A high score indicates a good clustering: points are more similar to other points from their cluster than from other clusters.
A score close to zero, as the ones we find in our experiments, indicates that there is little distinction beween the clusters.
On this very low baseline, Prokhorov still performs slightly better than the other distances.

Prokhorov     Bottleneck     1-Wasserstein     2-Wasserstein
test score     $0.7579$    $0.7603$    $0.803$    $0.7930$
computation time $[s]$     $ 85.51$    $163.2$    $102.3$    $258.6$

3DSegmentation, 5% training data, $K$Neighbors Classifier ($K=5$) scores, $f(t)=17t^2+3t$.

Prokhorov     Bottleneck     1-Wasserstein     2-Wasserstein
test score     $0.4733$    $0.3180$    $0.4063$    $0.2572$
computation time $[s]$     $ 2281 $    $ 2297 $    $ 1528 $    $ 2579 $

3DSegmentation, 5% training data, K-Medoids scores, $f(t)=5t^2+13t$.

Prokhorov     Bottleneck     1-Wasserstein     2-Wasserstein
test score     $ 0.004614 $    $-0.006915$    $-0.03000$    $-0.02572$
computation time $[s]$     $ 1864 $    $ 2055 $    $ 1629 $    $ 3126 $

3DSegmentation, 5% training data, silhouette scores, $f(t)=24 t^2$.

§.§ Synthetic Dataset

Finally, we consider the dataset introduced by <cit.>.
It contains six shape classes: A sphere, a torus, clusters, clusters within clusters, a circle and the unit cube.
From each class take 25 samples of 500 points.
Then add two levels of Gaussian noise ($\eta = 0.05, 0.1$) and the zeroth and first persistent homology of the Vietoris-Rips filtration are computed.
We compute the distance matrices and evaluate them based on the $K$-medoids and silhouette scores.
The results are displayed in Table <ref>.
We find that Prokhorov performs on par with Bottleneck and worse than Wasserstein.
Prokhorov takes similarly long as 1-Wasserstein; Bottleneck is faster and 2-Wasserstein is slower.

dim                      $\eta$                      $f(t)$                    Prokhorov     Bottleneck     1-Wasserstein     2-Wasserstein
KMedoids score             3*$0$     3*$0.05$     3*$7t$     $0.6083 $    $ 0.6 $    $ 0.975 $    $ 0.9417 $
Silhouette score                                                                                $ 0.2091 $    $ 0.1747 $    $ 0.6697 $    $ 0.5045 $
computation time $[s]$                                                                                $ 160.4 $    $ 52.63 $    $ 270.5 $    $ 1043 $
KMedoids score             3*$0$     3*$0.1$     3*$7t$      $0.6083 $    $ 0.6 $    $ 0.975 $    $ 0.9417 $
Silhouette score                                                                                $ 0.2091 $    $ 0.1747 $    $ 0.6697 $    $ 0.5045 $
computation time $[s]$                                                                                $ 143.4 $    $ 51.79 $    $ 265.2 $    $ 994.9 $
KMedoids score             3*$1$     3*$0.05$     3*$7t$    $0.825 $    $ 0.825 $    $ 1.0 $    $ 1.0 $
Silhouette score                                                                                $ 0.3785 $    $ 0.3785 $    $ 0.6365 $    $ 0.5998 $
computation time $[s]$                                                                                $ 23.02 $    $ 18.62 $    $ 20.74 $    $ 98.43 $
KMedoids score             3*$1$     3*$0.1$     3*$3t^2 + 5 t$     $0.6833 $    $ 0.6833 $    $ 0.9417 $    $ 0.95 $
Silhouette score                                                                                $ 0.2875 $    $ 0.2882 $    $ 0.4835 $    $ 0.4468 $
computation time $[s]$                                                                                $ 24.24 $    $ 19.84 $    $ 24.20 $    $ 127.6 $

Classification scores for the synthetic dataset.

§.§ Discussion

In order to explain the differences in the computation time, we note the size of the persistence diagrams in the various settings:

3D-Segmentation     Synthetic data      Synthetic data      Synthetic data     Synthetic data

$H_0$, $\eta = 0.05$     $H_0$, $\eta = 0.1$      $H_1$, $\eta = 0.05$     $H_1$, $\eta = 0.1$
Mean size              $ 11.84 $    $ 500 $    $ 500 $    $ 177.7 $    $ 189.9 $
standard deviation     $ 4.893 $     $ 0 $    $ 0 $    $ 40.53 $    $ 38.84 $

Cardinalities of the persistence diagrams for the considered experiments.

By inspecting Table <ref> wee see that the 3D segmentation dataset contains way smaller diagrams, on which the Prokhorov metric seems to perform well, both in terms of runtime and score.
On the bigger diagrams from the synthetic dataset, the Wasserstein metrics yield the highest scores.
Prokhorov barely outperforms Bottleneck in the scores, if at all, and at the cost of higher runtimes.
The difference in the computation time is caused by the evaluation of $f(t)$, which is the only difference between the Bottleneck and Prokhorov implementations.

Bottleneck and Prokhorov do not work well on zero-dimensional PDs.
There, every class is born at time zero, hence the PD is intrinsically one-dimensional and points are matched in linear order.
The bottleneck distance is less meaningful in this setting.
Moreover, the Prokhorov (and even more the Bottleneck) distance do not take points matched over a small distance into account.
This is a consequence of being designed to be robust against noise.
However, this data can actually contain meaningful information, which is picked up by the Wasserstein distances.
This is a possible explanation for the fact that Wasserstein yields better scores in the synthetic dataset.

Hence, the Prokhorov metric works best on rather small diagrams and runs fastest with simple (e. g. linear) parameter functions $f$.
Even then, one needs to take the additional time for tuning the parameter $f$ into account.

§ DISCUSSION AND OUTLOOK

Summarizing the results from the previous section, we find that the Prokhorov metric is well-suited for small persistence diagrams.
Large scale computations can be improved by the technique of entropic regularization from the theory of optimal transport <cit.>.
As the classical Prokhorov metric admits an optimal transport characterization, our discrete variant might be tractable using similar techniques.

A major aspect of the importance of the Bottleneck distance is its algebraic formulation in terms of interleavings.
We hope that this theory generalizes to incorporate the family of Prokhorov metrics.
An algebraic formulation would also provide a perspective on generalizations to multiparameter persistence.

Morally, stability theorems should involve related metrics on the input point cloud and on the persistence diagram side.
This motivates to investigate Prokhorov-type distances for point clouds in $\R^n$.
Such distances might be useful throughout data science.

§.§.§ Ackknowledgement

This work was in part supported by he Centre for Topological Data Analysis, EPSRC grant EP/R018472/1, and by the Dioscuri program initiated by the Max Planck Society, jointly managed with the National Science Centre (Poland), and mutually funded by the Polish Ministry of Science and Higher Education and the German Federal Ministry of Education and Research.
We thank Gesine Reinert for valuable discussions.

Missing 'biblatex' package
The bibliography requires the 'biblatex' package.

Microtome Publishing, Brookline, MA

issn1532-4435; 1533-7928/e
journaltitleJ. Mach. Learn. Res.
noteId/No 8
titlePersistence images: a stable vector representation of persistent homology

abstractWe study distributions of persistent homology barcodes associated to taking subsamples of a fixed size from metric measure spaces. We show that such distributions provide robust invariants of metric measure spaces, and illustrate their use in hypothesis testing and providing confidence intervals for topological data analysis.
titleRobust Statistics, Hypothesis Testing, and Confidence Intervals for Persistent Homology on Metric Measure Spaces
omputer Science - Computational Geometry

abstractPersistence modules are a central algebraic object arising in topological data analysis. The notion of interleaving provides a natural way to measure distances between persistence modules. We consider various classes of persistence modules, including many of those that have been previously studied, and describe the relationships between them. In the cases where these classes are sets, interleaving distance induces a topology. We undertake a systematic study the resulting topological spaces and their basic topological properties.
issn2367-1726, 2367-1734
journaltitleJ Appl. and Comput. Topology
titleTopological Spaces of Persistence Modules and Their Properties
55N99; 54A05; 18A25,Mathematics - Algebraic Topology,Mathematics - General Topology

title3D shape segmentation using TDA
typeIPython notebook

journaltitleACM Transactions on Graphics (Proc. SIGGRAPH)
titleA Benchmark for 3D Mesh Segmentation

Wiley Online Library

booktitleComputer Graphics Forum
titleGromov-Hausdorff stable signatures for shapes using persistence

abstractWe propose the use of persistent homology in a supervised way. We believe homological persistence is fundamentally not about decomposition theorems but a central role is played by a choice of metrics. Choosing a pseudometric between persistent vector spaces leads to a model. Fitting this model is what we believe supervised homological persistence is. We develop theory behind constructing such models and we give evidence of the usefulness of this approach in concrete data analysis tasks.
titleMetrics and Stabilization in One Parameter Persistence
hachólski and Riihimäki - 2020 - Metrics and stabilization in one parameter persist.pdf
ematics - Algebraic Topology

World Scientific, Singapore

issn0219-4988; 1793-6829/e
journaltitleJ. Algebra Appl.
noteId/No 1550066
titleDecomposition of pointwise finite-dimensional persistence modules.

GUDHI Editorial Board

booktitleGUDHI User and Reference Manual
titlePersistence representations

Providence, RI: American Mathematical Society (AMS)

titleComputational topology. An introduction.
pagesxii + 241

Springer US, New York, NY

issn0178-4617; 1432-0541/e
titleGeometry helps in bottleneck matching and related problems.

GUDHI Editorial Board

booktitleGUDHI User and Reference Manual
titleBottleneck distance

Springer Science & Business Media

titleMetric structures for Riemannian and non-Riemannian spaces

journaltitleArXiv preprint
titleOn Choosing and Bounding Probability Metrics

abstractThe present paper shows how to construct a maximum matching in a bipartite graph with n vertices and m edges in a number of computation steps proportional to $(m + n)\sqrt n $.
journaltitleSIAM J. Comput.
notePublisher: Society for Industrial and Applied Mathematics
titleAn $n^{5/2}$ Algorithm for Maximum Matchings in Bipartite Graphs

abstractMatplotlib is a 2D graphics package used for Python for application development, interactive scripting, and publication-quality image generation across user interfaces and operating systems.
journaltitleComputing in Science & Engineering
titleMatplotlib: A 2D graphics environment

Association for Computing Machinery (ACM), New York, NY

journaltitleACM J. Exp. Algorithm.
noteId/No 1.4
titleGeometry helps to compare persistence diagrams.

abstractPersistence diagrams (PDs) are now routinely used to summarize the underlying topology of complex data. Despite several appealing properties, incorporating PDs in learning pipelines can be challenging because their natural geometry is not Hilbertian. Indeed, this was recently exempliﬁed in a string of papers which show that the simple task of averaging a few PDs can be computationally prohibitive. We propose in this article a tractable framework to carry out standard tasks on PDs at scale, notably evaluating distances, estimating barycenters and performing clustering. This framework builds upon a reformulation of PD metrics as optimal transport (OT) problems. Doing so, we can exploit recent computational advances: the OT problem on a planar grid, when regularized with entropy, is convex can be solved in linear time using the Sinkhorn algorithm and convolutions. This results in scalable computations that can stream on GPUs. We demonstrate the efﬁciency of our approach by carrying out clustering with diagrams metrics on several thousands of PDs, a scale never seen before in the literature.
journaltitlearXiv:1805.08331 [cs, stat]
titleLarge Scale computation of Means and Clusters for Persistence Diagrams using Optimal Transport
arge Scale computation of Means and Clusters for .pdf:/home/nihell/snap/zotero-snap/common/Zotero/storage/RLR58JRU/Lacombe et al. - 2018 - Large Scale computation of Means and Clusters for .pdf:application/pdf
ine Learning,Statistics - Machine Learning,Computer Science - Computational Geometry

IOP Publishing, Bristol

journaltitleInverse Probl.
noteId/No 124007
titleProbability measures on the space of persistence diagrams

abstractOptimal transport (OT) theory can be informally described using the words of the French mathematician Gaspard Monge (1746-1818): A worker with a shovel in hand has to move a large pile of sand lying on a construction site. The goal of the worker is to erect with all that sand a target pile with a prescribed shape (for example, that of a giant sand castle). Naturally, the worker wishes to minimize her total effort, quantified for instance as the total distance or time spent carrying shovelfuls of sand. Mathematicians interested in OT cast that problem as that of comparing two probability distributions, two different piles of sand of the same volume. They consider all of the many possible ways to morph, transport or reshape the first pile into the second, and associate a "global" cost to every such transport, using the "local" consideration of how much it costs to move a grain of sand from one place to another. Recent years have witnessed the spread of OT in several fields, thanks to the emergence of approximate solvers that can scale to sizes and dimensions that are relevant to data sciences. Thanks to this newfound scalability, OT is being increasingly used to unlock various problems in imaging sciences (such as color or texture processing), computer vision and graphics (for shape manipulation) or machine learning (for regression, classification and density fitting). This short book reviews OT with a bias toward numerical methods and their applications in data sciences, and sheds lights on the theoretical properties of OT that make it particularly useful for some of these applications.
titleComputational Optimal Transport
uturi - 2020 - Computational Optimal Transport.pdf
ine Learning

journaltitleJournal of Machine Learning Research
titleScikit-learn: Machine Learning in Python

Chichester etc.: John Wiley &| Sons Ltd.

journaltitleWiley Ser. Probab. Math. Stat.
titleProbability metrics and the stability of stochastic models
pagesxiv + 494

abstractA new graphical display is proposed for partitioning techniques. Each cluster is represented by a so-called silhouette, which is based on the comparison of its tightness and separation. This silhouette shows which objects lie well within their cluster, and which ones are merely somewhere in between clusters. The entire clustering is displayed by combining the silhouettes into a single plot, allowing an appreciation of the relative quality of the clusters and an overview of the data configuration. The average silhouette width provides an evaluation of clustering validity, and might be used to select an ‘appropriate’ number of clusters.
journaltitleJournal of Computational and Applied Mathematics
titleSilhouettes: A graphical aid to the interpretation and validation of cluster analysis
Graphical display,cluster analysis,clustering validity,classification

shorttitlekmedoids 0.1.5-dev
titlekmedoids 0.1.5-dev : k-Medoids clustering with the FasterPAM algorithm
Scientific/Engineering - Artificial Intelligence,Software Development - Libraries

abstractClustering non-Euclidean data is diﬃcult, and one of the most used algorithms besides hierarchical clustering is the popular algorithm Partitioning Around Medoids (PAM), also simply referred to as k-medoids clustering.
journaltitlearXiv:2008.05171 [cs, stat]
shorttitleFast and Eager k-Medoids Clustering
titleFast and Eager k-Medoids Clustering: O(k) Runtime Improvement of the PAM, CLARA, and CLARANS Algorithms
chubert and Rousseeuw - 2020 - Fast and Eager k-Medoids Clustering O(k) Runtime .pdf:application/pdf
ine Learning,Statistics - Machine Learning

Providence, RI

American Mathematical Society (AMS)

annotationMSC2010: 90-02 = Research exposition (monographs, survey articles) pertaining to operations research and mathematical programming MSC2010: 28D05 = Measure-preserving transformations MSC2010: 35B65 = Smoothness and regularity of solutions to PDEs MSC2010: 35J60 = Nonlinear elliptic equations MSC2010: 49N90 = Applications of optimal control and differential games MSC2010: 49Q20 = Variational problems in a geometric measure-theoretic setting MSC2010: 90B20 = Traffic problems in operations research
journaltitleGraduate Studies in Mathematics
pagetotalxvi + 370 p
titleTopics in Optimal Transportation

Missing 'biblatex' package
The bibliography requires the 'biblatex' package.

Microtome Publishing, Brookline, MA

issn1532-4435; 1533-7928/e
journaltitleJ. Mach. Learn. Res.
noteId/No 8
titlePersistence images: a stable vector representation of persistent homology

abstractWe study distributions of persistent homology barcodes associated to taking subsamples of a fixed size from metric measure spaces. We show that such distributions provide robust invariants of metric measure spaces, and illustrate their use in hypothesis testing and providing confidence intervals for topological data analysis.
titleRobust Statistics, Hypothesis Testing, and Confidence Intervals for Persistent Homology on Metric Measure Spaces
omputer Science - Computational Geometry

abstractPersistence modules are a central algebraic object arising in topological data analysis. The notion of interleaving provides a natural way to measure distances between persistence modules. We consider various classes of persistence modules, including many of those that have been previously studied, and describe the relationships between them. In the cases where these classes are sets, interleaving distance induces a topology. We undertake a systematic study the resulting topological spaces and their basic topological properties.
issn2367-1726, 2367-1734
journaltitleJ Appl. and Comput. Topology
titleTopological Spaces of Persistence Modules and Their Properties
55N99; 54A05; 18A25,Mathematics - Algebraic Topology,Mathematics - General Topology

title3D shape segmentation using TDA
typeIPython notebook

journaltitleACM Transactions on Graphics (Proc. SIGGRAPH)
titleA Benchmark for 3D Mesh Segmentation

Wiley Online Library

booktitleComputer Graphics Forum
titleGromov-Hausdorff stable signatures for shapes using persistence

abstractWe propose the use of persistent homology in a supervised way. We believe homological persistence is fundamentally not about decomposition theorems but a central role is played by a choice of metrics. Choosing a pseudometric between persistent vector spaces leads to a model. Fitting this model is what we believe supervised homological persistence is. We develop theory behind constructing such models and we give evidence of the usefulness of this approach in concrete data analysis tasks.
titleMetrics and Stabilization in One Parameter Persistence
hachólski and Riihimäki - 2020 - Metrics and stabilization in one parameter persist.pdf
ematics - Algebraic Topology

World Scientific, Singapore

issn0219-4988; 1793-6829/e
journaltitleJ. Algebra Appl.
noteId/No 1550066
titleDecomposition of pointwise finite-dimensional persistence modules.

GUDHI Editorial Board

booktitleGUDHI User and Reference Manual
titlePersistence representations

Providence, RI: American Mathematical Society (AMS)

titleComputational topology. An introduction.
pagesxii + 241

Springer US, New York, NY

issn0178-4617; 1432-0541/e
titleGeometry helps in bottleneck matching and related problems.

GUDHI Editorial Board

booktitleGUDHI User and Reference Manual
titleBottleneck distance

Springer Science & Business Media

titleMetric structures for Riemannian and non-Riemannian spaces

journaltitleArXiv preprint
titleOn Choosing and Bounding Probability Metrics

abstractThe present paper shows how to construct a maximum matching in a bipartite graph with n vertices and m edges in a number of computation steps proportional to $(m + n)\sqrt n $.
journaltitleSIAM J. Comput.
notePublisher: Society for Industrial and Applied Mathematics
titleAn $n^{5/2}$ Algorithm for Maximum Matchings in Bipartite Graphs

abstractMatplotlib is a 2D graphics package used for Python for application development, interactive scripting, and publication-quality image generation across user interfaces and operating systems.
journaltitleComputing in Science & Engineering
titleMatplotlib: A 2D graphics environment

Association for Computing Machinery (ACM), New York, NY

journaltitleACM J. Exp. Algorithm.
noteId/No 1.4
titleGeometry helps to compare persistence diagrams.

abstractPersistence diagrams (PDs) are now routinely used to summarize the underlying topology of complex data. Despite several appealing properties, incorporating PDs in learning pipelines can be challenging because their natural geometry is not Hilbertian. Indeed, this was recently exempliﬁed in a string of papers which show that the simple task of averaging a few PDs can be computationally prohibitive. We propose in this article a tractable framework to carry out standard tasks on PDs at scale, notably evaluating distances, estimating barycenters and performing clustering. This framework builds upon a reformulation of PD metrics as optimal transport (OT) problems. Doing so, we can exploit recent computational advances: the OT problem on a planar grid, when regularized with entropy, is convex can be solved in linear time using the Sinkhorn algorithm and convolutions. This results in scalable computations that can stream on GPUs. We demonstrate the efﬁciency of our approach by carrying out clustering with diagrams metrics on several thousands of PDs, a scale never seen before in the literature.
journaltitlearXiv:1805.08331 [cs, stat]
titleLarge Scale computation of Means and Clusters for Persistence Diagrams using Optimal Transport
arge Scale computation of Means and Clusters for .pdf:/home/nihell/snap/zotero-snap/common/Zotero/storage/RLR58JRU/Lacombe et al. - 2018 - Large Scale computation of Means and Clusters for .pdf:application/pdf
ine Learning,Statistics - Machine Learning,Computer Science - Computational Geometry

IOP Publishing, Bristol

journaltitleInverse Probl.
noteId/No 124007
titleProbability measures on the space of persistence diagrams

abstractOptimal transport (OT) theory can be informally described using the words of the French mathematician Gaspard Monge (1746-1818): A worker with a shovel in hand has to move a large pile of sand lying on a construction site. The goal of the worker is to erect with all that sand a target pile with a prescribed shape (for example, that of a giant sand castle). Naturally, the worker wishes to minimize her total effort, quantified for instance as the total distance or time spent carrying shovelfuls of sand. Mathematicians interested in OT cast that problem as that of comparing two probability distributions, two different piles of sand of the same volume. They consider all of the many possible ways to morph, transport or reshape the first pile into the second, and associate a "global" cost to every such transport, using the "local" consideration of how much it costs to move a grain of sand from one place to another. Recent years have witnessed the spread of OT in several fields, thanks to the emergence of approximate solvers that can scale to sizes and dimensions that are relevant to data sciences. Thanks to this newfound scalability, OT is being increasingly used to unlock various problems in imaging sciences (such as color or texture processing), computer vision and graphics (for shape manipulation) or machine learning (for regression, classification and density fitting). This short book reviews OT with a bias toward numerical methods and their applications in data sciences, and sheds lights on the theoretical properties of OT that make it particularly useful for some of these applications.
titleComputational Optimal Transport
uturi - 2020 - Computational Optimal Transport.pdf
ine Learning

journaltitleJournal of Machine Learning Research
titleScikit-learn: Machine Learning in Python

Chichester etc.: John Wiley &| Sons Ltd.

journaltitleWiley Ser. Probab. Math. Stat.
titleProbability metrics and the stability of stochastic models
pagesxiv + 494

abstractA new graphical display is proposed for partitioning techniques. Each cluster is represented by a so-called silhouette, which is based on the comparison of its tightness and separation. This silhouette shows which objects lie well within their cluster, and which ones are merely somewhere in between clusters. The entire clustering is displayed by combining the silhouettes into a single plot, allowing an appreciation of the relative quality of the clusters and an overview of the data configuration. The average silhouette width provides an evaluation of clustering validity, and might be used to select an ‘appropriate’ number of clusters.
journaltitleJournal of Computational and Applied Mathematics
titleSilhouettes: A graphical aid to the interpretation and validation of cluster analysis
Graphical display,cluster analysis,clustering validity,classification

shorttitlekmedoids 0.1.5-dev
titlekmedoids 0.1.5-dev : k-Medoids clustering with the FasterPAM algorithm
Scientific/Engineering - Artificial Intelligence,Software Development - Libraries

abstractClustering non-Euclidean data is diﬃcult, and one of the most used algorithms besides hierarchical clustering is the popular algorithm Partitioning Around Medoids (PAM), also simply referred to as k-medoids clustering.
journaltitlearXiv:2008.05171 [cs, stat]
shorttitleFast and Eager k-Medoids Clustering
titleFast and Eager k-Medoids Clustering: O(k) Runtime Improvement of the PAM, CLARA, and CLARANS Algorithms
chubert and Rousseeuw - 2020 - Fast and Eager k-Medoids Clustering O(k) Runtime .pdf:application/pdf
ine Learning,Statistics - Machine Learning

Providence, RI

American Mathematical Society (AMS)

annotationMSC2010: 90-02 = Research exposition (monographs, survey articles) pertaining to operations research and mathematical programming MSC2010: 28D05 = Measure-preserving transformations MSC2010: 35B65 = Smoothness and regularity of solutions to PDEs MSC2010: 35J60 = Nonlinear elliptic equations MSC2010: 49N90 = Applications of optimal control and differential games MSC2010: 49Q20 = Variational problems in a geometric measure-theoretic setting MSC2010: 90B20 = Traffic problems in operations research
journaltitleGraduate Studies in Mathematics
pagetotalxvi + 370 p
titleTopics in Optimal Transportation
Our Implementation
We provide an implementation as a part of a custom gudhi fork at <https://github.com/nihell/gudhi-devel>.
It is a modification of the implementation of the Bottleneck distance <cit.>
Let us first illustrate how to use it before we come to runtime considerations.
The algorithm is implemented in C++ and comes with Python bindings.

It asks for three inputs: ,  and .
The two diagrams need to be presented as 2D numpy arrays.
The third parameter is a 1D numpy array representing the coefficients of a polynomial to be used as $f$.
Note that the zeroth entry needs to be zero in order to obtain a metric, compare Lemma <ref>.
However, setting the polynomial to be a constant integer one recovers the values of $D_{X,Y}$, which is a feature.
In the technical details, our approach follows <cit.>, which follows <cit.>.

In addition, we also add the Prokhorov metric to <cit.>, allowing for parallel computations of distance matrices and integration with .