Coresets for Minimum Enclosing Balls over Sliding Windows

TL;DR

This paper introduces efficient algorithms for maintaining coresets for the minimum enclosing ball over sliding windows in data streams, significantly improving speed while preserving approximation quality.

Contribution

The authors develop the first practical algorithms for sliding-window coresets for MEB, including SWMEB and SWMEB+, with extensions to reproducing kernel Hilbert spaces.

Findings

Achieve up to 10,000x speedup over batch algorithms.

Maintain small approximation errors for MEB coresets.

Support for kernels in reproducing kernel Hilbert spaces.

Abstract

\emph{Coresets} are important tools to generate concise summaries of massive datasets for approximate analysis. A coreset is a small subset of points extracted from the original point set such that certain geometric properties are preserved with provable guarantees. This paper investigates the problem of maintaining a coreset to preserve the minimum enclosing ball (MEB) for a sliding window of points that are continuously updated in a data stream. Although the problem has been extensively studied in batch and append-only streaming settings, no efficient sliding-window solution is available yet. In this work, we first introduce an algorithm, called AOMEB, to build a coreset for MEB in an append-only stream. AOMEB improves the practical performance of the state-of-the-art algorithm while having the same approximation ratio. Furthermore, using AOMEB as a building block, we propose two novel algorithms, namely SWMEB and SWMEB+, to maintain coresets for MEB over the sliding window with constant approximation ratios. The proposed algorithms also support coresets for MEB in a reproducing kernel Hilbert space (RKHS). Finally, extensive experiments on real-world and synthetic datasets demonstrate that SWMEB and SWMEB+ achieve speedups of up to four orders of magnitude over the state-of-the-art batch algorithm while providing coresets for MEB with rather small errors compared to the optimal ones.

Figures

Tables

Table 1. Table 1: Statistics of datasets

dataset	source	size	$m$	$\theta$	$\gamma$
Census	UCI	2,458,285	68	228.5	8984.58
CovType	LIBSVM	297,711	54	246.8	3.734
GIST	TEXMAX	1,000,000	960	45.77	4.0409
Gowalla	SNAP	6,442,892	2	$\approx$3000	7455.33
HIGGS	UCI	5,829,123	28	13.69	38.8496
SIFT	TEXMAX	1,000,000	128	12.14	298919.5
Synthetic	-	10,000,000	50	2.96	100.6

Table 2. Table 2: The average errors of different algorithms for Euclidean MEB

Algorithm	average error
	Census	CovType	GIST	Gowalla	HIGGS	SIFT	Synthetic
COMEB	0	0	0	0	0	0	0
CoreMEB	6.03e-04	3.66e-04	4.32e-04	2.01e-04	4.59e-04	4.97e-04	4.49e-04
SSMEB	5.90e-02	1.56e-01	1.10e-01	7.19e-03	1.59e-01	1.61e-01	9.49e-02
BBC	1.59e-03	8.05e-03	4.47e-04	1.53e-03	1.33e-02	5.83e-03	2.05e-02
DyMEB	9.14e-04	3.27e-03	4.51e-04	1.98e-04	2.19e-03	4.34e-03	5.96e-03
AOMEB	9.55e-04	4.55e-03	4.55e-04	2.23e-04	4.73e-03	2.86e-03	1.14e-02
SWMEB	6.12e-04	8.36e-03	2.31e-04	2.57e-04	3.87e-03	2.94e-03	1.11e-02
SWMEB+	8.98e-04	4.17e-03	1.81e-03	1.93e-04	9.67e-03	3.01e-03	1.52e-02

Table 3. Table 3: The average update time of different algorithms for Euclidean MEB

Algorithm	average update time (ms)
	Census	CovType	GIST	Gowalla	HIGGS	SIFT	Synthetic
COMEB	1150.7	1914.8	3813.7	13.52	146.67	12861.8	2793.7
CoreMEB	643.1	1275.1	2647.2	30.19	337.72	2789.1	1698.8
SSMEB	43.75	44.11	262.8	24.55	45.66	74.16	54.46
BBC	193.4	767.4	961.8	37.61	118.19	6807.2	1150.4
DyMEB	1639.4	2725.9	7498.8	80.32	3483.2	9195.1	2786.3
AOMEB	95.74	507.3	1755.5	17.82	125.48	1413.9	480.63
SWMEB	2.129	25.08	127.77	0.1593	4.861	57.61	19.9
SWMEB+	1.467	5.414	72.31	0.1887	2.592	14.37	5.679

Table 4. Table 4: The average errors of different algorithms for kernelized MEB

Algorithm	average error
	Census	CovType	GIST	Gowalla	HIGGS	SIFT	Synthetic
CoreMEB	9.31e-05	9.16e-05	9.76e-05	7.97e-05	9.64e-05	9.40e-05	9.53e-05
SSMEB	2.22e-01	1.73e-01	1.83e-01	2.26e-01	2.08e-01	1.67e-01	1.40e-01
BBC	1.21e-03	1.55e-03	1.78e-05	1.15e-02	9.02e-04	1.60e-03	1.48e-03
DyMEB	1.62e-04	6.07e-05	1.18e-04	2.07e-03	1.57e-04	9.91e-05	4.26e-05
AOMEB	5.63e-04	6.02e-04	1.62e-06	1.02e-04	6.08e-04	7.55e-04	4.59e-04
SWMEB	7.77e-04	4.71e-03	5.01e-03	1.19e-04	5.82e-03	3.12e-03	4.99e-03
SWMEB+	5.86e-04	1.01e-03	4.85e-04	1.67e-04	2.10e-03	6.78e-04	1.61e-03

Table 5. Table 5: The average update time of different algorithms for kernelized MEB

Algorithm	average update time (ms)
	Census	CovType	GIST	Gowalla	HIGGS	SIFT	Synthetic
CoreMEB	75041.6	91571.8	36522	5813.1	87875.1	149386	76223.1
SSMEB	208.67	251.12	2893.9	118.71	169.92	395.5	255.25
BBC	14914.5	34728.1	704506	410.30	59789.6	148539	44383.9
DyMEB	106093	256616	2484757	7142.1	345400	1333728	513587
AOMEB	2002.1	5716.6	42961.2	253.87	6471.1	19163.7	10556.2
SWMEB	146.09	545.9	4550.3	18.118	606.2	1705.7	1077.2
SWMEB+	26.847	52.06	372.4	4.769	55.85	130.46	87.21

Equations

$$r_t\ge d({\mathbf{c}}_t,{\mathbf{p}}_t)\ge d({\mathbf{p}}_t,{\mathbf{c}}_{t-1})-d({\mathbf{c}}_t,{\mathbf{c}}_{t-1})\ge (1+\frac{{\epsilon }_1}{2})r_{t-1}$$

$$r_t\ge d({\mathbf{c}}_t,\mathbf{q})\ge \sqrt{d^2({\mathbf{c}}_{t-1},\mathbf{q})+d^2({\mathbf{c}}_t,{\mathbf{c}}_{t-1})}>\sqrt{1+(\frac{{\epsilon }_1}{2}{)}^2}\cdot r_{t-1}$$

$$\begin{gathered} d({\mathbf{p}}_t,{\mathbf{c}}_n)\le d({\mathbf{c}}_t,{\mathbf{c}}_n)+(1+{\epsilon }_1)r_t\le (d({\mathbf{c}}_t,{\mathbf{c}}_n)+r_t)+{\epsilon }_1{r}_t \\ \le \sqrt{2}\cdot \sqrt{d^2({\mathbf{c}}_t,{\mathbf{c}}_n)+r_t^2}+{\epsilon }_1{r}_t\le \sqrt{2}\cdot d({\mathbf{c}}_n,\mathbf{q})+{\epsilon }_1{r}_t\le (\sqrt{2}+{\epsilon }_1)r_n \end{gathered}$$

$$d(\mathbf{c}_{1},\mathbf{p}_{t^{\prime}})\leq d(\mathbf{c}_{1},\mathbf{c}_{2})+d(\mathbf{c}_{2},\mathbf{p}_{t^{\prime}})\leq(1+\varepsilon_{1})r_{2}+\sqrt{r^{2}_{2}-r^{2}_{1}}\\ \leq(1+\varepsilon_{1})(1+\varepsilon_{2})r_{1}+\sqrt{(1+\varepsilon_{2})^{2}-1}\cdot r_{1}=\big{(}1+O(\varepsilon_{1}+\sqrt{\varepsilon_{2}})\big{)}r_{1}$$

$$d(\mathbf{p}_{t^{\prime}},\mathbf{c}_{t})\leq d(\mathbf{c}_{1},\mathbf{c}_{t})+(1+\varepsilon)r_{1}\leq\big{(}d(\mathbf{c}_{1},\mathbf{c}_{t})+r_{1}\big{)}+\varepsilon r_{1}\\ \leq\sqrt{2}\cdot\sqrt{d^{2}(\mathbf{c}_{1},\mathbf{c}_{t})+r_{1}^{2}}+\varepsilon r_{1}\leq\sqrt{2}\cdot d(\mathbf{c}_{t},\mathbf{q})+\varepsilon r_{1}\leq(\sqrt{2}+\varepsilon)r_{t}$$

$$d(\mathbf{c},{\mathbf{p}}^{\ast })\le d({\mathbf{c}}^{\ast },\mathbf{c})+d({\mathbf{c}}^{\ast },{\mathbf{p}}^{\ast })\le \sqrt{(r^{\ast }{)}^2-r^2}+r^{\ast }=(\mu +\sqrt{{\mu }^2-1})r$$

$$d({\mathbf{c}}_1,{\mathbf{p}}_2)\le d({\mathbf{c}}_2,{\mathbf{p}}_2)+d({\mathbf{c}}_1,{\mathbf{c}}_2)\le \mu r_2+(r_1-r_2)\le \mu r_1$$

$$\begin{split}r^{*}(P_{1}\cup P_{3})&\leq\frac{1}{2}\big{(}r^{*}(P_{1})+d(\mathbf{c}^{*}(P_{1}),\mathbf{c}^{*}(P_{3}))+r^{*}(P_{3})\big{)}\\ &\leq\frac{1}{2}\big{(}z+\sqrt{z^{2}-1}+z_{1}+z_{2}\big{)}\cdot r^{*}(P_{2})\end{split}$$

$$\frac{r^{\ast }(P_1\cup P_3)}{r^{\ast }(P_2\cup P_3)}\le \frac{1}{2}(z+\sqrt{z^2-1}+z_1+z_2)\le \frac{1}{2}(1+z+\sqrt{z^2-1})\le z+\frac{1}{2}$$

$$\frac{r^{\ast }(P_1\cup P_3)}{r^{\ast }(P_2\cup P_3)}\le \frac{z+\sqrt{z^2-1}+z_1+z_2}{2z_1}=f_1(z,z_1,z_2)$$

$$\max$$

$$1+z_2\le z_1,z\ge 1,z_2\ge 0$$

$$\sqrt{r^{\ast }(P_2\cup P_3{)}^2-r^{\ast }(P_2{)}^2}+\sqrt{r^{\ast }(P_2\cup P_3{)}^2-r^{\ast }(P_3{)}^2}\ge d({\mathbf{c}}^{\ast }(P_2),{\mathbf{c}}^{\ast }(P_3))$$

$$\frac{r^{*}(P_{1}\cup P_{3})}{r^{*}(P_{2}\cup P_{3})}\leq\frac{z_{2}(z+\sqrt{z^{2}-1}+z_{1}+z_{2})}{\sqrt{\big{(}(1-z_{1})^{2}+z_{2}^{2}\big{)}\big{(}(1+z_{1})^{2}+z_{2}^{2}\big{)}}}=f_{2}(z,z_{1},z_{2})$$

$$\max$$

$$1+z_2>z_1,z_1+z_2>1,z\ge 1,z_1>0,z_2>0$$

$$\textstyle\mathbf{p}\in\Big{(}\frac{\sqrt{2}+\varepsilon_{3}}{\mu}+\sqrt{\frac{2+\varepsilon_{3}}{\mu^{2}}-1}\Big{)}(\mu+\sqrt{\mu^{2}-1})\cdot B[x_{2},t]$$

$$\mathbf{p}\in\big{(}4+O(\sqrt{\varepsilon_{1}}+\sqrt{\varepsilon_{2}})\big{)}\cdot B^{*}[x_{2},t]$$

$$B^{*}[x_{2},t]\subset\big{(}\sqrt{2}+1+O(\sqrt{\varepsilon_{1}})\big{)}\cdot B[x_{2},t]$$

$$\min_{\mathbf{c},r}\ r^{2}\quad\text{s.t.}\quad\big{(}\mathbf{c}-\varphi(\mathbf{p}_{i})\big{)}^{\prime}\big{(}\mathbf{c}-\varphi(\mathbf{p}_{i})\big{)}\leq r^{2},\forall\mathbf{p}_{i}\in P$$

$$\underset{\mathbf{\alpha }}{\max }{\mathbf{\alpha }}^{\prime }\text{diag}(\mathbf{K})-{\mathbf{\alpha }}^{\prime }\mathbf{K}\mathbf{\alpha }\text{s.t.}\mathbf{\alpha }\ge 0,{\mathbf{\alpha }}^{\prime }1=1$$

$${\mathbf{c}}^{\ast }=\sum _{i=1}^n{\alpha }_i\cdot \phi ({\mathbf{p}}_i),(r^{\ast }{)}^2={\mathbf{\alpha }}^{\prime }\text{diag}(\mathbf{K})-{\mathbf{\alpha }}^{\prime }\mathbf{K}\mathbf{\alpha }$$

$$d({\mathbf{c}}^{\ast },\phi (\mathbf{q}){)}^2=\sum _{i,j=1}^n{\alpha }_i{\alpha }_{j}k({\mathbf{p}}_i,{\mathbf{p}}_j)+k(\mathbf{q},\mathbf{q})-2\sum _{i=1}^n{\alpha }_{i}k({\mathbf{p}}_i,\mathbf{q})$$

Full text

1 Introduction

Unprecedented growth of data poses significant challenges in designing algorithms that can scale to massive datasets. Algorithms with superlinear complexity often become infeasible on datasets with millions or billions of points. Coresets are effective approaches to tackling the challenges of big data analysis. A coreset is a small subset extracted from the original point set such that certain geometric properties are preserved with provable guarantees [3]. Instead of processing the original dataset, one can perform the computation on its coreset with little loss of accuracy. Various types of problems have been shown to be effective under coreset approximation, e.g., $k$-median and $k$-means clustering [20, 3, 9], non-negative matrix factorization (NMF) [17], kernel density estimation (KDE) [34, 26], and many others [6, 21, 10].

Coresets for minimum enclosing balls (MEB) [5, 33, 22, 4, 1, 32, 23] have received significant attention due to its wide applications in clustering [5, 22, 4], support vector machines [28], kernel regression [31], fuzzy inference [13], shape fitting [23], and approximate furthest neighbor search [25]. Given a set of points $P$, the minimum enclosing ball of $P$, denoted by $\mathsf{MEB}(P)$, is the smallest ball that contains all points in $P$. A subset $S\subset P$ is a $\mu$-coreset for $\mathsf{MEB}(P)$ if the distance between the center $\mathbf{c}^{*}(S)$ of $\mathsf{MEB}(S)$ and any point in $P$ is within $\mu\cdot r^{*}(S)$, where $r^{*}(S)$ is the radius of $\mathsf{MEB}(S)$ [4]. Existing studies [22, 4] show that there always exists a $(1+\varepsilon)$-coreset of size $O(\frac{1}{\varepsilon})$ ($0<\varepsilon<1$) for the MEB of any point set, which is independent of the dataset size and dimension.

Most existing methods [5, 22, 4, 32, 23] of coresets for MEB focus on the batch setting and must keep all points in memory when constructing the coresets. In many applications, such as network monitoring, financial analysis, and sensor data mining, one needs to process data in the streaming model [24] where the input points arrive one at a time and cannot be stored entirely. There have been several methods [2, 33, 1, 11] to maintain coresets for MEB in data streams. The state-of-the-art streaming algorithm [1] can maintain a $(\sqrt{2}+\varepsilon)$-coreset for MEB with a single pass through the input dataset. However, these algorithms only consider the append-only scenario where new points are continuously added to, but old ones are never deleted from, the stream. Hence, they fail to capture recency in time-sensitive applications since the computation may be performed on outdated data. To meet the recency requirement, the sliding window model [15, 8, 16, 9] that only considers the most recent $N$ points in the stream at any time is a popular approach for real-time analytics. One can trivially adapt append-only methods for the sliding window model but a complete coreset reconstruction is deemed inevitable once an expired point is deleted. To the best of our knowledge, there is no existing algorithm that can maintain coresets for MEB over the sliding window efficiently.

Our Results. In this paper, we investigate the problem of maintaining coresets for MEB in the sliding window model. In particular, our results are summarized as follows.

• •

In Section 3.1, we present the AOMEB algorithm to maintain a $(\sqrt{2}+\varepsilon)$-coreset of size $O(\frac{\log\theta}{\varepsilon^{2}})$ with $O(\frac{m\log\theta}{\varepsilon^{3}})$ computation time per point for the MEB of an append-only stream, where $m$ is the dimension of the points and $\theta$ is the ratio of the maximum and minimum distances between any two points in the input dataset. AOMEB shows better empirical performance than the algorithm in [1] while having the same approximation ratio.

• •

In Section 3.2, using AOMEB as a building block, we propose the SWMEB algorithm for coreset maintenance over a sliding window $W_{t}$ of the most recent $N$ points at time $t$. SWMEB divides $W_{t}$ into equal-length partitions. On each partition, it maintains a sequence of indices where each index corresponds to an instance of AOMEB. Theoretically, SWMEB can return a $(\sqrt{2}+\varepsilon)$-coreset for $\mathsf{MEB}(W_{t})$ with $O(\sqrt{N}\cdot\frac{m\log\theta}{\varepsilon^{3}})$ time and $O(\sqrt{N}\cdot\frac{\log\theta}{\varepsilon^{2}})$ space complexity, where $N$ is the window size.

• •

In Section 3.3, we propose the SWMEB+ algorithm to improve upon SWMEB. SWMEB+ only maintains one sequence of indices, as well as the corresponding AOMEB instances, over $W_{t}$. By keeping fewer indices, SWMEB+ is more efficient than SWMEB in terms of time and space. Specifically, it only stores $O(\frac{\log^{2}\theta}{\varepsilon^{3}})$ points with $O(\frac{m\log^{2}\theta}{\varepsilon^{4}})$ processing time per point, both of which are independent of $N$. At the same time, it can still return a $(9.66+\varepsilon)$-coreset for $\mathsf{MEB}(W_{t})$.

• •

In Section 3.4, we generalize our proposed algorithms to maintain coresets for MEB in a reproducing kernel Hilbert space (RKHS).

• •

In Section 4, we conduct extensive experiments on real-world and synthetic datasets to evaluate the performance of our proposed algorithms. The experimental results demonstrate that (1) AOMEB outperforms the state-of-the-art streaming algorithm [1] in terms of coreset quality and efficiency; (2) SWMEB and SWMEB+ can return coresets for MEB with rather small errors (mostly within $1\%$), which are competitive with AOMEB and other streaming algorithms; (3) SWMEB and SWMEB+ achieve 2 to 4 orders of magnitude speedups over batch algorithms while running between 10 and 150 times faster than AOMEB; (4) SWMEB+ further improves the efficiency of SWMEB by up to 14 times while providing coresets with similar or even better quality.

2 Preliminaries and Related Work

Coresets for MEB. For two $m$-dimensional points $\mathbf{p}=(p_{1},\ldots,p_{m})$, $\mathbf{q}=(q_{1},\ldots,q_{m})$, the Euclidean distance between $\mathbf{p}$ and $\mathbf{q}$ is denoted by $d(\mathbf{p},\mathbf{q})=\sqrt{\sum_{i=1}^{m}(p_{i}-q_{i})^{2}}$. An $m$-dimensional (closed) ball with center $\mathbf{c}$ and radius $r$ is defined as $B(\mathbf{c},r)=\{\mathbf{p}\in\mathbb{R}^{m}:d(\mathbf{c},\mathbf{p})\leq r\}$. We use $\mathbf{c}(B)$ and $r(B)$ to denote the center and radius of ball $B$. The $\mu$-expansion of ball $B(\mathbf{c},r)$, denoted as $\mu\cdot B$, is a ball centered at $\mathbf{c}$ with radius $\mu\cdot r$, i.e., $\mu\cdot B=B(\mathbf{c},\mu\cdot r)$.

Given a set of $n$ points $P=\{\mathbf{p}_{1},\ldots,\mathbf{p}_{n}\}\subset\mathbb{R}^{m}$, the minimum enclosing ball of $P$, denoted as $\mathsf{MEB}(P)$, is the smallest ball that contains all points in $P$. The center and radius of $\mathsf{MEB}(P)$ are represented by $\mathbf{c}^{*}(P)$ and $r^{*}(P)$. For a parameter $\mu>1$, a ball $B$ is a $\mu$-approximate MEB of $P$ if $P\subset B$ and $r(B)\leq\mu\cdot r^{*}(P)$. A subset $S\subset P$ is a $\mu$-coreset for $\mathsf{MEB}(P)$, or $\mu$-$\mathsf{Coreset}(P)$ for brevity, if $P\subset\mu\cdot\mathsf{MEB}(S)$. Since $S\subseteq P$ and $r^{*}(S)\leq r^{*}(P)$, $\mu\cdot\mathsf{MEB}(S)$ is always a $\mu$-approximate MEB of $P$.

Sliding Window Model. This work focuses on maintaining coresets for MEB in append-only streaming and sliding window settings. For a sequence of (possibly infinite) points $P=\langle\mathbf{p}_{1},\mathbf{p}_{2},\ldots\rangle$ arriving continuously as a data stream where $\mathbf{p}_{t}$ is the $t$-th point, we first consider the problem of maintaining a $\mu$-$\mathsf{Coreset}(P_{t})$ for $P_{t}=\{\mathbf{p}_{1},\ldots,\mathbf{p}_{t}\}$ at any time $t$.

Furthermore, we consider the count-based sliding window111In this paper, we focus on the count-based sliding window model. But our proposed approaches can be trivially extended to the time-based sliding window model [15]. on the stream $P$: given a window size $N$, the sliding window [15] $W_{t}$ at any time $t$ always contains the latest $N$ points, i.e., $W_{t}=\{\mathbf{p}_{t^{\prime}},\ldots,\mathbf{p}_{t}\}$ where $t^{\prime}=\max(1,t-N+1)$. We consider the problem of maintaining a $\mu$-$\mathsf{Coreset}(W_{t})$ for $W_{t}$ at any time $t$.

Related Work. We review the literature on MEB computation and coresets for MEB. Gärtner [19] and Fischer et al. [18] propose two pivoting algorithms that resemble the simplex method of linear programming for computing exact MEBs. Both algorithms have an exponential complexity w.r.t. the dimension $m$ and thus are not scalable for large datasets with high dimensions. Subsequently, a line of research work [5, 22, 4, 32, 23] studies the problem of building coresets to approximate MEBs. They propose efficient batch algorithms for constructing a $(1+\varepsilon)$-$\mathsf{Coreset}(P)$ of any point set $P$. The basic scheme used in these algorithms is presented in Algorithm 1. First of all, it selects the point $\mathbf{p}_{a}$ furthest from $\mathbf{p}_{1}$ and $\mathbf{p}_{b}$ furthest from $\mathbf{p}_{a}$ out of $P$, using $S=\{\mathbf{p}_{a},\mathbf{p}_{b}\}$ as the initial coreset (Lines 1 & 1). The center $\mathbf{c}$ and radius $r$ of $\mathsf{MEB}(S)$ can be computed from $\mathbf{p}_{a}$ and $\mathbf{p}_{b}$ directly (Line 1). Then, it iteratively picks the point $\mathbf{q}$ furthest from the current center $\mathbf{c}$, adds $\mathbf{q}$ to $S$, and updates $B(\mathbf{c},r)$ so that $B$ is $\mathsf{MEB}(S)$, until no point in $P$ is outside of the $(1+\varepsilon)$-expansion of $B$ (Lines 1–1). Finally, it returns $S$ as a coreset for $\mathsf{MEB}(P)$ (Line 1). Theoretically, Algorithm 1 terminates in $O(\frac{1}{\varepsilon})$ iterations and returns a $(1+\varepsilon)$-$\mathsf{Coreset}(P)$ of size $O(\frac{1}{\varepsilon})$ [4]. Compared with exact MEB solvers [19, 18], coreset-based approaches run in linear time w.r.t. the dataset size $n$ and dimension $m$, and achieve better performance on high-dimensional data. Nevertheless, they must store all points in memory and process them in multiple passes, which are not suitable for data stream applications.

Several methods are proposed to approximate MEBs or coresets for MEB in streaming and dynamic settings. Agarwal et al. [2] and Chan [11] propose algorithms to build $(1+\varepsilon)$-coresets for MEB in append-only streams. Though working well in low dimensions, both algorithms become impractical for higher dimensions (i.e., $m>10$) due to $O(1/\varepsilon^{O(m)})$ complexity. Zarrabi-Zadeh and Chan [33] propose a $1.5$-approximate algorithm to compute MEBs in append-only streams. Agarwal and Sharathkumar [1] design a data structure that can maintain $(\sqrt{2}+\varepsilon)$-coresets for MEB and $(1.37+\varepsilon)$-approximate MEBs over append-only streams. Chan and Pathak [12] propose a method for maintaining $(1.22+\varepsilon)$-approximate MEBs in the dynamic setting, which supports the insertions and deletions of random points. To the best of our knowledge, none of the existing methods can maintain coresets for MEB over the sliding window efficiently. All of them have to store the entire window of points and recompute from scratch for every window slide, which is expensive in terms of time and space.

3 Our Algorithms

In this section we present our algorithms to maintain coresets for MEB. We first introduce a $(\sqrt{2}+\varepsilon)$-approximate append-only streaming algorithm, called AOMEB, in Section 3.1. Using AOMEB as a building block, we propose the SWMEB algorithm with the same $(\sqrt{2}+\varepsilon)$-approximation ratio in Section 3.2. Furthermore, we propose a more efficient SWMEB+ algorithm that retains a constant approximation ratio in Section 3.3.

3.1 The AOMEB Algorithm

The AOMEB algorithm is inspired by CoreMEB [4] (see Algorithm 1) to work in the append-only streaming model. Compared with CoreMEB, which can access the entire dataset and optimally select the furthest point into the coreset at each iteration, AOMEB is restricted to process the dataset in a single pass and determine whether to include a point into the coreset or discard it immediately after seeing it. Therefore, AOMEB adopts a greedy strategy for coreset maintenance: adding a new point to the coreset once it is outside of the MEB w.r.t. the current coreset.

The pseudo code of AOMEB is presented in Algorithm 2. First of all, it takes $S_{1}=\{\mathbf{p}_{1}\}$ as the initial coreset with $B_{1}(\mathbf{p}_{1},0)$ as $\mathsf{MEB}(S_{1})$ (Line 2). Then, it performs a one-pass scan over the point set, using the procedure in Lines 2–2 for each point $\mathbf{p}_{t}$: It first computes the distance between $\mathbf{p}_{t}$ and $\mathbf{c}_{t-1}$. If $\mathbf{p}_{t}\in(1+\varepsilon_{1})\cdot B_{t-1}$, no update is needed; otherwise, it adds $\mathbf{p}_{t}$ to the coreset $S_{t-1}$ and updates $B_{t-1}$ to $\mathsf{MEB}(S_{t})$. Finally, after processing all points in $P$, it returns $S_{n}$ as the coreset $S$ for $\mathsf{MEB}(P)$ (Line 2).

Theoretical Analysis. Next, we provide an analysis of the approximation ratio and complexity of AOMEB. It is noted that the greedy strategy of AOMEB is also adopted by existing streaming algorithms, i.e., SSMEB [33] and blurred ball cover (BBC) [1]. Nevertheless, the update procedure is different: SSMEB uses a simple geometric method to enlarge the MEB such that both the previous MEB and the new point are contained while AOMEB and BBC recompute the MEB once the coreset is updated. As a result, AOMEB and BBC are less efficient than SSMEB but ensure a better approximation ratio. Compared with BBC, which keeps the “archives” of MEBs for previous coresets, AOMEB only maintains one MEB w.r.t. $S_{t}$ at time $t$. Therefore, AOMEB is more efficient than BBC in practice. Next, we will prove that AOMEB has the same $(\sqrt{2}+\varepsilon)$-approximation as BBC. First of all, we present the hemisphere property [5] that forms the basis of our analysis.

Lemma 1** (Hemisphere Property [5]).**

For a set of points $P\subset\mathbb{R}^{m}$, any closed half-space that contains $\mathbf{c}^{*}(P)$ must contain at least a point $\mathbf{p}\in P$ such that $d(\mathbf{c}^{*}(P),\mathbf{p})=r^{*}(P)$.

The proof of Lemma 1 can be found in Section 2 of [5]. Based on Lemma 1, we can analyze the complexity and approximation ratio of AOMEB theoretically.

Theorem 1**.**

For any $\mathbf{p}_{t}\in S_{n}$, it holds that $r_{t}\geq(1+\frac{\varepsilon_{1}^{2}}{8})r_{t-1}$.

Proof.

If $\mathbf{p}_{t}\in S_{n}$, then $d(\mathbf{p}_{t},\mathbf{c}_{t-1})>(1+\varepsilon_{1})r_{t-1}$. We discuss two cases of $d(\mathbf{c}_{t},\mathbf{c}_{t-1})$ separately. If $d(\mathbf{c}_{t},\mathbf{c}_{t-1})\leq\frac{\varepsilon_{1}}{2}r_{t-1}$, then

[TABLE]

If $d(\mathbf{c}_{t},\mathbf{c}_{t-1})>\frac{\varepsilon_{1}}{2}r_{t-1}$, then let $H$ be a hyperplane passing through $\mathbf{c}_{t-1}$ with $\mathbf{c}_{t}-\mathbf{c}_{t-1}$ as its normal. Let $H^{+}$ be the closed half-space, bounded by $H$, that does not contain $\mathbf{c}_{t}$. According to Lemma 1, there must exist a point $\mathbf{q}\in S_{t-1}\cap H^{+}$ such that $d(\mathbf{c}_{t-1},\mathbf{q})=r_{t-1}$. Thus,

[TABLE]

In addition, $1+\frac{\varepsilon_{1}}{2}>\sqrt{1+(\frac{\varepsilon_{1}}{2})^{2}}>1+\frac{\varepsilon_{1}^{2}}{8}$ as $\varepsilon_{1}\in(0,1)$. Therefore, we prove that $r_{t}\geq(1+\frac{\varepsilon_{1}^{2}}{8})r_{t-1}$ in both cases. ∎

Theorem 2**.**

For any $\mathbf{p}_{t}\in P$, it holds that $\mathbf{p}_{t}\in(\sqrt{2}+\varepsilon_{1})\cdot B_{n}$.

Proof.

For any $\mathbf{p}_{t}\in P$, we have either $\mathbf{p}_{t}\in S_{n}$ or $d(\mathbf{p}_{t},\mathbf{c}_{t})\leq(1+\varepsilon_{1})r_{t}$. If $\mathbf{p}_{t}\in S_{n}$, it is obvious that $\mathbf{p}_{t}\in(1+\varepsilon_{1})\cdot B_{n}$. If $\mathbf{p}_{t}\notin S_{n}$, we have $d(\mathbf{p}_{t},\mathbf{c}_{n})\leq d(\mathbf{c}_{t},\mathbf{c}_{n})+(1+\varepsilon_{1})r_{t}$. According to Lemma 1, there must exist a point $\mathbf{q}\in S_{t}$ such that $d^{2}(\mathbf{c}_{n},\mathbf{q})\geq d^{2}(\mathbf{c}_{t},\mathbf{c}_{n})+r_{t}^{2}$. Therefore,

[TABLE]

We conclude that $\mathbf{p}_{t}\in(\sqrt{2}+\varepsilon_{1})\cdot B_{n},\forall\mathbf{p}_{t}\in P$. ∎

Theorem 2 indicates that AOMEB returns a $(\sqrt{2}+\varepsilon)$-$\mathsf{Coreset}(P)$ where $\varepsilon=O(\varepsilon_{1})$ for an arbitrary point set $P$. According to Theorem 1, the radius of $\mathsf{MEB}(S_{t})$ increases by $1+O(\varepsilon^{2})$ times whenever a new point is added to $S_{t}$. After processing $\mathbf{p}_{1}$ and $\mathbf{p}_{2}$, the coreset $S_{2}$ contains both points with $r_{2}\geq d_{min}/2$ where $d_{min}=\min_{\mathbf{p},\mathbf{q}\in P\wedge\mathbf{p}\neq\mathbf{q}}d(\mathbf{p},\mathbf{q})$. In addition, the radius $r_{n}$ of $S_{n}$ is bounded by $r^{*}(P)<d_{max}$ where $d_{max}=\max_{\mathbf{p},\mathbf{q}\in P}d(\mathbf{p},\mathbf{q})$. Therefore, $r_{n}/r_{2}=O(\theta)$ where $\theta=d_{max}/d_{min}$ and the size of $S$ is $O(\frac{\log\theta}{\varepsilon^{2}})$. Finally, the update procedure for each point $\mathbf{p}_{t}$ spends $O(m)$ time to compute $d(\mathbf{c}_{t-1},\mathbf{p}_{t})$ and $O(\frac{m\log\theta}{\varepsilon^{3}})$ time to update $B_{t}$.

3.2 The SWMEB Algorithm

In this subsection, we present the SWMEB algorithm for coreset maintenance over the sliding window $W_{t}$. The basic idea is to adapt AOMEB for the sliding window model by keeping multiple AOMEB instances with different starting points over $W_{t}$. However, the key problem is to identify the appropriate indices, i.e., starting points, for these instances. A naive scheme, i.e., creating a set of indices that are evenly distributed over $W_{t}$, cannot give any approximation guarantee of coreset quality. Therefore, we design a partition-based scheme for index maintenance in SWMEB: dividing $W_{t}$ into equal-length partitions and keeping a sequence of indices on each partition such that at least one instance can provide an approximate coreset for $\mathsf{MEB}(W_{t})$ at any time $t$.

The procedure of SWMEB is illustrated in Figure 1. It divides $W_{t}$ into $l=N/L$ partitions $\{P_{1},\ldots,P_{l}\}$ of equal length $L$. It keeps a sequence of $s_{i}$ indices $\langle x_{i,1},\ldots,x_{i,s_{i}}\rangle$ from the end to the beginning of each partition $P_{i}$. As $W_{t}$ slides over time, old points in $P_{1}$ expire (colored in grey) while new points are temporarily stored in a buffer $Q$. The index $x_{1,s_{1}}$ on $P_{1}$, which is the closest to the beginning of $W_{t}$, will be deleted once it expires. When the size of $Q$ reaches $L$, it will delete $P_{1}$ and shift remaining partitions as all points in $P_{1}$ must have expired. Then, it creates a new partition $P_{l}$ for the points in $Q$ and the indices on $P_{l}$. Moreover, each index $x_{i,j}$ corresponds to an AOMEB instance $\mathcal{A}(x_{i,j})$ that processes $P[x_{i,j},t]=\{\mathbf{p}_{x_{i,j}},\ldots,\mathbf{p}_{t}\}$ at any time $t$. Specifically, $\mathcal{A}(x_{i,j})$ will process the points from the end of $P_{i}$ to $\mathbf{p}_{x_{i,j}}$ when $x_{i,j}$ is created and then update for each point till $\mathbf{p}_{t}$. Finally, the coreset is always provided by $\mathcal{A}(x_{1,s_{1}})$.

The pseudo code of SWMEB is presented in Algorithm 3. For initialization, the latest partition ID $l$ is set to [math] and the buffer $Q$ as well as the indices $X_{0}$ are set to $\varnothing$ (Line 3). Then, it processes all points in the stream $P$ one by one with the procedure of Lines 3–3, which can be separated into four phases as follows.

• •

Phase 1 (Lines 3–3): After adding a new point $\mathbf{p}_{t}$ to $Q$, it checks the size of $Q$. If $|Q|=L$, a new partition will be created for $Q$. When $t\leq N$, it increases $l$ by $1$ and creates a new partition $P_{l}$. Otherwise, $P_{1}$ must have expired and thus is dropped. Then, the partitions $\{P_{2},\ldots,P_{l}\}$ (and the indices on each partition) are shifted to $\{P_{1},\ldots,P_{l-1}\}$ and a new partition $P_{l}$ is created.

• •

Phase 2 (Lines 3–3): Next, it creates the indices and corresponding AOMEB instances on $P_{l}$. It runs an AOMEB instance $\mathcal{A}$ to process each point in $P_{l}$ inversely from $\mathbf{p}_{t}$ to $\mathbf{p}_{t-L+1}$. Initially, the number of indices $s_{l}$ on $P_{l}$ and the radius $r_{l}$ w.r.t. the latest index $x_{l,s_{l}}$ are [math]. We denote the coreset maintained by $\mathcal{A}$ after processing $\mathbf{p}_{t^{\prime}}$ as $S[t,t^{\prime}]$. Then, $\mathsf{MEB}(S[t,t^{\prime}])$ is represented by $B[t,t^{\prime}]$ with radius $r[t,t^{\prime}]$. If $r[t,t^{\prime}]\geq(1+\varepsilon_{2})r_{l}$, it will update $r_{l}$ to $r[t,t^{\prime}]$, add a new index $x_{l,s_{l}}$ to $X_{t}$, and use the snapshot of $\mathcal{A}$ after processing $\mathbf{p}_{t^{\prime}}$ as $\mathcal{A}(x_{l,s_{l}})$. After the indices on $P_{l}$ is created, $Q$ will be reset for new incoming points.

• •

Phase 3 (Lines 3–3): It checks whether $x_{1,s_{1}}$, i.e., the earliest index on $P_{1}$, has expired. If so, it will delete $x_{1,s_{1}}$ from $X_{t}$ and terminate $\mathcal{A}(x_{1,s_{1}})$ accordingly.

• •

Phase 4 (Lines 3–3): For each index $x_{i,j}\in X_{t}$ with $i\in[1,l]$ and $j\in[1,s_{i}]$ at time $t$, it updates the corresponding AOMEB instance $\mathcal{A}(x_{i,j})$ by processing $\mathbf{p}_{t}$.

Finally, it always returns $S[x_{1,s_{1}},t]$ from $\mathcal{A}(x_{1,s_{1}})$ as the coreset $S_{t}$ for $\mathsf{MEB}(W_{t})$ at time $t$ (Line 3).

Theoretical Analysis. In the following, we will first prove the approximation ratio of $S_{t}$ returned by SWMEB for $\mathsf{MEB}(W_{t})$. Then, we discuss the time and space complexity of SWMEB.

We first prove the following lemma that will be used in subsequent analyses.

Lemma 2**.**

For any two point sets $P_{1},P_{2}$ such that $P_{2}\subset P_{1}$, it must hold that $d^{2}(\mathbf{c}^{*}(P_{1}),\mathbf{c}^{*}(P_{2}))\leq(r^{*}(P_{1}))^{2}-(r^{*}(P_{2}))^{2}$.

Proof.

Obviously, Lemma 2 must hold when $\mathbf{c}^{*}(P_{1})=\mathbf{c}^{*}(P_{2})$. When $\mathbf{c}^{*}(P_{1})\neq\mathbf{c}^{*}(P_{2})$, we consider a hyperplane $H$ passing through $\mathbf{c}^{*}(P_{2})$ with $\mathbf{c}^{*}(P_{1})-\mathbf{c}^{*}(P_{2})$ as its normal. Let $H^{+}$ be the close half-space, bounded by $H$, that does not contain $\mathbf{c}^{*}(P_{1})$. According to Lemma 1, there must exist $\mathbf{q}\in P_{2}\cap H^{+}$ such that $d(\mathbf{c}^{*}(P_{2}),\mathbf{q})=r^{*}(P_{2})$. In addition, $d(\mathbf{c}^{*}(P_{1}),\mathbf{q})\leq r^{*}(P_{1})$ for $\mathbf{q}\in P_{1}$. Finally, as $\mathbf{c}^{*}(P_{1})\notin H^{+}$, $\mathbf{q}\in H^{+}$, and $\mathbf{c}^{*}(P_{1})-\mathbf{c}^{*}(P_{2})$ is the normal of $H$, we acquire $d^{2}(\mathbf{c}^{*}(P_{1}),\mathbf{q})\geq d^{2}(\mathbf{c}^{*}(P_{2}),\mathbf{q})+d^{2}(\mathbf{c}^{*}(P_{1}),\mathbf{c}^{*}(P_{2}))$. Thus, it holds that $d^{2}(\mathbf{c}^{*}(P_{1}),\mathbf{c}^{*}(P_{2}))\leq(r^{*}(P_{1}))^{2}-(r^{*}(P_{2}))^{2}$ and we conclude the proof. ∎

Theorem 3**.**

For any $\mathbf{p}\in W_{t}$, it holds that $\mathbf{p}\in(\sqrt{2}+\varepsilon)\cdot\mathsf{MEB}(S_{t})$ where $\varepsilon=O(\varepsilon_{1}+\sqrt{\varepsilon_{2}})$.

Proof.

According to Algorithm 3, the instance $\mathcal{A}(x_{1,s_{1}})$ is always used to return the coreset $S_{t}$. Since $\mathcal{A}(x_{1,s_{1}})$ has already processed the points from $\mathbf{p}_{x_{1,s_{1}}}$ to $\mathbf{p}_{t}$, these points must be contained in $(\sqrt{2}+\varepsilon_{1})\cdot\mathsf{MEB}(S_{t})$ according to Theorem 2. Thus, we only need to consider the points from $\mathbf{p}_{t-N+1}$ to $\mathbf{p}_{x_{1,s_{1}}-1}$. To create the indices on $P_{1}$, we process the points of $P_{1}$ with an AOMEB instance (see Line 3). Here we use $S_{1},B_{1}(\mathbf{c}_{1},r_{1})$ and $S_{2},B_{2}(\mathbf{c}_{2},r_{2})$ to denote the coresets of this instance and the corresponding MEBs after processing $\mathbf{p}_{x_{1,s_{1}}}$ and $\mathbf{p}_{t^{\prime}}$ ($t-N<t^{\prime}<x_{1,s_{1}}$) respectively. For each point $\mathbf{p}_{t^{\prime}}$, it holds that $\mathbf{p}_{t^{\prime}}\in(1+\varepsilon_{1})\cdot B_{2}$ and $r_{2}<(1+\varepsilon_{2})r_{1}$. In addition, as $S_{1}\subseteq S_{2}$, we have $d^{2}(\mathbf{c}_{1},\mathbf{c}_{2})\leq r^{2}_{2}-r^{2}_{1}$ from Lemma 2. Therefore,

[TABLE]

In addition, according to Lemma 1, there must exist a point $\mathbf{q}\in S_{1}$ such that $d^{2}(\mathbf{c}_{t},\mathbf{q})\geq d^{2}(\mathbf{c}_{1},\mathbf{c}_{t})+r_{1}^{2}$. Let $\varepsilon=O(\varepsilon_{1}+\sqrt{\varepsilon_{2}})$, we have

[TABLE]

We prove that $\mathbf{p}_{t^{\prime}}\in(\sqrt{2}+\varepsilon)\cdot B_{t}$ and thus conclude the proof. ∎

Theorem 3 shows that $S_{t}$ returned by SWMEB is a $(\sqrt{2}+\varepsilon)$-$\mathsf{Coreset}(W_{t})$ where $\varepsilon=O(\varepsilon_{1}+\sqrt{\varepsilon_{2}})$ at any time $t$. To analyze the complexity of SWMEB, we first consider the number of indices in $X_{t}$. For each partition, SWMEB maintains $O(\frac{\log\theta}{\varepsilon_{2}})$ indices where $\theta=d_{max}/d_{min}$. Thus, $X_{t}$ contains $O(\frac{l\log\theta}{\varepsilon_{2}})$ indices and the number of points stored by SWMEB is $O(\frac{l\log^{2}\theta}{\varepsilon_{1}^{2}\varepsilon_{2}}+L)$. Furthermore, the time of SWMEB to update a point $\mathbf{p}_{t}$ comprises (1) the time to maintain the instance w.r.t. each index in $X_{t}$ for $\mathbf{p}_{t}$ and (2) the amortized time to create the indices for each partition. Overall, the time complexity of SWMEB to update each point is $O(\frac{lm\log^{2}\theta}{\varepsilon_{1}^{3}\varepsilon_{2}})$. As $l=N/L$, the number of points maintained by SWMEB is minimal when $L=\frac{\log\theta}{\varepsilon_{1}}\cdot\sqrt{\frac{N}{\varepsilon_{2}}}$. In this case, the number of points stored by SWMEB is $O(\sqrt{N}\cdot\frac{\log\theta}{\varepsilon^{2}})$ and the time complexity of SWMEB to update one point is $O(\sqrt{N}\cdot\frac{m\log\theta}{\varepsilon^{3}})$.

3.3 The SWMEB+ Algorithm

In this subsection we present the SWMEB+ algorithm that improves upon SWMEB in terms of time and space while still achieving a constant approximation ratio. The basic idea of SWMEB+ is illustrated in Figure 2. Different from SWMEB, SWMEB+ only maintains a single sequence of $s$ indices $X_{t}=\{x_{1},\ldots,x_{s}\}$ over $W_{t}$. Then, each index $x_{i}$ also corresponds to an AOMEB instance $\mathcal{A}(x_{i})$ that processes a substream of points from $\mathbf{p}_{x_{i}}$ to $\mathbf{p}_{t}$. We use $S[x_{i},t]$ for the coreset returned by $\mathcal{A}(x_{i})$ at time $t$ and $B[x_{i},t]$ centered at $\mathbf{c}[x_{i},t]$ with radius $r[x_{i},t]$ for $\mathsf{MEB}(S[x_{i},t])$. Furthermore, SWMEB+ maintains the indices based on the radii of the MEBs. Specifically, given any $\varepsilon_{2}>0$, for three neighboring indices $x_{i},x_{i+1},x_{i+2}$, if $r[x_{i},t]\leq(1+\varepsilon_{2})r[x_{i+2},t]$, then $x_{i+2}$ is considered as a good approximation for $x_{i}$ and thus $x_{i+1}$ can be deleted. In this way, the radii of the MEBs gradually decreases from $x_{1}$ to $x_{s}$, with the ratios of any two neighboring indices close to $(1+\varepsilon_{2})$. Any window starting between $x_{i}$ and $x_{i+1}$ is approximated by $\mathcal{A}(x_{i+1})$. Finally, SWMEB+ keeps at most one expired index (and must be $x_{1}$) in $X_{t}$ to track the upper bound for the radius $r^{*}(W_{t})$ of $\mathsf{MEB}(W_{t})$. The AOMEB instance corresponding to the first non-expired index ($x_{1}$ or $x_{2}$) provides the coreset for $\mathsf{MEB}(W_{t})$.

The pseudo code of SWMEB+ is presented in Algorithm 4. In the initialization phase, $X_{0}$ and $s$ are set to $\varnothing$ and [math] respectively (Line 4). Then, all points in $P$ are processed one by one with the procedure of Lines 4–4, which includes four phases as follows.

• •

Phase 1 (Lines 4–4): Upon the arrival of $\mathbf{p}_{t}$ at time $t$, it creates a new index $x_{s}=t$ and adds $x_{s}$ to $X_{t}$; accordingly, an AOMEB instance $\mathcal{A}(x_{s})$ w.r.t. $x_{s}$ is initialized to process the substream beginning at $\mathbf{p}_{t}$.

• •

Phase 2 (Lines 4–4): When there exists more than one expired index (i.e., earlier than the beginning of $W_{t}$), it deletes the first index $x_{1}$ and terminates $\mathcal{A}(x_{1})$ until there is only one expired index. Note that it shifts the remaining indices after deletion to always guarantee $x_{i}$ is the $i$-th index of $X_{t}$.

• •

Phase 3 (Lines 4–4): For each $x_{i}\in X_{t}$, it updates the instance $\mathcal{A}(x_{i})$ for $\mathbf{p}_{t}$. The update procedure follows Line 2–2 of Algorithm 2. After the update, $\mathcal{A}(x_{i})$ maintains a coreset $S[x_{i},t]$ and its MEB $B[x_{i},t]$ by processing a stream $P[x_{i},t]=\langle\mathbf{p}_{x_{i}},\ldots,\mathbf{p}_{t}\rangle$.

• •

Phase 4 (Lines 4–4): It executes a scan of $X_{t}$ from $x_{1}$ to $x_{t}$ to delete the indices that can be approximated by their successors. For each $x_{i}\in X_{t}$ ($i\leq s-2$), it checks the radii $r[x_{i},t]$ and $r[x_{i+2},t]$ of $B[x_{i},t]$ and $B[x_{i+2},t]$. If $r[x_{i},t]\leq(1+\varepsilon_{2})r[x_{i+2},t]$, then it deletes the index $x_{i+1}$ from $X_{t}$, terminates $\mathcal{A}(x_{i+1})$, and shifts the remaining indices accordingly.

After performing the above procedure, it returns either $S[x_{1},t]$ (when $x_{1}$ has not expired) or $S[x_{2},t]$ (when $x_{1}$ has expired) as the coreset $S_{t}$ for $\mathsf{MEB}(W_{t})$ at time $t$.

Theoretical Analysis. The strategy of index maintenance based on the ratios of radii is inspired by Smooth Histograms [8] for estimating stream statistics over sliding windows. However, Smooth Histograms cannot be applied to our problem because it requires an oracle to provide a $(1+\varepsilon)$-approximate function value in any append-only stream [16] but any practical solution (i.e., [1] and AOMEB) only gives a $(\sqrt{2}+\varepsilon)$-approximation for $r^{*}(P)$ of an append-only stream $P$. In addition, Smooth Histograms are also used for submodular maximization in the sliding window model [16, 29, 30]. Nevertheless, such an extension is still not applicable for our problem because the radius function $r^{*}(\cdot)$ is not submodular in the view of set functions, which is shown by Example 1. In the following, we will prove that SWMEB+ still has a constant approximation ratio by an analysis that is different from [8, 16].

Example 1**.**

A function $f(\cdot)$ is submodular if $f(P\cup\{\mathbf{p}\})-f(P)\geq f(Q\cup\{\mathbf{p}\})-f(Q)$ for any set $P\subseteq Q$ and point $\mathbf{p}\notin Q$. In Figure 3, for $P_{1}=\{\mathbf{p}_{1},\mathbf{p}_{2},\mathbf{p}_{3}\}$ and $P_{2}=\{\mathbf{p}_{2},\mathbf{p}_{3}\}$, $r^{*}(P_{1}\cup\{\mathbf{p}_{4}\})-r^{*}(P_{1})>0$ but $r^{*}(P_{2}\cup\{\mathbf{p}_{4}\})-r^{*}(P_{2})=0$. As $P_{2}\subset P_{1}$, we show that $r^{*}(\cdot)$ is not a submodular function.

We first present two lemmas that will be used in our analysis of SWMEB+.

Lemma 3**.**

For a set of points $P$, let $B^{*}(\mathbf{c}^{*},r^{*})=\mathsf{MEB}(P)$ and $B(\mathbf{c},r)=\mathsf{MEB}(S)$ where $S\subseteq P$ is the coreset returned by Algorithm 2. If $\frac{r^{*}}{r}=\mu$ ($1\leq\mu\leq\sqrt{2}+\varepsilon_{1}$), then $B^{*}\subset(\mu+\sqrt{\mu^{2}-1})\cdot B$.

Proof.

According to Lemma 2, $d^{2}(\mathbf{c}^{*},\mathbf{c})\leq(r^{*})^{2}-r^{2}$ because of $S\subseteq P$.

Then $\forall\mathbf{p}^{*}\in B^{*}$, it holds that

[TABLE]

and we conclude the proof. ∎

Lemma 4**.**

For any two balls $B_{1}(\mathbf{c}_{1},r_{1})$ and $B_{2}(\mathbf{c}_{2},r_{2})$ with $B_{2}\subset B_{1}$, it holds that $\mu\cdot B_{2}\subset\mu\cdot B_{1},\forall\mu\geq 1$.

Proof.

First, we can prove $d(\mathbf{c}_{1},\mathbf{c}_{2})\leq r_{1}-r_{2}$ by contradiction. Then, based on the above result, we can see that, $\forall\mathbf{p}_{2}\in\mu\cdot B_{2}$,

[TABLE]

where the last inequality holds from $r_{1}\geq r_{2}$. Thus, $\forall\mathbf{p}_{2}\in\mu\cdot B_{2}$, $\mathbf{p}_{2}\in\mu\cdot B_{1}$. ∎

According to the above results, we prove the following lemma.

Lemma 5**.**

For any three point sets $P_{1},P_{2},P_{3}$ such that $P_{2}\subset P_{1}$ and $P_{3}\not\subset P_{1}$, if $\frac{r^{*}(P_{1})}{r^{*}(P_{2})}=z$, then it holds that $\frac{r^{*}(P_{1}\cup P_{3})}{r^{*}(P_{2}\cup P_{3})}\leq z+\frac{\sqrt{2}}{2}$.

Proof.

Let $z_{1}=\frac{r^{*}(P_{3})}{r^{*}(P_{2})}$, $z_{2}=\frac{d(\mathbf{c}^{*}(P_{2}),\mathbf{c}^{*}(P_{3}))}{r^{*}(P_{2})}$. First, we have

[TABLE]

Next, we consider three cases separately.

• •

Case 1 ($r^{*}(P_{3})+d(\mathbf{c}^{*}(P_{2}),\mathbf{c}^{*}(P_{3}))\leq r^{*}(P_{2})$): In this case, we have $r^{*}(P_{2}\cup P_{3})=r^{*}(P_{2})$ and $z_{1}+z_{2}\leq 1$. Thus, we can acquire

[TABLE]

• •

Case 2 ($r^{*}(P_{2})+d(\mathbf{c}^{*}(P_{2}),\mathbf{c}^{*}(P_{3}))\leq r^{*}(P_{3})$): In this case, we have $r^{*}(P_{2}\cup P_{3})=r^{*}(P_{3})$ and $1+z_{2}\leq z_{1}$. Thus, we can acquire

[TABLE]

We can transform it into the following problem:

[TABLE]

Solving Problem 2, we can get $\frac{r^{*}(P_{1}\cup P_{3})}{r^{*}(P_{2}\cup P_{3})}\leq z+\frac{1}{2}$.

• •

Case 3 (Otherwise): It holds that $z_{2}>0$, $1+z_{2}>z_{1}$, and $z_{1}+z_{2}>1$. In addition, we have $d(\mathbf{c}^{*}(P_{2}),\mathbf{c}^{*}(P_{3}))>0$ because the criteria of either Case 1 or 2 must be satisfied when $d(\mathbf{c}^{*}(P_{2}),\mathbf{c}^{*}(P_{3}))=0$. As $d(\mathbf{c}^{*}(P_{2}\cup P_{3}),\mathbf{c}^{*}(P_{2}))^{2}\leq r^{*}(P_{2}\cup P_{3})^{2}-r^{*}(P_{2})^{2}$ and $d(\mathbf{c}^{*}(P_{2}\cup P_{3}),\mathbf{c}^{*}(P_{3}))^{2}\leq r^{*}(P_{2}\cup P_{3})^{2}-r^{*}(P_{3})^{2}$ according to Lemma 2, we can get

[TABLE]

Considering Equations 1 and 3 collectively, we have

[TABLE]

We can transform it into the following problem:

[TABLE]

Solving Problem 4, we can get $\frac{r^{*}(P_{1}\cup P_{3})}{r^{*}(P_{2}\cup P_{3})}\leq z+\frac{\sqrt{2}}{2}$.

In all the three cases, we prove Lemma 5. ∎

Considering previous results collectively, we can prove the approximation ratio of SWMEB+ in the following theorem.

Theorem 4**.**

For any $\mathbf{p}\in W_{t}$, it holds that $\mathbf{p}\in(9.66+\varepsilon)\cdot\mathsf{MEB}(S_{t})$ where $\varepsilon=O(\sqrt{\varepsilon_{1}}+\sqrt{\varepsilon_{2}})$.

Proof.

We consider four different cases of $S_{t}$ returned by SWMEB+ at time $t$.

• •

Case 1 ($x_{1}\geq t-N+1$): In this case, we have $t<N$, $x_{1}=1$, and $W_{t}=P[x_{1},t]$. Therefore, it holds that $S_{t}=S[x_{1},t]$ and $\mathbf{p}\in(\sqrt{2}+\varepsilon_{1})\cdot\mathsf{MEB}(S_{t})$ for any $\mathbf{p}\in W_{t}$ according to Theorem 2.

• •

Case 2 ($x_{1}<t-N+1\wedge x_{2}=x_{1}+1$): In this case, we have $x_{2}=t-N+1$ and $W_{t}=P[x_{2},t]$. Similarly, $S_{t}=S[x_{2},t]$ and $\mathbf{p}\in(\sqrt{2}+\varepsilon_{1})\cdot\mathsf{MEB}(S_{t})$ for any $\mathbf{p}\in W_{t}$.

• •

Case 3 ($x_{1}<t-N+1\wedge\frac{r[x_{1},t]}{r[x_{2},t]}\leq 1+\varepsilon_{2}$): In this case, we have $W_{t}\subset P[x_{1},t]$ and $W_{t}\subset B^{*}[x_{1},t]=\mathsf{MEB}(P[x_{1},t])$. Then, we can acquire $1\leq\frac{r^{*}[x_{1},t]}{r[x_{1},t]},\frac{r^{*}[x_{2},t]}{r[x_{2},t]}\leq\sqrt{2}+\varepsilon_{1}$ from Theorem 2. Given $\frac{r^{*}[x_{2},t]}{r[x_{2},t]}=\mu$, it holds that \mathbf{p}\in\big{(}\frac{\sqrt{2}+\varepsilon_{3}}{\mu}+\sqrt{\frac{2+\varepsilon_{3}}{\mu^{2}}-1}\big{)}\cdot B^{*}[x_{2},t],\forall\mathbf{p}\in W_{t}, where $\varepsilon_{3}=O(\varepsilon_{1}+\varepsilon_{2})$. Additionally, we have $B^{*}[x_{2},t]\subset(\mu+\sqrt{\mu^{2}-1})\cdot B[x_{2},t]$ from Lemma 3. According to Lemma 4, it holds that $\forall\mathbf{p}\in W_{t}$,

[TABLE]

By finding the maximum of f(\mu)=\Big{(}\sqrt{\frac{2+\varepsilon_{3}}{\mu^{2}}-1}\Big{)}(\mu+\sqrt{\mu^{2}-1}) in Equation 5 on range $[1,\sqrt{2}+\varepsilon_{1}]$, we can acquire \mathbf{p}\in\big{(}3.36+\sqrt{\varepsilon_{3}}\big{)}\cdot B[x_{2},t] for any $\mathbf{p}\in W_{t}$.

• •

Case 4 ($x_{1}<t-N+1\wedge\frac{r[x_{1},t]}{r[x_{2},t]}>1+\varepsilon_{2}$): At time $t^{\prime}<t$ when $x_{1}$ and $x_{2}$ become neighboring indices, it holds that $\frac{r[x_{1},t^{\prime}]}{r[x_{2},t^{\prime}]}\leq 1+\varepsilon_{2}$. Then, $\frac{r^{*}[x_{1},t^{\prime}]}{r^{*}[x_{2},t^{\prime}]}\leq(1+\varepsilon_{2})(\sqrt{2}+\varepsilon_{1})$ where $r^{*}[x_{1},t^{\prime}]$ and $r^{*}[x_{2},t^{\prime}]$ are the radii of the MEBs for $P[x_{1},t^{\prime}]$ and $P[x_{2},t^{\prime}]$ respectively. Furthermore, according to Lemma 5, we have $\frac{r^{*}[x_{1},t]}{r^{*}[x_{2},t]}\leq(1+\varepsilon_{2})(\sqrt{2}+\varepsilon_{1})+\frac{\sqrt{2}}{2}$. According to Lemmas 3 and 4, we can acquire that $\forall\mathbf{p}\in W_{t}$,

[TABLE]

Moreover, we have $\frac{r^{*}[x_{2},t]}{r[x_{2},t]}\leq\sqrt{2}+\varepsilon_{1}$ and thus

[TABLE]

Combining Equations 6 and 7, we finally get \mathbf{p}\in\big{(}4\sqrt{2}+4+\varepsilon\big{)}\cdot B[x_{2},t]\approx(9.66+\varepsilon)\cdot B[x_{2},t],\forall\mathbf{p}\in W_{t}, where $\varepsilon=O(\sqrt{\varepsilon_{1}}+\sqrt{\varepsilon_{2}})$.

Considering all above cases collectively, we conclude the proof. ∎

According to Theorem 4, SWMEB+ can always return a $(9.66+\varepsilon)$-$\mathsf{Coreset}(W_{t})$ at any time $t$. In practice, since $r[x_{1},t]/r[x_{2},t]$ is bounded by $1+O(\varepsilon_{2})$ in almost all cases, the approximation ratio of SWMEB+ can be improved to $(3.36+\varepsilon)$ accordingly. Furthermore, for any $x_{i}\in X_{t}$ ($i\leq s-2$), either $r[x_{i+1},t]$ or $r[x_{i+2},t]$ is less than $(1+\varepsilon_{2})r[x_{i},t]$. In addition, it holds that $r[x_{1},t]\leq d_{max}$ and $r[x_{s-1},t]\geq 0.5\cdot d_{min}$. Therefore, the number of indices in $X_{t}$ is $O(\frac{\log\theta}{\varepsilon})$. Accordingly, the time complexity for SWMEB+ to update each point is $O(\frac{m\log^{2}\theta}{\varepsilon^{4}})$ while the number of points stored by SWMEB+ is $O(\frac{\log^{2}\theta}{\varepsilon^{3}})$, both of which are independent of $N$.

3.4 Generalization to Kernelized MEB

In real-world applications [28, 13, 31], it is required to compute the coreset for MEB in a reproducing kernel Hilbert space (RKHS) instead of Euclidean space. Given a symmetric positive definite kernel $k(\cdot,\cdot):\mathbb{R}^{m}\times\mathbb{R}^{m}\rightarrow\mathbb{R}$ and its associated feature mapping $\varphi(\cdot)$ where $k(\mathbf{p},\mathbf{q})=\langle\varphi(\mathbf{p}),\varphi(\mathbf{q})\rangle$ for any $\mathbf{p},\mathbf{q}\in\mathbb{R}^{m}$, the kernelized MEB of a set of points $P$ is the smallest ball $B^{*}(\mathbf{c}^{*},r^{*})$ in the RKHS such that the maximum distance from $\mathbf{c}^{*}$ to $\varphi(\mathbf{p})$ is no greater than $r^{*}$, which can be formulated as follows:

[TABLE]

However, it is impractical to solve Problem 8 directly in the primal form due to the infinite dimensionality of RKHS. We transform Problem 8 into the dual form as follows:

[TABLE]

where $\bm{\alpha}=[\alpha_{1},\ldots,\alpha_{n}]^{\prime}$ is the $n$-dimensional Lagrange multiplier vector, $\mathbf{0}=[0,\ldots,0]^{\prime}$ and $\mathbf{1}=[1,\ldots,1]^{\prime}$ are both $n$-dimensional vectors, $\mathbf{K}=[k(\mathbf{p}_{i},\mathbf{p}_{j})]_{i,j=1}^{n}$ is the $n\times n$ kernel matrix of $P$, and $\text{diag}(\mathbf{K})$ is the diagonal of $\mathbf{K}$. Problem 9 is known to be a quadratic programming [7] (QP) problem. According to the KKT conditions [7], the kernelized MEB $B^{*}$ of $P$, can be recovered from the Lagrange multiplier vector $\bm{\alpha}$ as follows:

[TABLE]

Here the center $\mathbf{c}^{*}$ is represented implicitly by each point $\mathbf{p}_{i}$ in $P$ and the corresponding $\alpha_{i}$’s. Then, the distance between $\mathbf{c}^{*}$ and $\varphi(\mathbf{q})$ for any $\mathbf{q}\in\mathbb{R}^{m}$ can be computed by:

[TABLE]

Next, we introduce how to generalize AOMEB in Algorithm 2 to maintain coresets for kernelized MEB. First, to represent the center $\mathbf{c}_{t}$ and compute the distance from $\varphi(\mathbf{p}_{t})$ to the center, it always keeps the Lagrange multiplier vector $\bm{\alpha}$. In Line 2, $\alpha_{1}$ is initialized to $1$ so that $\mathbf{c}_{1}=\varphi(\mathbf{p}_{1})$. Then, in Line 2, Equation 11 is used to compute the distance between $\mathbf{c}_{t-1}$ and $\varphi(\mathbf{p}_{t})$. If $\mathbf{p}_{t}$ is added to $S_{t}$, it will re-optimize Problem 9 on $S_{t}$ to adjust $\bm{\alpha}$ so that $B_{t}=\mathsf{MEB}(S_{t})$ in Line 1. Specifically, the Frank-Wolfe algorithm [32, 14] is used to solve Problem 9 as it is efficient for QP problems with unit simplex constraints. SWMEB and SWMEB+ can also maintain coresets for kernelized MEB by using the generalized AOMEB instance in Algorithms 3 and 4.

Theoretically, the generalized algorithms have the same approximation ratios and coreset sizes as the original ones. However, the time complexity will increase by a factor of $\frac{1}{\varepsilon}$ because the time to compute the distance from $\mathbf{c}_{t-1}$ to $\varphi(\mathbf{p}_{t})$ using Equation 11 is $O(\frac{m}{\varepsilon})$ instead of $O(m)$.

3.5 Discussion

For ease of presentation, we describe Algorithms 2–4 in the single-update-mode where the coreset is maintained for every new point. In practice, it is not required to update the coreset at such an intense rate. Here we discuss how to adapt these algorithms for the mini-batch-mode.

In the mini-batch-mode, given a batch size $b$, each update will add $b$ new points while deleting the earliest $b$ points at the same time. The adaptations of AOMEB in Algorithm 2 for the mini-batch-mode are as follows: (1) In Line 2, an initial coreset $S_{b}$ is built for the first $b$ points using CoreMEB in Algorithm 1. (2) In Lines 2–2, the coreset is updated for a batch of $b$ points collectively. Specifically, it adds the points not contained in $(1+\varepsilon_{1})\cdot B_{t-b}$ in the batch to $S_{t}$ and then updates $B_{t}$. To adapt SWMEB and SWMEB+ for the mini-batch-mode, each AOMEB instance should run in the mini-batch-mode as shown above. In addition, the indices and AOMEB instances in SWMEB and SWMEB+ are created and updated for batches instead of points. Note that the approximation ratio and coreset size will remain the same in the mini-batch-mode but the efficiency will be improved as fewer indices are created.

4 Experiments

In this section we evaluate the empirical performance of our proposed algorithms on real-world and synthetic datasets. First of all, we will introduce the experimental setup in Section 4.1. Then, we will present the experimental results on effectiveness and efficiency in Section 4.2. Finally, the experimental results on scalability are presented in Section 4.3.

4.1 Experimental Setup

Algorithms. We compare the following eight algorithms for computing MEBs or building coresets for MEB in our experiments.

• •

COMEB [18]: a combinatorial algorithm for computing the exact MEB of a set of points in the Euclidean space. It has an exponential time complexity w.r.t. the dimension $m$. Moreover, it is not applicable to kernelized MEB.

• •

CoreMEB [4]: a batch algorithm for constructing a $(1+\varepsilon)$-approximate coreset for the MEB of a set of points. The procedure is as described in Algorithm 1.

• •

SSMEB [33]: a $1.5$-approximation algorithm for computing a MEB in an append-only stream. We adopt the method described in [27] to compute a kernelized MEB by SSMEB.

• •

Blurred Ball Cover (BBC) [1]: an append-only streaming algorithm to maintain a $(\sqrt{2}+\varepsilon)$-approximate coreset for MEB.

• •

DyMEB [12]: a $1.22$-approximate dynamic algorithm for MEB computation. It keeps a data structure that permits to insert/delete random points without fully reconstructions.

• •

AOMEB: our append-only streaming algorithm presented in Section 3.1. It has the same $(\sqrt{2}+\varepsilon)$-approximation ratio as BBC.

• •

SWMEB: our first sliding-window algorithm presented in Section 3.2. It can maintain a $(\sqrt{2}+\varepsilon)$-coreset for MEB over the sliding window.

• •

SWMEB+: our second sliding-window algorithm presented in Section 3.3. It has higher efficiency than SWMEB at the expense of a worse approximation ratio.

We do not compare with the algorithms in [2, 11] because they cannot scale to the datasets with $m>10$. In our experiments, all algorithms run in the mini-batch-mode with batch size $b=100$. Furthermore, all algorithms, except for SWMEB and SWMEB+, cannot directly work in the sliding window model. The batch and append-only streaming algorithms store the entire sliding window and rerun from scratch for each update. DyMEB also stores the entire sliding window for tracking the expired point to delete. For each update, it must execute one deletion/insertion for every expired/arrival point in a mini-batch to maintain the coreset for MEB w.r.t. the up-to-date window. The implementations of our algorithms are available at https://github.com/yhwang1990/SW-MEB.

Datasets. The dataset statistics are listed in Table 1. We use 6 real-world datasets and 1 synthetic dataset for evaluation. All real-world datasets are downloaded from publicly available sources, e.g., UCI Machine Learning Repository222https://archive.ics.uci.edu/ml/index.php, LIBSVM333https://www.csie.ntu.edu.tw/~cjlin/libsvmtools/datasets, SNAP444http://snap.stanford.edu, and TEXMAX555http://corpus-texmex.irisa.fr. The generation procedure of the Synthetic dataset is as follows. We first decide the dimension $m$. By default, we set $m=50$. For testing the scalability of different algorithms w.r.t. $m$, we vary $m$ from 10 to 100 and from 1,000 to 10,000. Then, we generate a point by drawing the values of $m$ dimensions from a normal distribution $\mathcal{N}(0,1)$ independently. For kernelized MEB, we adopt the Gaussian kernel $k(\mathbf{p}_{i},\mathbf{p}_{j})=\exp(-d(\mathbf{p}_{i},\mathbf{p}_{j})^{2}/\gamma)$ where $\gamma=\frac{1}{n^{2}}\sum_{i,j=1}^{n}d(\mathbf{p}_{i},\mathbf{p}_{j})^{2}$ is the kernel width. In practice, we sample 10,000 points randomly from each dataset to compute $\gamma$. The results are also listed in Table 1. Note that the values of $\gamma$ vary with $m$ ($\gamma\approx 2m$) on the Synthetic dataset. More details on datasets are provided in Appendix A.3.

In an experiment, all points in a dataset are processed sequentially by each algorithm in the same order as a stream and the results are recorded for every batch of points.

Parameter Setting. The default window size $N$ is $10^{5}$ in all experiments except the ones for testing the scalability w.r.t. $N$, where we vary $N$ from $10^{5}$ to $10^{6}$. The parameter $\varepsilon_{1}$ in AOMEB, SWMEB, and SWMEB+ (as well as $\varepsilon$ in CoreMEB, BBC, and DyMEB) is $10^{-3}$ for Euclidean MEB and $10^{-4}$ for kernelized MEB. Then, we use $\varepsilon_{2}=0.1$ in SWMEB and $\varepsilon_{2}=\min(4^{i-1}\cdot\frac{\varepsilon_{1}}{10},0.1)$ for each $x_{i}\in X_{t}$ in SWMEB+. Finally, the partition size $L$ in SWMEB is $\frac{N}{10}$. The procedure of parameter tuning is shown in Appendix A.5.

Environment. All experiments are conducted on a server running Ubuntu 16.04 with a 1.9GHz processor and 128 GB memory. The detailed hardware and software configurations are in Appendices A.1 and A.2.

4.2 Effectiveness and Efficiency

The result quality of different algorithms is evaluated by average error computed as follows. For a sliding window $W_{t}$, each algorithm (except COMEB and SSMEB) returns a coreset $S_{t}$ for $\mathsf{MEB}(W_{t})$. The errors of coreset-based algorithms are acquired based on the definition of coresets: we compute the minimal $\lambda^{\prime}$ such that $W_{t}\subset\lambda^{\prime}\cdot\mathsf{MEB}(S_{t})$, and use $\varepsilon^{\prime}=\frac{\lambda^{\prime}\cdot r^{*}(S_{t})-r^{*}(W_{t})}{r^{*}(W_{t})}$ for the relative error. Since SSMEB directly returns MEBs, we compute its relative error according to the definition of approximate MEB, i.e., $\varepsilon^{\prime}=\frac{r(B^{\prime}_{t})-r^{*}(W_{t})}{r^{*}(W_{t})}$ where $B^{\prime}_{t}$ is the approximate MEB for $W_{t}$. It is noted that for kernelized MEB we use the radius of the MEB returned by CoreMEB when $\varepsilon=10^{-9}$ (i.e., the relative error is within $10^{-9}$) as $r^{*}(W_{t})$ since exact MEBs are intractable in a RKHS [28]. On each dataset, we take $100$ timestamps over the stream, obtain the result of each algorithm for the sliding window at each sampled timestamp, compute the relative errors of obtained results, and use the average error as the quality metric. The efficiency of different algorithms is evaluated by average update time. We record the CPU time of each algorithm to update every batch of $100$ points and compute the average time per update over the entire dataset.

The average errors and average update time of different algorithms for Euclidean and kernelized MEBs are presented in Tables 2–5 respectively. Note that COMEB cannot be applied to kernelized MEB and thus is not included in Tables 4 and 5.

In general, all algorithms except SSMEB can provide MEBs or coresets for MEB with rather small errors (at most $1.5\%$) in all cases. Furthermore, the update time for kernelized MEBs is one to two orders of magnitude longer than that for Euclidean MEB due to (1) higher time complexity for distance evaluation, (2) larger coresets caused by the infinite dimensionality of RKHS, and (3) smaller $\varepsilon_{1}$ (or $\varepsilon$) used. In terms of effectiveness, the error of COMEB for Euclidean MEB is always [math] because it can return exact results. Additionally, as CoreMEB guarantees a $(1+\varepsilon)$-approximation ratio theoretically, the errors of the coresets returned by CoreMEB are less than $10^{-3}$ for Euclidean MEB and $10^{-4}$ for kernelized MEB. With regard to the efficiency for Euclidean MEB, COMEB only runs faster than CoreMEB on low-dimensional datasets (i.e., Gowalla and HIGGS). CoreMEB can always outperform COMEB when $m\geq 30$. The reason is that COMEB has a higher time complexity than CoreMEB for its exponential dependency on $m$.

There are three append-only streaming algorithms in our experiments, namely BBC, SSMEB, and our proposed AOMEB. Among them, the result quality of SSMEB is not competitive with any other algorithms, though its update time is shorter than BBC and AOMEB. This is because the simple geometric method for update in SSMEB is efficient but largely inaccurate. The errors of BBC and AOMEB are slightly higher than those of CoreMEB while they are more efficient than CoreMEB. These experimental results are consistent with our theoretical analysis: AOMEB and BBC have a lower $(\sqrt{2}+\varepsilon)$-approximation ratio but only require a single-pass scan over the window. Finally, AOMEB can run up to $9$ times faster than BBC while having similar or better coreset quality because AOMEB maintains fewer MEBs than BBC, which leads to a more efficient update procedure.

The dynamic algorithm, i.e., DyMEB, shows slightly better coreset quality than AOMEB and BBC but runs even slower than CoreMEB. There are two reasons for such observations: first, the data structure maintained by DyMEB contains all points in $W_{t}$ for coreset construction, which naturally leads to good quality; second, the performance of DyMEB depends on the assumption that the probability of deleting any existing point is equal. However, the sliding window model always deletes the earliest point, which obviously violates this assumption. As a result, when expired points are deleted, DyMEB frequently calls for (partial) coreset reconstructions. In practice, the average update time of DyMEB even exceeds the time to build the coreset from scratch.

Finally, our sliding-window algorithms, namely SWMEB and SWMEB+, achieve two to four orders of magnitude speedups over CoreMEB across all datasets for both Euclidean and kernelized MEBs. In addition, they run $10$ to $150$ times faster than AOMEB. The reason for their superior efficiency is that they can maintain the coreset incrementally over the sliding window without rebuilding from scratch. In terms of effectiveness, the errors of SWMEB and SWMEB+ are slightly higher than those of AOMEB. Although both algorithms use AOMEB instances to provide the coresets, the index schemes inevitably cause quality losses since the points between the beginning of $W_{t}$ (i.e., $t-N+1$) and the first non-expired index (i.e., $x_{1,{s_{1}}}$ in SWMEB and $x_{2}$ in SWMEB+) are missing from the coresets. Lastly, SWMEB+ runs up to $14$ times faster than SWMEB owing to the fact that SWMEB+ maintains fewer indices than SWMEB. At the same time, SWMEB+ can return coresets with nearly equal quality to SWMEB, which means that the radius-based index scheme of SWMEB+, though having a worse theoretical guarantee, is competitive with the partition-based scheme of SWMEB empirically.

4.3 Scalability

We compare the scalability of different algorithms with varying the window size $N$ and dimension $m$. The performance for Euclidean and kernelized MEBs with varying $N$ is shown in Figures 4 and 5. Here we only present the results on the Census dataset. The experimental results on other datasets are shown in Appendix B.

First, we observe that the average errors of different algorithms basically remain stable w.r.t. $N$ (around $10^{-3}$ for SWMEB & SWMEB+). The update time of all algorithms, except SWMEB and SWMEB+, increases with $N$ because the number of points processed per update is equal to $N$. For SWMEB and SWMEB+, the number of indices in $X_{t}$ hardly changes $N$. Besides, the update frequency of any AOMEB instance decreases over time, when the coreset grows larger and fewer new points can be added. Since both algorithms contain “older” instances when $N$ is larger, they take less time for each update on average. In addition, the coreset sizes of our algorithms are $3$–$5$ times larger than that of CoreMEB because the greedy strategy used by AOMEB inevitably adds some redundant points to coresets. Despite this, the coreset sizes are at most $0.3\%\cdot N$. Moreover, the coreset size for kernelized MEBs is around $5$ times larger then that for Euclidean MEBs due to the infinite dimensionality of RKHS. In terms of space, SWMEB stores up to $100$ times more points than SWMEB+ since it not only maintains more AOMEB instances that SWMEB+ but also needs to keep a buffer $Q$, whose size is at most $10\%\cdot N$. Specifically, SWMEB stores at most $13\%\cdot N$ points while SWMEB+ only keeps up to $2,000$ points, which barely changes with $N$.

The update time and space usage of different algorithms with varying the dimension $m$ is shown in Figure 6. Here we only present the results for Euclidean MEBs as the trend is generally similar for kernelized MEBs. As plotted in log-log scale, we can observe that the update time of all algorithms except SSMEB increases almost linearly with $m$. Nevertheless, SWMEB and SWMEB+ still demonstrate their superiority in efficiency on high-dimensional datasets. Furthermore, both algorithms store more points when $m$ is larger because the coreset size grows with $m$. But even when $m=10,000$ the ratios of stored points are at most $25\%$ and $5\%$ for SWMEB and SWMEB+ respectively, whereas any other algorithms require to store the entire window of points.

5 Conclusion

We studied the problem of maintaining a coreset for the MEB of a sliding window of points in this paper. Firstly, we proposed the AOMEB algorithm to maintain a $(\sqrt{2}+\varepsilon)$-coreset for MEB in an append-only stream. Then, based on AOMEB, we proposed two sliding-window algorithms, namely SWMEB and SWMEB+, for coreset maintenance with constant approximation ratios. We further generalized our proposed algorithms for kernelized MEBs. Empirically, SWMEB and SWMEB+ improved the efficiency of the state-of-the-art batch algorithm by up to four orders of magnitude while providing coresets with rather small errors compared to the optimal ones. For future work, we plan to explore the applications of coresets to various data mining and machine learning problems in the streaming or sliding window model.

Appendix A Details on Experimental Setup

In this section, we present the details of our experimental setup for reproducibility.

A.1 Hardware Configuration

All the experiments are conducted on a server with the following specifications:

• •

CPU: Intel(R) Xeon(R) E7-4820 v3 @ 1.90GHz

• •

Memory: 128GB (8$\times$16GB) RAM 2133MHz DDR4 memory

• •

Hard Disk: a 480GB SATA-III solid-state drive

A.2 Software Environment

The server runs Ubuntu GNU/Linux 16.04.3 LTS 64-bit with kernel v4.11.0-rc2. All the code is written in Java 8 only using the standard libraries. No third-party software/library is required. The version of JVM for compilation is Java HotSpot(TM) 64-Bit Server VM 18.3 (build 10.0.2+13). In the experiments, each instance is limited to use a single thread for computation to ensure the fairness of comparison. In addition, we use a JVM option “-Xmx80000m” to restrict the maximum heap size used by each instance. Note that the memory usage is much less than 80GB and our purpose is to guarantee that the bottleneck is not in memory and I/O.

A.3 Datasets

The statistics of the datasets used in our experiments are listed in Table 1. Here we briefly describe the real-world datasets we use and the preprocessing procedures.

• •

Census contains a one percent sample of the Public Use Microdata Samples person records drawn from the full 1990 census sample. It is downloaded from UCI Machine Learning Repository.

• •

CovType is a dataset for predicting forest cover type from cartographic variables. It is available on UCI Machine Learning Repository and LIBSVM. We download the dataset from LIBSVM. In the preprocessing, we only retain the data with class label “0”.

• •

GIST is an image collection retrieved from TEXMAX.

• •

Gowalla is a collection of user check-ins over the period of February 2009 – October 2010 on gowalla.com. It is downloaded from SNAP. In the preprocessing, we extract the latitude and longitude of each check-in as a two-dimensional point and dispose other attributes.

• •

HIGGS is a dataset for distinguishing between a signal process which produces Higgs bosons and a background process which does not. It is downloaded from UCI Machine Learning Repository. In the preprocessing, we only retain the data with class label “1” (i.e., signal process).

• •

SIFT is an image collection retrieved from TEXMAX.

The procedure of generating Synthetic has been introduced in the main paper. We transform each dataset into a single file that stores the points in dense format: each point is represented by one line in the file and different dimensions of the point are split by a single space. After dataset preprocessing and format transformation, we shuffle all points of each dataset randomly. It is guaranteed that a dataset must be processed by all algorithms in the same order for the fairness of comparison.

A.4 Implementation Issues

The algorithms we compare in the experiments are as listed in Section 4.1. The implementations of our algorithms are available at https://github.com/yhwang1990/SW-MEB. Here we discuss the implementations of these algorithms.

We use the Java code published by the authors, which is available at https://github.com/hbf/miniball, for the implementation of COMEB [18]. All the other algorithms are implemented by ourselves. First of all, the basic scheme of CoreMEB [4] is presented in Algorithm 1. Our practical implementations are based on two improved versions of CoreMEB, i.e., Figure 2 in [23] for Euclidean MEB and Algorithm 4.1 in [32] for kernelized MEB. They use the same scheme as shown in Algorithm 1 but have lower computational costs and quicker convergence rates. The implementation of SSMEB is based on Section 2 of [33], which is extended to kernelized MEB according to Section 4.2 of [27]. We implement the Blurred Ball Cover (BBC) algorithm according to Section 2 of [1]. In addition, DyMEB is implemented based on Algorithms 1–3 in [12]. Finally, our proposed methods, i.e., AOMEB, SWMEB and SWMEB+, are implemented according to Algorithms 2–4 in this paper.

A.5 Parameter Tuning

The procedure of parameter tuning is as follows. Generally, there are three parameters used in our experiments: $\varepsilon_{1}$ (or $\varepsilon$) in all algorithms except COMEB and SSMEB, $\varepsilon_{2}$ in SWMEB and SWMEB+, the partition size $L$ in SWMEB.

First of all, we choose appropriate $\varepsilon_{1}$ for Euclidean and kernelized MEBs respectively. The parameter $\varepsilon_{1}$ determines the trade-off between coreset quality and efficiency. We test the effect of $\varepsilon_{1}$ by validating on $[10^{-1},10^{-2},\ldots,10^{-6}]$. For Euclidean MEB, the coreset quality hardly improves but the running time increases rapidly when $\varepsilon_{1}<10^{-3}$. Therefore, we set $\varepsilon_{1}=10^{-3}$ for Euclidean MEB. Using the same method, we set $\varepsilon_{1}=10^{-4}$ for kernelized MEB.

After choosing appropriate $\varepsilon_{1}$’s, we further test the remaining parameters in SWMEB and SWMEB+. Specifically, $\varepsilon_{2}$ is also selected from $[10^{-1},10^{-2},\ldots,10^{-6}]$. In SWMEB, $\varepsilon_{2}$ affects the index construction on each partition. When $\varepsilon_{2}\leq 0.01$, SWMEB suffers from high overhead because of containing too many indices. Thus, we use $\varepsilon_{2}=0.1$ for SWMEB. Additionally, to avoid the indices being too sparse, we restrict the maximal distance between any neighboring indices in the same partition to $\frac{L}{10}$. In SWMEB+, $\varepsilon_{2}$ adjusts the number of indices in $X_{t}$. Firstly, we observe that the coreset quality cannot be improved any more when $\varepsilon_{2}<\frac{\varepsilon_{1}}{10}$. Secondly, we scale $\varepsilon_{2}$ by a factor of $\lambda>1$ among indices, i.e., $\varepsilon_{2}=\lambda^{i-1}\cdot\frac{\varepsilon_{1}}{10}$ for $x_{i}\in X_{t}$, to reduce the index size without seriously affecting the quality. We select $\lambda$ from $\{2,4,8,16\}$ and use $\lambda=4$ for SWMEB+ since the quality seriously degrades when $\lambda>4$. Thirdly, we set $0.1$ as the upper bound of $\varepsilon_{2}$ to ensure the theoretical soundness. To sum up, we use $\varepsilon_{2}=\min(4^{i-1}\cdot\frac{\varepsilon_{1}}{10},0.1)$ for each $x_{i}\in X_{t}$ in SWMEB+. Finally, the partition size $L$ in SWMEB affects the balance between space and time complexity. If $L$ is smaller, $W_{t}$ will be divided into more partitions, which leads to more indices in $X_{t}$, but fewer points will be stored in the buffer $Q$, and vice versa. We try to select $L$ from range $[\frac{N}{5},\frac{N}{10},\ldots]$. The results show that SWMEB cannot scale to large datasets when $L<\frac{N}{10}$ for too many indices in $X_{t}$. Therefore, the partition size $L$ in SWMEB is set to $\frac{N}{10}$.

Appendix B Additional Experimental Results

We present additional experimental results for scalability on different datasets in Figures 7–11.