Oblivious subspace embeddings (OSEs) were introduced by Sarlos [Sar06] to solve linear algebra problems more quickly than traditional methods. An OSE is a distribution of matrices with
such that, for any d-dimensional subspace
, with “high” probability S preserves the norm of every vector in the subspace. OSEs are a generalization of the classic JohnsonLindenstrauss lemma from vectors to subspaces. Formally, we require that with probability
Because A is a “tall” matrix with more rows than columns, the system is overdetermined and there is likely no solution to Ax = b, but regression will find the closest point to b in the space spanned by A. The classic answer to regression is to use the Moore-Penrose pseudoinverse:
is the “pseudoinverse” of A (assuming A has full column rank, which we will typically do for simplicity). This classic solution takes time, where
is the matrix multiplication constant [CW90, Wil12, Gal14]:
time to compute
time to compute the inverse. OSEs speed up the process by replacing (2) with
for an OSE S on d + 1-dimensional spaces. This replaces the regression problem with an
problem, which can be solved more quickly since
lies in the d + 1-dimensional space spanned by b and the columns of A, with high probability S preserves the norm of
That is, S produces a solution which preserves the cost of the regression problem. The running time for this method depends on (1) the reduced dimension m and (2) the time it takes to multiply S by A. We can compute these for “standard” OSE types:
• If S has i.i.d. Gaussian entries, then is sufficient (and in fact,
is required [NN14]). However, computing
time, which is worse than solving the original regression problem (one can speed this up using fast matrix multiplication, though it is still worse than solving the original problem).
• If S is a subsampled randomized Hadamard transform (SRHT) matrix with random sign flips (see Theorem 2.4 in [Woo14] for a survey, and also see [CNW16] which gives a recent improvement) then m increases to we can compute SA using the fast Hadamard transform in O(nd log n) time. This makes the overall regression problem take
• If S is a random sparse matrix with random signs (the “Count-Sketch” matrix), then m = a decreasing function of the sparsity [CW13, MM13, NN13, BDN15, Coh16]. (The definition of a Count-Sketch matrix is, for any
and the column sparsity of matrix S is s. Independently in each column s positions are chosen uniformly at random without replacement, and each chosen position is set to
with probability
with probability 1/2.) Sparse OSEs can benefit from the sparsity of A, allowing for a running time of
nnz(A) denotes the number of non-zeros in A.
When n is large, the latter two algorithms are substantially faster than the naïve
1.1 Our Contributions
Despite the success of using subspace embeddings to speed up regression, often what practitioners are interested is not in preserving the cost of the regression problem, but rather in the generalization or prediction error provided by the vector . Ideally, we would like for any future (unseen) example
with high probability.
Ultimately one may want to use to do classification, such as regularized least squares classi-fication (RLSC) [RYP03], which has been found in cases to do as well as support vector machines but is much simpler [ZP04]. In this application, given a training set of examples with multiple (non-binary) labels identified with the rows of an
, one creates an
each column indicating the presence or absence of one of the r possible labels in each example. One then solves the multiple response regression problem
, and uses X to classify future examples. A commonly used method is for a future example a, to compute
where
are the columns of X. One then chooses the label
is maximum.
For this to work, we would like the inner products to be close to
is the solution to
is the solution to
-accurate OSE on d + r dimensional spaces [Sar06], which also satisfies so-called approximate matrix multiplication with error
, we get that
where is the spectral norm of
, which equals the reciprocal of the smallest singular value of A. To obtain a generalization error bound for an unseen example a, one has
which could be tight if given only the guarantee in (3). However, if the difference vector were distributed in a uniformly random direction subject to (3), then one would expect an
factor improvement in the bound. This is what our main theorem shows:
Theorem 1 (Main Theorem, informal)and matrix
are given. Let
be a subsampled randomized Hadamard transform matrix with
rows for an arbitrarily small constant
, and any fixed
with probability for an arbitrarily large constant C > 0. This implies that
with probability. If n > poly(d), then by first composing S with a Count-Sketch OSE with poly(d) rows, one can achieve the same guarantee.
(Here is a constant going to zero as n increases; see Theorem 10 for a formal statement of Theorem 1.)
Notice that Theorem 1 is considerably stronger than that of (4) provided by existing guarantees. Indeed, in order to achieve the guarantee (6) in Theorem 1, one would need to set existing OSEs, resulting in
rows. In contrast, we achieve only
rows. We can improve the bound in Theorem 1 to
is a matrix of i.i.d. Gaussians; however, as noted, computing
is slower in this case.
Note that Theorem 1 also makes no distributional assumptions on the data, and thus the data could be heavy-tailed or even adversarially corrupted. This implies that our bound is still useful when the rows of A are not sampled independently from a distribution with bounded variance.
The ) of Theorem 1 is achieved by applying (5) to the standard basis vectors
and applying a union bound. This
guarantee often has a more natural interpretation than the
guarantee—if we think of the regression as attributing the observable as a sum of various factors, (6) says that the contribution of each factor is estimated well. One may also see our contribution as giving a way for estimating the pseudoinverse
we get that
in the sense that each entry is within additive
a lot of work on computing entries of inverses of a matrix, see, e.g., [ADL
sparse (e.g. [Lee12]). In such cases, thresholding to the top k entries will yield an guarantee a factor
better than (3). One could ask if Theorem 1 also holds for sparse OSEs, such as the Count-Sketch. Surprisingly, we show that one cannot achieve the generalization error guarantee in Theorem 1 with high probability, say,
, using such embeddings, despite the fact that such embeddings do approximate the cost of the regression problem up to a
factor with high probability. This shows that the generalization error guarantee is achieved by some subspace embeddings but not all.
Theorem 2 (Not all subspace embeddings give the guarantee; informal version of Theorem 20). The Count-Sketch matrix with
rows and sparsity
—which is an OSE with exponentially small failure probability—with constant probability will have a result
that does not satisfy the
guarantee (6).
We can show that Theorem 1 holds for S based on the Count-Sketch OSE with
probability. We can thus compose the Count-Sketch OSE with the SRHT matrix and obtain an
time algorithm to compute
achieving (6). We can also compute
is a matrix of Gaussians, which is more efficient now that STA only has
rows; this will reduce the number of rows to
Another common method of dimensionality reduction for linear regression is leverage score sampling [DMIMW12, LMP13, PKB14, CMM15], which subsamples the rows of A by choosing each row with probability proportional to its “leverage scores”. With rows taken, the result
will satisfy the
) with probability
. However, it does not give a good
bound:
Theorem 3 (Leverage score sampling does not give the guarantee; informal version of Theorem 23)
small failure probability—with constant probability will have a result that does not satisfy the
guarantee (6).
Finally, we show that the rows that SRHT matrices use is roughly optimal:
Theorem 4 (Lower bounds for guarantees; informal versions of of Theorem 14 and Corollary 18). Any sketching matrix distribution over
matrices that satisfies either the
guarantee (3) or the
Notice that our result shows the necessity of the separation between the results originally defined in Equation (3) and (4) of Theorem 12 of [Sar06]. If we want to output some vector
that
, then it is known that
is necessary and sufficient. However, if we want to output a vector
, then we show that
is necessary and sufficient.
1.1.1 Comparison to Gradient Descent
While this work is primarily about sketching methods, one could instead apply iterative methods such as gradient descent, after appropriately preconditioning the matrix, see, e.g., [AMT10, ZF13, CW13]. That is, one can use an OSE with constant to construct a preconditioner for A and then run conjugate gradient using the preconditioner. This gives an overall dependence of
The main drawback of this approach is that one loses the ability to save on storage space or number of passes when A appears in a stream, or to save on communication or rounds when A is distributed. Given increasingly large data sets, such scenarios are now quite common, see, e.g., [CW09] for regression algorithms in the data stream model. In situations where the entries of A appear sequentially, for example, a row at a time, one does not need to store the full A but only the
Also, iterative methods can be less efficient when solving multiple response regression, where one wants to minimize . This is the case when
is constant and t is large, which can occur in some applications (though there are also other applications for which
is very small). For example, conjugate gradient with a preconditioner will take
time while using an OSE directly will take only
time (since one effectively replaces n with O (d) after computing
), separating t from d. Multiple response regression, arises, for example, in the RLSC application above.
1.1.2 Proof Techniques
Theorem 1. As noted in Theorem 2, there are some OSEs for which our generalization error bound does not hold. This hints that our analysis is non-standard and cannot use generic properties of OSEs as a black box. Indeed, in our analysis, we have to consider matrix products of the form for our random sketching matrix S and a fixed matrix U, where k is a positive integer. We stress that it is the same matrix S appearing multiple times in this expression, which considerably complicates the analysis, and does not allow us to appeal to standard results on approximate matrix product (see, e.g., [Woo14] for a survey). The key idea is to recursively reduce
using a property of S. We use properties that only hold for specifics OSEs S: first, that each column of S is unit vector; and second, that for all pairs
, the inner product between
is at most
with probability
Theorems 20 and 23. To show that Count-Sketch does not give the guarantee, we construct a matrix A and vector b as in Figure 1, which has optimal solution
with all coordinates
. We then show, for our setting of parameters, that there likely exists an index
ing the following property: the jth column of S has disjoint support from the kth column of S for all
except for a single k > d, for which
share exactly one common entry in their support. In such cases we can compute
explicitly, getting
. By choosing suitable parameters in our construction, this gives that
. The lower bound for leverage score sampling follows a similar construction.
Theorem 14 and Corollary 18. The lower bound proof for the guarantee uses Yao’s minimax principle. We are allowed to fix an
sketching matrix S and design a distribution over [A b]. We first write the sketching matrix
in its singular value decomposition (SVD). We choose the d + 1 columns of the adjoined matrix [A, b] to be random orthonormal vectors. Consider an
orthonormal matrix R which contains the columns of V as its first m columns, and is completed on its remaining
columns to an arbitrary orthonormal basis. Then
. Notice that
is equal in distribution to [A, b], since R is fixed and [A, b] is a random matrix with d + 1 orthonormal columns. Therefore,
is equal in distribution to
corresponds to the first m rows of an
uniformly random matrix with orthonormal columns.
A key idea is that if , then by a result of Jiang [J
submatrix of a random
orthonormal matrix has o(1) total variation distance to a
of i.i.d. N(0, 1/n) random variables, and so any events that would have occurred had G and h been independent i.i.d. Gaussians, occur with the same probability for our distribution up to an
factor, so we can assume G and h are independent i.i.d. Gaussians in the analysis.
The optimal solution in the sketch space equals
, and by using that SA has the form
, one can manipulate
to be of the form
, where the SVD of
We can upper bound
, since it is just the maximum singular value of a Gaussian matrix, where r is the rank of S, which allows us to lower bound
Then, since h is i.i.d. Gaussian, this quantity concentrates to
for a vector h of i.i.d. N(0, 1/n) random variables. Finally, we can lower bound
by the Pythagorean theorem, and now we have that
is the identity, and so this expression is just equal to the rank of
, which we prove is at least d. Noting that
for our instance, putting these bounds together gives
. The last ingredient is a way to ensure that the rank of S is at least d. Here we choose another distribution on inputs A and b for which it is trivial to show the rank of S is at least d with large probability. We require S be good on the mixture. Since S is fixed and good on the mixture, it is good for both distributions individually, which implies we can assume S has rank d in our analysis of the first distribution above.
1.2 Notation
For a positive integer, let [n] = {1, 2, . . . , n}. For a vector . For a matrix
to be the spectral norm of
to be the Frobenius norm of
to denote the Moore-Penrose pseudoinverse of
is its SVD (where
), is given by
In addition to notation, for two functions f, g, we use the shorthand
indicate that
) for an absolute constant
for constants c, C.
Figure 1: Our construction of A and b for the proof that Count-Sketch does not obey the guarantee.
Definition 5 (Subspace Embedding)-subspace embedding for the column space of an
is a matrix S for which for all
Definition 6 (Approximate Matrix Product)be a given approximation parameter. Given matrices A and B, where A and B each have n rows, the goal is to output a matrix C so that
. Typically C has the form
, for a random matrix S with a small number of rows. In particular, this guarantee holds for the subsampled randomized Hadamard transform
We first show that if S is a Gaussian random matrix, then it satisfies the generalization guarantee. This follows from the rotational invariance of the Gaussian distribution.
has full column rank. If the entries of
, then for any vectors
, we have, with probability
Because SA has full column rank with probability . Therefore
Thus it suffices to only consider vectors , or equivalently
. In such cases, SU will be independent of Sb, which will give the result. The proof is in Appendix A.
We first provide the definition of the subsampled randomized Hadamard transform(SRHT): let diagonal matrix with i.i.d. diagonal entries
, for which
in uniform on
. The matrix
is the Hadamard matrix of size
, and we assume n is a power of
coordinates of an n dimensional vector uniformly at random.
For other subspace embeddings, we no longer have that SU and Sb are independent. To analyze them, we start with a claim that allows us to relate the inverse of a matrix to a power series.
, and define
Suppose SA has linearly independent columns and
Proof.
We then bound the kth term of this sum:
be the subsampled randomized Hadamard transform, and let a be a unit vector. Then with probability
Hence, for , this is at most
with probability at least
We defer the proof of this lemma to the next section, and now show how the lemma lets us prove that SRHT matrices satisfy the generalization bound with high probability:
has full column rank with
subsampled randomized Hadamard transform with
vectors
with probability
Proof. Define
e.g., Theorem 2.4 of [Woo14] and references therein), so we condition on. Hence for
. In particular,
and we can apply Claim 8.
we may assume without loss of generality.
with probability. Suppose this event occurs, bounding the k = 0 term of (7). Hence it suffices to show that the
terms of (7) are bounded by
By approximate matrix product (Definition 6), we also have with
probability that
Combining with we have for any k that
Since this decays exponentially in k at a rate of , the sum of all terms greater than k is bounded by the kth term. As long as
On the other hand, by Lemma 9, increasing factor, we have for all k that
with probability at least , as long as
. Since the
can be expanded as a sum of
terms of this form, we get that
with probability at least , as long as
for a sufficiently large constant C. Combining with (9), the result holds as long as
for any
Combining Different Matrices. In some cases it can make sense to combine different matrices that satisfy the generalization bound.
be drawn from distributions of matrices that are
-approximate OSEs and satisfy the generalization bound (6). Then RS satisfies the generalization bound with a constant factor loss in failure probability and approximation factor.
We defer the details to Appendix B.
Proof. Each column of the subsampled randomized Hadamard transform has the same distribution as
is a random sign. It also has
with probability . See, e.g., [LDFU13]. By expanding the following product into a sum, and rearranging terms, we obtain
Note that is independent of
. We observe that in the above expression if
, then the sum over these indices equals
in this case for all c. Moreover, the sum over all indices conditioned on
equal to 0. Indeed, in this case, the expression can be factored into the form
random variable
Let W be a matrix with . We need Khintchine’s inequality:
Fact 12 (Khintchine’s Inequality)be i.i.d. sign random variables, and let
be real numbers. Then there are constants
We note that Khintchine’s inequality sometimes refers to bounds on the moment of though the above inequality follows readily by applying a Markov bound to the high moments.
We apply Fact 12 to each column of W, so that if -th column, we have by a union bound that with probability
simultaneously for all columns i. It follows that with the same probability,
We condition on this event in the remainder.
Thus, it remains to bound . By squaring
and using that
otherwise, we have,
for any Having computed the expectation of
, we now would like to show concentration. Consider a specific
Thus, we can apply Khintchine’s inequality recursively over all the from which it follows that with probability
, for each such
). We thus have with this probability, that
completing the proof.
We prove a lower bound for the guarantee, which immediately implies a lower bound for the
guarantee.
Definition 13. Given a matrix We say that an algorithm A(A, b, S) that outputs a vector
“succeeds” if the following property holds:
is a distribution over
with the property that for any
Proof. The proof uses Yao’s minimax principle. Let D be an arbitrary distribution over then E
averaging argument implies the existence of a fixed matrix
Thus, we must construct a distribution
cannot hold for any which does not satisfy
The proof can be split into three parts. First, we prove a useful property. Second, we prove a lower bound for the case
. Third, we show why
is necessary.
(I) We show that [SA, Sb] are independent Gaussian, if both [A, b] and S are orthonormal matrices. We can rewrite SA in the following sense,
where S is the complement of the orthonormal basis identity matrix, and
the left
submatrix of
. Thus, using [J
] as long as
(because of
the total variation distance between [SA, Sb] and a random Gaussian matrix is small, i.e.,
where each entry of H is i.i.d. Gaussian N(0, 1/n). (II) Here we prove the theorem in the case when (we will prove this is necessary in part III. Writing
in its SVD, we have
where . By a similar argument in Equation (12), as long as
we have that G
also can be approximated by a Gaussian matrix, where each entry is sampled from i.i.d. N(0, 1/n).
Similarly, also can be approximated by a Gaussian matrix, where each entry is
sampled from i.i.d. N(0, 1/n).
again, we have
where the the first equality follows by Equation (13), the second equality follows by the SVD of G, the third and fifth equality follow by properties of the pseudo-inversehas orthonormal rows and
is a diagonal matrix, and the fourth equality follows since
is an orthonormal basis.
Because each entry of is sampled from an i.i.d. Gaussian N(0, 1), using the result of [Ver10] we can give an upper bound for the maximum singular value of
with probability at least .99. Thus,
Because h is a random Gaussian vector which is independent of
where each entry of h is sampled from i.i.d. Gaussian N(0, 1/n). Then, using the Pythagorean Theorem,
from ; with probability 1/2, A, b are sampled from
In distribution
is a random
orthonormal basis and d is always orthogonal to A. In distribution identity matrix in the top-d rows and 0s elsewhere, while b is a random unit vector. Then, for any (A, b) sampled from
needs to work with probability at least 9/10. Also for any (A, b) sampled from
S needs to work with probability at least 9/10. The latter two statements follow since overall S succeeds on
with probability at least 19/20. Consider the case where A, b are sampled from distribution
Then consider which is the optimal solution to
where S can be decomposed into two matrices Plugging
into the original regression problem,
, which is at most
is a submatrix of S, the rank of S is also d.
It remains to define several tools which are used in the main proof of the lower bound.
Claim 15. For any matrix , if each entry of a vector
is chosen from an i.i.d Gaussian
Proof.
Let random variables. The random variables
t degree of freedom. Furthermore, the following tail bounds are known.
Fact 16 (Lemma 1 of [LM00])random variables. Then for any
and
Definition 17. Given a matrix We say that an algorithm B(A, b, S) that outputs a vector
“succeeds” if the following property holds:
is a distribution over
with the property that for any
[ADLPatrick Amestoy, Iain S. Duff, Jean-Yves L’Excellent, Yves Robert, François-Henry Rouet, and Bora Uçar. On computing inverse entries of a sparse matrix in an out-of-core environment. SIAM J. Scientific Computing, 34(4), 2012.
[AMT10] Haim Avron, Petar Maymounkov, and Sivan Toledo. Blendenpik: Supercharging lapack’s least-squares solver. SIAM J. Scientific Computing, 32(3):1217–1236, 2010.
[BDN15] Jean Bourgain, Sjoerd Dirksen, and Jelani Nelson. Toward a unified theory of sparse dimensionality reduction in euclidean space. In STOC, 2015.
[CMM15] Michael B Cohen, Cameron Musco, and Christopher Musco. Ridge leverage scores for low-rank approximation. arXiv preprint arXiv:1511.07263, 2015.
[CNW16] Michael B Cohen, Jelani Nelson, and David P. Woodruff. Optimal approximate matrix product in terms of stable rank. In Proceedings of the 43rd International Colloquium on Automata, Languages and Programming (ICALP 2016), Rome, Italy, July 12-15, arXiv preprint arXiv:1507.02268, 2016.
[Coh16] Michael B. Cohen. Nearly tight oblivious subspace embeddings by trace inequalities. In Proceedings of the Twenty-Seventh Annual ACM-SIAM Symposium on Discrete Algorithms, SODA 2016, Arlington, VA, USA, January 10-12, 2016, pages 278–287, 2016.
[CW90] Don Coppersmith and Shmuel Winograd. Matrix multiplication via arithmetic progressions. J. Symb. Comput., 9(3):251–280, 1990.
[CW09] Kenneth L. Clarkson and David P. Woodruff. Numerical linear algebra in the streaming model. In Proceedings of the forty-first annual ACM symposium on Theory of computing, pages 205–214. ACM, 2009.
[CW13] Kenneth L Clarkson and David P. Woodruff. Low rank approximation and regression in input sparsity time. In Proceedings of the forty-fifth annual ACM symposium on Theory of computing, pages 81–90. ACM, 2013.
[DMIMW12] Petros Drineas, Malik Magdon-Ismail, Michael W. Mahoney, and David P. Woodruff. Fast approximation of matrix coherence and statistical leverage. Journal of Machine Learning Research, 13:3475–3506, 2012.
[DMMS11] Petros Drineas, Michael W Mahoney, S Muthukrishnan, and Tamás Sarlós. Faster least squares approximation. Numerische Mathematik, 117(2):219–249, 2011.
[Gal14] François Le Gall. Powers of tensors and fast matrix multiplication. In International Symposium on Symbolic and Algebraic Computation, ISSAC ’14, Kobe, Japan, July 23-25, 2014, pages 296–303, 2014.
[JTiefeng Jiang et al. How many entries of a typical orthogonal matrix can be approximated by independent normals? The Annals of Probability, 34(4):1497–1529, 2006.
[LAKD08] Song Li, Shaikh S. Ahmed, Gerhard Klimeck, and Eric Darve. Computing entries of the inverse of a sparse matrix using the FIND algorithm. J. Comput. Physics, 227(22):9408–9427, 2008.
[LDFU13] Yichao Lu, Paramveer Dhillon, Dean Foster, and Lyle Ungar. Faster ridge regression via the subsampled randomized hadamard transform. In Proceedings of the Neural Information Processing Systems (NIPS) Conference, 2013.
[Lee12] Jeff Leek. Prediction: the lasso vs. just using the top 10 predictors. http://simplystatistics.org/2012/02/23/ prediction-the-lasso-vs-just-using-the-top-10/, 2012.
[LM00] Beatrice Laurent and Pascal Massart. Adaptive estimation of a quadratic functional by model selection. Annals of Statistics, pages 1302–1338, 2000.
[LMP13] Mu Li, Gary L Miller, and Richard Peng. Iterative row sampling. In Foundations of Computer Science (FOCS), 2013 IEEE 54th Annual Symposium on, pages 127–136. IEEE, 2013.
[MM13] Xiangrui Meng and Michael W Mahoney. Low-distortion subspace embeddings in input-sparsity time and applications to robust linear regression. In Proceedings of the forty-fifth annual ACM symposium on Theory of computing, pages 91–100. ACM, 2013.
[NN13] Jelani Nelson and Huy L Nguyên. Osnap: Faster numerical linear algebra algorithms via sparser subspace embeddings. In Foundations of Computer Science (FOCS), 2013 IEEE 54th Annual Symposium on, pages 117–126. IEEE, 2013.
[NN14] Jelani Nelson and Huy L. Nguyên. Lower bounds for oblivious subspace embeddings. In Automata, Languages, and Programming - 41st International Colloquium, ICALP 2014, Copenhagen, Denmark, July 8-11, 2014, Proceedings, Part I, pages 883–894, 2014.
[PKB14] Dimitris Papailiopoulos, Anastasios Kyrillidis, and Christos Boutsidis. Provable deterministic leverage score sampling. In Proceedings of the 20th ACM SIGKDD international conference on Knowledge discovery and data mining, pages 997–1006. ACM, 2014.
[RYP03] Ryan Rifkin, Gene Yeo, and Tomaso Poggio. Regularized least-squares classification. Nato Science Series Sub Series III Computer and Systems Sciences, 190:131–154, 2003.
[Sar06] Tamás Sarlós. Improved approximation algorithms for large matrices via random projections. In Symposium on, pages 143–152. IEEE, 2006.
[Ver10] Roman Vershynin. Introduction to the non-asymptotic analysis of random matrices. arXiv preprint arXiv:1011.3027, 2010.
[Wil12] Virginia Vassilevska Williams. Multiplying matrices faster than coppersmithwinograd. In Proceedings of the 44th Symposium on Theory of Computing Conference, STOC 2012, New York, NY, USA, May 19 - 22, 2012, pages 887–898, 2012.
[Woo14] David P. Woodruff. Sketching as a tool for numerical linear algebra. Foundations and Trends in Theoretical Computer Science, 10(1-2):1–157, 2014.
[ZF13] Anastasios Zouzias and Nikolaos M. Freris. Randomized extended kaczmarz for solving least squares. SIAM J. Matrix Analysis Applications, 34(2):773–793, 2013.
[ZP04] Peng Zhang and Jing Peng. Svm vs regularized least squares classification. In ICPR (1), pages 176–179, 2004.
Lemma 19. If the entries of
for any vectors a, b with probability
Proof. With probability 1, the matrix SA has linearly independent columns, and so
It is well-known (stated, for example, explicitly in Theorem 2.3 of [Woo14]) that with probability , the singular values of
. We condition on this event. It follows that
where the first equality uses that V is a rotation, the first inequality follows by sub-multiplicativity, and the second inequality uses that the singular values of SU are in the range with probability
The main observation is that since is statistically independent from Sb. Hence, Sb is distributed as
, conditioned on the vector
. It follows that conditioned on the value of
is distributed as
and so using (14) , with probability
In some cases it can make sense to combine different matrices that satisfy the generalization bound.
be drawn from distributions of matrices that are
-approximate OSEs and satisfy the generalization bound (6). Then RS satisfies the generalization bound with a constant factor loss in failure probability and approximation factor.
Proof. For any vectors we want to show
As before, it suffices to consider the case. We have with probability
suppose this happens. We also have by the properties of R, applied to SA and Sb, that
Because S is an OSE, we have . Therefore
We describe a few of the applications of combining sketches.
B.1 Removing dependence on n via Count-Sketch
One of the limitations of the previous section is that the choice of k depends on n. To prove that theorem, we have to assume that log d > log log n. Here, we show an approach to remove that assumption.
The main idea is instead of applying matrix directly, we pick two matrices
is FastJL matrix and C is Count-Sketch matrix with s = 1. We first compute
, then compute
. The benefit of these operations is S only needs to multiply with a matrix (CA) that has poly(d) rows, thus the assumption we need is log d > log log(poly(d)) which is always true. The reason for choosing C as a Count-Sketch matrix with
(2) The running time is
B.2 Combining Gaussians and SRHT
By combining Gaussians with SRHT matrices, we can embed into the optimal dimension with fast
embedding time.
B.3 Combining all three
By taking Gaussians times SRHT times Count-Sketch, we can embed into the optimal dimension embedding time.
Here we demonstrate an A and a b such that Count-Sketch will not satisfy the guarantee with constant probability, so such matrices cannot satisfy the generalization guarantee (6) with high probability.
be drawn as a Count-Sketch matrix with s nonzeros per column. There exists a matrix
such that, if
, then the “true” solution
and the approximation
have large
distance with constant probability:
Plugging in we find that
even though such a matrix is an OSE with probability exponential in s. Therefore there exists a constant c for which this matrix does not satisfy the generalization guarantee (6) with probability.
Proof. We choose the matrix A to be the identity on its top
set the value of the first d coordinates of vector and set the value to be
for the next
coordinates, with the remaining entries all zero. Note that
and
Let denote the kth column vector of matrix
. We define two events, Event I,
; Event II,
such that
, and all other
. Using Claim 21, with probability at least .99 there exists a k for which both events hold.
Given the constructions of A and b described early, it is obvious that
Conditioned on event I and II are holding, then denote . Consider the terms involving
in the quadratic form
it can be written as . Hence the optimal
, which is different from the desired
. Plugging in our requirement of
where the last inequality follows by . Thus, we get the result.
, with probability at least .99 there exists a
for which both event I and II hold.
Figure 2: Count Sketch matrix Event I, for any
. Event II, there exists a unique
intersect at exactly one location(row index).
, then for any
be an indicator that the entries of column i are disjoint from all
, so by Markov’s inequality, with probability .99, we have .99d columns having this property (indeed, the expected value of
is at most
). Define Event E to be that .99d columns of first d columns have the property that the entries of that column are disjoint from all the other
columns. Let S be the set of these .99d columns. Let N be the union of supports of columns in S.
Each column non-zero entries. Define event F( which is similar as event E) to be that
columns of the next
columns have the property that the entries of that column are disjoint from all the other
columns. By the same argument, since
with probability .99, we have
columns in
being disjoint from other columns in
. Condition on event F holding. Let L be the multiset union of supports of all columns in
be the union of supports of all columns in
, that is, the set union rather than the multiset union. Note that
because of
columns are disjoint from each other.
The intersection size x of N and M is hyper-geometrically distributed with expectation
By a lower tail bound for the hypergeometric distribution
where
where the last inequality follows by setting . Thus, we get with probability .99, the intersection size is at least
Now let W be the distinct elements in L\M, so necessarily . By an upper tail bound for the hypergeometric distribution, the intersection size y of W and N satisfies
where , we again get
If |W| = 0, then y = 0. Otherwise, we can set so that this probability is less than .01, and we get with probability .99, the intersection size y is at most
Note that we have that
are bounded by
suffices to ensure y is at most
, and earlier that x is at least
The probability one of the |S| blocks in N has two or more intersections with M is less than
times the probability two random distinct items in the intersection land in the block. This
probability is
So the expected number of such blocks is which is less than
, which we have. So, there are at least x/2 blocks which have intersection size exactly 1 with N. Note that the number of intersections of the |S| blocks with W is at most y, which is at most
, and therefore there exists a block, that is, a column among the first d columns, which intersects M in exactly one position and does not intersect W. This is our desired column. Thus, we complete the proof.
Not only does Count-Sketch fail, but so does leverage score sampling, which is a technique that takes a subsample of rows of A with rescaling. In this section we show an A and a b such that leverage score sampling will not satisfy the guarantee. We start with a formal definition of leverage scores.
Definition 22 (Leverage Scores). Given an arbitrary denote the
matrix consisting of the d left singular vectors of
denote the
th row of the matrix
is a row vector. Then the leverage scores of the rows of A are given by
The leverage score sampling matrix can be thought of as a square diagonal matrix with diagonal entries chosen from some distribution. If
, it means we do not choose the i-th row of matrix
, it means we choose that row of the matrix A and also rescale that row. We show that the leverage score sampling matrix cannot achieve
guarantee, nor can it achieve our notion of generalization error.
be a leverage score sampling matrix with m nonzeros on the diagonal. There exists a matrix
and a vector
such that, if
, then the “true” solution
and the approximation
have large
distance with constant probability:
Therefore there exists a constant c for which this matrix does not satisfy the generalization guarantee (6) with probability.
Proof. We choose the matrix A to be the identity on its top d rows, and L scaled identity matrices
for the next dL rows, where
(to normalize each column of A), which implies
. Choose some
. Set the value of the first d coordinates of vector b to be
and set the value to be
for the next
coordinates, with the remaining entries all zero. Note that
the following term,
where the optimal should be set to 1/d. The other involves the term:
where the optimal should be set to
. Because we are able to choose
Second, we compute . With high probability, there exists a j satisfying Equation (16), but after applying leverage score sampling, the middle term of Equation (16) is removed. Let
denote the leverage score of each of the top
denote the leverage score of each of the next Ld rows of A. We need to discuss the cases
separately. If m > d, then the following term involves
to be
Third, we need to compute . It is easy to see that
is an orthonormal matrix. The upper bound for
, and the lower bound is also a constant, which can be proved in the following way:
Before getting into the proof details, we define the key property of S being used in the rest of the proofs.
Definition 24 (All Inner Product Small(AIPS) Property)
we say that S satisfies the “AIPS” property.
is a subsampled Hadamard transform matrix, then the AIPS property holds with probability at least
Proof. From the structure of , we have with probability
. Applying a union bound over
pairs, we obtain that
The main idea for bounding is to rewrite it as
AIPS holds
AIPS does not hold ]. Because Pr[ AIPS does not hold] is at most 1/poly(n), the first term dominates the second term, which means we only need to pay attention to the first term. We repeatedly apply this idea until all the S are removed.
We start by boundinbg by squaring
and using that
otherwise. Then, we obtain,
where the first equality is from our conditioning, and the second equality is the definition of Hence,
where the first equality follows from the property just mentioned, and the inequality follows by including back the tuples of indices for which , using that each summand is non-negative.
Our next step will be to bound the term
where the first equality is the definition of B, the second equality follows by separating out the index , the third equality uses that
, that is, the squared norm of the
, the inequality follows since all rows of a projection matrix
norm at most 1, and the final equality uses that
no longer appears in the expression.
We now merge our bounds for the terms A and B in the following way:
where the first two inequalities and first equality are by definition of A and B above. The first inequality follows by induction, since at this point we have replaced , and can repeat the argument, incurring another multiplicative factor of O(d(log n)/r)+1. Repeating the induction in this way we arrive at the last inequality. Finally, the last inequality follows by plugging in the definition of
, using that
where the inequality follows since each row of has norm at most 1, and a is a unit vector. The final result is that