New Frameworks for Offline and Streaming Coreset Constructions

Introduction

Coresets are an important technique in machine learning, data sciences, and statistics for representing a large dataset with a much smaller amount of memory. Coresets are often used as a pre-processing dimensionality technique to improve the downstream efficiency of algorithms, both space and time. Informally speaking, a coreset $S$ of an input set $P$ of underlying points $p_{1},\ldots,p_{n}$ is a smaller number of weighted representatives of $P$ that can be used to approximate the cost of any query from a set of a given queries. For example, in the common $k$-means clustering problem, the coreset must approximate $\sum_{i=1}^{n}d(p_{i},C)^{2}$ for every query $C$, where $C$ is a set of $k$ points and $d(p_{i},C)$ is taken to be the smallest Euclidean distance from $p_{i}$ to any point in $C$. Thus to use a coreset $S$ to approximately solve the $k$-means clustering problem, it suffices to find the optimal clustering on $S$ rather than find the optimal clustering on $P$. Because the size of $S$ is much smaller than the size of $P$, i.e., $|S|\ll|P|$, then finding an optimal clustering on $S$ instead of $P$ will be much more efficient.

More generally, coreset is a set of points $P^{\prime}$ with corresponding weight function $w(\cdot)$ such that $\sum_{p^{\prime}_{i}\in P^{\prime}}w(p^{\prime})d(p^{\prime}_{i},C)^{2}$ is a $(1\pm\epsilon)$ approximation to $\sum_{i=1}^{n}d(p_{i},C)^{2}$. Coresets have been extensively studied in $k$-means clustering [BHI02, HM04, FS05, FS08, FL11, FS12, FSS13, BLUZ19, HV20, FSS20], subspace approximation [DRVW06, DV07, FL11, FMSW10a, FSS13, CW15, SW18], and a number of other geometric problems and applications [AHY06, FFS06, Cla08, DDH+08, AB09, PT18, HJLW18, ABB+19, MSSW18, BDM+18, MOB+20], due to the increasing availability of big data and the necessity for scalable methods to process this information.

The most common algorithmic procedure to designing a coreset is the following simple template. An algorithm first approximately evaluates the sensitivity of each point in the dataset. Informally, the sensitivity of a point quantities how important or distinct that point is, with respect to the given objective function on which we would like to optimize. Approximating the sensitivity of each point can often be done efficiently, so that the time to construct a coreset is often a lower order term compared to the runtime of the post-processing algorithm. The template then samples a fixed number of points, so that each point in the dataset with probability proportional to the sensitivity of the point. This approach is called sensitivity sampling and the fixed number of points is often a monotonically increasing function of the total sensitivity, defined to be the sum of the sensitivities of each point. Hence, if the numbered of sampled points is much smaller than the number of input points, this approach allows for compact dimensionality reduction, leading to improved performance of post-processing algorithms.

Since the total sensitivity is a central quantity to coreset techniques, the total sensitivity for various objective functions has been well-studied and completely characterized in some cases. However, it is not quite known what the optimal dependency between the total sensitivity and the size of the coreset should be; that is, what is optimal monotonically increasing function of the total sensitivity that governs the number of sampled points? Clearly smaller functions lead to smaller coresets, which lead to more efficient post-processing functions. Many recent coreset constructions in the past decade require constructing coresets whose size depends quadratically on the total sensitivity. In this paper, we show this dependency is not optimal; we introduce in Theorem 1.1 a generic construction whose dependency on the total sensitivity $t$ is only $O(t\log t)$ rather than $O(t^{2})$ [FL11]. Because the dependency is already black-boxed into the design of many coreset constructions, our results automatically improve many existing coreset algorithms simply by lowering the number of required samples, without modifying any other property of the algorithm; we are only showing that the worst-case theoretical guarantee of these algorithms is significantly and universally better than previously thought.

We show that the common sensitivity sampling framework only needs to sample $O(t\log t)$ points, where $t$ is the total sensitivity.

Let $d$ be the dimension of a query space $(P,w,Q,f)$. For each point $p$, let $m(p)$ be an upper bound on the sensitivity of point $p$. Let $t(p)=\sum_{p\in P}m(p)$, and $\varepsilon,\delta\in(0,1)$. Then by sampling $O\left(\frac{t}{\varepsilon^{2}}\left(d\log t+\log\left(\frac{1}{\delta}\right)\right)\right)$ i.i.d. points from $P$ and rescaling each sampled point by $\frac{1}{m(p)}$, the resulting sample is an $\epsilon$-coreset for $P$.

In contrast, previous analysis showed that the sensitivity sampling framework required $O(t^{2})$ points [FL11]. We emphasize that our results are purely theoretical; we show that any worst-case guarantee that could previously be achieved with $O(t^{2})$ samples can actually be achieved with only $O(t\log t)$ samples. Hence our results can be universally plugged into any existing coreset construction algorithm simply by requiring a lower number of samples. Moreover, our results are optimal, since it can be shown by standard coupon-collector arguments that $\Omega(t\log t)$ samples are necessary in some cases.

We show that the results of [LLS01] also succeed if the VC-dimension of $\mathcal{F}$ is $d$, rather than the pseudo-dimension. This implies that the algorithm outputs a $(p,\epsilon)$-approximation of all functions in $\mathcal{F}$, which means if the function is too small, then the resulting data structure can only provide an additive error guarantee rather than a multiplicative relative error guarantee, but if the function is adequately large, then the resulting data structure provides a multiplicative error guarantee.

Fortunately, we show this guarantee suffices to obtain an $\varepsilon$-coreset. We break the query space into partitions, based on how much a point contributes to a query, compared to the total contribution to a query across all the points. If the contribution of a partition is large, then our $(p,\epsilon)$-approximation guarantees a good approximation to this partition. Now if the contribution of a partition is small, then two things can happen. Either it is possible that the sum of the contributions of all of the “small” partitions is large, in which case our $(p,\epsilon)$-approximation again guarantees a good approximation, or the sum of the contributions remains insignificant. In this case, we only have an additive approximation of the contributions for these points, but because the sum of the contributions is insignificant, an additive error on these points translates to a small relative error on the entire objective function.

Theorem 1.1 has applications to many problems in machine learning. In Section 3, we describe applications to model fitting problems. Specifically, we consider the $(j,k)$-projective clustering problems such as $k$-median/$k$-means, $k$-line clustering, $j$-subspace approximation, and the integer $(j,k)$-projective clustering problem.

Informally, the goal is to find a model $F$ in a restricted family $\mathcal{F}$ of set of $k$-tuples of affine $j$-subspaces that minimizes $\sum_{p\in P}d(p,F)$, where $P$ is a set of input points. Here, $F$ is the union of $k$ $j$-flats so that if $j=0$, then each $j$-flat reduces to a point and the $(j,k)$-projective clustering problem becomes the $k$-median problem with the appropriate metric. On the other hand, with an alternate distance function (the squared Euclidean distance), $d(\cdot,\cdot)$, the $(j,k)$-projective clustering objective becomes the $k$-means problem. When $j=1$ and $k$ is fixed, the objective becomes the $k$-line clustering problem but if $j$ is fixed and $k=1$, then the objective instead becomes the subspace approximation problem. Finally, in the integer $(j,k)$-projective clustering problem, all points in $P$ are assumed to have integer coordinates from some predetermined range. Our results subsume earlier versions online that have not received independent verification and are summarized in Figure 1.

Although our primary contribution is theoretical, we complement our worst-case guarantees with empirical evaluations on both small and large-scale datasets, which we describe in Section 5.

2 Preliminaries

Let $X$ be called a ground set. Let $P\subseteq X$ be a (possibly ordered) multi-set and $w:P\to[0,\infty)$ be a function that maps every $p\in P$ to a weight $w(p)\geq 0$. The pair $(P,w)$ is called a weighted set in $X$. If $w(p)=1$ for every $p\in P$ then the (un)weighted set $(P,w)$ may be denoted by $P$ for short.

The order of the points in $P$ can be arbitrary in this paper. However, even if $P$ contains only a single copy of each point, the corresponding coreset may contain multiple instances of some points. Hence, we consider coresets as multi-sets, although duplicated points can usually be replaced by a single weighted point without changing the claimed results. The union and intersection are also implied to be over multi-sets in this paper.

Usually the loss function has specific properties such as being a pseudo distance function $D$, as will be defined later. However, it may also be more complicated such as a subtraction between pseudo distance functions, which will also be used in this paper. This is also why it may return a negative number. In general, we will be interested in approximating $\sum_{p\in P}w(p)f(p,q)$ for every query $q$ in the query space up to an additive error of $\varepsilon$.

For a set $X$, a query function $Q$, and a cost function $f$, we define the VC-dimension of the range space that it induced, as defined below. The classic VC-dimension was defined for sets and subset and here we generalize it to query spaces, following [FL11].

The following definition of sensitivity is central to our paper, as we shall show that the coreset size of our algorithm is proportional to the total sensitivity of the input set.

Let $(P,w,Q,f)$ be a query space over a ground set $X$, where $w:P\to[0,\infty)$ and $f:P\times Q(P)\to[0,\infty)$. Then we define the sensitivity of a point $p\in P$ by $s(p)=\sup_{C\in Q}w(p)|f(p,C)|$.

Then the total sensitivity of an input set is the natural definition:

Let $(P,w,Q,f)$ be a query space over a ground set $X$, where $w:P\to[0,\infty)$ and $f:P\times Q(P)\to[0,\infty)$. We define the total sensitivity of $P$ by $\sum_{p\in P}s(p)$, where $s(p)$ is the sensitivity of $p$.

We next define two related concepts, the $(\nu,\alpha)$-samples and relative $(p,\varepsilon)$-approximations.

[LLS01] Let $\alpha,\nu>0$. For every $a,b\geq 0$, we define the distance function $d_{\nu}(a,b)=\frac{|a-b|}{a+b+\nu}$. Let $(P,w,Q,f)$ be a query space over a ground set $X$, where $w:P\to[0,\infty)$ and $f:P\times Q(P)\to[0,\infty)$. Then the weighted set $(S,u)$ is called a $(\nu,\alpha)$-sample for $(P,w,Q,f)$ if $(S,u,Q,f)$ is a query space, and for every $q\subseteq Q(S)$, $d_{\nu}(\overline{f}(P,w,;q),\overline{f}(S,u,q))\leq\alpha$, where $\overline{f}(P,w,q)=\frac{\sum_{p\in P}|w(p)\cdot f(p,q)|}{\sum_{p\in P}w(p)}$.

[HS11] Let $0<p,\varepsilon<1$. Let $(P,w,Q,f)$ be a query space over a ground set $X$, where $w:P\to[0,\infty)$ and $f:P\times Q(P)\to[0,\infty)$. Then the weighted set $(S,u)$ is called a $(\nu,\alpha)$-sample for $(P,w,Q,f)$ if $(S,u,Q,f)$ is a query space, and for every $q\subseteq Q(S)$,

$(1-\varepsilon)\overline{f}(P,w,q)\leq\overline{f}(S,u,q)\leq(1+\varepsilon)\overline{f}(P,w,q)$, for $\overline{f}(P,w,q)\geq p$

$\overline{f}(P,w,q)-\varepsilon p\leq\overline{f}(S,u,q)\leq\overline{f}(P,w,q)+\varepsilon p$, for $\overline{f}(P,w,q)\leq p$,

where $\overline{f}(P,w,q):=\frac{\sum_{p\in P}|w(p)\cdot f(p,q)|}{\sum_{p\in P}w(p)}$.

We recall the equivalence between $(\nu,\alpha)$-samples and relative $(p,\varepsilon)$-approximations:

[HS11] Let $(X,\mathcal{R})$ be a range space. If $(Z,\mathcal{R})$ is a $(\nu,\alpha)$-sample for $(X,\mathcal{R})$ with $0<\alpha<\frac{1}{4}$ and $\nu>0$, then $Z$ is a relative $(\nu,4\alpha)$-approximation for $(X,\mathcal{R})$.

Let $P^{\prime}=(P,w)$ be a weighted set in $X$, and $\varepsilon>0$ be an approximation error. The weighted set $C^{\prime}=(C,u)$ in $X$ is an $\varepsilon$-coreset for a query space $(P^{\prime},Q,f)$ if for every $q\in Q(C)$ we have $(1-\varepsilon)\sum_{p\in C}u(p)f(p,q)\leq\sum_{p\in P}w(p)f(p,q)\\ \leq(1+\varepsilon)\sum_{p\in C}u(p)f(p,q)$.

Sensitivity Sampling

In this section, we show that provable worst-case guarantees for constant factor approximation can be achieved using the sensitivity sampling framework to construct coresets of size $O(t\log t)$, where $t$ is the total sensitivityWe also achieve optimal dependence on $\frac{1}{\epsilon}$ for a $(1+\epsilon)$-approximation, but we omit these factors for ease of discussion. This improves on previous analysis that the sensitivity sampling framework to sample $O(t^{2})$ points to construct coresets that guaranteed constant factor approximation [FL11]. We remark that our result is purely theoretical and does not require novel algorithmic implementation. Instead, our result shows that the parameters in existing coreset construction algorithms can be improved while still guaranteeing worst-case performance.

Recall that the sensitivity of a point $p\in P$ is defined by $s(p)=\sup_{C\in Q(P)}w(p)|f(p,C)|$, where $w:P\to[0,\infty)$ is a weight function $p$ and $f:P\times Q(P)\to[0,\infty)$ is a loss function between an input point and a query set. However, determining the exact sensitivity of a point can be time-consuming, so we instead define $m(p)\geq s(p)$ to be an upper bound on the sensitivity. It turns out that upper bounds $m(p)$ that are within a constant factor approximation of the exact sensitivity of a point are often efficiently computable. Then $t:=t(P)=\sum_{p\in P}m(p)$ is an upper bound on the total sensitivity, which is the sum of the sensitivities of all points in $P$.

We now formalize the sensitivity sampling framework broadly used in algorithmic design. We form a sample $S$ by picking the first point of $S$ to be $p\in P$ with probability $\frac{m(p)}{t(P)}$ and reweighting the sampled point with the inverse of the sampling probability. We show that repeatedly sampling points from $P$ with replacement until $S$ has $\frac{ct}{\varepsilon^{2}}\left(d\log t+\log\left(\frac{1}{\delta}\right)\right)$ points suffices to obtain an $\varepsilon$-coreset for $P$ with probability $1-\delta$ if the underlying query space has VC dimension $d$. The sensitivity sampling framework appears in full in Algorithm 1.

We first recall the following definition of pseudo-dimension:

[LLS01] show that if $\mathcal{F}$ has pseudo-dimension $d$, then $O\left(\frac{1}{\alpha^{2}\nu}\log\frac{1}{\nu}\right)$ samples suffices to simultaneously obtain an $(\nu,\alpha)$-sample to expectation of all functions in $\mathcal{F}$ with constant probability. Namely, [LLS01] show the following two lemmas:

Let $\mathcal{F}$ be a set of functions from $X$ to $$,$\mu$be a probability distribution over$X$and$\nu>0$,$0<\alpha<1$, and$N\geq\frac{2}{\alpha^{2}\nu}$. For any integer$N>0$, let$\Gamma_{N}$denote the set of all permutations of$\{1,\ldots,2N\}$so that for each$i\leq N$, either$i$and$N+i$are fixed, or$i$and$N+i$are swapped. Let$U$be the uniform distribution over$\Gamma_{m}$. Then

Let $d$ be the pseudo-dimension of $F\subseteq^{2N}$, where $N\geq\frac{125(2d+1)}{\alpha^{2}\nu}$ for any $\alpha,\nu>0$. Let $U$ be the uniform distribution over $\Gamma_{N}$. Then

Observe that combining Lemma 2.2 and Lemma 2.3 and solving for $N$ recovers the bound from [LLS01] of $O\left(\frac{1}{\alpha^{2}\nu}\left(d\log\frac{1}{\nu}+\log\frac{1}{\delta}\right)\right)$. We need an analog of their sampling result for VC-dimension rather than pseudo-dimension. As it turns out, the only place [LLS01] uses pseudo-dimension in Lemma 2.3 is a black-box reduction from the following lemma to bound the size of $\mathcal{F}$:

Moreover, [Hau95] proved the exact same statement when $d$ is the VC-dimension of $F$, rather than the pseudo-dimension.

Specifically, Lemma 2.5 follows from Corollary 1 in [Hau95] because $e(d+1)(2e/\varepsilon)^{d}\leq(2e^{3}/\varepsilon)^{d}<(41/\varepsilon)^{d}$.

Thus by using Lemma 2.5 rather than Lemma 2.4, we can recover Lemma 2.3 using VC-dimension rather than pseudo-dimension in the following formulation of Theorem 2.6. Hence, we can relate the sampling complexity of learning a class of functions to their VC-dimension:

2 Reduction to ε𝜀\varepsilon-Coresets

We now show that an $(\nu,\alpha)$-sample to a class of functions $\mathcal{F}$ suffices to achieve an $\varepsilon$-coreset under the appropriate parameters. The proof partitions the points in an input set $P$ by their contribution to $\sum_{p\in P}w(p)f(p,X)$ for some $X$ in the query space. A subset $S_{i}$ that contributes a large fraction towards $\sum_{p\in P}w(p)f(p,X)$ will be well-estimated by the $(\nu,\alpha)$-sample. On the other hand, if $S_{i}$ is not well-estimated by the $(\nu,\alpha)$-sample, then its contribution towards $\sum_{p\in P}w(p)f(p,X)$ must be small, so that intuitively, the additive error from the $(\nu,\alpha)$-sample is also small. Thus we can show that the sample is actually an $\varepsilon$-coreset.

Let $d$ be the dimension of a query space $(P,w,Q,f)$. Suppose that $m:P\to[0,\infty)$ such that $m(p)\geq\sup_{C\in Q(P)}w(p)|f(p,C)|$. Let $t\geq\sum_{p\in P}m(p)$, and $\varepsilon,\delta\in(0,1)$. Let $c\geq 1$ be a sufficiently large constant, and let $S$ be a sample of

Since $Z$ is $\left(\frac{1}{t},\varepsilon\right)$-approximation for $(X,\mathcal{R})$ and $A_{2}$ only consists of indices $i$ such that $S_{i}$ contributes at least $\frac{1}{t}$ fraction of the mass, then

Now if $\sum_{i=\in A_{1}}\mu_{i}\geq\frac{1}{t}$, then since $Z$ is a $\left(\frac{1}{t},\varepsilon\right)$-approximation for $(X,\mathcal{R})$, then

Hence combining with (1) and noting that $P=P_{1}\cup P_{2}$, then we have

On the other hand, if $\sum_{i=\in A_{1}}\mu_{i}<\frac{1}{t}$, then since $Z$ is a $\left(\frac{1}{t},\varepsilon\right)$-approximation for $(X,\mathcal{R})$, it follows that

Since we have $t\sum_{p\in P}m(p)\geq\sum_{p\in P}s(p)=\sum_{p\in P}\sup_{C\in Q(P)}w(p)|f(p,C)|\geq 0$, then $\frac{\varepsilon}{t}\leq\frac{\varepsilon}{\sum_{p\in P}s(p)}$ and thus

Moreover, $\sum_{i=\in A_{1}}\mu_{i}<\frac{1}{t}$ implies that $\frac{\sum_{p\in S_{i}}w(p)f(p,X)}{\sum_{p\in P}s(p)}<\frac{\varepsilon}{t}\leq\frac{\varepsilon}{\sum_{p\in P}s(p)}$. Therefore, we have $\sum_{p\in S_{i}}w(p)f(p,X)<1$, so that

Again combining with (1) and noting that $P=P_{1}\cup P_{2}$, then

Note that $P_{2}$ contains at most $\log t$ sets $S_{i}$, it suffices to obtain $\left(\frac{1}{t},\varepsilon\right)$-sample for $O(\log t)$ sets, each with failure probability $\Theta\left(\frac{\delta}{\log t}\right)$, using Theorem 2.6. By a union bound, the total failure probability is at most $\delta$. $\Box$

The intuition is that the class of functions $\mathcal{F}$ represents the objective in the query space, so that a particular $f\in\mathcal{F}$ represents the objective for a particular query in the query space. For objectives like $k$-means or $k$-median clustering, each $f$ represents the objective for a separate set of $k$ centers. Then the goal is to learn $\mathcal{F}$ simultaneously with a small number of samples.

The domain $X$ for the class of functions $\mathcal{F}$ translates exactly to the ground set $X$, which is the input points for objectives like $k$-means or $k$-medians. We first note that sampling a point of $X$ and then rescaling by the (inverse of the) sampling probability provides an unbiased estimator to the objective. [LLS01] then states that if we sample uniformly over $X$, we can obtain a $(1+\varepsilon)$-approximation to the objective by bounding the variance through a small number of samples. Then the idea of sensitivity sampling is that instead of uniformly sampling points from $X$, we sample each point of $X$ according to its sensitivity, but still rescale by the (inverse of the) sampling probability. Now the expectation of the samples is still the objective, but the variance is much smaller and so we require a smaller number of samples.

The real workhorse in this bound is Lemma 2.2 by [LLS01], which uses the chaining technique of Kolmogorov and refined by Talagrand [Tal94]. Crucially, the usage of chaining by [LLS01] manages to simultaneously learn a large number of functions in a class $\mathcal{F}$ without needing to union bound over a net over the functions in $\mathcal{F}$. It is precisely this technique that avoids a quadratic dependency on $t$ from the union bound.

Applications

Our theoretical worst-case guarantee has a wide range of applications due to the prevalence of the coreset technique and how well-studied the total sensitivity is of various optimization problems. Note that for any problem whose sensitivity is known to be $t$, Theorem 1.1 gives an improvement on dependency of $t$ from $O(t^{2})$ to $O(t\log t)$ for the number of sampled points. [VX12b] considers sensitivity sampling for shape fitting problems, focusing on $(j,k)$-projective clustering problems, such as $k$-median/$k$-means, $k$-line clustering, $j$-subspace approximation, and the integer $(j,k)$-projective clustering problem. We show that our results imply more efficient coreset constructions using the total sensitivity bounds on these problems obtained by [VX12b].

For the $k$-line center problem and the integer general $(j,k$)-projective clustering problem, [VX12a] showed the following upper bounds on the total sensitivity:

Then Theorem 3.2 and Theorem 1.1 together imply efficient coresets for both the $k$-line center problem and the integer $(j,k)$-projective clustering problem.

We remark that Theorem 3.4 is subsumed by Theorem 3.6 below, as [VX12b] tighten the total sensitivity upper bound for the $k$-line center problem by showing that $f(d,k)$ in Theorem 3.2 is independent of $d$.

From Theorem 3.5 and Theorem 1.1, we have

[VX12b] bounded the total sensitivity for the $(0,k$)-projective clustering problem, which includes $k$-median and $k$-means.

For either the $k$-median problem or the $k$-means problem, there exists an algorithm that outputs a set of $O\left(\frac{d}{\varepsilon^{2}}k\log k\right)$ weighted points that is an $\epsilon$-coreset, with probability at least $\frac{2}{3}$.

[VX12b] also bounded the total sensitivity for the $j$-subspace fitting problem.

By Theorem 3.9 and Theorem 1.1, we conclude that

Coreset for k𝑘k-clustering

This section considers tighter bounds for $k$-clustering and may be skipped for general applications. We introduce $\rho$-pseudo distances and define the importance of a point as a generalization of the sensitivity. Using the importance, we then give an analog of Theorem 1.1 with sharper bounds for $\rho$-pseudo distances. As a result, we obtain stronger bounds for coreset constructions for $k$-clustering.

We reiterate that our results in this section mirror those of Section 2. We only provide theoretical guarantees on the number of samples required by the sensitivity sampling framework for $k$-clustering. We specify the framework in full in Algorithm 2.

We first require the following definition of $(\alpha,\beta)$-assignment for bicriteria algorithms.

Let $X$ be a ground set and $(P,w,Q,g)$ be a query space where $g:X^{2}\to[0,\infty)$. Let $Q^{*}\in Q(P)$ be a query that minimizes the loss of $P$ over every query,

Let $\alpha,\beta>0$, and $B\subseteq X$ such that $|B|\leq\beta|Q^{*}|$. A function $\mathcal{B}:P\to B$ is an $(\alpha,\beta)$-assignment for $(P,w,Q,g)$ if

Every $b\in B$ is called a center and its cluster in $P$ is $\mathcal{B}^{-1}(b)=\left\{p\in P\mid\mathcal{B}(p)=b\right\}$.

Intuitively, an $(\alpha,\beta)$-assignment is just a bicriteria clustering with approximation factor $\alpha$ and an overselection of centers by a factor of $\beta$. Thus for $k$-means clustering, any $\alpha$-approximation algorithm that chooses $\beta k$ centers can be used to determine an $(\alpha,\beta)$-assignment for each of the points.

The following definition is especially useful for $k$-means clustering.

Let $X$ be a ground set and $\rho\geq 1$. A symmetric function $g:X^{2}\to[0,\infty)$ is $\rho$-pseudo distance over $X$ if for every $(p,q,x)\in X^{3}$

For a finite set $Q\subseteq X$, we denote $g(p,Q):=\min_{x\in Q}g(p,x)$.

Note that the above definition does not assume that $g(p,p)=0$ for every $p\in P$. The inequality in Definition 4.2 is sometimes called “weak triangle inequality”.

Let $g$ be a $\rho$-pseudo distance over a ground set $X$. For every pair of points $p,q\in X$ and a finite set $Q\subseteq X$,

Proof : For every $p,q\in X$, and a center $x_{p}\in Q$ that is closest to $p$, i.e. $g(p,Q)=g(p,x_{p})$, we have

We now give a generalization of the notion of sensitivity in the form of importance.

Let $(P,w,Q,g)$ be a query space where $w:P\to[0,\infty)$ and $g$ is a $\rho$-pseudo distance. Let $\mathcal{B}:P\to B$ be an $(\alpha,\beta)$-assignment for $(P,w,Q,g)$. For every $p\in P$ and $Q\in Q(P)$, we define

For every center $b\in B$ and a point $p\in\mathcal{B}^{-1}(b)$ in its cluster we define

We now show that the importance $m$ satisfies a similar function as the notion of sensitivity.

Let $\mathcal{B}:P\to B$, $m$ and $f$ be defined as in Definition 4.4. Then for every $p\in P$ we have

Proof : For simplicity, we denote $p^{\prime}=\mathcal{B}(p)$ and $q^{\prime}=\mathcal{B}(q)$ for every $p,q\in P$, and $P_{b}=\mathcal{B}^{-1}(p)$ for every $b\in B$. We prove the claim for a point $p\in P_{b}$ in the cluster of some center $b\in B$, as in Definition 4.4. Let $Q\in Q(P)$ and assume $w(p)f(p,Q)>0$, otherwise the lemma trivially holds. We need to upper bound

where the first inequality holds by Lemma 4.3, and the second inequality holds since $\mathcal{B}$ is an $(\alpha,\beta)$-assignment. To bound the last term, note that

where the first inequality is by Lemma 4.3, (4) holds since $p^{\prime}=b$ and since $g$ is symmetric by definition, and (5) holds since $\mathcal{B}$ is an $(\alpha,\beta)$-assignment.

Dividing by $\sum_{q\in P_{b}}w(q)\cdot\sum_{q\in P}w(q)g(q,Q)$ yields

Substituting this in (3) yields the desired result

As a warm-up, we now prove an analog of Theorem 1.1 for $\rho$-pseudo distances that provides tighter bounds, due to the tighter setting of $t$ from the $(\alpha,\beta)$ assignment.

$(P,w,Q,g)$ be a query space, where $g$ is a $\rho$-pseudo distance and $w:P\to[0,\infty)$.

$\mathcal{B}:P\to B$ be an $(\alpha,\beta)$ assignment for $(P,w,Q,g)$.

$d$ be the VC-dimension of $(P,w,Q,f)$, where $f$ was defined in (2).

$c\geq 1$ be a sufficiently large constant, $\varepsilon,\delta\in(0,1)$, and

$(C,u)$ be the output of a call to $\textsc{Coreset}(P,w,\mathcal{B},s)$; see Algorithm 2.

Then, $C\subseteq P$, $u:C\to[0,\infty)$ and, with probability at least $1-\delta$, $(C,u)$ is an $\varepsilon$-coreset of size $|C|=s$ for $(P,w,Q,g)$.

Proof : Let $p\in P$. By Lemma 4.5, for every $Q\in Q(P)$

The probability of choosing $p$ to be, say, the first point in $C$, is

Using the last inequality and (6), we apply Theorem 2.7 to obtain that, with probability at least $1-\delta$, we have that for all $Q\in Q(C)$,

implies that $(C,u)$ is an $\varepsilon$-coreset as desired. $\Box$

In this section, we define a $(\rho,\psi,\phi)$-pseudo distance function, which serves as a generalization of $\rho$-pseudo distances, and give smaller coreset constructions for $(\rho,\psi,\phi)$-pseudo distance functions. Our algorithms appear in Algorithm 3 and Algorithm 4.

We define the following generalization of distance to handle $k$-clustering.

Let $g$ be a $\rho$-pseudo distance over $X$ as in Definition 4.2. For $\phi>0$ and $\psi\geq 0$, $g$ is also a $(\rho,\psi,\phi)$-pseudo distance function if for every $(p,q,x)\in X^{3}$ we have

Intuitively, the $\rho$-pseudo distance handles loss functions such as the squared Euclidean distance that do not satisfy the triangle inequality but rather a generalized version of the triangle inequality.

Let $g:X^{2}\to[0,\infty)$ be a $(\rho,\psi,\phi)$-pseudo distance function. Then for every finite set $M\subseteq X$ and $p,q\in M$ we have

where the first inequality is by the definition of $g(p,Q)=\min_{x\in Q}g(p,Q)$, and the second inequality is by Definition 4.7.

where the last inequality is by the assumption of this case. Combining (7) and (8) yields that (in both cases)

Let $(P,w,Q,g)$ be a query space where $w:P\to[0,\infty)$, and $g:X^{2}\to[0,\infty)$ is a $(\rho,\psi,\phi)$-pseudo distance function. Let $\mathcal{B}:P\to B$ be an $(\alpha,\beta)$-assignment for $(P,w,Q,g)$. Let $\varepsilon\in(0,1)$ such that $\psi<\varepsilon/(4\rho(\alpha+1))$, and

Consider the variables in Definition 4.9. Let

and $d$ be the VC-dimension of $(P,w,Q,h)$. Let $c^{\prime}$ be a sufficiently large constant,

and $(C\cup B,u)$ be the output of a call to algorithm $\textsc{Smaller-Coreset}(P,w,\mathcal{B},s)$; see Algorithm 3. Then, $C\subseteq P$, $u:C\to[0,\infty)$, and with probability at least $1-\delta$, we have that for all $X\in Q(C)$,

Proof : Let $Q\in Q(C)$ and extend the function $u$ as defined in Algorithm 3 to be $u(p)=0$ for every $p\in P\setminus C$. Also define $v(p)=w(p)-u(p)$ and $p^{\prime}=\mathcal{B}(p)$ for every $p\in P$. The difference in the loss between taking the original points or its coreset $C\cup B$ is

where (12) is by the triangle inequality, and (13) holds since

where the last equality is by the definition of $u$ in Line 3 of Algorithm 3.

By letting $F=F(Q)$ as defined in (10), (13) is bounded by $\left|\sum_{p\in P}v(p)(g(p,Q)-g(p^{\prime},Q))\right|$, which equals

where the last inequality is by the triangle inequality. We now bound each of the last terms.

where (16) holds since $\mathcal{B}$ is an $(\alpha,\beta)$-assignment, i.e.,

(17) is by Lemma 4.8, and (18) is by (10) and the assumption $p\in P\setminus F$.

$t\geq 2$, and for every constant $c>0$ there is a sufficiently large constant $c^{\prime}$ such that

Plugging these bound in Theorem 2.7 with $\varepsilon/4$, and the query space $(P,w,Q,h)$ yields that, with probability at least $1-\delta/2$, we have that for all $Q\in Q(C)$,

Assume that the last equation indeed holds (which happens with probability at least $1-\delta/2$). By this and the definition of $g$, for every $Q\in Q(C)$, (14) is bounded by

Bound on (15): Since $|v(p)|=|w(p)-u(p)|\leq w(p)+u(p)$, and using the triangle inequality

where the first inequality is by Lemma 4.8, the second holds since $p\in F$, and the last inequality is by (9).

Our $(\alpha,\beta)$-assignment approximates the sum of distances to a query up to an additive error as follows.

where (23) is by Lemma 4.3, and (24) is by (19). Hence,

where (25) is by (21), (26) holds by (23), (27) holds by (22), and (28) holds since $F\subseteq P$.

It is left to bound the rightmost term in (28). Let $z:P\times B\to[0,\infty)$ such that for every $b\in B$ and $p\in P_{b}$

where (30) holds since $P=\bigcup_{b\in B}P_{b}$, in (31) we simply multiplied and divided by $\sum_{q\in P_{b}}w(q)g(q^{\prime},Q)$, and (32) holds since $p^{\prime}=q^{\prime}$ for every $p,q\in P_{b}$.

where (33) holds since $v(p)=w(p)-u(p)$, and (34) is by definition (29) of $z$.

Let $b\in B$, $t^{\prime}=2|B|$, and for every $p\in P$, let $m(p)=w(p)z(p,b)$. Hence, for every $p\in P$,

and for every constant $c\geq 1$ there is a sufficiently large $c^{\prime}$ such that

Substituting the query space $(P,w,\left\{b\right\},z)$, $\varepsilon=1/2$, $d=1$, and $\delta/|B|$ instead of $\delta$ in Theorem 2.7, yields that with probability at least $1-\delta/(2|B|)$, we have

Assume the event that (35) holds for every $b\in B$ occurs, which happens with probability at least $\delta/2$, by the union boundInstead of using the union bound, we could simply choose $B$ as the set of queries, $\delta$ instead of $\delta/(2|B|)$ and $d=\log|B|$. However, in this would introduce a term of $d\log t=O(\log^{2}|B|)$ in the coreset size compared to the current $\log|B|$ term.. Plugging (35) in (34) yields

Combining the last inequalities bounds (15) with probability at least $1-\delta/2$, as

where (37) holds by (28), (38) by (30), and (39) by (36).

Finally, replacing (14) and (15) with (20) and (39) respectively, proves that, with probability at least $1-\delta/2-\delta/2=1-\delta$ we have

By this and (13), it follows that $(C,u)$ approximates $X$ as desired. $\Box$

We now handle the specific case where $g:X^{2}\to[0,\infty)$ is a $(\rho,\psi,\phi)$-pseudo distance function.

Consider the variables in Theorem 4.10, where $s$ is replaced by

Let $(C,u)$ be the output of a call to algorithm $\textsc{Coreset}(P,w,\mathcal{B},s)$; see Algorithm 3.

Then, $C\subseteq P$, $u:C\to[0,\infty)$, and with probability at least $1-\delta$, $(C,u)$ is an $\varepsilon$-coreset of size $s$ for $(P,w,Q,g)$.

Proof : Let $\varepsilon^{\prime}=\varepsilon/(1+\rho(\alpha+1))$. After replacing $\delta$ with $\delta/2$ and $\varepsilon$ with $\varepsilon^{\prime}$ in Theorem 4.10, we obtain that with probability at least $1-\delta/2$,

Assume that this event indeed occurs and the inequality holds, and let $Q\in Q(C)$.

We will bound the error by excluding $B$ from this coreset, i.e.,

where (41) is by the triangle inequality, and (42) is by (40). The rightmost term is

where (43) is by the definition of $u$ in Line 3 of Algorithm 3.

The bound on the rightmost term is similar to (35), after replacing the bound $1/2$ with $\varepsilon^{\prime}$, which is the reason for the largest size $s$ of the coreset. Specifically, let $b\in B$, $t^{\prime}=2|B|$, and for every $p\in P$, let $m(p)=w(p)z(p,b)$, where $z$ is defined in (29). Hence, for every $p\in P$,

and for every constant $c\geq 1$ there is a sufficiently large $c^{\prime}$ such that

Substituting the query space $(P,w,\left\{b\right\},z)$, $d=1$, and $\delta/|B|$ instead of $\delta$ in Corollary 2.7, yields that with probability at least $1-\delta/(2|B|)$, we have

Assume the event that (35) holds for every $b\in B$ indeed occurs, which happens with probability at least $\delta/2$. Substituting the value of $z(p,b)$ from (35) and multiplying by $\sum_{q\in P_{b}}w(q)$ yields

where (46) is by (45), and (47) is by the property of $(\alpha,\beta)$-assignment in (24).

Combining the previous inequalities all together yields the desired result

where (48) is by (42), and (49) is by (47).

Using the union bound on previous assumptions, this holds with probability at least $1-\delta/2-\delta/2=1-\delta$. $\Box$

2 Positively weighted coresets

In this section, we give a construction for a coreset that is guaranteed to output positive weights associated with each sampled point.

Consider the variables in Theorem 4.10. Let

and let $(C\cup B,u^{\prime})$ be the output of a call to algorithm $\textsc{Conditional-Coreset}(P,w,\mathcal{B},s,\varepsilon^{\prime})$; see Algorithm 4. Then, $C\subseteq P$, $u:C\times Q(C)\to[0,\infty)$, and with probability at least $1-\delta$, we have that for all $X\in Q(C)$,

Proof : The proof is the same as the proof of Theorem 4.10 except for replacing $u(p)$ with $u^{\prime}(p,Q)$ everywhere, and replacing the bound on (25) by

where the first inequality is similar to (27), and the equality is since $u^{\prime}(p,Q)=0$ for every

Empirical Evaluations

This concludes our discussion of the general sensitivity sampling framework. Although our contribution is primarily theoretical, we nevertheless performed empirical evaluations in Python 3.6 via the Numpy and Scipy.sparse libraries on a desktop machine with an Intel i7-6850K CPU @ 3.60GHZ, 64GB RAM. We consider coreset constructions based on sensitivity sampling for bicriteria algorithms (Algorithm 1), general loss functions that satisfy the weak triangle inequality (Algorithm 2), and the conditional normalized distance (Algorithm 3). We compared Algorithms 1-3 to uniform sampling on $k$-means clustering on both relatively small offline data and large-scale streaming data that cannot fit into memory. Algorithms 1-3 each require a bicriteria algorithm to approximate the importance of each point; we use kmeans++ with $\alpha=O(\log k)$ and $\beta=1$ to approximate the importances, so that the runtime is linear.

For experiments on small offline datasets, we compared our coreset constructions for $k$-means clustering in Algorithms 1-3 vs. uniform sampling on the datasets: (i) Gyroscope data and (ii) Accelerometer data. Collected by [AGO+13b], and can be found on [AGO+13a], the experiments have been carried out with a group of 30 volunteers within an age bracket of 19-48 years. Each person performed six activities (walking, walking upstairs, walking downstairs, sitting, standing, laying) while wearing a Samsung Galaxy S II smartphone on the waist. Using its embedded gyroscope (resp. accelerometer), 3-axial angular velocity (resp. linear acceleration) were captured at a constant rate of 50Hz. The experiments have been video-recorded to label the data manually. Data was collected from $n=7352$ measurements; each instance consists of measurements from $d=3$ dimensions: $x$, $y$, $z$, each in a size of 128.

We ran Algorithms 1-3 and uniform sampling on the above six datasets with different sample/coreset size, between 1000 to 7000, with $k=100$ and $k=200$. The multiplicative approximation error (empirical $\varepsilon$) was calculated by $\varepsilon:=\frac{\text{\footnotesize{COST}}(A,Q_{C})-\text{\footnotesize{COST}}(A,Q_{A})}{\text{\footnotesize{COST}}(A,Q_{A})}$, where $A$ is the matrix whose rows correspond to the $n$ input points, $Q_{A}$ corresponds to the $k$ centers of the whole data (two Lloyd’s iterations after kmeans++ initialization) and $Q_{C}$ is the clustering of the coreset. Our results show a significant improvement of our algorithms over uniform sampling; we present the gyroscope data evaluations in Figure 2 and the accelerometer data evaluations in Figure 3.

2 Evaluations on Streaming Data

To handle large-scale streaming data that cannot fit into memory, our system separates the $n$ points of the data into chunks of a desired size of coreset, called $m$. We use a merge-and-reduce framework on a binary tree, e.g. [FMSW10b], where each node is a coreset of the union of the data represented by its children nodes and the bottom layer of the tree consists of consecutive chunks of the data of size $4516$. Thus the root of the tree is a coreset of the whole data. We build a tree of height $10$ for our data, dividing the $n=4624611$ input points across $1024$ chunks of size $4516$.

We compared uniform sampling to Algorithms 1 and 3 for $k$-means clustering on a created document-term matrix of Wikipedia (parsed enwiki-latest-pages-articles.xml.bz2-rss.xml from [wic19]), i.e. sparse matrix with 4624611 rows and 100k columns where each cell $(i,j)$ equals the value of how many appearances the word number $j$ has in article number $i$. We use a standard dictionary of the 100k most common words in Wikipedia [Dic12]. We concatenated the coreset received in each floor and compared the received error in each floor. The error we determined was calculated by the formula $\frac{\text{\footnotesize{COST}}(A,Q_{C})-\text{\footnotesize{COST}}(A,Q_{A})}{\text{\footnotesize{COST}}(A,Q_{A})}$, where $A$ is the original data matrix, $Q_{A}$ is the clustering of the whole data (Lloyd’s iterations until 1% convergence, after ++ initialization) and $Q_{C}$ is the clustering of the coreset. We used two values of $k$, 100 and 200. We present our results in Figure 4. We concatenated the coreset received in each floor and compared the received error in each floor. The error we determined was calculated by the formula $\frac{\text{\footnotesize{COST}}(A,Q_{C})-\text{\footnotesize{COST}}(A,Q_{A})}{\text{\footnotesize{COST}}(A,Q_{A})}$, where $A$ is the original data matrix, $Q_{A}$ is the clustering of the whole data (Lloyd’s iterations until 1% convergence, after ++ initialization) and $Q_{C}$ is the clustering of the coreset. We present our results in Figure 4 for $k=100$ and $k=200$. Similar to the offline evaluations, we obtain better results for our algorithms than uniform sampling. However, unlike than the offline data, here the conditional normalized algorithm gets much better results than the general sensitivity sampling algorithm.

Algorithms.

The algorithms we compared are uniform sampling, and our Algorithm 2 and 4.

Results.

We concatenated the coreset received in each floor and compared the received error in each floor. The error we determined was calculated by the formula $\frac{\text{\footnotesize{COST}}(A,Q_{A})-\text{\footnotesize{COST}}(A,Q_{C})}{\text{\footnotesize{COST}}(A,Q_{A})}$, where $A$ is the original data matrix, $Q_{A}$ is the clustering of the whole data (Lloyd’s iterations until 1% convergence, after ++ initialization) and $Q_{C}$ is the clustering of the coreset. We used two values of $k$, 100 and 200. We present our results in Figure 4.

Discussion.

Indeed also for this dataset we got better results for our algorithm than uniform sampling. However, unlike than in Section 5, here Algorithm 2 gets much better results than Algorithm 1.

References