#  $k$-Median clustering under discrete Fréchet and Hausdorff distances
††thanks: The authors would like to thank Pankaj Agarwal, Kamesh Munagala and
anonymous reviewers for helpful discussions and feedback. A shorter version of
this paper will appear in SoCG 2020.

Abhinandan Nath111Part of the work was done when the author was a graduate
student at Duke University
Mentor Graphics, Fremont, USA
<EMAIL_ADDRESS>Erin Taylor
Duke University, Durham, USA
<EMAIL_ADDRESS>
###### Abstract

We give the first near-linear time $(1+\varepsilon)$-approximation algorithm
for $k$-median clustering of polygonal trajectories under the discrete Fréchet
distance, and the first polynomial time $(1+\varepsilon)$-approximation
algorithm for $k$-median clustering of finite point sets under the Hausdorff
distance, provided the cluster centers, ambient dimension, and $k$ are bounded
by a constant. The main technique is a general framework for solving
clustering problems where the cluster centers are restricted to come from a
_simpler_ metric space. We precisely characterize conditions on the simpler
metric space of the cluster centers that allow faster
$(1+\varepsilon)$-approximations for the $k$-median problem. We also show that
the $k$-median problem under Hausdorff distance is NP-Hard.

##  1 Introduction

We study the $k$-median problem for an arbitrary metric space
${\mathcal{X}}=(X,\mathtt{d})$, where the cluster centers are restricted to
come from a (possibly infinite) subset $C\subseteq X$. We call it the
$(k,C)$-median problem. We prove general conditions on the structure of $C$
that allow us to get efficient $(1+\varepsilon)$-approximation algorithms for
the $(k,C)$-median problem for any $\varepsilon>0$. As applications of our
framework, we give $(1+\varepsilon)$-approximation algorithms for the metric
space defined over polygonal trajectories and finite point sets in
$\mathbb{R}^{d}$ under the discrete Fréchet and Hausdorff distance
respectively, where the cluster centers have bounded complexity. For
trajectories, our algorithm runs in near-linear time in the number of input
points (Theorem 13) and is exponentially faster than the previous best
algorithm ([10], Theorem 11). For point sets, ours is the first
$(1+\varepsilon)$-approximation algorithm that runs in time polynomial in the
number of input points (Theorem 15) for bounded dimensions and cluster
complexity. Our results are summarized in Table 1. We also show that the
$k$-median problem under Hausdorff distance problem is NP-Hard.

Table 1: Our results for $(1+\varepsilon)$-approximate $k$-median for $n$ input trajectories/point sets in $\mathbb{R}^{d}$, each having at most $m$ points. Each cluster center can have at most $l$ points. While stating running times, we assume $k,l$ and $d$ are constants independent of $n,m$, and $\tilde{O}$ hides logarithmic factors in $n,m$. Metric space |  Our result |  Previous best
---|---|---
Polygonal traj., discrete Fréchet |  $\tilde{O}\left(nm\right)$ |  $\tilde{O}\left(n^{dkl+1}m\right)$ [10]
Point sets, Hausdorff |  $nm^{O(dl)}$ |  –

The $k$-median problem has been very widely studied. We are given a set $P$ of
$n$ elements from a metric space. The goal is now to select $k$ _centers_ so
that the sum of the distances of each point to the nearest cluster center is
minimized. In the simplest setting, $P\subseteq\mathbb{R}^{d}$ under the
Euclidean metric. In this paper, each individual element of $P$ is itself a
collection of points in $\mathbb{R}^{d}$, e.g., a curve traced by a moving
object, or a point cloud. Since the objective of clustering is to group
_similar_ objects into the same cluster and to summarize each cluster using
its cluster center, it is important that a meaningful distance function is
used to compare two input elements. In our work, we look at the widely used
discrete Fréchet and Hausdorff distances for trajectories and point sets
respectively.

Trajectories can model a variety of systems that change with time. As such,
trajectory data is being collected at enormous scales. As a first step,
clustering is hugely important in understanding and summarizing the data. It
involves partitioning a set of trajectories into clusters of similar
trajectories, and computing a representative trajectory (a center) per
cluster. It can be viewed as a compression scheme for large trajectory
datasets, effectively performing non-linear dimension reduction. If the
centers have low complexity, this representation can reduce uncertainty and
noise found in individual trajectories. Information provided by the set of
centers is useful for trajectory analysis applications such as similarity
search and anomaly detection [31]. The Hausdorff distance is another widely
used shape-based distance [8, 23]. In shape matching applications, we may want
to cluster similar shapes into one group (where a shape is represented by a
point cloud).

Figure 1: The red, green and blue trajectories are similar to each other and
form a single cluster. If the cluster center trajectory is restricted to have
four vertices, it will look like the black trajectory at the bottom. However
if the center is unrestricted, it will contain a lot of vertices inside the
four noticeable noisy _clumps_ of vertices of the input trajectories, thereby
overfitting to the input.

The $k$-median problem is hard to solve exactly, even in Euclidean space [16,
29]. There is a long line of work on both constant factor and
$(1+\varepsilon)$-approximations, with varying running time dependence on $k$
and the ambient dimension (if applicable). Many of these algorithms require
the underlying metric space to have bounded doubling dimension (e.g., see
[1]). However, it can be shown that both the discrete Fréchet and Hausdorff
distances do not have doubling dimension bounded by a constant (Appendix A).
We circumvent this problem by considering cluster centers from a somewhat
_simpler_ metric space compared to the input metric space. For trajectories
and point sets, we restrict the centers to have low complexity, i.e., a
bounded number of points. This approach has been used before ([9, 10, 13]). It
has the added benefit of preventing the cluster center from overfitting to the
elements of its cluster. This is crucial, since real-life measurements are
noisy and error-prone, and without any restrictions the cluster center can
inherit noise and high complexity from the input (see Fig. 1). As another
example, in clustering financial time-series data using Hausdorff distance
[6], frequent intra-day fluctuations may not be useful in capturing long-term
trends, and we want to avoid them by retricting the cluster centers’
complexity. However, our work differs from previous approaches in that we
precisely characterize general conditions on the simpler metric space for the
cluster centers, which leads to faster $(1+\varepsilon)$-approximation
algorithms for the $k$-median problem for the discrete Fréchet and Hausdorff
distances.

#### Problem definition.

A _metric space_ ${\mathcal{X}}=(X,\mathtt{d})$ consists of a set $X$ and a
distance function $\mathtt{d}:X\times X\rightarrow\mathbb{R}_{\geq 0}$ that
satisfies the following properties : (i) $\mathtt{d}(x,x)=0$ for all $x\in X$;
(ii) $\mathtt{d}(x,y)=\mathtt{d}(y,x)$ for all $x,y\in X$; and (iii)
$\mathtt{d}(x,z)\leq\mathtt{d}(x,y)+\mathtt{d}(y,z)$ for all $x,y,z\in X$.

Given subsets $P,C\subseteq X$, the _$(k,C)$ -median_ problem is to compute a
set $C^{\prime}\subseteq C$ of $k$ _center_ points that minimizes

$\displaystyle\sum_{p\in P}\mathtt{d}(p,C^{\prime}),$

where $\mathtt{d}(p,C^{\prime})=\min_{c\in C^{\prime}}\mathtt{d}(p,c)$. Here
$P$ is finite, but $C$ need not be.

Let $T^{l}$ (resp. $U^{l}$) be the set of all trajectories (resp. point sets)
in $\mathbb{R}^{d}$, where each trajectory (resp. point set) has at most $l$
points. Thus, $T=\bigcup_{l>0}T^{l}$ and $U=\bigcup_{l>0}U^{l}$ are the set of
all trajectories and finite point sets in $\mathbb{R}^{d}$ respectively. As
special cases of the $(k,C)$-median problem, we discuss the $(k,T^{l})$-median
and the $(k,U^{l})$-median problems for the metric spaces
${\mathcal{T}}=(T,\mathtt{d}_{F})$ and ${\mathcal{U}}=(U,\mathtt{d}_{H})$
respectively, i.e., each center trajectory or point set can have at most $l$
points. Here, $\mathtt{d}_{F}$ and $\mathtt{d}_{H}$ denote the discrete
Fréchet and Hausdorff distances respectively.

#### Challenges and ideas.

As mentioned before, both the discrete Fréchet and Hausdorff metrics do not
have low doubling dimension, so a number of previous techniques do not
directly apply to our setting to yield efficient algorithms. One approach
would be to embed these metrics into other metric spaces. Backurs and
Sidiropoulos [4] give an embedding of the Hausdorff metric over point sets of
size $s$ in $d$-dimensional Euclidean space, into $l_{\infty}^{s^{O(s+d)}}$
with distortion $s^{O(s+d)}$. However, both the distortion and the resultant
dimension are too high for many applications. It is not known if the Fréchet
distance can be embedded into an $l_{p}$ space using finite dimension.

We circumvent the problem by restricting the cluster centers to come from a
subset $C$ of the original space $X$, namely trajectories and point sets
defined by a bounded number of points each. We show that if every metric ball
in $C$ can be _covered_ by a small number of metric balls of a fixed smaller
radius, then we can use a sampling-based algorithm similar to the one in
Ackermann et al. [1]; we make this precise by introducing the notion of
_coverability_ of $C$. We crucially show that the centers of the balls in the
_cover_ can be arbitrary, and need not come from $C$. This is more general
than $C$ having bounded doubling dimension. This allows us to approximate the
optimal $(1,C)$-median of $P$ using the $(1,C)$-median of a constant sized
random sample of $P$, allowing us to use the framework of [1].

It is not known how to efficiently compute the optimal $(1,C)$-median under
the discrete Fréchet and Hausdorff distances. However, we show that all we
need is to compute a constant number of candidate centers in time independent
of the size of the input, at least one of which is a good approximation to the
optimal $(1,C)$-median of the input. We show how to compute these candidates
for a _coverable_ set $C$. Then, we can apply the sampling technique of [1]
and recursively use this property to find $k$ centers that approximate the
cost of the $(k,C)$-median optimal solution. Although our work heavily relies
on the framework of Ackermann et al. [1], it is a significant improvement from
existing work on clustering under the discrete Fréchet distance, and the first
such result for clustering under the Hausdorff distance.

#### Previous work.

Trajectory clustering has a lot of practical applications, e.g., discovering
frequent movement patterns in trajectory data. As such, there has been work on
trajectory clustering [17, 22, 32], and possibly computing a representative
trajectory for each cluster. Many proposed algorithms and models have no
provable performance gaurantees and are experimental in nature.

Driemel et al. [13] started the rigorous study of clustering trajectories
under the continuous Fréchet distance under the classic $k$-clustering
objectives. However, they only deal with 1D-trajectories. They introduce the
$(k,l)$-clustering problem, i.e., clustering trajectories with $k$ cluster
centers such that each center can have at most $l$ points. For
$1D$-trajectories, they give $(1+\varepsilon)$-approximation algorithms for
both the $(k,l)$-center and $(k,l)$-median problem that run in near-linear
time for constant $\varepsilon,k,l$. Buchin et al. [9] study the
$(k,l)$-center clustering problem for trajectories in $\mathbb{R}^{d}$ under
the discrete and continuous Fréchet distances, and give both upper and lower
bounds. The most closely related work to ours is the one by Buchin, Driemel
and Struijs [10], where they give algorithms for the $(k,l)$-median problem
under discrete Fréchet; however their running times are much slower (see Table
1). They also show that the $1$-median problem under discrete Fréchet distance
is NP-Hard, and W[1]-Hard in the number of input trajectories.

On the other hand, clustering under Hausdorff distance has received much less
attention. There is work on hierarchical clustering of financial time series
data using Hausdorff distance [6]. Chen et al. [11] use the DBSCAN algorithm
[15] while Qu et al. [30] use spectral clustering for trajectories and using
the Hausdorff distance. We are not aware of any theoretical analysis for
$k$-median clustering of point sets under the Hausdorff distance.

In general metric spaces, a polynomial time
$(1+\sqrt{3}+\varepsilon)$-approximation algorithm to the $k$-median problem
exists [28], whereas no polynomial time algorithm can achieve an approximation
ratio less than $(1+2/e)$ unless NP $\subseteq$ DTIME $[n^{O(\log\log n)}]$
[24]. For the Euclidean $k$-median problem in $d$ dimensions, Guruswami and
Indyk [18] showed that there is no PTAS for $k$-median if both $k$ and $d$ are
part of the input. Arora et al. [3] gave the first PTAS when $d$ is fixed,
whose running time was subsequently improved in [25]. Kumar et al. [27] gave a
$(1+\varepsilon)$-approximate algorithm with runtime
$2^{(k/\varepsilon)^{O(1)}}dn$, this was extended by Ackermann et al. [1] to
those metric spaces for which the optimal $1$-median can be approximated using
a constant-sized random sample; this holds true for doubling metric spaces.
There are coreset-based approaches with running times linear in $n$ and either
exponential in $d$ and polynomial in $k$ [20, 19], or vice versa [5].
Recently, Cohen-Addad et al. [12] gave a PTAS for $k$-median in low-
dimensional Euclidean and minor-free metrics using local search.

##  2 Preliminaries, definitions and an overview

We formally define the discrete Fréchet and Hausdorff distances. We also
define two properties on $C$ and $\mathtt{d}$ which allow us to design
efficient clustering algorithms. Finally we give an overview of our algorithm.

#### Discrete Fréchet and Hausdorff distance.

Consider two finite sets $\zeta_{1}$ and $\zeta_{2}$. A _correspondence_
${\mathcal{C}}$ between $\zeta_{1}$ and $\zeta_{2}$ is a subset of
$\zeta_{1}\times\zeta_{2}$ such that every element of $\zeta_{1}$ and
$\zeta_{2}$ appears in at least one pair in ${\mathcal{C}}$. For
$\zeta_{1},\zeta_{2}\subseteq\mathbb{R}^{d}$, the _Hausdorff distance_ [21] is
defined as

$\displaystyle\mathtt{d}_{H}(\zeta_{1},\zeta_{2})=\min_{{\mathcal{C}}\in\Xi(\zeta_{1},\zeta_{2})}\max_{(p,q)\in{\mathcal{C}}}\left\lVert
p-q\right\rVert,$

where $\left\lVert.\right\rVert$ denotes the $l_{2}$ norm, and
$\Xi(\zeta_{1},\zeta_{2})$ is the set of all possible correspondences between
$\zeta_{1}$ and $\zeta_{2}$.

A trajectory $\gamma$ is a finite sequence of points $\langle
p_{1},p_{2}\ldots\rangle$ in $\mathbb{R}^{d}$. A correspodence can be defined
for a pair of trajectories $\gamma_{1}=\langle p_{1},p_{2},\ldots\rangle$ and
$\gamma_{2}=\langle q_{1},q_{2},\ldots\rangle$ by treating each trajectory as
a point sequence; such a correspondence is said to be _monotonic_ if it also
respects the ordering of points in the trajectories, i.e., if
$(p_{i_{1}},q_{j_{1}}),(p_{i_{2}},q_{j_{2}})\in{\mathcal{C}}$, then $i_{2}\geq
i_{1}\Rightarrow j_{2}\geq j_{1}$. The _discrete Fréchet distance_ [14]
between $\gamma_{1}$ and $\gamma_{2}$ is defined as

$\displaystyle\mathtt{d}_{F}(\gamma_{1},\gamma_{2})=\min_{{\mathcal{C}}\in\Xi_{M}(\gamma_{1},\gamma_{2})}\max_{(p,q)\in{\mathcal{C}}}\left\lVert
p-q\right\rVert,$

where $\Xi_{M}(\gamma_{1},\gamma_{2})$ is the set of all monotone
correspondences between $\gamma_{1}$ and $\gamma_{2}$. Our algorithms cluster
in the metric space defined by the discrete Fréchet and Hausdorff distance.

#### Strong and weak sampling properties.

We define two properties that make efficient clustering algorithms possible.
These are generalizations of the _strong_ and _weak_ sampling properties
defined by Ackermann et al.(see Theorem 1.1 and Property 4.1 in [1]). The
major difference is that the cluster centers are restricted to a subset $C$ of
the metric space $X$. These properties allow fast approximation of the
$(1,C)$-median using only a constant sized random sample of the input. We
later show how to get an efficient $(k,C)$-median algorithm using the fast
$(1,C)$-median algorithm as a subroutine. We denote the optimal $(1,C)$ median
of any set $P\subseteq X$ by $c_{P}$.

######  Definition 1 (Strong sampling property).

Let $0<\varepsilon,\delta<1$ be arbitrary. $(X,C,\mathtt{d})$ is said to
satisfy the strong sampling property for $\varepsilon,\delta$ iff

(i) For any finite $P\subseteq X$, $c_{P}$ can be computed in time depending
only on $|P|$.

(ii) There exists a positive integer $m_{\delta,\varepsilon}$ depending on
$\delta,\varepsilon$ such that for any $P\subseteq X$, the optimal
$(1,C)$-median $c_{S}$ of a uniform random multiset $S\subseteq P$ of size
$m_{\delta,\varepsilon}$ satisfies

$\displaystyle\Pr\left[\sum_{p\in
P}\mathtt{d}(p,c_{S})\leq(1+\varepsilon)\sum_{p\in
P}\mathtt{d}(p,c_{P})\right]\geq 1-\delta.$

The strong sampling property characterizes those instances in which the
optimal $(1,C)$-median of a constant-sized random sample is a good
approximation to the optimal $(1,C)$-median of the whole set. However, in many
cases it is impossible to efficiently solve the $(1,C)$-median exactly (e.g.,
when $C,X=\mathbb{R}^{d}$ and $\mathtt{d}$ is the Euclidean metric). This is
also true for the discrete Fréchet and Hausdorff distances for polygonal
trajectories and finite point sets respectively. The following definition
becomes helpful then.

######  Definition 2 (Weak sampling property).

Let $0<\varepsilon,\delta<1$ be arbitrary. $(X,C,\mathtt{d})$ is said to
satisfy the weak sampling property for $\varepsilon,\delta$ iff there exist
positive integers $m_{\delta,\varepsilon}$ and $t_{\delta,\varepsilon}$
depending on $\delta,\varepsilon$ such that for any $P\subseteq X$ and a
uniform random multiset $S\subseteq X$ of size $m_{\delta,\varepsilon}$, there
exists a set $\Gamma(S)\subseteq C$ of size $t_{\delta,\varepsilon}$ that
satisfies

$\displaystyle\Pr\left[\exists c\in\Gamma(S)\mid\sum_{p\in
P}\mathtt{d}(p,c)\leq(1+\varepsilon)\sum_{p\in
P}\mathtt{d}(p,c_{P})\right]\geq 1-\delta.$

Furthermore, $\Gamma(S)$ can be computed in time depending on
$\delta,\varepsilon,|S|$ but independent of $|P|$.

The weak sampling property characterizes those instances in which it is
possible to generate a constant number of candidate cluster centers in time
independent of the size of the input set, and at least one of which is
guaranteed to be a good approximation to the optimal $1$-center of the input
set. We later show that the discrete Fréchet and Hausdorff distances satisfy
the weak sampling property. We use $m_{\delta,\varepsilon}$ to denote the size
of the random sample for both the strong and weak sampling properties.

#### Algorithm overview.

We give an overview of our algorithm, denoted Cluster, in Fig. 2. It is
similar to the algorithm Cluster from [1], with the small but crucial
difference being that the set of candidate cluster centers $\overline{C}_{S}$
comes from $C$; this also changes how the candidates are generated. We show
that with a careful choice of $C$, our instance satisfies one of the sampling
properties (Definitions 1, 2), and we can apply the framework of Ackermann et
al. [1].

The algorithm takes as input the set of points $\overline{P}\subseteq P$ that
are yet to be assigned cluster centers, the number of cluster centers
$\overline{k}$ still to be computed, and the centers $\overline{C}\subseteq C$
already computed. It returns the final set of cluster centers. To solve the
$(k,C)$-median problem for $P$, we call Cluster$(P,k,\\{\\})$ with values of
$\alpha,m_{\delta,\varepsilon},F$ that we will specify later.

Briefly, the algorithm has two phases. In the pruning phase no new centers are
added. Rather, the set $N$ containing half of the points of $\overline{P}$
closest to $\overline{C}$ are removed from $\overline{P}$, and the algorithm
is called recursively on $\overline{P}\setminus N$. In the sampling phase, new
centers are added. The algorithm first samples a uniformly random multiset $S$
of $\overline{P}$ of size $\frac{2}{\alpha}m_{\delta,\varepsilon}$ for some
constant $\alpha$ and $m_{\delta,\varepsilon}$ to be defined later. Then for
each subset $S^{\prime}\subset S$ of size $m_{\delta,\varepsilon}$, a set of
candidate centers $F(S^{\prime})$ is generated; the function $F$ varies
depending on certain conditions satisfied by $(X,C,\mathtt{d})$; in particular
$F(S^{\prime})=\\{c_{S^{\prime}}\\}$ for the strong sampling property, and
$F(S^{\prime})=\\{\Gamma(S^{\prime})\\}$ for the weak sampling property. Each
candidate center is in turn added to $\overline{C}$, and the algorithm is run
recursively. Finally, the solution with the lowest cost is returned.

|  $\textsc{Cluster}(\overline{P},\overline{k},\overline{C})$:
---
input: Point set $\overline{P}$, remaining number of centers $\overline{k}$,
computed centers $\overline{C}$
if $\overline{k}=0$: return $\overline{C}$
else:
if $\overline{k}\geq|\overline{P}|$: return $\overline{C}\cup\overline{P}$
else:
/* Pruning phase */
$N\leftarrow$ set of $\frac{1}{2}|\overline{P}|$ minimal points
$p\in\overline{P}$ w.r.t $\mathtt{d}(p,\overline{C})$
$C^{*}\leftarrow\textsc{Cluster}(\overline{P}\setminus
N,\overline{k},\overline{C})$
/* Sampling phase */
$S\leftarrow$ uniform random multisubset of $\overline{P}$ of size
$\frac{2}{\alpha}m_{\delta,\varepsilon}$
$\overline{C}_{S}\leftarrow\bigcup_{S^{\prime}\subset
S,|S^{\prime}|=m_{\delta,\varepsilon}}F(S^{\prime})$
for all $\overline{c}\in\overline{C}_{S}$ :
$C^{\overline{c}}\leftarrow\textsc{Cluster}(\overline{P},\overline{k}-1,\overline{C}\cup\\{\overline{c}\\})$
return $C^{\overline{c}}$ or $C^{*}$ with the lowest cost
Figure 2: Algorithm Cluster computes a $(k,C)$-median clustering.

The rest of the paper is organized as follows. In Section 3, we show that for
$(X,C,\mathtt{d})$ satisfying the strong and weak sampling properties, the
algorithm Cluster (using the appropriate $\alpha,m_{\delta,\varepsilon},F$)
computes a $(1+\varepsilon)$-approximation to the optimal $(k,C)$-median. In
Section 4, we prove sufficient conditions on $C$ for the sampling properties
to hold. In Section 5, we give clustering algorithms for the discrete Fréchet
and Hausdorff distances using the framework developed in previous sections. In
Section 6, we show that the $k$-median problem for Hausdorff distance is NP-
Hard.

##  3 Clustering via sampling

We show that the Cluster algorithm (Fig. 2) computes an approximate
$(k,C)$-median for instances satisfying the sampling properties and for
appropriate $\alpha,m_{\delta,\varepsilon}$ and $F$, i.e., if we use
$F(S^{\prime})=\\{c_{S^{\prime}}\\}$ for the strong sampling property and
$F(S^{\prime})=\Gamma(S^{\prime})$ for the weak sampling property. The
analysis closely follows that of Ackermann et al. [1], so we leave detailed
proofs to Appendix B.

The _superset sampling lemma_ (Appendix B.1, Lemma 17) shows how to draw a
uniform random multiset from $P^{\prime}\subseteq P$ while only knowing $P$
(without explicitly knowing $P^{\prime}$), provided $P^{\prime}$ contains a
constant fraction of the points of $P$. Using this lemma and the strong and
weak sampling properties, we have the following. We leave the proof of the
_superset sampling lemma_ and the following lemma in Appendices B.1 and B.2
respectively.

######  Lemma 3.

Let $\alpha<\frac{1}{4k}$ be an arbitrary positive constant. Suppose
$(X,C,\mathtt{d})$ satisfies the strong or weak sampling property (Definitions
1, 2) for some $\varepsilon,\delta\in(0,1)$. Given $P\subseteq X$, algorithm
Cluster run with input $(P,k,\\{\\})$ and appropriate $F$ computes a set
$\tilde{C}\subseteq C$ of size $k$ such that

$\displaystyle\Pr\left[\sum_{p\in P}\mathtt{d}(p,\tilde{C})\leq(1+8\alpha
k^{2})(1+\varepsilon)\sum_{p\in
P}\mathtt{d}(p,C^{*})\right]\geq\left(\frac{1-\delta}{5}\right)^{k},$

where $C^{*}$ is an optimal solution to the $(k,C)$-median problem for $P$.

#### Running time.

We characterize the running time of Cluster in terms of $\mathtt{d}$, $C$ and
$F$. Let $h(C)$ denote the maximum time required to compute
$\mathop{\mathrm{arg\,min}}_{y\in C}\mathtt{d}(x,y)$ for any $x\in X$, and let
$t(C)$ denote the maximum time required to compute $\mathtt{d}(x,y)$ for any
$x\in X,y\in C$. Finally, let $w(m)=\max_{S\subseteq X,|S|=m}|F(S)|$, and let
$f(m)$ be the maximum number of operations needed to compute $F(S)$ for any
$S\subseteq X$ of size $m$, where computing $\mathtt{d}$ between points in $X$
and $C$, and computing the closest point in $C$ to any point in $X$ count as
one operation each. The proof is similar to the running time analysis from
[1], and is given in Appendix B.3.

######  Lemma 4.

Suppose $(X,C,\mathtt{d})$ satisfies the strong or weak sampling properties
(Definitions 1 and 2) for some $\epsilon,\delta\in(0,1)$. Given $P\subseteq X$
containing $n$ points, algorithm Cluster runs in time

$\displaystyle n\cdot
2^{O\left(km_{\delta,\varepsilon}\log\left(\frac{1}{\alpha}m_{\delta,\varepsilon}\right)\right)}\cdot\left(w(m_{\delta,\varepsilon})\cdot
f(m_{\delta,\varepsilon})\right)^{O(k)}\cdot(h(C)+t(C)).$

By setting $\alpha=\frac{\varepsilon}{8k^{2}}$, the approximation factor in
Lemma 3 becomes $(1+3\varepsilon)$. Moreover, the error probability can be
made arbitrarily small by running the Cluster algorithm $2^{\Theta(k)}$ times
and taking the minimum cost solution, without changing the asymptotic running
time. We thus get the following. Note that $w(m_{\delta,\varepsilon})$ takes
on values $1$ and $t_{\delta,\varepsilon}$ for the strong and weak sampling
properties respectively.

######  Theorem 5.

Suppose $(X,C,\mathtt{d})$ satisfies the strong sampling property (Definition
1) for some $\varepsilon,\delta\in(0,1)$. Further, suppose $c_{S}$ can be
computed in $a(m_{\delta,\varepsilon})$ operations, where computing
$\mathtt{d}$ between points in $X$ and $C$, and computing the closest point in
$C$ to any point in $X$ count as one operation each. Given $P\subseteq X$
having $n$ points and $k\in\mathbb{N}$, with probability $\geq 1-\delta$, a
$(1+3\varepsilon)$-approximate solution to the $(k,C)$-median problem for $P$
can be computed in time

$\displaystyle n\cdot
2^{O\left(km_{\delta,\varepsilon}\log\left(\frac{k}{\varepsilon}m_{\delta,\varepsilon}\right)\right)}\cdot
a(m_{\delta,\varepsilon})^{O(k)}\cdot(h(C)+t(C)).$

######  Theorem 6.

Suppose $(X,C,\mathtt{d})$ satisfies the weak sampling property (Definition 2)
for some $\varepsilon,\delta\in(0,1)$. Further, suppose $\Gamma(S)$ can be
computed in $b(m_{\delta,\varepsilon})$ operations, where computing
$\mathtt{d}$ between points in $X$ and $C$, and computing the closest point in
$C$ to any point in $X$ count as one operation each. Given $P\subseteq X$
having $n$ points and $k\in\mathbb{N}$, with probability $\geq 1-\delta$, a
$(1+3\varepsilon)$-approximate solution to the $(k,C)$-median problem for $P$
can be computed in time

$\displaystyle n\cdot
2^{O\left(km_{\delta,\varepsilon}\log\left(\frac{k}{\varepsilon}m_{\delta,\varepsilon}\right)\right)}\cdot\left(t_{\delta,\varepsilon}\cdot
b(m_{\delta,\varepsilon})\right)^{O(k)}\cdot(h(C)+t(C)).$

##  4 Covering metric spaces

We specify sufficient conditions on $C$ for the strong and weak sampling
properties to hold. These conditions characterize how well can certain subsets
of $C$ be _covered_ using a small number of sets.

Let ${\mathcal{X}}=(X,\mathtt{d})$ be a metric space. Given $x\in X$, let
${\mathcal{B}}_{\mathtt{d}}(x,r)=\\{x^{\prime}\in
X\mid\mathtt{d}(x,x^{\prime})\leq r\\}$ denote the ball of radius $r$ (under
$\mathtt{d}$) centered at $x$; we will drop the subscript $\mathtt{d}$ if it
is clear from the context. An _r-cover_ of a subset $X^{\prime}\subseteq X$
for some $r>0$ is a set $Y\subseteq X$ such that
$X^{\prime}\subseteq\bigcup_{y\in Y}{\mathcal{B}}(y,r).$ Note that the
elements of $Y$ need not be in $X^{\prime}$. Also note that if $Y$ is an
$r$-cover for $X^{\prime}$, it is also an $r$-cover for any subset of
$X^{\prime}$.

A subset $Y\subseteq X$ is said to be _g-coverable_ for some non-decreasing
function $g:\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}_{\geq 0}$ iff for all
$y\in Y$ and $r>r^{\prime}>0$, there exists an $r^{\prime}$-cover of
${\mathcal{B}}(y,r)\cap Y$ of size at most $g(r/r^{\prime})$. Note that if $Y$
is $g$-coverable then any subset $Y^{\prime}\subseteq Y$ is also
$g$-coverable.

Intuitively, if $C$ has a small cover, then for any metric ball in $C$ there
exists a small set of points (not necessarily from $C$), termed the cover,
such that the distance from any point in the ball to a point in the cover is
smaller than the radius of the ball.

#### Sufficient conditions for the strong sampling property.

The following theorem gives sufficient conditions for the strong sampling
property to hold in terms of coverability of $C$. The proof is similar to
Lemma 3.4 of [1] but has been adapted to our setting, see Appendix C for more
details.

######  Theorem 7.

If $C$ is $g$-coverable, and for any $P\subseteq X$, $c_{P}$ can be computed
in time depending only on $|P|$, then $(X,C,\mathtt{d})$ satisfies the strong
sampling property (Definition 1) for any $\varepsilon,\delta\in(0,1)$. Here,
the constant $m_{\delta,\varepsilon}=m_{\delta,\varepsilon,g}$ also depends on
$g$.

From Theorems 5 and 7, we get the following.

######  Corollary 8.

Suppose $C$ is $g$-coverable and the optimal $(1,C)$-median of any subset of
$X$ can be computed in time depending on the size of the subset. Let
$\varepsilon,\delta\in(0,1)$. Given $P\subseteq X$ having $n$ points and
$k\in\mathbb{N}$, with probability $\geq 1-\delta$, a
$(1+3\varepsilon)$-approximate solution to the $(k,C)$-median problem for $P$
can be computed in time

$\displaystyle n\cdot
2^{O\left(km_{\delta,\varepsilon,g}\log\left(\frac{k}{\varepsilon}m_{\delta,\varepsilon}\right)\right)}\cdot
a(m_{\delta,\varepsilon,g})^{O(k)}\cdot(h(C)+t(C)),$

where $m_{\delta,\varepsilon,g}$ is a constant depending only on
$\varepsilon,\delta,g$, and $a(m)$ is the number of operations needed to
compute the optimal $(1,C)$-median of $m$ points in $X$, where computing
$\mathtt{d}$ between points in $X$ and $C$, and computing the closest point in
$C$ to any point in $X$ count as one operation each.

#### Sufficient conditions for the weak sampling property.

We give sufficient conditions for the weak sampling property to hold in terms
of coverability of $C$.

For $Y\subseteq X$, let $\theta_{Y}\left(\frac{r}{r^{\prime}}\right)$ be the
number of operations required to compute an $r^{\prime}$-cover of
${\mathcal{B}}(y,r)\cap Y$ for any $y\in Y$ (if such a cover exists) where
computing $\mathtt{d}$ between points in $X$ and $C$, and computing the
closest point in $C$ to any point in $X$ count as one operation each (we
assume that the number of operations can be expressed in terms of
$\frac{r}{r^{\prime}}$).

The following lemma will be helpful, and states that if $C$ has a small cover
and if we have a good estimate of the _cost_ of the optimal $(1,C)$-median,
then we can construct a small set of points in $C$ such that at least one of
them is a good approximation to the optimal $(1,C)$-median. Further, this can
be done in time independent of $|P|$. Both of these properties are necessary
for the weak sampling property (Definition 2).

######  Lemma 9.

Let $P\subseteq X$. Suppose $C$ is $g$-coverable. Then given $a,b$ such that
$a\leq\frac{1}{|P|}\sum_{p\in P}\mathtt{d}(c_{P},p)\leq b$, we can compute a
set $Q\subseteq C$ of size $O(g\left(\frac{4b}{\varepsilon\delta a}\right))$
such that

$\displaystyle Pr\left[\exists q\in Q\mid\sum_{p\in
P}\mathtt{d}(p,q)\leq(1+\varepsilon)\sum_{p\in
P}\mathtt{d}(p,c_{P})\right]\geq 1-\delta.$

Further, $Q$ can be computed in
$O\left(\theta_{C}\left(\frac{4b}{\varepsilon\delta
a}\right)+g\left(\frac{4b}{\varepsilon\delta a}\right)\right)$ operations,
where computing $\mathtt{d}$ between points in $X$ and $C$, and computing the
closest point in $C$ to any point in $X$ count as one operation each.

###### Proof.

For any $x\in X$, we define $x^{\prime}=\mathop{\mathrm{arg\,min}}_{y\in
C}\mathtt{d}(x,y)$. Consider a point $q\in P$ chosen uniformly at random. By
Markov’s inequality, $\mathtt{d}(q,c_{P})\leq\frac{1}{\delta|P|}\sum_{p\in
P}\mathtt{d}(p,c_{P})$ with probability $\geq 1-\delta$. In such a case,

$\displaystyle\mathtt{d}(q^{\prime},c_{P})\leq\mathtt{d}(q^{\prime},q)+\mathtt{d}(q,c_{p})\leq
2\mathtt{d}(q,c_{P})\leq\frac{2}{\delta|P|}\sum_{p\in P}\mathtt{d}(p,c_{P}).$

Thus, $c_{P}\in{\mathcal{B}}(q^{\prime},\frac{2b}{\delta})$ with probability
at least $1-\delta$.

Let $C^{\prime}$ be an $\left(\frac{\varepsilon a}{2}\right)$-cover of
${\mathcal{B}}(q^{\prime},\frac{2b}{\delta})\cap C$. Since $C$ is
$g$-coverable, $|C^{\prime}|=g\left(\frac{4b}{\varepsilon\delta a}\right)$. We
will argue that $Q=\\{q^{\prime}\\}\cup\\{c^{\prime}\mid c\in C^{\prime}\\}$
is the required solution. Let $x=\mathop{\mathrm{arg\,min}}_{y\in
C^{\prime}}\mathtt{d}(y,c_{P})$. Since $C^{\prime}$ is an
$\left(\frac{\varepsilon a}{2}\right)$-cover,
$\mathtt{d}(x,c_{P})\leq\frac{\varepsilon a}{2}$. Also,
$\mathtt{d}(x,x^{\prime})\leq\mathtt{d}(x,c_{P})\leq\frac{\varepsilon a}{2}$.
Thus,
$\mathtt{d}(x^{\prime},c_{P})\leq\mathtt{d}(x,x^{\prime})+\mathtt{d}(x,c_{P})\leq\varepsilon
a$. We then have

$\displaystyle\sum_{p\in P}\mathtt{d}(p,x^{\prime})$
$\displaystyle\leq\sum_{p\in
P}\left(\mathtt{d}(p,c_{P})+\mathtt{d}(x^{\prime},c_{P})\right)\leq\left(\sum_{p\in
P}\mathtt{d}(p,c_{P})\right)+\varepsilon a|P|$
$\displaystyle\leq(1+\varepsilon)\sum_{p\in P}\mathtt{d}(p,c_{P}).$

Computing $C^{\prime}$ takes $\theta_{C}\left(\frac{4b}{\varepsilon\delta
a}\right)$ operations. Computing the output set takes
$|C^{\prime}|+1=O\left(g\left(\frac{4b}{\varepsilon\delta a}\right)\right)$
operations. ∎

The next theorem shows that with a small random sample of $P$, one of two
things can happen. Either one of the samples is close to an approximate
$(1,C)$-median, or we can approximate the cost of the optimal $(1,C)$-median
in time independent of $|P|$. This along with Lemma 9 shows that the weak
sampling property holds if $C$ is $g$-coverable. The proof is inspired by the
proof of Theorem 1 in [26], and it also shows how to compute $\Gamma$ for the
weak sampling property (Definition 2).

######  Theorem 10.

If $C$ is $g$-coverable, then $(X,C,\mathtt{d})$ satisfies the weak sampling
property (Definition 2) for $0<\varepsilon<\frac{4}{9}$ and
$\left(1-\frac{5}{18}\varepsilon\right)<\delta<1$. Further, the constants
$m_{\delta,\varepsilon}=1+\frac{4}{\varepsilon}$ and
$t_{\delta,\varepsilon}=O\left(g\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)\right)$,
and the number of operations needed to compute $\Gamma(S)$ is
$O\left(\frac{1}{\varepsilon}+\theta_{C}\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)+g\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)\right)$,
where $\delta_{1}=\frac{\varepsilon}{2}-\frac{9}{5}(1-\delta)$; and computing
$\mathtt{d}$ between points in $X$ and $C$, and computing the closest point in
$C$ to any point in $X$ count as one operation each.

###### Proof.

Let $\varepsilon_{1}=\frac{\varepsilon}{4}$ and
$\bar{r}=\frac{1}{|P|}\sum_{p\in P}\mathtt{d}(p,c_{P})$. Also, let
$x^{\prime}=\mathop{\mathrm{arg\,min}}_{y\in C}\mathtt{d}(x,y)$ for any $x\in
X$.

Let $Q\subseteq P$ be a uniform random multiset of size
$\frac{1}{\varepsilon_{1}}$, and $q\in P$ be another point chosen uniformly at
random. We will show that $Q\cup\\{q\\}$ plays the role of $S$ in Definition
2.

Using Markov’s inequality and union bound, we have

$\displaystyle\Pr\left[\mathtt{d}(q,c_{P})>\frac{\bar{r}}{2\varepsilon_{1}^{2}}\right]<2\varepsilon_{1}^{2}\text{
and }\Pr\left[\exists p\in
Q\mid\mathtt{d}(p,c_{P})>\frac{\bar{r}}{2\varepsilon_{1}^{2}}\right]<\left(\tfrac{1}{\varepsilon_{1}}\right)2\varepsilon_{1}^{2}=2\varepsilon_{1}.$

Thus with probability $\geq 1-2\varepsilon_{1}-2\varepsilon_{1}^{2}$, $q$ and
$Q$ are in
${\mathcal{B}}\left(c_{P},\frac{\bar{r}}{2\varepsilon_{1}^{2}}\right)$. We
assume that this event happens. Now, by definition of $q^{\prime}$, we have
$\mathtt{d}(q,q^{\prime})\leq\mathtt{d}(q,c_{P})$. Hence,

$\displaystyle\mathtt{d}(q^{\prime},c_{P})\leq\mathtt{d}(q,q^{\prime})+\mathtt{d}(q,c_{P})\leq
2\mathtt{d}(q,c_{P})\leq\frac{\bar{r}}{\varepsilon_{1}^{2}}.$

Let
${\mathcal{B}}_{1}={\mathcal{B}}\left(c_{P},\frac{\bar{r}}{\varepsilon_{1}^{2}}\right)$,
${\mathcal{B}}_{2}={\mathcal{B}}\left(q^{\prime},\varepsilon_{1}\bar{r}\right)$
and $P^{\prime}=P\cap{\mathcal{B}}_{1}$. Then,
$q^{\prime}\in{\mathcal{B}}_{1}$ and $Q\subseteq P^{\prime}$. We consider two
cases now.

Case 1: $P^{\prime}$ has at least $2\varepsilon_{1}|P^{\prime}|$ points
outside ${\mathcal{B}}_{2}$. For any $p\in Q$, the probability $p$ is outside
${\mathcal{B}}_{2}$ is $2\varepsilon_{1}$. Thus, with probability at least
$2\varepsilon_{1}$, there exists $p\in Q$ such that
$\mathtt{d}(p,q^{\prime})\geq\varepsilon_{1}\bar{r}$ and hence $\sum_{p\in
Q}\mathtt{d}(p,q^{\prime})\geq\varepsilon_{1}\bar{r}$. Also,
$\mathtt{d}(p,q^{\prime})\leq\frac{2\bar{r}}{\varepsilon_{1}^{2}}$ for any
$p\in Q$. Hence, $\sum_{p\in
Q}\mathtt{d}(p,q^{\prime})\leq\frac{2\bar{r}}{\varepsilon_{1}^{3}}$.

Let $\delta_{1}=2\varepsilon_{1}-\frac{9}{5}(1-\delta)$. We can now use Lemma
9 with $a=\frac{\varepsilon_{1}^{3}}{2}\sum_{p\in Q}\mathtt{d}(p,q^{\prime})$
and $b=\frac{1}{\varepsilon_{1}}\sum_{p\in Q}\mathtt{d}(p,q^{\prime})$ to
compute a set $Q_{1}\subseteq C$ of
$O\left(g\left(\frac{4b}{\varepsilon\delta_{1}a}\right)\right)=O\left(g\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)\right)$
candidate centers, one of which is a $(1+\varepsilon)$-approximate center with
probability at least $1-\delta_{1}$. The total probability of getting a good
set of candidate centers is
$(2\varepsilon_{1}-\delta_{1})(1-2\varepsilon_{1}-2\varepsilon_{1}^{2})>1-\delta$.

Case 2: $P^{\prime}$ has at most $2\varepsilon_{1}|P^{\prime}|$ points putside
${\mathcal{B}}_{2}$. We further consider two cases.

Case 2(a): $\mathtt{d}(q^{\prime},c_{P})\leq 4\varepsilon_{1}\bar{r}.$ Then

$\displaystyle\sum_{p\in P}\mathtt{d}(p,q^{\prime})\leq\sum_{p\in
P}\left(\mathtt{d}(p,c_{p})+\mathtt{d}(q^{\prime},c_{p})\right)\leq(1+4\varepsilon_{1})\sum_{p\in
P}\mathtt{d}(p,c_{P})\leq(1+\varepsilon)\sum_{p\in P}\mathtt{d}(p,c_{P}).$

Case 2(b): $\mathtt{d}(q^{\prime},c_{P})>4\varepsilon_{1}\bar{r}.$ Suppose we
assign all points from $c_{P}$ to $q^{\prime}$. By an averaging argument, we
have $|P^{\prime}|\geq(1-\varepsilon_{1}^{2})|P|$. Then, the number of points
of $P$ that are outside ${\mathcal{B}}_{2}$ is at most

$\displaystyle|P\setminus P^{\prime}|+2\varepsilon_{1}|P^{\prime}|=$
$\displaystyle|P|-(1-2\varepsilon_{1})|P^{\prime}|$ $\displaystyle\leq$
$\displaystyle|P|-(1-2\varepsilon_{1})(1-\varepsilon_{1}^{2})|P|$
$\displaystyle\leq$
$\displaystyle(\varepsilon_{1}^{2}+2\varepsilon_{1}(1-\varepsilon_{1}^{2}))|P|.$

Thus,
$|P\cap{\mathcal{B}}_{2}|\geq(1-\varepsilon_{1}^{2}-2\varepsilon_{1}(1-\varepsilon_{1}^{2}))|P|$.
Now, for $p\in P\cap{\mathcal{B}}_{2}$, the decrease in cost on switching from
$c_{P}$ to $q^{\prime}$ is at least

$\displaystyle\mathtt{d}(p,c_{P})-\mathtt{d}(p,q^{\prime})\geq\mathtt{d}(p,c_{P})-\varepsilon_{1}\bar{r}\geq\mathtt{d}(q^{\prime},c_{P})-2\varepsilon_{1}\bar{r}.$

For $p\in P\setminus{\mathcal{B}}_{2}$, the increase in cost on switching from
$c_{P}$ to $q^{\prime}$ is at most

$\displaystyle\mathtt{d}(p,q^{\prime})-\mathtt{d}(p,c_{P})\leq\mathtt{d}(q^{\prime},c_{P}).$

The overall decrease in cost is

$\displaystyle|P\cap{\mathcal{B}}_{2}|(\mathtt{d}(q^{\prime},c_{P})-2\varepsilon_{1}\bar{r})-|P\setminus{\mathcal{B}}_{2}|\mathtt{d}(q^{\prime},c_{P})>0$

for our choice of $\varepsilon_{1}$ and
$\mathtt{d}(q^{\prime},c_{P})>4\varepsilon_{1}\bar{r}$. But $c_{P}$ is the
optimal $(1,C)$-median of $P$, a contradiction. Hence case 2(b) cannot occur.

From the two cases above, we can see that $Q_{1}\cup\\{q^{\prime}\\}$ plays
the role of $\Gamma(S)$ in Definition 2.

Sampling $Q,q$ takes time $O(\frac{1}{\varepsilon})$. Computing $\sum_{p\in
P}\mathtt{d}(p,q^{\prime})$ involves $O(\frac{1}{\varepsilon})$ computations
of $\mathtt{d}$ between a pair of points from $X$, at least one of which comes
from $C$. Computing the set of candidates takes
$O\left(\theta_{C}\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)+g\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)\right)$
operations. Thus, total number of operations needed is
$O\left(\frac{1}{\varepsilon}+\theta_{C}\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)+g\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)\right)$.
Moreover,
$m_{\delta,\varepsilon}=|Q\cup\\{q\\}|=1+\frac{1}{\varepsilon_{1}}=1+\frac{4}{\varepsilon}$,
and
$t_{\delta,\varepsilon}=|Q_{1}\cup\\{q^{\prime}\\}|=O\left(g\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)\right)$.
∎

From Theorems 6 and 10, we get the following.

######  Corollary 11.

Suppose $C$ is $g$-coverable. Let
$\varepsilon\in(0,\frac{4}{9}),\delta\in(1-\frac{5}{18}\varepsilon,1)$. Given
$P\subseteq X$ of size $n$ and $k\in\mathbb{N}$, with probability $\geq
1-\delta$, a $(1+3\varepsilon)$-approximate solution to the $(k,C)$-median
problem for $P$ can be computed in time

$\displaystyle n\cdot
2^{O\left(\frac{k}{\varepsilon}\log\left(\frac{k}{\varepsilon}\right)\right)}\cdot\left(b(m_{\delta,\varepsilon})\cdot
t_{\delta,\varepsilon}\right)^{O(k)}\cdot(h(C)+t(C)),$

where
$b(m_{\delta,\varepsilon})=O\left(\frac{1}{\varepsilon}+\theta_{C}\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)+g\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)\right)$,
$t_{\delta,\varepsilon}=O\left(g\left(\frac{2048}{\delta_{1}\varepsilon^{5}}\right)\right)$,
and $\delta_{1}=\frac{\varepsilon}{2}-\frac{9}{5}(1-\delta)$.

##  5 Clustering discrete Fréchet and Hausdorff distances

In this section, we show how the results from the previous section can be used
to cluster trajectories and point sets under the discrete Fréchet and
Hausdorff distances respectively.

#### Clustering under discrete Fréchet distance.

Recall that $T^{l}$ is the set of all trajectories in $\mathbb{R}^{d}$ having
at most $l$ points each; thus $T=\bigcup_{l>0}T^{l}$ is the set of all
trajectories in $\mathbb{R}^{d}$. Given ${\mathcal{T}}=(T,\mathtt{d}_{F})$ and
trajectories $P\subseteq T$, the $(k,l)$-median problem [13, 10] is equivalent
to the $(k,T^{l})$-median problem in our setting, i.e., the center
trajectories contain at most $l$ points. We show that $T^{l}$ is $g$-coverable
for some $g$ that depends on $l$.

######  Lemma 12.

$T^{l}$ is $g$-coverable under $\mathtt{d}_{F}$ for $g(x)=l^{2l}\cdot
x^{O(dl)}$. Further, an $r^{\prime}$-cover of
${\mathcal{B}}_{\mathtt{d}_{F}}(\gamma,r)\cap T^{l}$ for $r>r^{\prime}>0$ and
$\gamma\in T^{l}$ can be computed in
$l^{2l}\cdot\left(\frac{r}{r^{\prime}}\right)^{O(dl)}$ time.

###### Proof.

Let $r>r^{\prime}>0$ be arbitrary. Let $\gamma=\langle
p_{1},\ldots,p_{l^{\prime}}\rangle\in T^{l}$ for some $l^{\prime}\leq l$.
Since the Euclidean metric in $\mathbb{R}^{d}$ has doubling dimension $O(d)$,
for any $p\in\mathbb{R}^{d}$ there exist
$\left(\frac{r}{r^{\prime}}\right)^{O(d)}$ points in the Euclidean ball
${\mathcal{B}}_{E}(p,r)$ centered at $p$ such that any point in
${\mathcal{B}}_{E}(p,r)$ is at most $r^{\prime}$ distance away from one of
these points; denote these points by $B_{p}(r,r^{\prime})$. Consider the set
of points $\bigcup_{i=1}^{l^{\prime}}B_{p_{i}}(r,r^{\prime})$; this set has
cardinality $l^{\prime}\cdot\left(\frac{r}{r^{\prime}}\right)^{O(d)}$.

Next, consider the set of all trajectories $T^{\prime}$ defined by at most
$2l$ points from $\bigcup_{i=1}^{l^{\prime}}B_{p_{i}}(r,r^{\prime})$ and
containing at least one point from $B_{p_{i}}(r,r^{\prime})$ for every $i$;
further these points respect the ordering of the sets that they belong to,
i.e., if $p\in B_{p_{i}}(r,r^{\prime})$, $q\in B_{p_{j}}(r,r^{\prime})$, and
$i<j$, then $p$ appears before $q$ in the trajectory (for points coming from
the same set $B_{p_{i}}(r,r^{\prime})$ all possible orderings are considered).
Note that
$|T^{\prime}|\leq\left(l^{\prime}\cdot\left(\frac{r}{r^{\prime}}\right)^{O(d)}\right)^{2l}$.
Further, $\mathtt{d}_{F}(\gamma,\gamma^{\prime})\leq r$ for all
$\gamma^{\prime}\in T^{\prime}$.

We will show that for any
$\gamma^{\prime\prime}\in{\mathcal{B}}_{\mathtt{d}_{F}}(\gamma,r)\cap T^{l}$,
there exists $\gamma^{\prime}\in T^{\prime}$ such that
$\mathtt{d}_{F}(\gamma^{\prime},\gamma^{\prime\prime})\leq r^{\prime}$. Thus
$T^{\prime}$ is the desired cover, completing the first part of our proof. Let
$\gamma^{\prime\prime}=\langle q_{1},\ldots,q_{l^{\prime\prime}}\rangle$ for
some $l^{\prime\prime}\leq l$. By definition of $\mathtt{d}_{F}$ and the fact
that $\mathtt{d}_{F}(\gamma,\gamma^{\prime\prime})\leq r$, each $q_{i}$ has a
corresponding sequence of points $\langle
p_{j_{i}},p_{j_{i}+1},\ldots,p_{j^{\prime}_{i}}\rangle$ each of which is at
most $r$ distance away from $q_{i}$. Moreover, for all $1\leq
i<l^{\prime\prime}$ we have $j^{\prime}_{i}\leq j_{i+1}\leq j^{\prime}_{i}+1$,
and $j_{1}=1,j^{\prime}_{l^{\prime\prime}}=l^{\prime}$.

For each $q_{i}$ and $j\in\\{j_{i},j_{i}+1,\ldots,j^{\prime}_{i}\\}$, let
$u_{j}$ denote the point in $B_{p_{j}}(r,r^{\prime})$ that is closest to
$q_{i}$. Note that $q_{i}\in{\mathcal{B}}_{E}(p_{j},r)$ and $\left\lVert
q_{i}-u_{j}\right\rVert\leq r^{\prime}$. Consider the sequence of points
$\gamma(q_{i})=\langle
u_{j_{i}},u_{j_{i}+1},\ldots,u_{j^{\prime}_{i}}\rangle$. Then,
$\mathtt{d}_{F}(\langle q_{i}\rangle,\gamma(q_{i}))\leq r^{\prime}$. Let
$\gamma^{\prime}$ be the trajectory obtained by concatenating
$\gamma(q_{1}),\ldots,\gamma(q_{l^{\prime\prime}})$ in order. Then we get
$\mathtt{d}_{F}(\gamma^{\prime},\gamma^{\prime\prime})\leq r^{\prime}$.
Further, by construction $\gamma^{\prime\prime}\in T^{\prime}$.

As far as running time is concerned, computing the set
$\bigcup_{i=1}^{l^{\prime}}B_{p_{i}}(r,r^{\prime})$ takes time
$l^{\prime}\cdot\left(\frac{r}{r^{\prime}}\right)^{O(d)}$. From this set,
computing $T^{\prime}$ takes time
$l^{2l}\cdot\left(\frac{r}{r^{\prime}}\right)^{O(dl)}$. ∎

For a trajectory having $m$ points, computing $\mathtt{d}_{F}$ to any
trajectory in $T^{l}$ takes time $O(ml)$ using the standard dynamic
programming algorithm, whereas computing the closest trajectory in $T^{l}$
under $\mathtt{d}_{F}$ takes time $O\left(lm\log
m\log\left(\frac{m}{l}\right)\right)$ time (see [7], Theorem 3). This, along
with Corollary 11 and Lemma 12 give the following.

######  Theorem 13.

Let $\varepsilon\in(0,\frac{4}{9}),\delta\in(1-\frac{5}{18}\varepsilon,1)$.
Given a set of $n$ trajectories $P\subseteq T$ each having at most $m$ points
and $k\in\mathbb{N}$, with probability $\geq 1-\delta$, the algorithm Cluster
(Fig. 2) computes a $(1+3\varepsilon)$-approximate solution to the
$(k,T^{l})$-median problem for $P$ under the discrete Fréchet distance in time

$\displaystyle nm\log m\log\left(\frac{m}{l}\right)\cdot
2^{O\left(\frac{k}{\varepsilon}\log\left(\frac{k}{\varepsilon}\right)\right)}\cdot\left(\frac{l}{\delta_{1}\varepsilon}\right)^{O(kdl)},$

where $\delta_{1}=\frac{\varepsilon}{2}-\frac{9}{5}(1-\delta)$.

#### Clustering under Hausdorff distance.

Recall that $U^{l}$ is the set of all point sets in $\mathbb{R}^{d}$
containing at most $l$ points each. Thus, $U=\bigcup_{l>0}U^{l}$ is the set of
all finite point sets of $\mathbb{R}^{d}$. Given
${\mathcal{U}}=(U,\mathtt{d}_{H})$ and subsets $P\subseteq U$, we show how to
approximately solve the $(k,U^{l})$-clustering problem for $P$ and $k,l>0$.

######  Lemma 14.

$U^{l}$ is $g$-coverable under $\mathtt{d}_{H}$ for $g(x)=l^{l}\cdot
x^{O(dl)}$. Further, an $r^{\prime}$-cover of
${\mathcal{B}}_{\mathtt{d}_{H}}(\zeta,r)\cap U^{l}$ for $r>r^{\prime}>0$ and
$\zeta\in U^{l}$ can be computed in
$l^{l}\cdot\left(\frac{r}{r^{\prime}}\right)^{O(dl)}$ time.

###### Proof.

Let $r>r^{\prime}>0$ be arbitrary, and let
$\zeta=\\{p_{1},p_{2},\ldots,p_{l^{\prime}}\\}$ for some $l^{\prime}\leq l$.
Since the Euclidean metric has doubling dimension $O(d)$, there exist
$\left(\frac{r}{r^{\prime}}\right)^{O(d)}$ points in the Euclidean ball
${\mathcal{B}}_{E}(p,r)$ such that any point in ${\mathcal{B}}_{E}(p,r)$ is at
most $r^{\prime}$ distance away from one of these points; denote these points
by $B_{p}(r,r^{\prime})$.

Consider all subsets of size $l$ of the set
$\bigcup_{i=1}^{l^{\prime}}B_{p_{i}}(r,r^{\prime})$, denote this set by
$U^{\prime}$. Note that
$|U^{\prime}|=\left(l^{\prime}\cdot\left(\frac{r}{r^{\prime}}\right)^{O(d)}\right)^{l}$.

Next, consider any point set
$\zeta^{\prime}\in{\mathcal{B}}_{\mathtt{d}_{H}}(\zeta,r)\cap U^{l}$. By
definition of Hausdorff distance, the points of $\zeta^{\prime}$ (there are at
most $l$ of them) must lie in
$\bigcup_{i=1}^{l^{\prime}}{\mathcal{B}}_{E}(p_{i},r)$. Thus, for each
$p\in\zeta^{\prime}$, there is some $i\in\\{1,\ldots,l^{\prime}\\}$ such that
$\left\lVert p-q\right\rVert\leq r^{\prime}$ for some $q\in
B_{p_{i}}(r,r^{\prime})$. Thus, there exists $\zeta^{\prime\prime}\in
U^{\prime}$ such that $\mathtt{d}_{H}(\zeta^{\prime},\zeta^{\prime\prime})\leq
r^{\prime}$. Then, $U^{\prime}$ is the desired cover, completing the first
part of our proof.

Computing $\bigcup_{i=1}^{l^{\prime}}B_{p_{i}}(r,r^{\prime})$ takes time
$l^{\prime}\cdot\left(\frac{r}{r^{\prime}}\right)^{O(d)}$. Computing
$U^{\prime}$ from it takes time
$\left(l^{\prime}\cdot\left(\frac{r}{r^{\prime}}\right)^{O(d)}\right)^{l}$. ∎

Computing $\mathtt{d}_{H}$ between two point sets of size $m_{1}$ and $m_{2}$
in $\mathbb{R}^{d}$ takes time $O(m_{1}m_{2})$. Given $\zeta\in U$ of size
$m$, computing the closest point to it in $U^{l}$ (under $\mathtt{d}_{H}$)
boils down to finding $l$ disks in $\mathbb{R}^{d}$ of minimum radius such
that all the points of $\zeta$ lie inside the union of these disks; the
centers of these disks give the desired set in $U^{l}$. This is the $l$-center
problem in $\mathbb{R}^{d}$ for the Euclidean metric. Solving this problem is
NP-Hard when $l$ is part of the input. Note that the total number of subsets
of $\zeta$ induced by disks in $\mathbb{R}^{d}$ is $m^{O(d)}$. By looking at
$l$ such subsets at a time, we can pick the one that covers $\zeta$ and
minimizes the radius of the largest disk; this takes total time $m^{O(dl)}$.
This along with Corollary 11 and Lemma 14 give the following.

######  Theorem 15.

Let $\varepsilon\in(0,\frac{4}{9}),\delta\in(1-\frac{5}{18}\varepsilon,1)$.
Given a set of $n$ point sets $P\subseteq U$ each having at most $m$ points
and $k\in\mathbb{N}$, with probability $\geq 1-\delta$, the algorithm Cluster
(Fig. 2) computes a $(1+3\varepsilon)$-approximate solution to the
$(k,U^{l})$-median problem for $P$ under the Hausdorff distance can be
computed in time

$\displaystyle nm^{O(dl)}\cdot
2^{O\left(\frac{k}{\varepsilon}\log\left(\frac{k}{\varepsilon}\right)\right)}\cdot\left(\frac{l}{\delta_{1}\varepsilon}\right)^{O(kdl)},$

where $\delta_{1}=\frac{\varepsilon}{2}-\frac{9}{5}(1-\delta)$.

#### Remark.

We believe the running time of Theorem 15 can be improved to be similar to
that of Theorem 13, since a $\left(1+O(\varepsilon)\right)$-approximate
nearest neighbor (instead of the actual nearest neighbor) in $C$ for any $x\in
X$ should suffice. For the Hausdorff case, this involves using a fast
approximation for the $l$-center problem [2].

##  6 Hardness of $k$-median clustering under Hausdorff distance

We prove the following hardness result for $k$-median clustering under the
Hausdorff distance.

######  Theorem 16.

The $k$-median clustering problem for finite point sets under the Hausdorff
distance is NP-Hard.

###### Proof.

We reduce the Euclidean $k$-median problem, which is known to be NP-Hard [29].
The reduction is fairly straightforward – for each input point $p$ of an
instance of the Euclidean $k$-median problem, we have a singleton set
$\\{p\\}$ as input to the Haudorff $k$-median problem. Any solution to the
Euclidean $k$-median problem is also a solution to the Hausdorff $k$-median
problem of the same cost – we replace cluster center $c$ in the Euclidean
version by the cluster center $\\{c\\}$ for the Haudorff version, and for each
$p$ assigned to $c$, we assign $\\{p\\}$ to $\\{c\\}$.

On the other hand, consider a solution to the instance of Hausdorff $k$-median
problem. In particular, let $S$ be a cluster center that is assigned the sets
$\\{\\{p_{1}\\},\ldots,\\{p_{n}\\}\\}$. The cost of this single cluster is

$\displaystyle\sum_{i=1}^{n}\mathtt{d}_{H}(\\{p_{i}\\},S)=\sum_{i=1}^{n}\max_{s\in
S}\left\lVert s-p_{i}\right\rVert,$

by the definition of $\mathtt{d}_{H}$. Thus, replacing $S$ by a singleton set
$\\{s\\}$ for any $s\in S$ does not increase the cost of clustering. Hence we
can assume that all cluster centers are singleton sets. We can then construct
a solution for the Euclidean $k$-median problem by assigning $p$ to $s$, where
$\\{s\\}$ is the cluster center that $\\{p\\}$ was assigned to in the
Hausdorff clustering solution. This does not increase the cost of the
clustering as well. ∎

##  7 Conclusion

We have given a framework for clustering where the cluster centers are
restricted to belong to a simpler metric space. We characterized general
conditions on this simpler space that allow us to obtain efficient
$(1+\varepsilon)$-approximation algorithms for the $k$-median problem. As
special cases, we gave efficient algorithms for clustering trajectories and
point sets under the discrete Fréchet and Hausdorff distances respectively.

We believe the general framework can be extended to other metric spaces as
well. The next step would be to see if it can be applied to the continuous
Fréchet distance, and to non-metric distance measures such as dynamic time
warping. It would be interesting to provide other characterizations on the
metric space for the cluster centers (as alternatives to the notion of
covering discussed in this paper) that are amenable to efficient clustering
algorithms.

## References

  * [1] Marcel R. Ackermann, Johannes Blömer, and Christian Sohler.  Clustering for metric and nonmetric distance measures.  ACM Trans. Alg., 6(4):59:1–59:26, 2010.
  * [2] Pankaj K Agarwal and Cecilia M Procopiuc.  Exact and approximation algorithms for clustering.  Algorithmica, 33(2):201–226, 2002.
  * [3] Sanjeev Arora, Prabhakar Raghavan, and Satish Rao.  Approximation schemes for Euclidean k-medians and related problems.  In Proc. ACM Symp. Th. Comput., volume 98, pages 106–113, 1998\.
  * [4] Arturs Backurs and Anastasios Sidiropoulos.  Constant-distortion embeddings of hausdorff metrics into constant-dimensional l_p spaces.  In APPROX/RANDOM. Schloss Dagstuhl-Leibniz-Zentrum fuer Informatik, 2016.
  * [5] Mihai Bādoiu, Sariel Har-Peled, and Piotr Indyk.  Approximate clustering via core-sets.  In Proc. ACM Symp. Th. Computing, pages 250–257. ACM, 2002.
  * [6] Nicolas Basalto, Roberto Bellotti, Francesco De Carlo, Paolo Facchi, Ester Pantaleo, and Saverio Pascazio.  Hausdorff clustering of financial time series.  Physica A: Statistical Mechanics Applications, 379(2):635–644, 2007\.
  * [7] Sergey Bereg, Minghui Jiang, Wencheng Wang, Boting Yang, and Binhai Zhu.  Simplifying 3D polygonal chains under the discrete Fréchet distance.  In Lat. Amer. Symp. Theoret. Informatics, pages 630–641. Springer, 2008.
  * [8] P. C. Besse, B. Guillouet, J. Loubes, and F. Royer.  Review and perspective for distance-based clustering of vehicle trajectories.  IEEE Transactions on Intelligent Transportation Systems, 17(11):3306–3317, 2016.
  * [9] Kevin Buchin, Anne Driemel, Joachim Gudmundsson, Michael Horton, Irina Kostitsyna, Maarten Löffler, and Martijn Struijs.  Approximating (k,l)-center clustering for curves.  In Proc. ACM-SIAM Symp. Disc. Alg., pages 2922–2938. SIAM, 2019\.
  * [10] Kevin Buchin, Anne Driemel, and Martijn Struijs.  On the hardness of computing an average curve.  arXiv preprint arXiv:1902.08053, 2019.
  * [11] Jinyang Chen, Rangding Wang, Liangxu Liu, and Jiatao Song.  Clustering of trajectories based on hausdorff distance.  In Proc. Int. Conf. Electronics Comm. Control, pages 1940–1944. IEEE, 2011.
  * [12] Vincent Cohen-Addad, Philip N Klein, and Claire Mathieu.  Local search yields approximation schemes for k-means and k-median in euclidean and minor-free metrics.  SIAM J. Comput., 48(2):644–667, 2019.
  * [13] Anne Driemel, Amer Krivošija, and Christian Sohler.  Clustering time series under the Fréchet distance.  In Proc. ACM-SIAM Symp. Disc. Alg., pages 766–785. SIAM, 2016\.
  * [14] Thomas Eiter and Heikki Mannila.  Computing discrete Fréchet distance.  Technical report, Information Systems Dept., Technical University of Vienna, 1994.
  * [15] Martin Ester, Hans-Peter Kriegel, Jörg Sander, and Xiaowei Xu.  A density-based algorithm for discovering clusters in large spatial databases with noise.  In Proc. ACM Int. Conf. Know. Disc. Data Mining, volume 96, pages 226–231, 1996.
  * [16] Tomás Feder and Daniel Greene.  Optimal algorithms for approximate clustering.  In Proceedings of the twentieth annual ACM symposium on Theory of computing, pages 434–444. ACM, 1988.
  * [17] Scott Gaffney and Padhraic Smyth.  Trajectory clustering with mixtures of regression models.  In Proc. ACM Int. Conf. Know. Disc. Data Mining, volume 99, pages 63–72, 1999.
  * [18] Venkatesan Guruswami and Piotr Indyk.  Embeddings and non-approximability of geometric problems.  In Proc. ACM-SIAM Symp. Disc. Alg., volume 3, pages 537–538, 2003\.
  * [19] Sariel Har-Peled and Akash Kushal.  Smaller coresets for k-median and k-means clustering.  Disc. Computat. Geom., 37(1):3–19, 2007.
  * [20] Sariel Har-Peled and Soham Mazumdar.  On coresets for k-means and k-median clustering.  In Proc. ACM Symp. Th. Computing, pages 291–300. ACM, 2004.
  * [21] Felix Hausdorff.  Grundzuge der mengenlehre, volume 61.  American Mathematical Soc., 1978.
  * [22] Chih-Chieh Hung, Wen-Chih Peng, and Wang-Chien Lee.  Clustering and aggregating clues of trajectories for mining trajectory patterns and routes.  Int. J. Very Large Databases, 24(2):169–192, 2015.
  * [23] Daniel P Huttenlocher, Gregory A Klanderman, and William J Rucklidge.  Comparing images using the hausdorff distance.  IEEE Transactions on pattern analysis and machine intelligence, 15(9):850–863, 1993.
  * [24] Kamal Jain, Mohammad Mahdian, and Amin Saberi.  A new greedy approach for facility location problems.  In Proc. ACM Symp. Th. Computing, pages 731–740. ACM, 2002.
  * [25] Stavros G Kolliopoulos and Satish Rao.  A nearly linear-time approximation scheme for the euclidean k-median problem.  SIAM J. Comput., 37(3):757–782, 2007.
  * [26] Amit Kumar, Yogish Sabharwal, and Sandeep Sen.  Linear time algorithms for clustering problems in any dimensions.  In Int. Coll. Automata Lang. Programming, pages 1374–1385. Springer, 2005.
  * [27] Amit Kumar, Yogish Sabharwal, and Sandeep Sen.  Linear-time approximation schemes for clustering problems in any dimensions.  J. ACM, 57(2):5, 2010.
  * [28] Shi Li and Ola Svensson.  Approximating k-median via pseudo-approximation.  SIAM J. Computing, 45(2):530–547, 2016.
  * [29] Nimrod Megiddo and Kenneth J Supowit.  On the complexity of some common geometric location problems.  SIAM J. Comput., 13(1):182–196, 1984.
  * [30] Lin Qu, Fan Zhou, and YW Chen.  Trajectory classification based on hausdorff distance for visual surveillance system.  J. Jilin University, 6:1618–1624, 2009.
  * [31] Cynthia Sung, Dan Feldman, and Daniela Rus.  Trajectory clustering for motion prediction.  In 2012 IEEE/RSJ International Conference on Intelligent Robots and Systems, pages 1547–1552. IEEE, 2012.
  * [32] Hongteng Xu, Yang Zhou, Weiyao Lin, and Hongyuan Zha.  Unsupervised trajectory clustering via adaptive multi-kernel-based shrinkage.  In Proc. IEEE Int. Conf. Comp. Vision, pages 4328–4336, 2015.

##  Appendix A Doubling dimension of Fréchet and Hausdorff distance

We show that the discrete Fréchet distance does not have a doubling dimension
bounded by a constant. The Hausdorff case can be shown similarly.

Suppose ${\mathcal{T}}=(T,\mathtt{d}_{F})$ has a constant doubling dimension
$D$. Then by definition, for any trajectory $\gamma\in T$ and any $r>0$, there
exist $2^{D}$ trajectories $T^{\prime}\subset T$ such that any trajectory
$\gamma^{\prime}\in{\mathcal{B}}_{\mathtt{d}_{F}}(\gamma,r)$ is at most at
distance $\frac{r}{2}$ from a trajectory in $T^{\prime}$. By the pigeon hole
principle, any set of $2^{D}+1$ or more trajectories in
${\mathcal{B}}_{\mathtt{d}_{F}}(\gamma,r)$ will have at least two trajectories
that have the same closest trajectory in $T^{\prime}$, and are therefore at
most $r$ distance apart (by the triangle inequality).

Consider a trajectory $\gamma^{\prime\prime}$ consisting of a sequence of
points $\langle p_{1},p_{2},\ldots,p_{m}\rangle$ in a straight line in
$\mathbb{R}^{d}$, with a distance of at least $3r$ between every consecutive
pair of points, for some $m,r>0$. For each $p_{i}$, let
$P_{i}=\\{a_{i},b_{i}\\}$ be such that $\left\lVert a_{i}-b_{i}\right\rVert>r$
and $\left\lVert a_{i}-p_{i}\right\rVert,\left\lVert
b_{i}-p_{i}\right\rVert\leq r$ (e.g., $a_{i}$ and $b_{i}$ can be the
diametrically opposite points of the Euclidean ball of radius $r$ centered at
$p_{i}$). Consider the set of trajectories $T^{\prime\prime}=\\{\langle
q_{1},q_{2},\ldots,q_{m}\rangle\mid q_{i}\in P_{i}\\}$. Then by construction,
$T^{\prime\prime}\subset{\mathcal{B}}_{\mathtt{d}_{F}}(\gamma^{\prime\prime},r)$,
and $\mathtt{d}_{F}(\gamma_{1},\gamma_{2})>r$ for any
$\gamma_{1},\gamma_{2}\in T^{\prime\prime}$. Further,
$|T^{\prime\prime}|=2^{m}\geq 2^{D}+1$ for sufficiently large $m$,
contradicting the fact that ${\mathcal{T}}=(T,\mathtt{d}_{F})$ has doubling
dimension $D$.

##  Appendix B Proofs of Section 3

###  B.1 Proof of the superset sampling lemma

We state the superset sampling lemma for our setting. As in Section 3, we
define $F(S^{\prime})$ to be $\\{c_{S^{\prime}}\\}$ or $\Gamma(S^{\prime})$
for instances satisfying the strong or weak sampling property respectively.
The lemma states that we can draw a uniform sample from $P^{\prime}\subseteq
P$ without explicitly knowing $P^{\prime}$, provided $P^{\prime}$ contains a
sufficient fraction of points of $P$. In our algorithm in Fig. 2, we use this
property in the _sampling phase_ in order to sample points from a single
cluster, without explicitly knowing the points in the said cluster.

######  Lemma 17.

Suppose $0<\alpha<\frac{1}{4}$. Suppose $(X,C,\mathtt{d})$ satisfies either
the strong or weak sampling property. Let $P\subseteq X$ of size $n$ and
$P^{\prime}\subseteq P$ with $|P^{\prime}|\geq\alpha n$. Let $S\subseteq P$ be
a uniform sample multiset of size at least
$\frac{2}{\alpha}m_{\delta,\varepsilon}$. Then with probability at least
$\frac{1-\delta}{5}$, there exists a subset $S^{\prime}\subseteq S$ with
$|S^{\prime}|=m_{\delta,\varepsilon}$ satisfying the following.

$\displaystyle\exists c\in F(S^{\prime})\text{ s.t. }\sum_{p\in
P}\mathtt{d}(p,c)\leq(1+\varepsilon)\sum_{p\in P}\mathtt{d}(p,c_{P}).$

Here, $F(S^{\prime})$ is either $\\{c_{S^{\prime}}\\}$ or $\Gamma(S^{\prime})$
for $(X,C,\mathtt{d})$ satisfying the strong or weak sampling property
respectively.

###### Proof.

This proof is similar to that of Lemma 2.1 in [1]. Define $Y$ to be a random
variable denoting the number of points from $P^{\prime}$ contained in $S$.
Note that $\mathbb{E}[Y]\geq 2m_{\delta,\varepsilon}$, since $S$ is sampled
with replacement. By applying a Chernoff bound, we obtain

$\displaystyle\Pr[Y<m_{\delta,\varepsilon}]\leq\Pr\left[Y<\frac{\mathbb{E}[Y]}{2}\right]\leq
e^{-m_{\delta,\varepsilon}/4}\leq e^{-1/4}<\frac{4}{5},$

since $m_{\delta,\varepsilon}\geq 1$. Thus, with probability at least
$\frac{1}{5}$, $S$ containts at least $m_{\delta,\varepsilon}$ points of
$P^{\prime}$. This along with the strong or weak sampling property finishes
the proof. ∎

###  B.2 Proof of Lemma 3

The following is similar to the proof of Theorem 2.2 [1], as the correctness
does not depend on how candidate centers are found. We include it for
completeness.

###### Proof.

We prove the lemma for $k=2$; it generalizes in a straightforward manner for
$k>2$. We assume $n=|P|$ is a power of $2$ for simplicity. Suppose
$P_{1},P_{2}$ are the clusters with centers $C^{*}=\\{c^{*}_{1},c^{*}_{2}\\}$
corresponding to the optimal $(2,C)$-median for $P$. Assume
$|P_{1}|\geq\frac{1}{2}P_{2}$. Let $opt_{1}(Q)$ and $opt_{2}(Q)$ denote the
_value_ of the optimal $(1,C)$-median and $(2,C)$-median respectively.

By Lemma 17, during the sampling phase of the initial call to Cluster,
$\overline{C_{S}}$ contains a $\overline{c_{1}}$ such that $\sum_{p\in
P_{1}}\mathtt{d}(p,\overline{c_{1}})\leq(1+\varepsilon)opt_{1}$. We consider
two cases: when the algorithm selects $\overline{c_{1}}$ and recurses with
$(\overline{P},1,\\{\overline{c_{1}}\\})$ where $\overline{P}$ contains a
suitably large number of points from $P_{2}$, and when this does not occur.

Case 1: Suppose there exists a call with
$(\overline{P},1,\\{\overline{c_{1}}\\})$ such that
$|P_{2}\cap\overline{P}|\geq\alpha|\overline{P}|$. Then, by Lemma 17, during
this call $\overline{C_{S}}$ contains a $\overline{c_{2}}$ such that
$\sum_{p\in
P_{2}\cap\overline{P}}\mathtt{d}(p,\overline{c_{2}})\leq(1+\varepsilon)opt_{1}(P_{2}\cap\overline{P})$.
In this case, we upper bound the cost of
$\overline{C}=\\{\overline{c_{1}},\overline{c_{2}}\\}$. Let
$N=P\setminus\overline{P}$ be the points removed by the pruning phase between
the sampling of $\overline{c_{1}}$ and $\overline{c_{2}}$. $P_{1},P_{2}\cap
N$, and $P_{2}\cap\overline{P}$ form a partition of $P$. We have that:

$\sum_{p\in P}\mathtt{d}(p,\overline{C})\leq\sum_{p\in
P_{1}}\mathtt{d}(p,\overline{c_{1}})+\sum_{p\in P_{2}\cap
N}\mathtt{d}(p,\overline{c_{1}})+\sum_{p\in
P_{2}\cap\overline{P}}\mathtt{d}(p,\overline{c_{2}}).$ (1)

By the definition of $\overline{c_{1}}$, we can bound the first term of
Equation 1:

$\sum_{p\in P_{1}}\mathtt{d}(p,\overline{c_{1}})\leq(1+\varepsilon)\sum_{p\in
P_{1}}\mathtt{d}(p,c^{*}_{1}).$

Next, we bound the last term of Equation 1. By the selection of
$\overline{c_{2}}$:

$\sum_{p\in
P_{2}\cap\overline{P}}\mathtt{d}(p,\overline{c_{2}})\leq(1+\varepsilon)\sum_{p\in
P_{2}\cap\overline{P}}\mathtt{d}(p,c_{P_{2}\cap\overline{P}})\leq(1+\varepsilon)\sum_{p\in
P_{2}}d(p,c^{*}_{2}).$

Now, we only need to bound the middle term of Equation 1. We first assume
$N\neq\emptyset$, otherwise we are done. Suppose there are $t$ recursive calls
(and thus pruning phases) between sampling $\overline{c_{1}}$ and
$\overline{c_{2}}$. Then, $N=N^{(1)}\cup\dots\cup N^{(t)}$ where
$|N^{(i)}|=\frac{n}{2^{i}}$, since in each pruning phase
$\frac{1}{2}|\overline{P}|$ points are removed. Intuitively, each $N^{(i)}$
contains only a few points from $P_{2}$. Let $\overline{P}^{(0)}=P$ and
$\overline{P}^{(i)}=\overline{P}^{(i-1)}\setminus N^{(i)}$. Recall that when
$|P_{2}\cap\overline{P}|\geq\alpha|\overline{P}|$, we assign these points to
$\overline{c_{2}}$. Thus,
$|P_{2}\cap\overline{P}^{(i)}|<\alpha|\overline{P}^{(i)}|$ for all $i<t$.
Thus, we bound the number of points of $P_{2}$ in each $N^{(i)}$ for $i\leq
t$:

$|P_{2}\cap
N^{(i)}|\leq|P_{2}\cap\overline{P}^{(i-1)}|<\alpha|\overline{P}^{(i-1)}|=2\alpha\frac{n}{2^{i}}.$
(2)

because $N^{(i)}\subset\overline{P}^{(i-1)}$. We can also bound the number of
points of $P_{1}$ in each $N^{(i)}$:

$|P_{1}\cap N^{(i)}|\geq|N^{(i)}|-|P_{2}\cap
N^{(i)}|\geq(1-2\alpha)\frac{n}{2^{i}}.$ (3)

We first show that assigning $P_{2}\cap N$ to $\overline{c_{1}}$ has small
cost. If $p\in N^{(i)}$ and $p^{\prime}\in N^{(i+1)}$, it must be that
$\mathtt{d}(p,\overline{c_{1}})\leq\mathtt{d}(p^{\prime},\overline{c_{1}})$
since minimal points with respect to $\overline{C}$ are chosen at each step to
be removed. Thus we can sum over such $p,p^{\prime}$, and for all $i<t$

$\frac{1}{|P_{2}\cap N^{(i)}|}\sum_{p\in P_{2}\cap
N^{(i)}}\mathtt{d}(p,\overline{c_{1}})\leq\frac{1}{|P_{1}\cap
N^{(i+1)}|}\sum_{p\in P_{1}\cap N^{(i+1)}}\mathtt{d}(p,\overline{c_{1}}).$

Combining this with Equations 2 and 3, for all $i<t$ we obtain

$\sum_{p\in P_{2}\cap
N^{(i)}}\mathtt{d}(p,\overline{c_{1}})\leq\frac{4\alpha}{1-2\alpha}\sum_{p\in
P_{1}\cap N^{(i+1)}}\mathtt{d}(p,\overline{c_{1}}).$ (4)

Finally, we bound the cost of assigning points in $P_{2}\cap N^{(t)}$ to
$\overline{c_{1}}$. We have
$|P_{1}\cap\overline{P}^{(t)}|=|\overline{P}^{(t)}|-|P_{2}\cap\overline{P}^{(t)}|\geq|\overline{P}^{(t)}|-|P_{2}\cap\overline{P}^{(t-1)}|>(1-2\alpha)\frac{n}{2^{t}}$.
We combine this with our previous observations:

$\sum_{p\in P_{2}\cap
N^{(t)}}\mathtt{d}(p,\overline{c_{1}})\leq\frac{2\alpha}{1-2\alpha}\sum_{p\in
P_{1}\cap\overline{P}^{(t)}}\mathtt{d}(p,\overline{c_{1}}).$ (5)

Now, we have computed a bound for assigning points of $P_{2}$ to
$\overline{c_{1}}$ in each $N^{(i)}$. Combining these:

$\displaystyle\sum_{p\in P_{2}\cap N}\mathtt{d}(p,\overline{c_{1}})$
$\displaystyle=\sum_{i=1}^{t}\sum_{p\in P_{2}\cap
N^{(i)}}\mathtt{d}(p,\overline{c_{1}})$
$\displaystyle\leq\frac{4\alpha}{1-2\alpha}\sum_{i=1}^{t-1}\sum_{p\in
P_{1}\cap
N^{(i+1)}}\mathtt{d}(p,\overline{c_{1}})+\frac{2\alpha}{1-2\alpha}\sum_{p\in
P_{1}\cap\overline{P}^{(t)}}\mathtt{d}(p,\overline{c_{1}})$ $\displaystyle\leq
8\alpha\sum_{i=1}^{t-1}\sum_{p\in P_{1}\cap
N^{(i+1)}}\mathtt{d}(p,\overline{c_{1}})+8\alpha\sum_{p\in
P_{1}\cap\overline{P}^{(t)}}\mathtt{d}(p,\overline{c_{1}})$ $\displaystyle\leq
8\alpha\sum_{p\in P_{1}}\mathtt{d}(p,\overline{c_{1}})$ $\displaystyle\leq
8\alpha(1+\varepsilon)\sum_{p\in P_{1}}\mathtt{d}(p,c^{*}_{1}),$

for $\alpha\leq\frac{1}{4}$.

Thus for Case 1, we can bound the algorithm cost from Equation 1:

$\displaystyle\sum_{p\in P}\mathtt{d}(p,\overline{C})$
$\displaystyle\leq\sum_{p\in P_{1}}\mathtt{d}(p,\overline{c}_{1})+\sum_{p\in
P_{2}\cap N}\mathtt{d}(p,\overline{c_{1}})+\sum_{p\in
P_{2}\cap\overline{P}}\mathtt{d}(p,\overline{c_{2}})$
$\displaystyle\leq(1+\varepsilon)\sum_{p\in
P_{1}}\mathtt{d}(p,c^{*}_{1})+8\alpha(1+\varepsilon)\sum_{p\in
P_{1}}\mathtt{d}(p,c^{*}_{1})+(1+\varepsilon)\sum_{p\in
P_{2}}\mathtt{d}(p,c^{*}_{2})$
$\displaystyle\leq(1+8\alpha)(1+\varepsilon)opt_{2}(P)$

Case 2: If there is no recursive call with $|P_{2}\cap R|\geq\alpha|R|$ the
pruning phase will be called recursively $\lceil\log n\rceil$ times until
there is a single point $q\in\overline{P}$. $q$ can be assigned to a cluster
by itself at a cost of $0$. Then, $N=P\setminus\\{q\\}$ and the proof of the
previous case also bounds the cost of assigning $P\cap N$ to
$\overline{c_{1}}$.

For $k>2$, consider the algorithm as each center is added. Let
$\overline{C_{i}}=\\{\overline{c_{1}},\dots\overline{c_{i}}\\}$ be the $i$
medians already approximated, corresponding to supercluster $P^{\prime}_{1}$,
with $P^{\prime}_{2}$ consisting of clusters whose medians are yet to be
found. Similar analysis as above shows that the cost of points incorrectly
assigned to centers in $\overline{C_{i}}$ in the pruning phases can be bounded
by $8\alpha k\sum_{p\in P^{\prime}_{1}}\mathtt{d}(p,\overline{C_{i}})$. See
the proof of Theorem 2.5 of [1] for full details. ∎

###  B.3 Proof of Lemma 4

###### Proof.

We will first count the number of operations required, where computing
$\mathtt{d}$ between points in $X$ and $C$, and computing the closest point in
$C$ for any point in $X$ both count as one operation each.

Let $T(n,k)$ denote the number of operations required by Cluster with $n$
input points and $k$ medians to be found. For $k=0$, we clearly have
$T(n,0)=O(1)$. For $n\leq k$, we can put each input point in its own cluster,
and return the points in $C$ closest to each input point as the cluster
medians. Thus, $T(n,k)\leq O(n)$.

Let us consider the case $n>k\geq 1$. In the sampling phase, the number of
candidate centers generated is
$2^{O\left(m_{\delta,\varepsilon}\log\left(\frac{1}{\alpha}m_{\delta,\varepsilon}\right)\right)}\cdot
w(m_{\delta,\varepsilon})$, each taking $f(m_{\delta,\varepsilon})$
operations. Each of the candidate centers is then tried recursively, each
taking $T(n,k-1)$ operations. The pruning phase takes $O(n)$ operations. After
pruning, the algorithm is called once for the remaining point set, requiring
$T(n/2,k)$ operations. We thus have

$\displaystyle T(n,k)\leq$ $\displaystyle
2^{O\left(m_{\delta,\varepsilon}\log\left(\frac{1}{\alpha}m_{\delta,\varepsilon}\right)\right)}\cdot
w(m_{\delta,\varepsilon})\cdot\left(f(m_{\delta,\varepsilon})+T(n,k-1)\right)+T(n/2,k)+O(n)$
$\displaystyle\leq$ $\displaystyle
2^{O\left(m_{\delta,\varepsilon}\log\left(\frac{1}{\alpha}m_{\delta,\varepsilon}\right)\right)}\cdot
w(m_{\delta,\varepsilon})\cdot f(m_{\delta,\varepsilon})\cdot
T(n,k-1)+T(n/2,k)+O(n).$

Solving the recurrence yields

$\displaystyle T(n,k)=n\cdot
2^{O\left(km_{\delta,\varepsilon}\log\left(\frac{1}{\alpha}m_{\delta,\varepsilon}\right)\right)}\cdot\left(w(m_{\delta,\varepsilon})\cdot
f(m_{\delta,\varepsilon})\right)^{O(k)}.$

Taking into account the time to compute $\mathtt{d}$ and the closest point in
$C$, we get the total running time to be $n\cdot
2^{O\left(km_{\delta,\varepsilon}\log\left(\frac{1}{\alpha}m_{\delta,\varepsilon}\right)\right)}\cdot\left(w(m_{\delta,\varepsilon})\cdot
f(m_{\delta,\varepsilon})\right)^{O(k)}\cdot(h(C)+t(C))$. ∎

##  Appendix C Proof of Theorem 7

The following two useful lemmas hold for any metric space
${\mathcal{X}}=(X,\mathtt{d})$.

######  Lemma 18 ([1], Lemma 3.2).

Let $c\in X$, $P\subseteq X$ of size $n$, and $\delta>0$. A uniform sample
multiset $S\subseteq P$ of size $m$ satisfies

$\displaystyle\Pr\left[\exists q\in S\mid\mathtt{d}(q,c)\geq\frac{1}{\delta
n}\sum_{p\in P}\mathtt{d}(p,c)\right]\leq m\delta.$

######  Lemma 19 ([1], Lemma 3.3).

Let $\varepsilon\in(0,1]$, $P\subseteq X$ of size $n$, and $b,c\in X$ be such
that $\sum_{p\in P}\mathtt{d}(p,b)>(1+\tfrac{4}{5}\varepsilon)\sum_{p\in
P}\mathtt{d}(p,c)$. A uniform sample multiset $S\subseteq P$ of size $m$
satisfies

$\displaystyle\Pr\left[\sum_{s\in S}\mathtt{d}(s,b)\leq\sum_{s\in
S}\mathtt{d}(s,c)+\frac{\varepsilon m}{5n}\sum_{p\in
P}\mathtt{d}(p,c)\right]<\exp\left(-\frac{\varepsilon^{2}m}{144}\right).$

Intuitively, Lemma 18 says that if we are sampling points uniformly, the
probability of sampling a point far away from a fixed point increases with the
size of the sample. Lemma 19 states that if the average distance of a point
set to a point $c$ is smaller than to a point $b$, then the average distance
of a uniform sample of the point set to $c$ relative to $b$ is also small,
with some probability that can be lower-bounded.

We now prove Theorem 7.

###### Proof of Theorem 7.

Let $P,S,\varepsilon,\delta$ be as specified in Definition 1. Let $|P|=n$ and
$|S|=m$; the value of $m$ will be set later in the proof.

Let $r=\frac{6m}{\delta n}\sum_{p\in P}\mathtt{d}(p,c_{P})$. Let
$r^{\prime}=\frac{r}{3},U={\mathcal{B}}(c_{P},r),$ and
$U^{\prime}={\mathcal{B}}(c_{P},r^{\prime})$. By Lemma 18, for all $s\in S$,
$\mathtt{d}(c_{P},s)\leq r^{\prime}$ with probability at least $1-\delta/2$.
Thus, $S\subseteq U^{\prime}$ with probability at least $1-\delta/2$. If
$c_{S}\notin U$, then $\sum_{s\in S}\mathtt{d}(s,c_{S})>2r^{\prime}m$ (by
triangle inequality). However since $S\subseteq U^{\prime}$, we have
$\sum_{s\in S}\mathtt{d}(s,c_{P})\leq r^{\prime}m$, contradicting the claim
that $c_{S}$ is the optimal $(1,C)$-median of S. Thus, with probability at
least $1-\delta/2$, $c_{S}\in U$.

Since $C$ is $g$-coverable, there exists a
$\left(\frac{\varepsilon}{5n}\sum_{p\in P}\mathtt{d}(p,c_{P})\right)$-cover of
$U\cap C$ of size $g\left(\frac{30m}{\delta\varepsilon}\right)$. Let
$C^{\prime}$ be such a cover. Define

$\displaystyle C^{\prime}_{bad}=\\{c\in C^{\prime}\mid\sum_{p\in
P}\mathtt{d}(p,c)>(1+\tfrac{4}{5}\varepsilon)\sum_{p\in
P}\mathtt{d}(p,c_{P})\\}.$

Setting $m$ to be the sufficiently large constant $m_{\delta,\varepsilon,g}$
and using the union bound and Lemma 19 we get 222We assume that $g$ is such
that for large enough $m$,
$g\left(\frac{30m}{\delta\varepsilon}\right)\exp\left(\frac{-\varepsilon^{2}m}{144}\right)<\delta/2$.

$\displaystyle\Pr\left[\exists c\in C^{\prime}_{bad}\mid\sum_{s\in
S}\mathtt{d}(s,c)\leq\sum_{s\in S}\mathtt{d}(s,c_{P})+\frac{\varepsilon
m}{5n}\sum_{p\in P}\mathtt{d}(p,c_{P})\right]$ $\displaystyle<$ $\displaystyle
g\left(\frac{30m}{\delta\varepsilon}\right)\exp\left(\frac{-\varepsilon^{2}m}{144}\right)$
$\displaystyle<$ $\displaystyle\delta/2.$

Thus, with probability at least $1-\delta/2$, for all $c\in C^{\prime}_{bad}$
we have

$\displaystyle\sum_{s\in S}\mathtt{d}(s,c)>\sum_{s\in
S}\mathtt{d}(s,c_{P})+\frac{\varepsilon m}{5n}\sum_{p\in
P}\mathtt{d}(p,c_{P}).$

Let $c^{\prime}$ be the closest point in $C^{\prime}$ to $c_{S}$. By
definition of $C^{\prime}$ and the fact that $c_{S}\in U\cup C$, we have
$\mathtt{d}(c_{S},c^{\prime})\leq\frac{\varepsilon}{5n}\sum_{p\in
P}\mathtt{d}(p,c_{P})$. We then have

$\displaystyle\sum_{s\in S}\mathtt{d}(s,c^{\prime})$
$\displaystyle\leq\sum_{s\in
S}\left(\mathtt{d}(s,c_{S})+\mathtt{d}(c_{S},c^{\prime})\right)\leq\sum_{s\in
S}\mathtt{d}(s,c_{S})+\frac{\varepsilon m}{5n}\sum_{p\in
P}\mathtt{d}(p,c_{P})$ $\displaystyle\leq\sum_{s\in
S}\mathtt{d}(s,c_{P})+\frac{\varepsilon m}{5n}\sum_{p\in
P}\mathtt{d}(p,c_{P})<\sum_{s\in S}\mathtt{d}(s,c)$

for all $c\in C^{\prime}_{bad}$ from previous inequality. Thus
$c^{\prime}\notin C^{\prime}_{bad}$ and hence

$\displaystyle\sum_{p\in
P}\mathtt{d}(p,c^{\prime})\leq(1+\tfrac{4}{5}\varepsilon)\sum_{p\in
P}\mathtt{d}(p,c_{P}).$

We then conclude

$\displaystyle\sum_{p\in P}\mathtt{d}(p,c_{S})\leq\sum_{p\in
P}\mathtt{d}(p,c^{\prime})+n\mathtt{d}(c^{\prime},c_{S})\leq(1+\varepsilon)\sum_{p\in
P}\mathtt{d}(p,c_{P}).$

This event holds with probability at least $(1-\delta/2)^{2}>1-\delta$. ∎