OSNAP: Faster Numerical Linear Algebra Algorithms via Sparser Subspace Embeddings

An oblivious subspace embedding (OSE) given some parameters $\epsilon$, d is a distribution $\mathcal{D}$ over matrices $\Pi \in \mathbb{R}^{m×n}$ such that for any linear subspace $W \subseteq \mathbb{R}^n$ with dim(W) = d, $\mathbb{P}_{\Pi \sim \mathcal{D}}(\forall x \in W ||\Pi x||_2 \in (1 \pm \epsilon)||x||_2) > 2/3$. We show that a certain class of distributions, Oblivious Sparse Norm-Approximating Projections (OSNAPs), provides OSE's with $m = O(d^{1+\gamma}/\epsilon^2)$, and where every matrix $\Pi$ in the support of the OSE has only $s = O_{\gamma}(1/\epsilon)$ non-zero entries per column, for $\gamma > 0$ any desired constant. Plugging OSNAPs into known algorithms for approximate least squares regression, $\ell_p$ regression, low rank approximation, and approximating leverage scores implies faster algorithms for all these problems. Our main result is essentially a Bai-Yin type theorem in random matrix theory and is likely to be of independent interest: we show that for any fixed $U \in \mathbb{R}^{n×d}$ with orthonormal columns and random sparse $\Pi$, all singular values of $\Pi U$ lie in $[1 - \epsilon, 1 + \epsilon]$ with good probability. This can be seen as a generalization of the sparse Johnson-Lindenstrauss lemma, which was concerned with d = 1. Our methods also recover a slightly sharper version of a main result of [Clarkson-Woodruff, STOC 2013], with a much simpler proof. That is, we show that OSNAPs give an OSE with $m = O(d^2/\epsilon^2)$, $s = 1$.Engineering and Applied Science

Sketching via hashing: from heavy hitters to compressed sensing to sparse fourier transform

Sketching via hashing is a popular and useful method for processing large data sets. Its basic idea is as follows. Suppose that we have a large multi-set of elements S=[formula], and we would like to identify the elements that occur “frequently" in S. The algorithm starts by selecting a hash function h that maps the elements into an array c[1…m]. The array entries are initialized to 0. Then, for each element a ∈ S, the algorithm increments c[h(a)]. At the end of the process, each array entry c[j] contains the count of all data elements a ∈ S mapped to j

New Constructions of RIP Matrices with Fast Multiplication and Fewer Rows

In this paper, we present novel constructions of matrices with the restricted isometry property (RIP) that support fast matrix-vector multiplication. Our guarantees are the best known, and can also be used to obtain the best known guarantees for fast Johnson Lindenstrauss transforms. In compressed sensing, the restricted isometry property is a sufficient condition for the efficient reconstruction of a nearly k-sparse vector $x \in \mathbb{C}^d$ from m linear measurements $\Phi x$. It is desirable for m to be small, and further it is desirable for $\Phi$ to support fast matrix-vector multiplication. Among other applications, fast multiplication improves the runtime of iterative recovery algorithms which repeatedly multiply by $\Phi$ or $\Phi^*$. The main contribution of this work is a novel randomized construction of RIP matrices $\Phi \in \mathbb{C}^{m×d}$, preserving the $\ell_2$ norms of all k-sparse vectors with distortion $1 + \epsilon$, where the matrix-vector multiply $\Phi x$ can be computed in nearly linear time. The number of rows m is on the order of $\epsilon^{−2}klogd log^2(klogd)$, an improvement on previous analyses by a logarithmic factor. Our construction, together with a connection between RIP matrices and the Johnson-Lindenstrauss lemma in [Krahmer-Ward, SIAM. J. Math. Anal. 2011], also implies fast Johnson-Lindenstrauss embeddings with asymptotically fewer rows than previously known. Our construction is actually a recipe for improving any existing family of RIP matrices. Briefly, we apply an appropriate sparse hash matrix with sign flips to any suitable family of RIP matrices. We show that the embedding properties of the original family are maintained, while at the same time improving the number of rows. The main tool in our analysis is a recent bound for the supremum of certain types of Rademacher chaos processes in [Krahmer-Mendelson-Rauhut, Comm. Pure Appl. Math. to appear].Engineering and Applied Science

Toward a unified theory of sparse dimensionality reduction in Euclidean space

Let $\Phi\in\mathbb{R}^{m\times n}$ be a sparse Johnson-Lindenstrauss transform [KN14] with $s$ non-zeroes per column. For a subset $T$ of the unit sphere, $\varepsilon\in(0,1/2)$ given, we study settings for $m,s$ required to ensure $$\mathop{\mathbb{E}}_\Phi \sup_{x\in T} \left|\|\Phi x\|_2^2 - 1 \right| < \varepsilon ,$$ i.e. so that $\Phi$ preserves the norm of every $x\in T$ simultaneously and multiplicatively up to $1+\varepsilon$. We introduce a new complexity parameter, which depends on the geometry of $T$, and show that it suffices to choose $s$ and $m$ such that this parameter is small. Our result is a sparse analog of Gordon's theorem, which was concerned with a dense $\Phi$ having i.i.d. Gaussian entries. We qualitatively unify several results related to the Johnson-Lindenstrauss lemma, subspace embeddings, and Fourier-based restricted isometries. Our work also implies new results in using the sparse Johnson-Lindenstrauss transform in numerical linear algebra, classical and model-based compressed sensing, manifold learning, and constrained least squares problems such as the Lasso

Sparsity lower bounds for dimensionality reducing maps

We give near-tight lower bounds for the sparsity required in several dimensionality reducing linear maps. First, consider the Johnson-Lindenstrauss (JL) lemma which states that for any set of n vectors in Rd there is an A∈Rm x d with m = O(ε-2log n) such that mapping by A preserves the pairwise Euclidean distances up to a 1 pm ε factor. We show there exists a set of n vectors such that any such A with at most s non-zero entries per column must have s = Ω(ε-1log n/log(1/ε)) if m < O(n/log(1/ε)). This improves the lower bound of Ω(min{ε-2, ε-1√(logm d)) by [Dasgupta-Kumar-Sarlos, STOC 2010], which only held against the stronger property of distributional JL, and only against a certain restricted class of distributions. Meanwhile our lower bound is against the JL lemma itself, with no restrictions. Our lower bound matches the sparse JL upper bound of [Kane-Nelson, SODA 2012] up to an O(log(1/ε)) factor. Next, we show that any m x n matrix with the k-restricted isometry property (RIP) with constant distortion must have Ω(k log(n/k)) non-zeroes per column if m=O(k log (n/k)), the optimal number of rows for RIP, and k < n/polylog n. This improves the previous lower bound of Ω(min{k, n/m}) by [Chandar, 2010] and shows that for most k it is impossible to have a sparse RIP matrix with an optimal number of rows. Both lower bounds above also offer a tradeoff between sparsity and the number of rows. Lastly, we show that any oblivious distribution over subspace embedding matrices with 1 non-zero per column and preserving distances in a d dimensional-subspace up to a constant factor must have at least Ω(d2) rows. This matches an upper bound in [Nelson-Nguyên, arXiv abs/1211.1002] and shows the impossibility of obtaining the best of both of constructions in that work, namely 1 non-zero per column and d ⋅ polylog d rows.Engineering and Applied Science

Recovering the Optimal Solution by Dual Random Projection

Random projection has been widely used in data classification. It maps high-dimensional data into a low-dimensional subspace in order to reduce the computational cost in solving the related optimization problem. While previous studies are focused on analyzing the classification performance of using random projection, in this work, we consider the recovery problem, i.e., how to accurately recover the optimal solution to the original optimization problem in the high-dimensional space based on the solution learned from the subspace spanned by random projections. We present a simple algorithm, termed Dual Random Projection, that uses the dual solution of the low-dimensional optimization problem to recover the optimal solution to the original problem. Our theoretical analysis shows that with a high probability, the proposed algorithm is able to accurately recover the optimal solution to the original problem, provided that the data matrix is of low rank or can be well approximated by a low rank matrix.Comment: The 26th Annual Conference on Learning Theory (COLT 2013

Toward a unified theory of sparse dimensionality reduction in Euclidean space

Let $\Phi \in \mathbb{R}^{m×n}$ be a sparse Johnson-Lindenstrauss transform [KN14] with s non-zeroes per column. For a subset T of the unit sphere, $\epsilon \in (0,1/2)$ given, we study settings for m,s required to ensure $\underset {\Phi_{x \in T}} {\mathbb{E} sup} \mid || \Phi x ||^2_2 - 1 \mid < \epsilon$, i.e.\ so that $\Phi$ preserves the norm of every $x \in T$ simultaneously and multiplicatively up to $1+\epsilon$. We introduce a new complexity parameter, which depends on the geometry of T, and show that it suffices to choose s and m such that this parameter is small. Our result is a sparse analog of Gordon's theorem, which was concerned with a dense $\Phi$ having i.i.d. gaussian entries. We qualitatively unify several results related to the Johnson-Lindenstrauss lemma, subspace embeddings, and Fourier-based restricted isometries. Our work also implies new results in using the sparse Johnson-Lindenstrauss transform in numerical linear algebra, classical and model-based compressed sensing, manifold learning, and constrained least squares problems such as the Lasso.Engineering and Applied Science