Multipass Linear Sketches for Geometric LP-Type Problems

In this paper, we explore using linear sketches to solve approximate LP-type problems in high accuracy regimes, in which , where $d$ is the dimensionality of the problem and $\epsilon$ is the approximation accuracy. This regime is applicable to many real world big data applications where $1/\epsilon$ might be several orders of magnitude larger than $d$, such as machine learning tasks of regression or classification. Take for example [38], where $d=2$, $1/\epsilon \in [{10}^2,{10}^5]$, and $n={10}^5$, or [45] where for some data sets $1/\epsilon$ is taken to be much larger than $d$ yet small compared to $n$.

We present an algorithm in the multiple linear sketch model for solving certain classes of LP-type problems, which can be utilized in several prominent big data models. The multiple linear sketch model, which we define further in Section 1.2, can be thought of as a multipass extension of the linear sketch model, where one is allowed to choose a fresh linear sketch in each pass of a multipass streaming algorithm in order to have a better representation of the data. There are two significant benefits of our sketching approach. Firstly, our algorithm works in the presence of input point deletions, such as in the turnstile model. Secondly, it allows our algorithm to work in the presence of duplicated points in the input, which generally can prove problematic for sampling based algorithms.

We are given a universe $S={\{-\mathrm{\Delta },-\mathrm{\Delta }+1,\mathrm{\dots },\mathrm{\Delta }\}}^d$ of $d$-dimensional points, where each point can be represented by a word of size $d\log \left(\mathrm{\Delta }n\right)$. This defines our cost model, and we give our space complexity in the number of points/words. We use a linear sketch to reduce an input set $P\subseteq S$ of $n$ input points into a point set $Q=\{{(1+\epsilon )}^l\cdot 𝒆\mid 𝒆\in E_{{𝒑}_c},l\in \left[{\log }_{1+\epsilon }\delta \right]\}$, where $E_{{𝒑}_c}$ is a metric $\epsilon$-net centered at a point ${𝒑}_c\in P$ and $\delta$ is an $O(1)$ approximation for the distance of the furthest point in $P$ to ${𝒑}_c$. Our algorithm is also applicable when $S={\{-1,-1+\frac{1}{\mathrm{\Delta }},\mathrm{\dots },1\}}^d$. Now, we reduce $P$ directly into the point set $Q=E$ where $E$ is the metric $\epsilon$-net centered at the origin. Intuitively, this has the effect of snapping each input point to the closest of a small number of net points, thereby reducing our input space. This linear sketch can be though of as an $N\times \left|S\right|$ matrix $A$ for $N=l\left|E\right|\ll n$ where each row corresponds to one possible point ${(1+\epsilon )}^l\cdot 𝒆$, containing 1’s for all the points $𝒑\in S$ with direction closest to $𝒆$ and norm closest to ${(1+\epsilon )}^l$. The major benefit of using this linear sketch is that we don’t actually store our metric $\epsilon$-net. We can calculate where each point is snapped to quickly on the fly, so we do not suffer the exponential in $d$ space complexity of $\epsilon$-kernel based formulations.

This linear sketch of the data can now be used to run our algorithm on, treating it as a smaller virtual input compared to the original input. Our general algorithm for solving LP-type problems follows the sampling idea presented by [8], which in turn is based on Clarkson’s algorithm for linear programs in low dimensions [21], where input points are assigned weights, and a weighted sampling procedure is used in order to solve the LP-type problem on smaller samples from the input. By modifying the weights throughout multiple iterations, the algorithm quickly converges on samples that contain a basis for the LP-type problem, for which the solution is the solution of the entire input set.

As our algorithm is based on linear sketches, we can utilize it in many big data models. We can solve problems on input streams that we can make multiple linear scans over in the Multipass Streaming model. Additionally, our algorithm works for streams with deletions, such as in the Multipass Strict Turnstile model. Finally, we can also formulate our algorithm for distributed inputs, minimizing communication load in the Coordinator and Parallel Computation models. A usage of these models is seeing the symmetric difference of data over time across multiple servers. Over two servers with large datasets, communication can be thought of as passes over a stream with deletions. Further, we can find the symmetric difference of server datasets in the Parallel Computation model as well. This is helpful in capturing change in data to predict current trends, see e.g. [24].

We give the following result for LP-type problems that are bounded in combinatorial and VC-dimension. More detailed bounds and time bounds are given in Sections A.1, A.2, A.3, and A.4.

For most problems, we have that $N\in {\left(1/\epsilon \right)}^{O(d)}$. With this, we get the following corollary.

Our algorithm provides powerful results when we set the parameter $s$ equal to $\log N$. As we have that $\log N\in \mathrm{poly}(d,\log \left(1/\epsilon \right))$, this allows us to achieve algorithms with pass and space complexities polynomial in $d$ and polylogarithmic in $1/\epsilon$.

We can use this general algorithm to find approximate solutions for many LP-type problems which satisfy some elementary properties, as described in Section 3. We emphasize that most LP-type problems fulfill these properties. We present applications of our general algorithm for the following problems. Note that the Linear SVM, Bounded LP, and Bounded SDP problems do not have the $\log \log (\mathrm{\Delta }d)$ pass increase in the strict turnstile model.

• $\mathrm{■}$

  MEB: Given $n$ points $P\subseteq {\{-\mathrm{\Delta },-\mathrm{\Delta }+1,\mathrm{\dots },\mathrm{\Delta }\}}^d$, our objective is to find a minimal ball enclosing the input points. We present an algorithm to output a ball that encloses $P$ with radius at most $(1+O(\epsilon ))$ times optimal. For the MEB problem, we have $N\in {\left(1/\epsilon \right)}^{O(d)}$.

• $\mathrm{■}$

  Linear SVM: Given $n$ labeled points $P\subseteq {\{-1.-1+\frac{1}{\mathrm{\Delta }},\mathrm{\dots },1\}}^d\times \{\pm 1\}$ and some $\gamma >0$ for which the points are either linearly inseparable or $\gamma$-separable, our objective is to find a separating hyperplane of $P$ with the largest margin. We present an algorithm to output either that the points are linearly inseparable, or a separating hyperplane $(𝒖,b)$ where ${\Vert 𝒖\Vert }^2$ is at most $(1+O(\epsilon ))$ times optimal. For the linear SVM problem, we have $N\in {\left(1/\epsilon \gamma \right)}^{O(d)}$.

• $\mathrm{■}$

  Bounded LP: Given an objective vector $𝒄\in {\{-1,-1+\frac{1}{\mathrm{\Delta }},\mathrm{\dots },1\}}^d$ such that $\Vert 𝒄\Vert \in O(1)$ and $n$ constraints $P\subseteq {\{-1.-1+\frac{1}{\mathrm{\Delta }},\mathrm{\dots },1\}}^{d+1}$, where for each constraint $({𝒂}_i,b_i)\in P$, $\Vert {𝒂}_i\Vert ,\left|b_i\right|\in O(1)$, our objective is to find a solution $𝒙$ with $\Vert 𝒙\Vert \in O(1)$ that maximizes ${𝒄}^T𝒙$ while satisfying the constraints. We present an algorithm to output a solution $O(\epsilon )$ within each constraint, and for which the objective is at most $O(\epsilon )$ smaller than optimal. For bounded LP problems, we have $N\in {\left(1/\epsilon \right)}^{O(d)}$.

• $\mathrm{■}$

  Bounded SDP: Given an objective matrix $C\in {\{-1.-1+\frac{1}{\mathrm{\Delta }},\mathrm{\dots },1\}}^{d\times d}$ such that ${\Vert C\Vert }_F\le 1$ and $n$ constraints $P\subseteq {\{-1.-1+\frac{1}{\mathrm{\Delta }},\mathrm{\dots },1\}}^{d\times d}\times \{-1,-1+\frac{1}{\mathrm{\Delta }},\mathrm{\dots },1\}$, where for each constraint $(A^{(i)},b^{(i)})\in P$, ${\Vert A^{(i)}\Vert }_2,\left|b^{(i)}\right|\le 1$, ${\Vert A^{(i)}\Vert }_F\le F$, and the number of nonzero entries of $A^{(i)}$ is bounded by $S$, our objective is to find a positive semidefinite solution $X$ of unit trace that maximizes the objective ${⟨C,X⟩}_F$ while satisfying the constraints. We present an algorithm to output a solution that is $O(\epsilon )$ within each constraint, and for which the objective is at most $O(\epsilon )$ smaller than optimal. As we are working with $d\times d$ matrices, the combinatorial and VC dimensions are $O(d^2)$, and we have $N\in \left(\genfrac{}{}{0pt}{}{d^2}{S}\right)\cdot {\left(1/\epsilon \right)}^{O(S)}+{\left(1/\epsilon \right)}^{O(d)}$. $\parallel \cdot {\parallel }_2$, $\parallel \cdot {\parallel }_F$, and ${⟨\cdot ,\cdot ⟩}_F$ refer to the spectral and Frobenius norms and Frobenius inner product respectively, explained further in Appendix A.8.

Table 1 compares our results to those of prior work, including the exact LP-type problem algorithm of [8], and various $\epsilon$-approximation LP-type problems.

For the approximate MEB problem, there are many results, as discussed in Section 1.2. One such is the Zarrabi-Zadeh result, which is among many that utilize an $\epsilon$-kernel [53]. While these algorithms are single pass, they require space complexity of ${\left(1/\epsilon \right)}^{O(d)}$. Thus, we present a large increase in space efficiency over them. We also present a lower bound on the space complexity for any one-pass algorithm for the MEB problem of ${\left(1/\epsilon \right)}^{\mathrm{\Omega }(d)}$ in Section 4, further motivating our multi-pass approach. While the result by Bâdoiu and Clarkson has space complexity polynomial in $d$, it is also polynomial in $1/\epsilon$ in passes and space [9]. Thus, in the high accuracy regime, we present an exponential improvement over their result.

Another approach to the $(1+\epsilon )$-MEB problem is by using the ellipsoid algorithm with a similar metric $\epsilon$-net construction. While this algorithm also can achieve $\mathrm{poly}(d,\log \left(1/\epsilon \right))$ pass and space complexity, we present two main improvements over using this approach. Most importantly, as our algorithm is in the multiple linear sketch model, we are able to apply our result to streams with addition and removal of points, and to models where the input is distributed. As the ellipsoid algorithm is not in the multiple linear sketch model, it cannot, for example, find the MEB of the symmetric difference of pointsets from two machines. Further, the bulk of our algorithm utilizes ${\mathrm{\ell }}_0$ samplers and estimators and simple arithmetic, which can allow our algorithm to be practical to implement.

A strong result to linear programs is by [34]. They use the interior point method for solving linear programs in order to get their result. There are two significant differences of our result to theirs. Firstly, as with the ellipsoid method, while this result might present an improvement in the multipass streaming model, our algorithm is able to handle point additions and deletions, as well as distributed inputs. Secondly, their tilde notation hides logarithmic dependence on $n$. Our result completely removes this dependence in both the pass and space complexities, showing that $n$ is not an inherent parameter in this problem.

Our improvements over the work by [8] is similar. By using linear sketches we are able to extend their results to models with point additions and deletions such as the turnstile model, and additionally deal with duplicate points due to our usage of ${\mathrm{\ell }}_0$ samplers and estimators.

As we present in Section A.7.1, our algorithm can be used to solve the linear classification problem. Note that our result finds an exact classifier. We compare our result to the result by [23]. In the high accuracy regime, this represents an exponential improvement over their result in pass and space complexity.

Another application of our framework is for additive $\epsilon$-approximation bounded SDP problems. As we present in Section A.8.1, our formulation can be used to solve for the saddle point problem within the unit simplex. For this problem, there is a result by [26] in the number of arithmetic operations. We achieve an exponential improvement in $1/\epsilon$ dependence over their result, significant in the high accuracy regime. Another result is by [43] who achieve a significant result for SDPs with high dimensionality. In the high accuracy regime with a large number of constraints $n\gg 1/\epsilon$, this presents a loss in number of passes polynomial in $d$, while achieving a significant gain in space complexity, as our result does not depend on $n$ at all.

Approximation algorithms for LP-type problems in big data models have been well-studied, where the trade-off between approximation ratio, number of passes, and space complexity has given rise to different techniques being used for different regimes.

For the MEB problem, there is the folklore result of a $1$-pass algorithm, where one stores the first point and a running furthest point from it, giving a $2$-approximation. There is also a classic $1$-pass $3/2$-approximation algorithm by Zarrabi-Zadeh and Chan, that uses $O(d)$ space complexity [54]. Allowing $O(d)$ passes, an algorithm by Agarwal and Sharathkumar, with further analysis by Chan and Pathak, achieves a $1.22$-approximation with $O\left(\left(d/{\epsilon }^3\log \left(1/\epsilon \right)\right)\right)$ space complexity [2, 19]. For the $(1+\epsilon )$-approximation problem, much work is done in the single-pass model. These algorithms are based on storing extremal points, called an $\epsilon$-kernel, and achieve a space complexity of ${\left(1/\epsilon \right)}^{O(d)}$ [1, 18, 3, 53]. There is also a body of work on bounding the size of these $\epsilon$-kernels, as they are useful to achieve low-space representations of large data [32, 10, 52, 22]. [7] presents a lower bound of ${\left(1/\epsilon \right)}^{\mathrm{\Omega }(d)}$ on the size of an $\epsilon$-kernel. In the multi-pass model, there is the classic result by Bâdoiu and Clarkson [9]. This algorithm is widely used in big data models, for its linear dependence on $d/\epsilon$ for space and $1/\epsilon$ for pass complexity, and simplicity [38, 45].

Linear programs are frequently studied in big data models, with a major focus on exact multi-pass LP results. A significant recent result is by [8], who in the multi-pass model achieve a result with $O(d\log n)$ passes and $\mathrm{poly}(d,\log n)$ space usage. Their result is also applicable to other LP-type problems as well as parallel computation models. Their algorithm builds on Clarkson’s sampling based method [21], and utilizes reservoir sampling to sample a $\mu$-net [20]. We build on the ideas presented by [8], to solve LP-type problem approximations. There is also recent work on multipass approximations for LPs, with a result by [34] who achieve an $\tilde{O}\left(\sqrt{d}\log \left(1/\epsilon \right)\right)$ pass $\tilde{O}\left(d^2\right)$ space algorithm for solving $(1-\epsilon )$-approximation LPs. For packing LPs [4] gives a $(1+\epsilon )$-approximation algorithm, and for covering LPs [28] gives a $(1+\epsilon )$-approximation algorithm. The linear classification problem is also well studied, as a fundamental machine learning problem. Its classic result is the mistake bound of the perceptron algorithm [39]. More recently, a result by [23] achieves $O\left(1/{\epsilon }^2\log n\right)$ iterations to find a linear classifier.

There is a lot of study on approximate semidefinite programs, focusing on different regimes. [26] give a result on solving bounded SDP saddle point problems with certain assumptions about the problem input. We use similar assumptions to compare our results to theirs. For the regime with a low number of constraints and a large dimension, [43] achieves an $O\left(\sqrt{d}\log d/\epsilon \right)$ pass algorithm with $\tilde{O}\left(n^2+d^2\right)$ space complexity.

As our algorithm handles models with point deletions, we must make sure to construct our metric $\epsilon$-net around a non-deleted point and find an approximation for the distance of the furthest non-deleted point from it, to define the size of our net. To do this, we use ${\mathrm{\ell }}_0$ samplers, which work in the presence of deletions. We first sample any non-deleted point ${𝒑}_c$, which becomes the center of our net. Now, by construction of $P$, we have an upper and lower bound on the distance of any point from ${𝒑}_c$. We simply binary search over the powers of $2$ in this range and sample a non-deleted point with norm larger than the current power of $2$ in each iteration. By the end of the binary search, which takes $\log \log \left(\mathrm{\Delta }d\right)$ iterations, we end up with a $2$-approximation for $\delta$, considering only non-deleted points. Alternatively, one can find this using $2$ iterations. First, we similarly sample a non-deleted input point to become ${𝒑}_c$. Second, we use an ${\mathrm{\ell }}_0$ sampler for each power of $2$, and in one pass over the points add them to the appropriate sampler. This way, instead of doing the steps of a binary search, we perform them all at once. While this negates the extra pass complexity, it adds a space complexity in words logarithmic in $\mathrm{\Delta }$, since we need to store that many samplers. We will choose the $\log \log \left(\mathrm{\Delta }d\right)$ iteration solution. We explain this procedure in detail in Section A.2.

When adapting our sampling procedure to distributed models, one cannot take a sample of $m$ points by asking machines to sample $m$ points each, as this would lead to a load of $m\cdot k$ points. Further, since we are doing weighted sampling, machines need to sample in accordance with other machines’ weights. In order to fix these problems, we employ the following protocol, which we describe for the coordinator model, noting that it can be implemented in the parallel computation model by having one machine assume the coordinator role. Each machine sends the coordinator the total weight of their subset. The coordinator then calculates how many of the $m$ sample points $B$ each machine shall sample, and sends them this number. Now, the machines can sample and send only their part of the sample, meaning $m$ points are sent in total now, instead of per machine.

In these big data models, we are able to solve many useful low dimensional LP-type problems. We show how to apply our algorithm in three separate ways for the problems. For the MEB problem, we give a multiplicative $(1+O(\epsilon ))$-approximation to the optimal ball. For the Linear SVM problem, we show how our algorithm can be presented as a Monte Carlo randomized algorithm. We achieve a one-sided error where if a point set is inseparable, we always output as such, and if a point set is separable, we output a $(1+O(\epsilon ))$-approximation for the optimal separating hyperplane with high probability. For Linear Programming problems, instead of a multiplicative approximation, we give an additive $\epsilon$-approximation to bounded LP problems. Full complexity results for the models are presented in their respective subsections of Section A. and an overview is shown in Theorem 1.

We can solve SDP problems up to additive $\epsilon$-approximation. For SDP problems, notice that by taking the vectorization $𝒙$ of a solution $X$, we can rewrite a given SDP problem as an LP in $d^2$ dimensions. We must also define the positive semidefinite constraint $X⪰0$. This is equivalent to constraining ${𝒚}^{T}X𝒚\ge 0$ for all vectors $𝒚$ where $\Vert y\Vert =1$. Notice that this requires an infinite number of constraints. To solve this, we create another lattice along with our original metric $\epsilon$-net. Using this lattice, we can instead use the constraints ${𝒛}^{T}X𝒛\ge 0$ for vectors $𝒛$ on the lattice. Just like property 3, we show that satisfying these constraints satisfies the original semidefinite constraints up to $\epsilon$. Thus, we are able to successfully reduce the SDP problem into an LP, which we can solve up to additive $\epsilon$-approximation, while retaining a small number of constraints.

Being able to solve bounded LP problems up to additive $\epsilon$-approximation allows us to use our algorithm for many important related problems. For example, we are able to solve the linear classification problem by finding a classifier within additive $\epsilon$ of an optimal classifier for a set of points. Since in the problem, the optimal classifier has a separation of at least $\epsilon$, this means that our found classifier is guaranteed to be a separator for the original points.