Distributed Least Squares in Small Space via Sketching and Bias Reduction
Authors: Sachin Garg
Kevin Tan
Michal Derezinski
Venue: NeurIPS 2024 poster (2024)
Keywords: Matrix Sketching
Least squares
Randomized Linear Algebra
Random Matrix Theory
Primary area: learning_theory
TL;DR: We give a sparse sketching method running in optimal space and current matrix multiplication time, recovering a nearly-unbiased least squares estimator using two passes over the data.

Abstract: Matrix sketching is a powerful tool for reducing the size of large data matrices. Yet there are fundamental limitations to this size reduction when we want to recover an accurate estimator for a task such as least square regression. We show that these limitations can be circumvented in the distributed setting by designing sketching methods that minimize the bias of the estimator, rather than its error. In particular, we give a sparse sketching method running in optimal space and current matrix multiplication time, which recovers a nearly-unbiased least squares estimator using two passes over the data. This leads to new communication-efficient distributed averaging algorithms for least squares and related tasks, which directly improve on several prior approaches. Our key novelty is a new bias analysis for sketched least squares, giving a sharp characterization of its dependence on the sketch sparsity. The techniques include new higher moment restricted Bai-Silverstein inequalities, which are of independent interest to the non-asymptotic analysis of deterministic equivalents for random matrices that arise from sketching.

Distributed Least Squares in Small Space via Sketching and Bias Reduction
Sachin Garg
Kevin Tan Department of Statistics University of Pennsylvania kevtan@uni ch. edu
Computer Science & Engineering University of Michigan sachg@umich. ch. edu
Michal Dereziński
Computer Science & Engineering
University of Michigan derezinQumichi edu
Abstract
Matrix sketching is a powerful tool for reducing the size of large data matrices. Yet there are fundamental limitations this size reduction when we want recover an accurate estimator for a task such as least square regression. We show that these limitations can be circumvented in the distributed setting by designing sketching methods that minimize the bias of the estimator, rather than its error. In particular, we give a sparse sketching method running in optimal space and current matrix multiplication time, which recovers nearly-unbiased least squares estimator using two passes over the data. This leads to new communication--fficee averaging algorithms for least squares and related tasks, which directly improve several prior approaches. Our key novelty is a new bias analysis for sketched least squares, giving sharp characterization of its dependence on the sketch sparsity. The techniques include new higher-moment restricted Sai-Silverstein inequalities, which of independent interest to the non-asymptotic analysis of deterministic equivalents for random matrices that arise from sketching.
distributed
Introduction
Matrix sketching is a powerful collection of randomized techniques for compressing large data matrices, developed over long line of works as part of Randomized Numerical Linear Algebra [RandNLA, c.g., 45, 26, 37, 39, 21]. Sketching can be used to reduce the large dimension of data matrix A Ra by applying a random sketching matrix (operator)
m to obtain the sketch A SA mxd where m < n. For example, sketching can be used to approximate the solution to the least squares problem, x" = argmin, (xx where L(x) ||Ax - bll², by using sketched estimator x argming [xx blla. where A SA and b Sb.
Perhaps the simplest form of sketching is subsampling, where the sketching operator S selects random sample of the rows of matrix A. However, the real advantage of sketching as framework emerges as we consider more complex operators S. such as sub-Gaussian matrices [1]. randomized Hadamard transforms [3], and sparse random matrices [12]. These approaches ensure higher quality and more robust compression of the data matrix, e-g., leading to provable e-approximation guarantees for the estimate X in the least squares task, i.e., [(x) < (1 + <x". Nevertheless, there are fundamental limitations to how far we can compress data matrix using sketching while ensuring an c-approximation. These limitations pose challenge particularly in space- limited computing
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

Subsample
(d/e)) x
Sketch
((d) d
((1/e) rows
((1/e)
(1/e)) rows
Figure 1: Illustration of the leverage score sparsification algorithm used in Theorem 1. Each row of the sketch mixes ((1/e) leverage score samples from A. Remarkably, the e-error guarantee of the subsampled estimator is retained e-bias of the sketched estimator.
environments, such as for streaming algorithms where observe the matrix A, say, one at time, and we have limited space for storing the sketch [11].
One strategy for overcoming the fundamental limitations of sketching as a compression tool is to look beyond the single approximation guarantee provided by a sketching-based estimator x, and consider how its broader statistical properties can be leveraged in a given computing environment. To that end, many recent works have demonstrated both theoretically and empirically that sketching- estimators often exhibit not only approximation robustness but also statistical robustness, for instance enjoying sharp confidence intervals, effectiveness of statistical inference tools such as bootstrap and cross- validatio,
well accuracy boosting techniques such as distributed averaging [e.g.. 35, 25, 32, 33). Yet, these results have had limited impact on the traditional computational complexity analysis in RandNLA and sketching literature, as many of them either impose additional assumptions, focus on sharpening the constant factors, or require using more expensive sketching techniques. In this work, demonstrate that the statistical properties of sketching-based estimators can in fact have a substantial impact on the computational trade-offs that arise in RandNLA.
Our key motivating example is the above mentioned least squares regression task. It is well understood that for an n x d least squares task, to recover an e-approximate solution out of an m X sketch, we need sketch size at least m ((d/e). e). This has been formalized in the streaming setting with lower bound of (2 -d22 d2 log(nd)) bits of space required, when all of the input numbers use (log(nd) bits of precision [11]. One setting where this be circumvented is in the distributed computing model where the bits can be spread out across many machines, so that the per-machine space can be smaller. Here, one could for instance hope that
maintain small (dd) d sketches in O(1/e) machines and then combine their estimates to recover an e-approximate solution. A simple and attractive approach to average the estimates x; produced by the individual machines, returning = 9 E X; as this only requires each machine to communicate (d (nd) bits of information about its sketch. This approach requires the sketching-based estimates x; to have sufficiently small bias for the averaging scheme to be effective. While this has been demonstrated empirically in many cases, existing theoretical results still require relatively expensive sketching methods to recover low-bias estimators, leading an unfortunate trade-off in the distributed averaging scheme between the time and space complexity required.
In this work, address the time-space trade-off in distributed averaging of sketching- based estima- tors, by giving a sharp characterization of how their bias depends on the sparsity of the sketching matrix. Remarkably, data down to the minimum size of Odd log(nd)) bits at no extra computational cost, while still being able to recover an e-approximate solution for least squares and related problems. Importantly, our results require the sketching matrix to be slightly denser than is necessary for obtaining approximation guarantees on single estimate, and thus, cannot be recovered by standard RandNLA sampling methods such as approximate leverage score sampling [27].
show that in the distributed streaming environment one can compress the
Before state the full result in the distributed setting, we give main technical contribution. which is the following efficient construction of a low-bias least squares estimator in single pass using only (d (n)) bits of space, assuming all numbers use Ollog(nd) bits of precision. Below, 7 > 0 denotes an arbitrarily small constant.

---

Reference
Method
Total runtime
Parallel passes
Folklore
Gaussian sketch
O(nd
[43]
Leverage Score Sampling
nnz( A) + O vvo)
1
[5]
Determinantal Point Process
nnz (B) + (dd))
log"(n// 1/C
[9]
Weighted Mb-SGD (sequential)
nnz A) + O(d2/e
This work (Thm. 1)
Leverage Score Sparsification
nnz A) + O(d2/e)
Table 1: Comparison of time complexities and parallel passes over the data required for differ-
ent methods to obtain a (1 <)-approximation in Odd log(nd)) bits of space for an x d least squares problem (A, b), given a preconditioner such that k(AP) = O(1) (see Section 3 for our computational model). We include the fully sequential Weighted Mb-SGD reference.
Theorem 1. Given streaming access to A E Raxd and b E R" and direct access to pre- conditioner matrix P Idddd such that s(AP) < a, within a single pass over (A, b), in O(7-"nnz(A) E - polylooo) time and O(d log(nd)) bits of space, we can construct randomized estimator x for the least squares solution x" argmin ||Ax blI such that:
(Bias)
||A - -bll S (1 + o) Axx - -bl?,
(Variance)
B|IAR b1l S |Ax - bll?.
The above construction assumes access to preconditioner matrix P with s(AP) S (where is denotes the condition number). Such matrix can be obtained efficiently with a = 0(1) ) separate single pass, leading to two-pass algorithm described later in Theorem 2.
Remark 1.
Our estimator & is constructed at the end of the data pass from a sketch (A, b) where À = SA and
b = Sb, by minimizing | b1|² using preconditioned conjugate gradient. Here, S is carefully constructed sparse sketching matrix which is inspired by the so-called leverage score sparsified (LESS) embeddings [18]. Leverage scores represent the relative importances of the rows of A which are commonly used for subsampling in least squares (see Definition 1), and their estimates can be easily obtained in single pass by using the preconditioner matrix P.
Our time complexity bound of Ó(nnz(A) + dz (f/c) matches the time it would take (for single machine) to subsample Ö(d/c) rows of A according to the approximate leverage scores and produce an estimator that achieves the e-crror bound ||Ax - bil² I (1 + e) - bll2. However, this strategy requires either maintaining (dd bits of space for the sketch, or computing x directly along the way, blowing up the runtime to (dd² /e). Since approximate leverage score sampling leads to significant least squares bias, averaging can only improve this to 0(d² /E) (see Table 1). An alternate strategy would be to combine leverage score sampling with a preconditioned mini-batch stochastic gradient descent (Weighted Mb-SGD), with mini-batches chosen so that they fit Ódd) space. This achieves the same time and space complexity as our method, but due to the streaming access A and the sequential nature of SGD, requires O(1/e) data passes.
Instead, our algorithm essentially mixes an Ö(d/e) size leverage score sample into an O(d) size sketch, merging Ö(L/e) rows of A into a single row of the sketch (see Figure 1). This results in better data compression compared to direct leverage score sampling, with only (dd) bits of space, while retaining the same Ö(nnz(A) A) + d/ec) runtime complexity as the above approaches and requiring only single data pass. The resulting estimator & can no longer recover the c-error bound, but remarkably, its expectation E[ still does. To turn this into an improved estimator in a distributed model, we can simply average
1/c such estimators, i.c., & = , D x, obtaining Ax - b1l2 < (1 + 2e) Axx - 2 As shown in Table 1, ours is the first result in this model to achieve Ō(d)) space in a single pass and faster than current matrix multiplication time (dd").
Finally, by incorporating a preconditioning scheme, we illustrate how our construction can be used to design the first algorithm that solves least squares in current matrix multiplication time, constant parallel passes and Odd log(nd)) bits of space. We note that the O(d²) cost comes only from the worst-case complexity of constructing the preconditioner P. which can often be accelerated in practice. The computational model used in Theorem 2 is described in detail in Section 3.
3

---

Theorem 2. Given A E Raxd and b E R in the parallel computing model. using two parallel passes with q machines, we can compute such that with probability 0.9
-Ax-bli bl1 S (1+++00//0))
-b]
in ' A) + d' + + bits of communication. In particular, choosing q 1/e we recover an -appoppprxinnation.
time, O log(nd)) bits of space and O(dlog(nd)
Remark 2. While the parallel computing model in Theorem 2 assunes that all machines have streaming access to the entire data matrix, this result can be easily extended to the setting where A has been randomly down-sampled or partitioned into separate size chunks A A,,... Ag. and each machine constructs an estimate &; based on a sketch of its own chunk. Then, with the same computational guarantees as Theorem 2, the averaged estimator & - X; with probability 0.9 enjoys a guarantee of:
A - 51 (1+ 1/9 Bebunk)) )1Axx) b[1
where Behunk is the bias that would be incurred if we solved each chunk exactly and averaged those solutions. Using existing guarantees for uniform down- down-sampling of least squares [Theorem 20, 43), one can bound this bias as Behunk = ( + (H)),, where p is the coherence of the data matrix and n is the chunk size. For sufficiently large chunks, this bias is negligible compared to the error e. We prove the above high probability guarantee Theorem in Appendix E.
Further applications. Our least squares analysis can be extended to other settings where prior works [e.g.. 17, 18, 16. 22] have analyzed randomized estimators based on sparse sketching via techniques from asymptotic random matrix theory. The primary and most direct application involves correcting inversion bias in the so-called sketched inverse covariance estimate (À "A)'', which was the motivating task of [18]. with applications including distributed second- optimization and statistical uncertainty quantification, where quantities like (À A) 'x are approximated.
Theorem 3 (informal Theorem 5). Given A xd and its LESS embedding S with sketch size m Cd and non-zeros per row, the inverse covariance sketch ( (m A S SA) 1 is (c, 6)- unbiased estimator of (A A) (see Definition 3) for c Ō(1 vd/s) vd and 5 1/poly(d). m
This result should be compared with c Ö( d/s) vi obtained by [18]. Thus, we get direct improvement for very sparse sketches, i.e., s= o(d). This can be immediately translated into an improved local convergence guarantee for Distributed Newton Sketch which a second-order convex minimization algorithm used in settings where the Hessian matrix can be expressed as A "A for tall matrix A, c.g., in generalized lincar models like logistic regression. In this method, following Corollary 16 of [18] as well related results [16, 44, 19], we use sketching and averaging to estimate a Newton step:
X
1 q
Σ
H
g
where g: is the gradient at x, and H, H., are Hessian sketches constructed by independent machines. The following corollary, which is a direct improvement over Corollary 16 of [18], shows that Distributed Newton Sketch on generalized linear model task can achieve fast local convergence rate of the form f(x==1)-((")) o(1) with (dd)) time per iteration (see Appendix B.1 for details). [xx)
Corollary (Distributed Newton Sketch). Consider fxx) = E" 4(xx Φ.) + |xx]] where l convex twice continuously differentiable functions, such that f has a Lipschitz Hessian, À >> O. and o; is the ith row of an n. x d data matrix &. Given >> O, there is a neighborhood U. around the minimizer x* = argming [x) such that, for any Xs E Uc, using two parallel passes with Õ(1/e) machines, we can compute Distributed Newton Sketch update Xe+1 such that
fxx(11) f(x*) S e [f(x)) - f(x")].
4) + + 2 + -1 polyloog) time, Odd log(nd)) bits of space and (dlo((nd)
bits of communication.

---

Our Techniques. At the core of our analysis are techniques inspired by asymptotic random matrix theory (RMT) in the proportional limit [e.g,, see 6]. Here, in order to establish the limiting spectral distribution (such as the Marchenko-Pastiee law) of a random matrix À Á whose dimensions diverge to infinity, one aims to show the convergence of the Stieltjes transform of its resolvent matrix (A*A - 2I)=1 Recently, [18] showed that these techniques can be adapted to sparse sketching matrices (via leverage score sparsification) in order to characterize the bias of the sketched inverse covariance (À A)- where À SA.
Our main contribution is two-fold. First, we show that similar argument can also be applied to analyze the bias of the least squares estimator, x = (À Ā) A "b Unlike the inverse covariance, this estimator no longer takes the form of resolvent matrix, but its bias is also associated with the inverse, which means that we can use a leave-one-out argument to characterize the effect of removing a single row of the sketch on the estimation bias. Our second main contribution is to improve the sharpness of the bounds relative to the sparsity of the sketching matrix by combining a careful application of Hölder's inequality with a higher moments analysis of the restricted Bai- Silverstein inequality for quadratic forms. Those improvements are not only applicable to the least squares analysis, but also to all existing RMT-style results for LESS embeddings. including the aforementioned inverse covariance estimation, as well as applications in stochastic optimization, resulting in the sketching cost of LESS embeddings dropping below matrix multiplication time.
2 Related Work
Randomized numerical linear algebra. RandNLA sketching techniques have been developed over long line of works, starting from fast least squares approximations of [42]; for an overview, see [45, 26, 37, 39, 21] among others. Since then, these methods have been used in designing fast algorithms not only for least squares but also many other fundamental problems in numerical linear algebra and optimization including low-rank approximation [13, 34], 1p regression [14], solving linear systems [28, 24] and more. Using sparse random matrices for matrix sketching also has long history, including cata-oblivious sketching methods such as CountSketch [12]. OSNAP [40], and more [38]. Leverage score sparsification (LESS) was introduced by [18] data-dependent sparse sketching method to enable RMT-style analysis for sketching (see below).
Unbiased estimators for least squares. To put our results in a proper context, let us consider other approaches for producing near-unbiased estimators for least squares, see also Table First, a well known folklore result states that the least squares estimator computed from dense Gaussian sketching matrix is unbiased. The bias of other sketching methods, including leverage score sampling and OSNAP, has been studied by [43], showing that these methods need a Ve-error guarantee to achieve an e-bias which leads to little improvement unless e is extremely small and the sketch size is sufficiently large. Another approach of constructing unbiased estimators for least squares, first proposed by [23], is based on subsampling with non-i.i.d. importance sampling distribution based on Determinantal Point Processes [DPPs, 31, 20]. However, despite significant efforts [30, 7, 51, sampling from DPPs remains quite expensive: the fastest known algorithm requires running a Markov chain for polylog(n/e many steps, each of which requires separate data pass and takes ((d)) time. Other approaches have also been considered which provide partial bias reduction for i.i.d. RandNLA subsampling schemes in various regimes that are either much more expensive or not directly comparable to ours [2, 44].
Statistical and RMT analysis of sketching. Recently, there has been significant interest in statisti- cal and random matrix theory (RMT) analysis of matrix sketching: see [21] for an overview. These approaches include both asymptotic analysis via limiting spectral distributions and deterministic cquivalents [35, 25, 32, 33], and non-asymptotic analysis under statistical assumptions [36, 41, 4]. A number of works have shown that the RMT-style techniques based on deterministic equivalents can be made rigorously non non-asyme asymptotic for certain sketching methods such as dense sub-Gaussian [17], LESS matrices [16, 18], and other sparse matrices [9]. which has been applied to low-rank approximation, fast subspace embeddings and stochastic optimization. Our new analysis can be viewed as general strategy for directly improving the sparsity required by LESS embeddings (and thereby, the sketching time complexity) in many of these applications, specifically those that rely on analysis inspired by the calculus of deterministic equivalents via generalized Stieltjes transforms.
5

---

3 Preliminaries
Notations. In all our results, we use lowercase letters to denote scalars, lowercase boldface for vectors, and uppercase boldface for matrices. The norm || I| denotes the spectral norm for matrices and the Euclidean norm for vectors, whereas || || denotes the Frobenius norm for matrices. We use z to denote the p.s.d. ordering of matrices.
Computational model. We next clarify the computational model that is used in Theorem 2. We consider central data server storing (A, b), and q machines. The jth machine has handle Stream(j), which can be used to open a stream and to read the next row/label pair (a,, b.) in the stream. After full pass, the machine The machines can operate their streams entirely asynchronously. and each has its own limited local storage space, e.g., in Theorem 2 can communicate some information back to the server, e.g., in Theorem 2, they communicate their final estimate vectors Xj. using (dd (on))) bits of communication. Then, the server computes the final estimate, in our case via averaging. X Ei= X,, which can be done either directly or via map-reduce type architecture.
re-open the handle and begin another pass over the data.
Odd log(nd)) bits of space per machine. At the end, they
We define the parallel passes required by such an algorithm as the maximum number of times the stream is opened by any single machine. We analogously define fine/spacccommuiiationn costs by taking a maximum over the costs required by any single machine (for communication, this refers only to the number of bits sent from the machine back to the server).
Definitions and useful lemmas. In our framework, we construct a sparse sketching matrix S where sparsification is achieved using a probability distribution over rows of data matrix A, that is proportional to the leverage scores of A. The next definition [following, e.g., 9] provides the explicit definition of exact and approximate leverage scores for our setting.
Definition 1 ((31,,2)-apppe appprxmaae leverage scores). Fix matrix A E Baaxd and consider matrix U E Raxd with orthonormal columns spanning the column space of A. Then, the leverage scores 1,, sii n are defined as the row norms squared U,, i.e., ||u I2, where u, the ith row of U. Furthermore, consider fixed 31, B2 > 1. Then I are called (allaeddii -aprrxxixi leverage scores for A if the following holds for all i
4 < i, and 54532dd
3
i=1
The approximate leverage scores can computed by first constructing a preconditioner matrix P E Raxd such that x s(AP) O(1), which takes O(nnz(A) A) d=") in single pass, and then relying on the following norm approximation scheme.
Lemma (Based Lemma 7.2 from [10]). Given A E Ruad and P E Rdxd, using single pass over A in time ((n=z((nz(() A) d')) for small constant 0, compute estimates 11....i I.
such that with probability 2 0.95:
n AP < S Ollog(n))lle? API2 Vi and
In the next definition, give the sparse sketching strategy used in analysis. This approach is similar to the original leverage score sparsification proposed by [18], except: 1) we adapted so that it can be implemented effectively in single pass, and 2) use it in much sparser regime (fewer non-zeros per row).
Definition ((s, (1, , (2) B embedding). Fix matrix E Raxd and some s 2 0. Let the tuple 20,
(I,,... In) denote (31, approxxxxee leverage scores for A. Let p: min
d J. We define
a (s, B1, 32) Bap-approx leverage score sparsifier as follows.
{-(/ VPI
b1
by
E
where
by ~ Bernoull(py).
Pr
Moreover, we define the (s, B1, -eeeeraaee score sparsified (LESS) embedding of size m as matrix Idaxd with i.i.d. rows Tm (81, B2) capapoxixie leverage score sparsifier and y, 12" consist of random +1 entries.
txx, X; such that x; = diag(E.)y; where & denotes a randomly generated

---

A key property of a sketching matrix is the subspace embedding property, defined below. It was recently shown by [9] that LESS embeddings require only polylogarithmicaly many non-zeros per row of S to prove that S is a subspace embedding for the data matrix A with the optimal m = (d) sketching dimension. For our main results, it is sufficient to use oj = 0(1) below.
Lemma 2 (Subspace embedding for LESS, Theorem 1.3, [9]). Fix n, 8 > O. Consider 31, 32 1 and a full rank matrix A E R=xd Then for a (31, 32) -leverage score sparsified embedding S E Rm=x with s 2 (d/6)/" and m = ((dd + log 1/6)/n²), with probability 1 - 6 we have
1
= (1 ++1) ATA.
1+n
4 Least Squares Bias Analysis
In this section we provide an outline of the bias analysis for the sketched least squares estimator constructed using a LESS embedding, leading to the proofs of our main results, Theorems and 2. In particular, we prove the following main technical result (detailed proof Appendix C).
Theorem 4 (Bias of LESS- least squares). Fix A RRd and let $ be (s, 31, 3p (,,),,8)-LESS LLESS embedding of size m for A. Let S satisfy (1) with 17 = and probability - â where ó < m Then there exists an event & with probability at least - ô such that
d
( ) (n//)) L(x").
[(Eg(x)) L(xx) = 0
m2
Remark 3. Thus, the bias [(Eg(x)) L(xx) of the LESS estimator using ((3,325) (s) non-
zeros per row of S is of the order Ô 1ix + A -((x). By comparison, the standard expected loss bound which holds for sketched least squares (including this estimator) [[[[(x)] L(x*)= S
), and the best known bound on the bias of most standard sketched estimators (e.g..
leverage score sampling) is Õ d² m²) L(xx), given by [43]. So, our result recovers the standard bias bound for = and improves on for s >> 1 by a factor of min(s, (s, d). At the end of the section, we discuss how to deal with the lower order term Ö (,At)) ) to reduce the bias further.
Proof sketch. Using a standard argument, we can replace the matrix A with the matrix U E R2ax consisting of orthonormal columns spanning the column space of A, and assume that = poly(d) Let S be an (s, 31, (,, L-ESS embedding for U. Also, let b E I be a vector of responses/labels corresponding to n rows in U. Let x = argmin ||SUx - Sb|l2. Furthermore for any x E RR we find the loss b- Uxx We also define Q = (~U S SU) as the sketched inverse covariance matrix with scaling 7 = representing the standard correction accounting for inversion bias. We condition on the high probability event & guaranteed in Lemma 2 and consider LEE ( -x(x") By Pythagorean theorem, we have L(Es - L(x*) = [|U [x)) - Ux*]| Note that by the normal equations we have x SU SU) S Sb "QUTS S Sb, and also S S Σm x.x_. These two facts lead to writing the bias as follows:
as [(x)
Ux bll Additionally, we use to denote the residual
L(Eg[x)) L(x") = 11 [EX
Using leave-one-out technique, we replace Q with Q matrix S without the ith row, by noting that Q the Sherman-Momee formula. This leads to the following relation:
where S_, denotes
2 U x U)- and applying
m
[(Ex(x))- L(x*)= < 2 Ee IQ_ U r1| +2
Eg
[
I %
1
1Zor
|Zar
where 7/ = 1 UU_ x;. Due to the subspace embedding assumption and assuming
m large enough, we have |Q|| = O(1) and also ||Q = O(1). The first term |Zr]|| is quite straightforward to bound since. if not for the conditioning the high probability event E, we would
have QQQ_.U x,X; r] = E[Q [Q_ = O, which follows from (b- - (b-Ux*) = O. We get an
upper bound on || as O d2 log(d/6)
|||r|22, which is sufficient for us.
smm
m
7

---

The central novelty of our analysis lies in bounding ||zr|| for (s, 81. 32)-LESS embeddings, which the dominant term. Our key observation is that, when examining random variable of the form x v for some vector v, the dependence on the sparsity of row x, only arises when considering moments higher than ([lo'(n]]' because otherwise we can simply rely on the fact that [F|x;x;] Thus,
when decomposing ||Zr]|l' we must carefully separate the contribution of near-second moments vs the contribution of higher moments to the overall bound.
To obtain this separation, we start by applying Hölder's inequality on ||Zr|| with p (log(n))
and = 1 Ollog(n) to get
|Zzr| [BE112 2=112-1107) 111) 1/8
1/9 sup Eg [v Q-;U X0,," r|] ||v||==
7/
Furthermore applying Cauchy-Schward inequality on the second term leads to
|Zzr|| < E31 7 111)) 1/P 7
1/24 (Es ||x] ||x/||| r1) 1/29
(Es 1|Q U'-,U" x,,
Unlike [18], we exploit the fact that ||x;|| 1/0(lo(((m)) = O(1) and get a constant upper bound on (Es ||-,U*
U x, ( However, this results in much more careful argument, requiring now an upper bound on (Esll 2 11')) /p for p = (llog(n)). First, we observe that
1/P
Bell 7
111)
- h
h
(Eg (6% 7)
where 7 = 1 + Eg (x, UQ_;U x.). In particular. for the second term, we have
(Es 1(x- - 7) 1/0 < 7 (IEs [(tr(Q-.) x; UQ_ U'"'"" + (EEs [tr(Q-i) 10
Egtr((Q -E))(QQ.))))/001
(3)
To bound the first of these two terms, we prove new version of the Restricted Bai-Silverstein
inequality (Lemma 3) for [s, (2),,,,))-La LLESS embeddings. Unlike [18], we provide a proof with any and any (81, (2) values. Furthermore, utilizing the subspace embedding guarantee from Lemma we prove a much more general result where the number of non-zeros in the approximate leverage score sparsifier can be much smaller than d (proof Appendix D).
Lemma 3 (Restricted Bai- Rai-Silverstein for (s, B1, -LESS embeddings). Let E N be fixed and U E Ruxd be such that U = Let x, diag(E)y: where y: E R" has independent #1 entries and & is an (s, B1, (a))approximate leverage score sparsifier for U. Then for any matrix with 0 3 C 3 (11) and any > 0 we have for an absolute constant e >> 0.
<. vdp² (1+
dplog(d/8)
(B (tr) - x UCU xxp)) (CCTT,P)"
8
Using Lemma 3, we upper bound the first term squared in (3) as Ô m (1 9) Moreover, also using Lemma we get a matching upper bound on 71 We then design martingale concentration argument to prove high probability upper bound on the last remaining term, |tr(Q-) Extr(Q-i).
which implies the desired moment bound (proof Appendix B), concluding the proof of Theorem
Lemma 4. For given 8 0 and matrix Q have with probability 1 ód:
d
(r(QQ=) Etr(E s c'y.
Vm log (m/8).
Completing the proof of Theorem 1. First, suppose that 2 O(polylog(d/// so that the bias bound can be achieved from Theorem 4. Our implementation mainly based on the online construction of approximate leverage scores, given the preconditioner P, using Lemma 1. Briefly, this construction proceeds by first sketching P using d x O(1/r) Gaussian matrix G to produce the matrix P = PG, and then, for each observed row a; of A, we compute = ||a P1 Assuming without loss of generality that d (nly(n) and adjusting 7, the estimates satisfy 3, O(ad)).
8

---

Next, we sample the non-zero entries of S corresponding to the observed row a,. i.e., the i-th column of S. Note that for this we only need to know the single leverage score estimate I,. Crucially for our analysis, the entries of this column need to be sampled i.i.d., which can be done in time proportional to the number of non-zeros in that column by first sampling corresponding Binomial distribution to determine how many non-zeros we need, then picking a random subset that size, and then sampling the random +1 values. Altogether, the cost of constructing the sketch is O(> "(y-"nnz(A () + (133dd))) Ory--unz(A A) +ace =f polylog(d)) by setting = (opolylog(d))))))) Finally, once we construct the sketch, at the end of the pass can run conjugate gradient preconditioned with P on the sketched problem, which takes (aad)).
We note that in the (somewhat artificial) regime where we require extremely small bias, i.e., e (ppolylog(d) /d), the bound claimed in Theorem 1 can still be obtained, since in this case for small enough y we have d = (dd' /V() with w 2.5, so we can rely on direct leverage score sampling (which corresponds to s = 1), and instead of maintaining the sketch, we compute the estimator x = (À A b directly along the way. This involves performing separate d x d matrix multiplication after collecting each d leverage score samples, to gradually compute Á A, and then inverting the matrix at the end. From Theorem 4, we see that it suffices to set sketch size m Ödd/ Vi), which leads to the desired runtime.
Completing the proof of Theorem 2. Here, we use a slightly modified variant of Lemma 2, given as Theorem 1.4 in [9], which shows that using a single pass we can compute sketch Å in time ((nnz(A) (A) + d²),, which satisfies the subspace embedding property (1) with r = 2/2 Then, we can perform the QR decomposition À QR and set P R-1 in additional time (dd) to obtain the desired preconditioner. Next, we use Theorem I to construct q i.i.d. estimators x; in second parallel pass, and finally, the estimators are aggregated compute = 1
q E x,, which satisfies E| - b < (1 ++ ((1/))) - bbl?. Applying Markov's 's inequality concludes the proof.
Experiments
In this section, we illustrate empirically how our results point to practical free lunch phenomenon in distributed averaging of sketching- based estimators. As mentioned in Section our construction from Theorem essentially works by taking subsample of the data and then mixing groups of those rows together to produce an even smaller sketch (see Figure 1). According to our theory, while the small sketch does not recover the same e-small error as the larger subsample, it does recover an e-small bias. Moreover, this happens without incurring any additional computational cost, as the cost of the sketching is proportional to the cost of simply reading the subsampled rows. This suggests that we can use sparse sketching to compress a data subsample down to a small size while retaining the least squares performance in a distributed averaging environment.
To verify this, we evaluate the effectiveness of distributed averaging of sketched least squares estimators on several benchmark datasets. Specifically, we visualize the relative error of the averaged sketch-and-solyy estimator L(x*)
against the number of machines g used to generate the estimate * q Ei= x; Each estimate x; is constructed with the same sparsification strategy used by LESS, except that instead of sparsifying the sketch with leverage scores, we instead sparsify them with uniform probabilities (which often sufficient in practice). Following [16), we call the resulting method LESSUniform. Within each dataset, we perform four simulations, designed so that the total sketching cost stays the same for all four test cases, by simultancously changing sketch size and sparsity. Concretely, vary these so that the product (sketch size nnz per row) stays the same in each case, so as to ensure that the total cost of sketching is fixed in each plot.
In Figure 2, on the X-axis we plot the number g of estimators being averaged, so that the bias of a single estimator appears on the right-hand side of the plot (large q), whereas the variance (error) appears on the left-hand side (q 1). In each plot, the line with nnz per row = 1 (and large sketch size) corresponds to uniform subsampling, whereas the remaining ones are sketches produced by compressing that subsample. The plot shows that decreasing the sketch size (i.e., compressing the sample) does increase the error of single estimator (as expected), however it also shows that the bias of these estimators remains essentially unchanged regardless of the sketch size (since all lines

---

MSD
Abalone
sinem 2000, mne per mome
Sketch aizem126, nez per rowml
Skatc! sicow1000. rin per rowe2
Skatch sizom6é miz per mow
Loxs/Least Squares Loss -1
Skotch sizow400, 1i2 5
Skatchh dize-32, niz per rewed
Sketch 200,
Less/Least Squares Loss
Sketch size-116 16, per rovn
20
10--
10
10
10
101
10
10
10
109
10
10
10a
10
Number of Machines
Number of Machines
Figure 2: Distributed averaging experiment on YearPrediction/See and Abalone datasets [8], showing that sparse sketching can be used to compress the data while preserving near-unbiasedness without increasing the estimation cost (see Appendix for results on the Boston dataset).
meet as q + oo), confirming that suitable sparse sketches that "compress" rows of data can preserve near-unbiasedness without increasing the cost.
This phenomenon may not occur for all sketching methods. Figure 4 within Appendix F showcases a few more interesting results. We first further demonstrate that suitable sketches that "compress" rows of data A into a single row of can preserve near-unbiasedness without increasing the cost. In particular, this desirable phenomenon that LESSUniform enjoys also extends to LESS proper, as well as the Gaussian and Subgaussian (Rademacher) sketches. In fact, we observe that LESS enjoys similar desirable performance as the Gaussian and Subgaussian sketches and virtually no least squares bias, while retaining the computational speedups of sparse sketching, suggesting that it attains the best of both worlds.
However, when we decrease the number of subsamples within leverage score subsampling, the bias introduced by subsampling increases expected. This happens as the number of subsamples is reduced without increasing the amount of "compression" as one would with LESS or LESSUniform. We also show that the subsampled randomized Hadamard transform (SRHT) can exhibit some amount of least squares bias as the sketch size decreases. The numerical results shown for the bias introduced by leverage score subsampling and the SRHT complement the lower bounds established in [18].
6 Conclusions
We gave new sparse sketching method that, using two passes over the data, produces a nearly- unbiased least squares estimator, which can be used to improve upon the space-time trade-offs of solving least squares in parallel distributed environments via simple averaging. In particular, our algorithm is the first to require only (dd log(nd)) bits of space and current matrix multiplication time Odd*) while obtaining = o(1) approximation in few passes. Our techniques are of broader interest to sketching- bresehing-based optimization algorithms, including Distributed Newton Sketch.
Acknowledgments
This work was partially supported by NSF CAREER CCF-2338655.
'One can think of LESS and LESSUniform generalizations of leverage score subsampling and uniform subsampling respectively, where we mix Ó(l/c) subsamples into a single row of the sketch. This suggests that sparse sketching, as an extension of subsampling, yields desirable bias reduction without significantly increasing computational costs.
10

---

References
[1] Dimitris Achlioptas. Database-
catabase-friendly friendly random projections: Johnson-Lindenstrane with binary coins. Journal of computer and System Sciences, 66(4):671-687, 2003.
[2] Naman Agarwal, Brian Bullins, and Elad Hazan. Second-order stochastic optimization for machine learning in linear time. The Journal of Machine Learning Research, 18(1)4448-4188, 2017.
[3] Nir Ailon and Bernard Chazelle. The fast Johnson--Lindenstrann transform and approximate nearest neighbors. SIAM Journal computing, 39(1 (9(1):302-322, 2009.
[4] Ahmed El Alaoui and Michael W. Mahoney. Fast randomized kernel ridge regression with sta- tistical guarantees. In Proceedings of the 28th International Conference on Neural Information Processing Systems, pages 775-783, 2015.
[5] Nima Anari, Yang Liu, and Thuy-Duong Vuong. Optimal sublinear sampling of spanning trees and determinantal point processes via average-case entropic independence. In 2022 IEEE 63rd Annual Symposium on Foundations of Computer Science (FOCS), pages 123-134. IEEE, 2022.
[6] Zhidong Bai and Jack W Silverstein. Spectral analysis of large dimensional random matrices, volume 20. Springer, 2010.
[7] Daniele Calandriello, Michal Derezinski, and Michal Valko. Sampling from a k-dpp without looking at all items. Advances in Neural Information Processing Systems, 33:6889-6888 2020.
[8] Chih- Chung Chang and Chih-Jen Lin. LIBSVM: library for support vector machines. ACM Transactions on Intelligent Systems and Technology, 2:27:1-27:27, 2011.
[9] Shabarish Chenakkod, Michat Dereziński, Xiaoyu Dong, and Mark Rudelson. Optimal embed- ding dimension for sparse subspace embeddings. 56th Annual ACM Symposium on Theory of Computing, 2024.
[10] Nadiaa Chepurko, Kenneth L Clarkson, Praneeth Kacham, and David P Woodruff. Near-optimal algorithms for linear algebra in the current matrix multiplication time. In Proceedings of the 2022 Annual ACM-SIAM Symposium Discrete Algorithms (SODA), pages 3043-3068. SIAM, 2022.
[11] 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, 2009.
[12] Kenneth L Clarkson and David P Woodruff. Low rank approximation and regression in input sparsity time. In Proceedings of the pages 81 81-90, 2013.
forty forty-fifth annual ACM symposium on Theory of Computing,
[13] Michael B Cohen, Cameron Musco, and Christopher Musco. Input sparsity time low-rank approximation via ridge leverage score sampling. In Proceedings of the Twenty-Eighth Annual ACM-SIAM Symposium on Discrete Algorithms, pages 1758 1777. SIAM, 2017.
[14] Michael B Cohen and Richard Peng. Lp row sampling by lewis weights. Proceedings of the symposium on Theory of computing, pages 183-192, 2015.
[15] Michal Dereziński. Stochastic variance-reduced newton: Accelerating finite-sum minimization with large batches. arXiv preprint arXiv:2206.02700 2022.
[16] Michat Dereziński, Jonathan Lacotte, Mert Pilanci, and Michael W Mahoney, Newton-LESS: Sparsification without trade-off's for the sketched newton update. Advances in Neural Informa- tion Processing Systems, 34:2835-2847, 2021.
[17] Michal Dereziński, Feynman T Liang, Zhenyu Liao, and Michael W Mahoney. Precise expres- sions for random projections: Low-rank approximation and randomized newton. Advances in Neural Information Processing Systems, 33, 2020.

---

[18] Michat Dereziński, Zhenyu Liao, Edgar Dobriban, and Michael Mahoney, Sparse sketches with small inversion bias. In Conference on Learning Theory, pages 1467-1510. PMLR, 2021.
[19] Michal Dereziáski and Michael W Mahoney. Distributed estimation of the inverse Hessian by determinantal averaging. In Advances Neural Information Processing Systems 32, pages
1 11401-11411 2019.
[20] Michat Dereziński and Michael Mahoney. Determinantal point processes in randomized numerical linear algebra. Notices the American Mathematical Society, 68(1):34-45, 2021. 2021.
[21] Michat Derezinski and Michael W Mahoney. Recent and upcoming developments randomized numerical linear algebra for machine learning. In Proceedings of the 30th ACM SIGKDD Conference on Knowledge Discovery and Data Mining, pages 6470-6479, 2024.
[22] Michal Dereziński and Elizaveta Rebrova. Sharp analysis of sketch-and -ketch-and-projeett methods via a connection to randomized singular value decomposition. SIAM Journal on Mathematics of Data Science, 6(1):127-153, 2024.
[23] Michat Dereziński and Manfred K. Warmuth. Unbiased estimates for linear regression via volume sampling. In Advances in Neural Information Processing Systems 30, pages 3087-3096, 2017.
[24] Michat Derezińisi and Jiaming Yang. Solving dense linear systems faster than via precondition- ing. In 56th Annual ACM Symposium on Theory of Computing, 2024.
[25] Edgar Dobriban and Sifan Liu. Asymptotics for sketching in least squares regression. Advances in Neural Information Processing Systems, 32, 2019.
[26] Petros Drineas and Michael Mahoney. Randnla: randomized numerical linear algebra. Communications of the ACM, 59(6):80-90 2016.
[27] Petros Drineas, Michael Mahoney, and Muthukrishnan. Sampling algorithms for & regression and applications. In Proceedings of the seventeenth annual ACM-SIAM symposium on Discrete algorithm, pages 1127-1136, 2006.
[28] Robert M Gower and Peter Richtárik. Randomized iterative methods for linear systems. SIAM Journal on Matrix Analysis and Applications, 36(4):1660-1690 2015.
[29] William B Johnson, Gideon Schechtman, and Joel Zinn. Best constants in moment incqualities for linear combinations of independent and exchangeable random variables. The Annals of Probability, pages
1985.
[30] Alex Kulesza and Ben Taskar. k-DPPs: Fixed-Size Determinantal Point Processes. In Pro- ceedings of the 28th International Conference Machine Learning, pages 1193-1200, June 2011.
[31] Alex Kulesza and Ben Taskar. Determinantal Point Processes for Machine Learning. Now Publishers Inc., Hanover, MA, USA, 2012.
[32] Jonathan Lacotte, Sifan Liu, Edgar Dobriban, and Mert Pilanci. Optimal iterative sketching methods with the subsampled randomized Hadamard transform. Advances in Neural Information Processing Systems, 33:9725-975 10119725-775..20 2020.
[33] Daniel LeJeune, Pratik Patil, Hamid Javadi, Richard G Baraniuk, and Ryan J Tibshirani. Asymptotics of the sketched pseudoinverse. arXiv preprint arXiv:2211.03751, 2022.
[34] Yi Li and David Woodruff. Input-sparsity low rank approximation in schatten norm. In International Conference on Machine Learning, pages 6001-6009. PMIR, 2020.
[35] Miles E Lopes, Shusen Wang, and Michael W Mahoney. A bootstrap method for error estimation in randomized matrix multiplication. The Journal of Machine Learning Research, 20(1):1434- 1473, 2019.
12

---

[36] Ping Ma, Yongkai Chen, Xinlian Zhang, Xin Xing, Jingyi Ma, and Michael W Mahoney. Asymptotic analysis of sampling estimators for randomized numerical linear algebra algorithms. The Journal of Machine Learning Research, 23(1):9770-0014 2022.
[37] Per-Gunnar Martinsson and Joel A Tropp. Randomized numerical linear algebra: Foundations and algorithms. Acta Numerica, 29:403-5777 572,, 2020.
[38] Xiangrui Meng and Michael W. Mahoney. Low-distortion subspace embeddings in input- sparsity time and applications robust linear regression. In Proceedings of the Symposium Theory of Computing, STOC 13, pages 91-100, 2013.
[39] R. Murray, J. Demmel, M. W. Mahoney, N. B. Erichson, M. Melnichenko, O. A. Malik, L. Grigori, M. Dereziński, M. E. Lopes, T. Liang, and H. Luo. Randomized Numerical Linear Algebra perspective on the field with an eye to software. Technical Report arXiv preprint arXiv:2302.1 1474, 2023.
[40] Jelani Nelson and Huy Nguyên. Osnap: Faster numerical linear algebra algorithms via sparser subspace embeddings. In 2013 ieee 54th annual symposium on foundations of computer science, pages 117-126. IEEE, 2013.
[41] G. Raskutti and M. W. Mahoney. A statistical perspective on randomized sketching for ordinary least-squares. Journal of Machine Learning Research, 17(214):1-31, 2016.
[42] Tamas Sarlos. Improved approximation algorithms for large matrices via random projections. In 2006 47th annual IEEE symposium on foundations of computer science (FOCS 06), pages 143-152. IEEE, 2006.
[43] S. Wang, A. Gittens, and M. W. Mahoney. Sketched ridge regression: Optimization perspective, statistical perspective, and model averaging. Journal of Machine Learning Research, 18(218):1 50, 2018.
[44] Shusen Wang, Fred Roosta, Peng Xu, and Michael W Mahoney, GIANT: globally improved approximate newton method for distributed optimization. Advances in Neural Information Processing Systems, 31:2332-223 22222222222222 2018.
[45] David P Woodruff. Sketching as tool for numerical linear algebra. Foundations and Trends@ in Theoretical Computer Science, 10(1-2):1-157, 2014.
13

---

A Detailed preliminaries
We start by providing several classical results, used in our analysis. The following formula provides way compute the inverse of matrix A after a rank-1 update, given the inverse before the update.
Lemma (Sherman-Morrison formula). For an invertible matrix A Raxd and vector u, R A + uv is invertible if and only if 1 + v Auuu O. If this holds then,
A A-'u uv A-1 1 11
(A (A+ A-
In particular,
A-'u u 1+v
(A + uv -
A 'u
The following inequality provides a crucial tool for writing expectation of the product of two random variables as the product of higher individual moments.
Lemma 6 (Hölder's inequality). For real-valued random variables X and Y,
E
XY|
(E[XP))) /p (E((Y1)))
where p. q 0 are Hölder's conjugates, i.e. ! p + 1 g 1.
The following technical lemmas provide concentration results for the sum of random quantities. We collect these results here and then refer to them while using in our analysis.
Lemma 7 (Matrix Chernoff Inequality). For i 1, 2, ... consider sequence Zi of d d positive semi-definite random matrices such that E[ Σ, Z,] Ia and |Za||
s R. Then for any > 0, we
have
Pr
Amax
E 1=1
2
21++0)
S d. exp
<deexp
(2+c)R) (x
( ne
where Xmaz denotes taking the maximum eigenvalue.
Lemma (Azuma's inequality). If (Yo, Y1, Y2,...) martingale with IY, - Yj-1l S oj Jhenfor
any m, À > 0 we have
X (- "5 &
Pr ( Yol 2 A) 2- exp
Lemma 9 (Rosenthal's inequality ([29), Theorem 2.5 and Corollary 2.6)). Let 1 S < so and X2, Xn, are nonnegative, independent random variables with finite th moments then,
1/P
E
Σ
x
2p
max
{
Σ ELXXI.
E EEX?
log(p)
Furthermore, for mean-zero independent and symmetric random variables have
1/P
{ ( EXX))
1/2
1/p
E
Σ
X,
2p
max
E
E[X?!
log(p)
Lemma 10 (Bai-Silverstein's Inequality Lemma B.26 from [6]). Let B be d d be fixed matrix
and x be a random vector of independent entries. Let [[ax]] = 0 and [[a"]] = l,and E| < vn.
Then for any p Z 1,
E x Bx Lr(B)) S (2p)². ((u)tr(BB"))/ + u2 (2 (BB"))))) )02
14

---

B Inversion bias analysis
We the notion of an (e unbiased estimator of A defined in [18].
Definition 3 ((e, -unbbiade estimator). For e, 0, a random positive definite matrix B gdxd is called an (e, 8) unbiased estimator of A if there exists an event & with Pr(E) 2 1 - ô such that,
A
Ee B 3 (1+ and, B 3 0() A,
when conditioned on the event E.
In this section, we give formal statement and proof for Theorem 3, which then used in the proof of Theorem 4. We replace A with U such that U consists of d orthonormal columns spanning the column space of A. Here S.. E Rm*x denotes LESS sketching matrix with independent rows , = y" diag(E) where consists of +1 Rademacher entries and & is an (s, 31. 52)-
approximate leverage score sparsifier. Note that E[xx* | I, We assume that the sketching matrix S. consists of m 2 10d i.i.d rows and 3 divides m. Also, we assume that the Sn satisfies the subspace embedding condition for (Theorem 2) with v = Let Q = (7U* S,,, U)- where
m-dd
Theorem (Small inversion bias for 11,81,30 (2,) -LESS embeddings). Let å > 0 satisfy ma and m 2 Od)). Let S., E Rm an (1. - embedding for data matrix E Rexd such such that,
that U = I. Then there exists event with Pr(E)
1 I Q 21 when conditioned on & and
I+c I3 EE (1++))
2
where = 0 (0llo""00m/1
1
Proof. Let S-; denote Sm without the 1th row, and S denote S. with the 1th and jth rows removed. Let Q-; ; = (TUS $_,S_;U) -1 and Q-ú = (+U -h We proceed with the same proof strategy adopted in [18]. We define the events & as follows:
E U
Σ
x
U T - 1,2,3 .E= 1U>
==
Note that event Ej means that the sketching matrix with just (1/3)rd rows (scaled to maintain unbiasedness of the sketch) from Sm, satisfies a lower spectral approximation U I. Also we notice that events & &2. &g are independent, and for any pair (i,j) there exists at least one event Ex, k E (1, 2, 3) such that & is independent of both x; and x,. Furthermore conditioned on &g we have
Q- and Q-ij 6 Ia.
Note that as guaranteed in Theorem 2, we have Pr(Ex) 2 1 for all k and therefore Pr(E) 2 d with 8' 6/3. Let Eg denote the expectation conditioned the event E.
I I E-E [Q] = -E [Q] + (E S.
-Eg Q EEg QU
U]
Q U
-E&(QQ + TEE
x UQ U 'x;
m
Es(Q-, U (xxx²" - DJU
+ Eg QQ_
Q Es [ (2
1)Q_,U x U
%
Zo
Z.
22
where 74 = 1 2XX The second equality follows by noting that ST Sm m Σ X and using linearity of expectation. The third equality holds due the application of
Sherman Morrison's 's (Lemma 5) formula on Q = SAUTST,S_.U UT upper bounding Zo and use the following from [18].
_
XU)- We start by
15

---

Lemma 11 (Upper bound on Zool) For any lo > 0 we have,
DJ 12 k8' + J Pr(x UU x, 2 x)dx
| (x
Using Chebyshev's inequality we have, Pr(x, UU x; 2 x) I Varlx; UU x:) Using Restricted x²
ani-Silverstein inquuality for (s, 31, appprxximat LESS embeddings i.e. Lemma (proved in
Lemma 16) with p 2, we have Var(x; UUTx) xx) S cd dlog(d/6)
for some absolute constant
Let m² and < mi we get,
1
d
df log(d/0)
(a
dl log(d/0)
Zoll O
-0
(4)
m
m
smn²
m
sm2
We use the bound on the term
|| directly from [18], provided below as a Lemma.
Lemma 12 (Upper bound on |/.] [18).
| [Q-1- QII O(1/m).
(5)
remains upper bound ||Z2|l.
|Z2| = || [(2 2 1)Q- x $
sup
Ee[ - 11.1vQQ-UUU x Uz|]
7
||v||=|,
||||==1
Applying Hölder's incquality with p (log(n)) and g - 1 1+ for A [[lox(n]] get,
1/p
Eg || Y -11)
(Eg ||v Q- U X Uzl 1/4
/Zll S
sup
|w||=1, ||
74
7 - -1p]]
(Re I Q-;Uxxx x; 201) 1/20
(Re [|x] Uz| (0]] 1/2g
< Re
sup || ||
sup
%
|
0(1)
0(1)
(6)
where we used Cauchy-Schwarz inequality on Es ||v Q-;U' -' x,x; Uz Now note that x; diag(E) Ys, we have |x,|| poly(n) and therefore ||x, = 0(1). We now show the terms involving exponents depending on q are ((1) highlighted in the incquality (6).
1/20
1/20
sup Ee ||v Q-UU'x; 1210)) ||v|==
sup Es || Q- x 41 Q- x, 2011
||v||=1
1/29
sup Ee [ x Q- ||x: ||x:||20|
||v|==
0(1) (Es ||Q x || 1QQ 1/24
(()) (2 Eg [IQ
U x; | Q- UT (2t1] 1/24
Here E' an event independent of x,. Without loss of generality. we can assume that E = Eı ^ &g. We first condition Q-, and take expectation over x;, Also note that Q and x; are independent and event E' is independent of x,, and furthermore Eg [x] xx]
F[x x;] 1. We get,
1/2g
sup E. I U x; [29] ||v||=1
< 0(1) (E. ||| - 1/24
Now we use that conditioned on E", 1 + & and |
= 1, we get,
1/2g
sup Ee [|v U Uxx1111 x3120] ||v|=1
=(1) 0(1)
16

---

Similarly, using |xx||> = (11) and g [x] x;} = E[x; x;} = 1,
sup
(Ee
||]] Un 20] 1/24 = O())
|| =1
Es -1)")
1/p
Now we prove an upper bound p is even. We have,
Without loss of generality, we assume that
(7 -11 -1)7 s 2 Eg
( -1)
Eg
where E' is event independent of xs. The above can be upper bounded as,
Bg (7: 7/
(2-1))
1/p
$ (Egg [[> - -b))o))
- - -7))))
< 7 - 71 + (Bg' [(7: 17x-9)99 -))))/))
where 7 = 1+ E. (x] UQ- x,). As E' is independent x, and Q_ is independent of x,, get Ee' (x, UQ_ x1) = (r(trQQQ)) Therefore, 7 + 2 E (Q-+). We now aim to
upper bound (Bg' - -)) as.
(Be' 1(% - 7)
[T [(tr(Q-)) x, UQ- x;) 1/8 (Eg' [tr(Q-;) Eg (QQ
Using our new Restricted Bai-Silverstein inequality from Lemma 3 (restated as Lemma 16 and proven in Appendix D), we have
UQQ xx)x))))) <c. <c. (14 dp Aiplog(d//) (d/5)
(Eg' [(tr(Q-)) x UQ_
We now consider (Eg' Q ( show that (r(QQ.)) -e. tr ( I on this high-probability event we have,
In Lemma 4 (restated below as Lemma 13), -dlog" (m3 (m/8) with probability at least - Conditioned
c'z vm
dy Vm
(Eg' [tr(Q-<) Ee to (r(2tt))))) <
dlog 4.5 (m//8)
for absolute constant c' > 0. Therefore we get with probability at least 1 5.
(Eg' [(7 - >> < 2 c-piva. [
dplog(d/8)
d7 Vm
d dogg dlog (m//))
As m d, we get (Eg' [(% - 2/1 (idlog++3 (Viloog m
1+ V
Also using the analysis
in [18] for upper bounding | 71. get matching upper bound 17 7/ follows:
ly 71 = = [Vdllg = ((dlog".((//0)
d (1+ (1+]
m
Substituting these bounds in (7) and then in (6) we get,
Z O Vdllog V
1 d
(8)
Combining the upper bounds for Zo, Zi and Ze using relations (4,5,8) we conclude our proof.
We now provide the proof of Lemma 4, which restate in the following Lemma.
17

---

Lemma 13. For given ā > 0 and matrix Q-s we have with probability 1 8:
cy
[tr(Q-:) g. tr tr()):)l
dllog (m/5) (log)(((mm/))
mm
for an absolute constant e' > 0.
Proof. Writing tr (r(Q-)) - ggg tr (Q-1) as finite sum, we have,
tr (0r(((1) -)= = Σ Ee tr((Q-i) (r((Q-))-BB - tr( (r(Q-;) E..trrQQ
1=1
Denoting X, ,tr(Qii)) with Xo = 18 tr tr(Q-i), we have the following formulation
m
tr(Q - tr (Q-1) 2 X, X,-1
with f Xj-1- The random sequence Xj forms martingale and tr(Q-:) Eg tr Q-1) Xm= Xo. We find an upper bound on XXy Xy-1l. To achieve that we note,
X, X--1 , (t((Q-)) --
-E Eg )-1)(tr(Q--)) Q-i) tr tr(Q-i))).
Therefore with w/ = (Bg's . E j-)r-))t - -Q-1) and X, = -(58's EE (e))-))t((Qs)))) have,
|X Xj-il V b, + x²l < |ر + Ixyl.
From [18], we have 1x5l S m
We now prove an upper bound on 1/
2 QQ by Ux UQ m
0 < tr(QQ -)) - tr(Q = tr
1+ UQ UTxj
;j
UQ UX,
U'x,
2 x UQ- U
m
Now look the term x, UQ2 U Uxx For any 0 and any k 0, by Markov's inequality have,
E UQ 2 by
U 1"1
Pr (x, UQ X'x 2 <
2 E
UQ is -
2*-1 (((
at
a
Let E, be an event independent of both x, and xj and have probability at least 1 - 8'. Therefore we have E|.] < 2 [.. We get,
2* Es, ||x UQ", U'x, tr(Q2,)) 2* Ee, (tt (Q (1
Pr (xj UQ U'x 2 2) <
at
a
We now upper bound both terms the right-hand side separately. Considering the term Ex, [|x] UQ U"x X'x;- - tr(Q2,))) and using Lemma 3 we get,
Es, [x] U'x tr(Q2 ct 3.3. (dkklog(a/o) + +d)
k/2
(10)
8
Now considering the second term in (9), i.e., Es, (trr(22,)))) We use that conditioned on E,, we have Q_ij 61.- Therefore,
Es, [[tr(( < 6* d
(11)
18

---

Substituting (10) and (11) in (9),
2k c 1.3. (dklogld/8)
+ +d)
6/2
2 6 dh
Pr (x] Uxx 2 20 S
+
ak
a
2 c h3k dh (k loggd//)))/1
Pr Px U a)
a
for some potentially different constant c. Consider log(m/5) log(2)
and 4.c- k2 d v /lgg(//6))
and have with probability at least 1 86/m,
x ij X'x S 4 1. d. Vlog(kd/o).
This implies that for an absolute constant è' we have,
[tr(Q-n) - tr(Q s d 2 d log (m/0).
Therefore we now have an upper bound for |//l
My = g'
Q-i))) s 2d' 7
d 1og3-0mm///m (m//))
m
This means for all j. we have,
[X - X 46 2 d log 3.5 (m//8)
m
3.5 (m/ Then Σ of m with probability least 5. Consider e Ad' 2 dd 3.5 dl log"(m//// (m/8). Applying Azuma's inequality (Lemma 8) with X Vm 4d d log (m//). We get with probability at least and for potentially different absolute constant d'
x2 16c2.
c¹y
|Xm Xoo S vm
dloo' (m/0))
This concludes our proof.
B.1 Proof of Corollary 1
We consider the following Distributed Newton Sketch method:
X Xx
Pto
g+
where g: fx;) and H;
ATSTTTTSSA is a sketched estimate of the Hessian H, m-
m
2 fxx,,, when expressed as the matrix product H, = A, A, for the appropriately chosen x d matrix A, (as can be done for any generalized linear model). This update can viewed as an approximate Newton step with the Hessian inverse estimate H, = + >. H. As long as each H. is an (Ve,6/2g) 2//2//2uubbaseed estimator of H then using Lemma from [18] we get that the averaged Hessian inverse with q = ((1ee) satisfies:
I+C H. + 3 + H (=V(+00//))-
Using standard approximate Newton analysis (e.g., see Lemmas and 3, and the related discussion in [15)), this implies that when x, is in sufficiently small neighborhood around x" (determined solely by the strong convexity and Lipschitz constants of the Hessian of f), we have:
f(x++1) f(x++1)-((x)) ( = O(e). f(x,) - f(x*)
Adjusting the constants appropriately, relying on Theorem for constructing H. and following the complexity analysis from Theorem 2, recover the claim.
19

---

C Least squares bias analysis: Proof of Theorem 4
In this section, we aim to prove Theorem Let data matrix containing n data points and b E R" is a vector containing labels corresponding to n data points. We adopt the same notations as used in the proof of Theorem 5. Let x argmin || Ux S..bl1 Furthermore for any we use r to denote the residual b - Uxx. We aim to provide an upper bound on the bias introduced due this sketch and solve paradigm, i.e. (((])))) I(x*). Similar to Theorem we condition on the high probability event & and consider ((Ee(x)) - L(x*). By Pythagorean theorem, we have
argminx ||Ux bll where U R is the
R''
find the loss at as L(x) ||Xx bll?. Additionally,
[[Eg[x]) - L(x") = [Exx)) - Ux 11 Also,
L [ - L(x') [ { - Ux 1x = |B - x
=I|Ex(( (U S Sm U) U Sm b²bl" 2
|E S. Sm U- U S [b UU b)]
Note that b - :U = = Ux = r. We get,
T.(Bx(x)) I(xx) | (U S Sm UU S. rrl
Consider Q = (TUTS --
- (7U (00:==(UUSSSSSS... "T
m and
m-d
Eg (e/h(7 Sm U) U ST Smr)
=>BE QIEEIQUTTT
-UXx,
where we used linearity of expectation in the last line combined with S m = 1 Σ i= X.X.
Using Sherman-Morrison formula (Lemmaa 5) we have QU',, = I+ UQ_.U
Q_ x.
Denote
7/ = 1+ x UQ- U"x,, x; and substitute we get,
BE((US U) Smr] -Ee ( 2 Q- Uxx:TT
Es [(2-1) Q- Q-,U'xxx;1
(SUUSS U0- S S Be
So we get the following decomposition:
[(E=(x)) L(x*) = [E[(UTST] ST U) Sm
Rgg (Q_, x Eg
2 -
%/
(2--) %
< 2 |E (Q- U x
+2
Ee
Q Q-;U U x,
(12)
Zerl
IZar 1
Note that a similar decomposition was considered the proof of Theorem (see Appendix B) with slightly different Zo and Zz. We first bound |/rr|² the following argument. Without loss of generality, we assume that events E and Eg are independent of x, and E' E1 A Ez.
Zor E x, r] UUQ= 1e]
Pr(E)
E
U x; 1g, 1e (1 1 1e))
Pr(E)
EQQ U xxx, r 18, 1e,]
E U 18, 11, 1-8,)]
Pr(E)
Pr(E)
[E xxxxxxx1
Eg
QQ x, r 1-5))))
20

---

Now note that in the first term Q- is independent of x, and x; is also independent of the event E Using this with the fact that E X U = O. Therefore,
I we get Eg [Q- (QQUUXxx Eg QQ
U
since
1 12 1 E [Q-UU x r 3, 11|
- ?2 X r 1-
The last inequality holds because conditioned on E", know that Q|Q=il =ll S 6. Using Cauchy- Schwarz inequality we have,
|Zor| 2 S 72 ViBg U 1-c,] v B [r [TTx,x]11 rl
Note that E [r x,,, r] = ||-|| Also, can bound g [x] UU x; 1 Egl following:
Eg' |x| UU x; 1..1. 1-8a] = Pr(XTUU X; 1-80 > x)de
%/a UUX UU X, -15s > xddr + / Pr(x] x, 1-50 > x)der
S y0 + S Pr(x, TU > x)dx.
Using Chebyshev's inequality, Pr(x] UU x; 2 x) v Var(x] UU x.)
By Restricted Bai-Silverstein, for some absolute constant
Lemma 3 with p 2, have Var(x; UU x,) d d (log(d/0)
Let y 2 and 8' m
we get,
1 d
de log log(d/0)
-0(4
de log(d/0)
E E UU x 1-5, =
m m²
sm2
m2
sm?
This finishes upper bounding |Zo|| as:
d log(d/8) +
|Zar|1 = 0 m=
|||||
(13)
smn
Now we proceed with ||,r||
|Zzr| = Ee ( -1 Q- U xxXxxxx)
Applying Hölder's inequality with p = (log(()) and g = + A where A
[oog(]]] we get,
I
1/P
1/a
|Zzr|| <
Eg %
7 -
11))
sup Eg |v Q- xxx/1|² ||v|=1
7 - 1
1/p
E QQ-UUXx112) 1/2g (Iee ||x] r1120) 1/24
Egg
7/6
Eg
7
118)
1/P
(2. Eg [||Q_, x 11 1| 1Q xx)1 1/2g (22 Eg ||x]] r|| ||xx r1| =1) 1/2g
(Ux,11)))))))
(2-Eg-
%
Since
||x;||
(n), have |x,||2 22 O(1) and therefore we have IQ_.U x, 24
(()) Q|Q UT 12!2 Also Eg ||Q_, UT x1 Eg | 2| Using that conditioned E
1/24
have
||Q=>|| < 6, get Eg |Q, UT x;| {1(²))
(L). Similarly using ||x;||2< م ((1)
Eg |xx r (20) 1/20
get
(1) ||]]| This gives us:
1/p
Zzr|| < 0())
Eg
1l))
|||||
7/4
21

---

Now using (8) from the proof of Theorem 5, we get,
LE 7 - 11²) 1/p
VVlogt
1+ / )
=0
7/
m
Finally the bound for |zrr|² follows as:
|Zar|2 O (dlo0((n100 m²
(0(1los(110) (1 (1 ) |||||/
(14)
Combining (13) and (14) we conclude our proof.
D Higher-Moment Restricted Bai-Silvestein (Lemma 3)
In this section, we prove Lemma We need the following two auxiliary lemmas derive the main theorem of this section. The first lemma uses Matrix-Chernoff concentration inequality to upper bound the spectral norm of U Ue where & is a (s, 31, (32 appprxxmaeee leverage score sparsifier.
Lemma 14 (Spectral norm bound with leverage score sparsifier). Let U Ruxd has orthonormal columns. Let & be a (s, $1, 52) (a)-approximate leverage score sparsifier for U, and denote Ug = Ui Then for any > we have,
Pr 1E 2 (1+ 3d log(d//)) (
8 )
<6 if < d,
Pr U Ue 2 (1+ (l+3log(d/0)) Só ifs2d.
Proof. Writing Ug Ue as sum of matrices have,
UNU = Σ u
n
Σ Z,
P:
=11
1=1
and E [Z] = u.u. Also Σ EZ: p. u.u Note that Z/s independent random variables
where P; min(1, (1, +3,1, and Z:
U = la. If p; = then Z, = i and therefore d
|Z || ||:1| < 1. If < 1, we have | + ||,11 a3, 4. As 4 S 814,, we get
|Z, < x- Therefore |Z,|| < max d) for all Denote R max(1, d1. We use Matrix Chernoff (Lemma 7) to upper bound the largest eigenvalue of U Ug. For any e > 0, we have,
Pr
Amax
Σ Zi
n
2 21+
S exp (-12+
22
With R max(1, 1, و and depending on the case whether s < do
we get
3d log(////))
Pr
Amax
Σ Zi
> 1+
Só if < d,
1=1
Pr
Amax
Σ Zi
2 +1+ 3lgg(d/8)
<o d.
i=1
Let Ag denote the event
U
Ue
< 1+ logl//6) holding with probability at least 1 - ō, for small
8 > In the next result we upper bound the higher moments of the trace of U.CUg for any matrix
C 3 O(1) We first prove the upper bound in the case when the high probability event Ag does not occur.
22

---

Lemma 15 (Trace moment bound over small probability event), Let k E N be fixed. Let U E Ruxd have orthonormal columns. Let & be a (s, 11. (2)-appiappoo leverage score sparsifier for U. Ler Ue = Ui. Also let Pr(Ac) 2 (12djar and event E' be independent of the sparsifier &. Then we have,
Bg ( - (4k)ah
for any fixed matrix C such that 0 3 C 3 61.
Proof.
Eg [(HUUCCRR.CCO)) = ! P((((((_CUU 1-Ag 2 xlE')dr )dr
Pr(tr(U,cu(U,cCu)) 1-de 2 z)da.
The last equality holds because E' is independent of &. Consider some fixed y > 0.
y Pr((t((UCCU)))** 2k 1 Ag 2 dr
Eg [[r[U=CC <
Pr((tr(U&CUT)) 2h -14 2 x)
Sy 6 + B((tr( x² Zdv.
(15)
The last inequality holds because by Chebyshev's inequality Pr((I((UCCU)))) 2 x)
C(U Also note that,
4k
(HT(
Σ u Cu;
by Σ
Cu;
P;
For 1 i n, let R2 be random variables denoting in u Cu,, Then note that R, are independent p. random variables with E[R,] u Cu,, Also R, are non-negative random variables with finite (4k)th moment. Using Rosenthal's inequality (Lemma 9) we get,
4k
4k
n Σ
n
u Cu;
S 24h (4k)ak
Σ
E
Σ
E[R:]
Ps
Now E [R, tr(UCTT) and E[ can be found as follows: if p; 1 then R; u Cu,
and therefore R. = (n Cu,) 1k < ||u 221--1) u C S u C.ku,,
u, if p; < 1, we have,
1
E| E[R
= p:
(n Cu,)
Lkk
s4k-1 (31 (n. Cu,)
prth
dak-1
|||. 2(46 1) u Caklu
-
(Bi1))kk11
-1
-k u th
s4k-1
(311))4h-1
Now using 4 S ByT, we get,
-1 BiR I s u 1 u.
Therefore,
467
b u Cu,
max 1, ) tx((C*t) + (tr(C)))
IE
< 2 (4k)
Pi
23

---

Using the above inequality along with using C = 61 we get,
4k
E
Σ
u. Cui
2 (4k)a 6th dth = =(1)) th. (4k) ih dtt
spi.
Substituting the above bound in (15), it follows that:
s sy.d 8 (12)ah th (4k)ah dth 1 y
E
For > (12d)th dk and d < ,, we get the desired result.
We are now ready to prove the main result this section. The following result, which is central to analysis, upper bounds the high moments of deviation of quadratic form from its mean.
Lemma 16 (Restricted Bai-Siiverstein for (s, B1, -LESS embedding). Let p E N be fixed and U E Raxd have orthonormal columns. Let x; diag(E)y; where y, E R" has independent #1 entries and & is a (s, 31, 32) appprxxmaae leverage score sparsifier for U. Let Ue diag(E)U. Then for any matrix C 61 and any > we have,
(B [tr(CC) UCU x < <c p.v. (1+V
dp log(d//)
for an absolute constant > O.
Proof. Let x, diag(E)y; where y; is vector of Rademacher #1 entries. Denote Ug diag(E)U..
E. [[r((C) x. UCU [ E [tr(C)) y, U
E [tr(CC tr(U,CU:) (r(UCCU)) y, U
E ((U((C) +E (tr(UcCU)) y U. S 2p-1
(16)
First, consider Ti, substitute c h. and assume exponent p to be even.
E [t((C) - tr(U E (r(CUC")) (UC
E
Σ (5 1)u7 Cu;
by
E
1 1 Cu:
Pt
For 1 < i < consider random variables R₁ where R
in u Cu,. Furthermore let Y; R₂ uu Cu,, Then Σ Y, Σ C - 1)u7 Cu, and [[] = Y, are independent mean
random variables with finite th moments and therefore symmetric random variables, Lemma 9) to get,
Rosenthal's inequality (for
p/2
E
Y
< A(p)
Σ BYAP Σ
(17)
T
T?
where A(p) is a constant depending on We bound Ti and I separately, starting with T as,
T Σ E.p..

---

Recall that Y; = R. Cu, R₁ is always non non-negative because C is positive semi-definite matrix. Therefore we have,
B(Y))* S B(R.)).
We find the th moment of R2, E(R)) (u Cu,) )" if p; = 1f P. 1 then,
d
SP-1 (311))-1 u Cu,) dP-1 1
B(R,)* pf (u,
Cu.)
||u112(0-1) 1) u Ceu,
SP-1 (311))=1
d= sp u u
where in the last incquality we use | 1u:1|20001 -000011111001 = 14- & (311)p-1 get upper bound for T the following:
Summing over i from 1 to n we
E(Ra)* & max 1, ) tr tr(UCCU').
Σ E
(18)
Now we upper bound T.. Note that B(Y,)² = E|R = (u if p: = and if pj < 1 we have E[R.]² p. (. Cu,)² s d
814 ||u u CZu, S u C²u,. Summing over from 1 to we
get an upper bound for 7 as,
Σ B(Y))
D/2
I max (I. d
tr (UC
p/2
(19)
Substituting (18) and (19) in (17) we have,
E [[r((C) - trU)) < A(p) (max (1, (11 tr(UCCU)) + max (1. tu(UC*U)) p/2 )
(20)
Now aim to upper bound the term T3 in (16) i..., E ( - P condition over & and take expectation over yı. This requires using standard Bai-Silverstein inequality (Lemma 10), we get,
First,
E [(r(UeCU)) y. < B(p) E
+ U,
where B(p) is a constant depending on p. Since consists of +1 entries we have 14, 12n
S 1. Also
using tr(AB) & tr())tr(B) and considering the high probability event Ag capturing U U S
(1 max
3d (dd//6)) llog(d/))) We get the following,
E (t((UCC))) y, y 2B(p) E (U(U)CU(UCCU)))
-2B(p) (M(UUCCU LAe) 1/2 + +BB(p) [E(H(U,CUU,,CU) 1=Ae)
<2B(p) (1+ m 3dlog(d/0) , 3log(d/o)) ((I(U + 2B(p) (2p)
P/2
(21)
The first term the last inequality follows from the Matrix-Chermoff (Lemma 7) and the second term follows from Lemma 15 by considering p/2 (assuming that 5 is small enough so that Lemma 15 is satisfied). We now upper-bound bound E (t=(U)C² using Rosenthal's 's inequality for uncentered
p/2
(non-symmetric) random variables,
b Σ
= - =
u
C
1m p.
25

---

For S i S n consider independent random variables R;
uu u² We have,
P4
Σ R = t
Here R are positive random variables with finite (p/2)th moment. Using Rosenthal's inequality we get,
p/2
p/2
E
R2
S C(//2)
B( R2)p/2 +
E.R.
(22)
T2
72
It is straightforward to upper bound TE R2) as,
B(R))//= = (n, / if pı 1,
d 4 - u C"u p: 1.
Summing over from 1 to n we get upper bound for Ti as,
di d B(R S max ( (1.4-11 tu(U')"
(23)
Now we consider T2. It is simply given as,
9/2
Σ
ER'
=(t(UCC477))12
(24)
Combining (23) and (24) and substituting in (22) we get,
E (tr(U (tr(U, 2/2 < C(p/2) (max 1, (1 (T(UC*U*++++ (M(CCTT))))))))))))) (25)
di-1
Substituting (25) in (21) and let w max(L.
E CU) y ye < 2B(p) (2p)2p
+2B(p)C(p/2) (1 3w log( log(d/6))**² - ("CUC )/
(26)
Combining the bounds for Ti (20) and T2 (26) substituting in (16), and noting that tr tr(UCC*T)> = ur(C*) for any k, we get,
E [tr(CC) x, UCU UCUTx;)) < A(p) (WP-1 tr(CP) + (w tr ()2)
22 B(')(() (1 3w loggo(//))) w ( -
(tr (tr(C²))
(2)
Now we specify the various constants depending on p. We have A(p) < (2p) B(p) < (2p)",C(pp/) C(p/2) s pp/2, Also we use tr(Ck) < 6t since lute constant we have,
6I. This implies for an abso-
dp dplog(d/8)
E [itr(C) - UCU UCUTU()) 1/1 < c VVd. 1
where 8 > 0 is now arbitrary.
26

---

E Application to distributed settings with partitioned data
In this section, we extend our main result settings where the dataset is partitioned across multiple machines. We illustrate this in the distributed setting considered in [43]. Let A E IN xd and b E R with rank(A) d. In practice, we randomly shuffle the rows of A and b, and partition the data uniformly across q machines, with every machine getting chunks of size n denoted by (A,, b;). In this distributed setup, A, can be considered a uniformly subsampled sketch of A with sketch size n. However, similar to the analysis in [43], in the following analysis we assume that each machine constructs its sketch (A,, by) by uniformly sampling rows from A and b with replacement and independently from other machines. This setting is not ideally partition but we believe similar results can also be shown for the setting where we partition the dataset after reshuffling the rows. Let N max (A) denote the matrix coherence of A. Let S E R be a uniform sampling matrix with size where will assume that S also scales the samples by VN/n, which implies that E|S S] = I. Then the following holds:
Lemma 17 (Subspace embedding using uniform subsampling; based on Theorem 12 in [43]). Let S E R=*N be uniform sampling matrix of size n where n 2 O(plog(d/8)). Then with
probability 5 we have,
I+ A' 3 (1+ ATA
for 1) - O / /log(d/))/n).
Let x = A b denote the minimizer of least squares loss ||Ax - bbl² and x; = Alb, denote the solution to the sketched least squares problem at jth machine. Here A, = S,A and b; = S,b, where S; denotes the uniform subsampling matrix for 7th machine. We can use the model averaging result from [43] to prove the following result.
Lemma 18 (Model averaging result; adapted from Theorem 20 in [43]). Let there be q machines and each machine constructs its sketch (A,, b,) independently by uniformly sampling n rows from A and b, where n Oµ log(dg) 6)). Then with probability 0.99,
Ax AA bl( -blp s 1 (1+c( (1+ ((loog(d//))
1e² log' (d/2((/)))) | ||Ax -b12.
(27)
qn
n²
where x = -/E x and >> O is an absolute constant".
Since computing x, at any machine could be potentially expensive due to large n. we instead find x, at 8th machine using the sketching procedure from Theorem 2. Note that this requires applying our debiasing techniques to the least squares problem min AA, - b, Il2 at every machine. To that end, >i X; and consider Ep Ax bll², where Ep means conditional expectation given the let x = 7 sketches (A,, b;) at every machine. The following theorem states that if the bias from partitioning is much smaller than the sketching error in each machine, then our approach can be successfully used reduce the bias of the sketching estimate so that averaging will work as desired.
Theorem 6. Let there be q machines and each machine constructs its sketch (A,, b,) independently by uniformly sampling n. rows from A and b. where n 2 O(u log(dq/o)), Also, let ith machine construct &; via locally sketching A; and b, such that x; satisfies Theorem (with A; and ba) at ath machine. Then, the estimator & "9 Xe; averaged across the machines satisfies
- -pllAx-bil2 ( 1 1s (log(d/5) +
||AX - -/12.
+
qn
n²
with probability at least 0.9 and an absolute constant d' >
Proof. By the Pythagorean theorem we have:
Ep Ax b1l2 Ax" b112 (E/eII/AAI"""))
(28)
'The result Theorem 20 [43] also requires that || < for all i, which holds trivially if S, is a sampling matrix without replacement. However. in our setting, we are implicitly assuming that with high probability [|ST 11 = o( and absorb the extra resulting constant in c.
27

---

We proceed to upper bound Ep || A (x x²11 as,
Epl| - x )| < 2 (Ep]IA(x x²)1 + 11(1 - -xx((=-x))))))
We use Lemma 18 upper bound the second term and get,
E- 2 [BPlA([X- + of p log(d/00)
1² logggdd//// (d/ō1) ||Axx - -b112)
gn
n²
(29)
On the other hand,
Ep||A(X =EplA=E Ep + Ep x112
= -EplIA(x - Ep X| + | (Bp[)(]- -
(30)
We upper bound both terms in (30) separately. Proceeding with the first term we have,
| - Ep[x])|
Ep
A
x - Ex(xx)))
1
Ep
A
(x: - E>(
92
1
Σ
Epl|A(E: - Ep/2111
(31)
92
$
+ "
>
Epl|A, - Ex(&1111
(32)
q²
2(1 + n)
Ep (IA:( + -
-x;)))
q²
2(1 + n)
|A,x; - -b1l2
(33)
9²
2(1 + 1)
< q²
Σ
|A - -b11
(34)
The equality (31) holds as all X/s are independent. The inequality (32) holds due to the subspace
embedding property (Lemma 17). The inequality (33) holds due to the variance bound in Theorem 1.
The last inequality holds because x, minimizes the least squares loss || b, 12 at 7th machine.
We now upper bound the second term in (30)
A - A
(E [2]] -xxx x1
|| A Ep(x - x!I|
i=1
VI+ 7
(B.p - -x)il
(35)
Ve(1 +
g
|| b.ll
(36)
9
Ve(1 + m)) |A,' -bill
9
The incquality (35) holds due to the subspace embedding property (Lemma 17) and inequality (36) holds due to the bias bound in Theorem 1. Therefore we get,
(1 + m)
Σ A - -b:l)) 2
[Epl x3)
(37)
q2
28

---

Combining (28, 29, 34, 37) and assuming
7 < we get,
y (log(d/6))
r² log (d) ||Ax - -b|l2+
8
9
Ep||Ax- - bl12 v 1 2c
+
qp ΣΣ
| bull?
qn
2
4c
+
A,x* - -b.ll)
(38)
q²
==1
A
E
ΣE
A
bill²)
||Ax* - -bll?. we just use Markov's inequality to upperbound the |Ax* - bll? with probability 0.95. For the last term, have
second term in (38) by 160
Σ
A
bill) = 9
A bill2 (q² 1
Σ |A,x Dil1-11A,x byll $.J.
Again 1 E ||A,x* - 1 20| - bll with probability 0.95. The last term is an average
q of gf q random variables and we can provide uniform bound for the expectation of of them using the independence among different machines. We have E|[] A,x* - b;ll ||A,x** bylll < /B||| 1 , by can upper bound the average of these q² q random variables with 20 times their expectation with probability 0.95, finishing the proof of the theorem.
||Ax bbl Again using Markov's inequality we
F Further experimental details
Here, we provide small set of numerical experiments to empirically examine the relative error of distributed averaging estimates from individual machines to return an estimator x 1 E/=1 x; First, we show that when the sestching-based estimates x; have a non-negligible bias, so that the distributed averaging estimator X remains inconsistent in the number of machines even with an unlimited number of machines, as long as the space on each machine is limited, the averaged estimator's 's performance will be limited by the bias of the individual estimates. Second, we show that one can use sketching to compress a data subsample at no extra computational cost, without increasing its bias, which we refer as the free lunch in distributed averaging via sketching.
We examine three benchmark regression datasets, Abalone (4177 rows, features), Boston (506 rows, 13 features), and YearPredi ct onMDD (truncated to the first 2500 rows. with 90 features), from the 1 libsvm repository from [8]. All the experiments were run on an 19-13900k processor with 1 128GB of RAM and an NVIDIA RTX 3090 GPU. We repeat the experiment 100 times, with the shaded region representing the standard error. We visualize the relative error of the averaged sketch- and-solve estimator 4(x) M(x")
against the number of machines q used to generate the estimate
[(x')
Each estimate 5; was constructed with the same sparsification strategy used by LESS, except that instead of sparsifying the sketch with leverage scores, we instead sparsify them with uniform probabilities. Following [16], we call the resulting method LESSUniform. Within each dataset, we perform four simulations, each with different sketch sizes and different numbers of nonzero entries per row. We vary these so that the product (sketch size x nnz per row) stays the same, so as to ensure that the total cost of sketching is fixed in each plot.
As expected, decreasing the sketch size while increasing the number of nonzeros per row (effectively increasing the amount of "compression" occurring here by sparse sketching) increases the error in all three datasets. However, remarkably, increasing the amount of "compression" does not seem to increase the amount of bias. We can therefore conclude that sparse sketches preserve near- unbiasedness, while enabling us reduce the sketch size from subsampling without incurring any additional computational cost. The increase in error can be mitigated in distributed setting by increasing the number of estimates/machinss.
E1 Comparison with other sketching estimators
To further illustrate the phenomenon that suitable sketched least squares estimators enjoy small bias, provide further set of numerical experiments below in Figure 4 for the Boston and Abalone
29

---

Boston
10
in other experiments, but still potentially much Sketch per
Sketch size-170, per row-2
Loss/Least Squares Loss 1
Sketch sire=70, row=
Sketch size=30, row=13
10-1
10-2
10
10
10a
102
103
10
Number of Machines
Figure Comparison of the relative error of the distributed averaging estimator of sketch-and- solve least squares estimates where the sketches are constructed with sparse sketching matrices with uniform probabilities (LFSSUniform) on libsvm dataset Boston (see Figure for results on YearPredi celicciionnn ctionMSD and Abalone). For each dataset, the computational cost of sketching is the same all four parameter settings. Remarkably, sketching to smaller size appears to preserve near-unbiasedness without incurring any additional computational cost.
datasets. Recall that in Figures 2 and 3, we show that increasing the amount of "compression" of the sample, by reducing the sketch size when using sparse sketching matrices with uniform probabilities, does not increase the bias of the sketched estimator. Figure 4 further demonstrates that Gaussian and Subgaussian sketches, which are computationally much more expensive than our proposed sparse sketches, exhibit virtually no bias in sketched least squares (as expected). A similar conclusion can be made for the proper LESS method, which uses leverage score estimates, and is therefore more expensive than the LESSUniform method cheaper than Subgaussian sketches.
However, the unbiasedness does not hold for all sparse sketches, we demonstrate below in Figure 4. At a high level, we show that sketches constructed using leverage score subsampling and the SRHT can exhibit a non-negligible level of least squares bias.
Gaussian and Subgaussian sketches, we see in the figure, enjoy near-unbiasedness while not introducing additional bias as the sketch size decreases. This does not hold for leverage score subsampling. When we decrease the number of subsamples within leverage score subsampling, the bias introduced by subsampling increases as suggested by the lower bound in Theorem 10 of [18]. Intuitively, this happens as the number of subsamples is reduced without increasing the amount of "compression" as one would with LESS or LESSUniform
We also show that the subsampled randomized Hadamard transform (SRHT) also introduces increased bias as the sketch size decreases. This complements the lower bound established for the failure of SRHT and other data-oblivious sparse embeddings like CountSketch to satisfy the restricted Bai-Silverstein property that is structural condition for near-umbiasedness sketches [18].
For completeness, also demonstrate the desirable performance of LESS proper. Figure 4 demon- strates that the desirable phenomenon that LESSUniform enjoys also extends to LESS proper. In fact, observe that LESS enjoys similar desirable performance the Gaussian and Subgaussian sketches and does not increase (in fact, it has minimal) least squares bias, while enjoying the computational speedups that sparse sketches also enjoy, suggesting a best-of-both-wool property.
30

---

Boston Leverage Score Subsampling
Abalone Leverage Score Subsampling
105
102
LosslLase Squares Loss
10
10
10
10
10-1
10
10-4
10-
10
10'
10"
10
30
10
107
101
30
10
10
101
Number of Machines
Number of Machines
Baston Dataset- SRHT Sketch
Abalone Dataset- SRHT Sketch
Sistal 540-240
Skoul 525-125
Stater -2--138
Seatch
Lossteast Squares Lncs-1
Sieeet size
Losst Losst Laast Squares LoSS-1
sketch
10
10-4
10-4.
10
10'
10
10*
10'
Number of machines
Number of machines
Boston Dataset- Gaussian Sketch
Abalone Dataset- Gaussian Sketch
10-1
10-4
State? 118
Starch Eow
Less/ Squares Lazs1
Soilt? S1e=12
Loskteast Squares Loss-1
Shotch sizc
Sbrich Sue-16
20
10-3
10
10-1
30
10"
10
10-
10'
10"
26"
10*
Number of machines
Number of machines
Boston Dataset- Subgaussian Sketch
Abalone Dataset- Subgaussian Sketch
Siatal S4c-240
Shouh 320-
Statir
Statch Eoe
Losgleastt Squares Lass-1
Stetct size
Loss Lastt Squares LoSS-1
sletch sloc-
20-5
10 Number of machines
10'
10
10'
Number of machines
Boston LESS
Abalone LESS
105
Skatch ebre 310. mnt per ron
Losst Squ Squares Loss -1
10"
Soatch dite 135, nts periow 1
Skatch so0- mi 2
Soatch stee- od. miz per 2ow
Lost/ Squares Loss
Skatch si nns per 1
30
Soatch stca- 30. me per roumd
Skutch 30 nre per 12
Soatch slo0- ma per row-d
10
20
2011
10-2
10:5
10
30
Number Machines
Number Machines
Figure 4: Distributed averaging experiment for the Boston and Abalone datasets using five different sketching techniques: Leverage Score Subsampling, Subsampled Randomized Hadamard Transform (SRHT), Gaussian sketches, Subgaussian sketches, and LESS. Both Gaussian and Subgaussian sketches exhibit no observable least squares bias, whereas Leverage Score Subsampling and SRHT exhibit small but measurable level of bias. LESS enjoys the similar high performance as the Gaussian and Subgaussian sketches and does not increase (in fact, it has minimal) Icast squares bias, but as sparse sketch also allows for the significant improvements in runtime that sparse sketches enjoy, demonstrating that it enjoys the -oost-of-othhworldd
31

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
The paper provides concise mathematical proofs for the claims made in the abstract.
Justification:
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
It is fine to include aspirational goals as motivation long is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [NA]
Justification: The results proved in the paper hold in very general setting, without any limitations or requiring strong assumptions.
Guidelines:
The answer NA means that the paper has limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create a separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (c.g.. independence assumptions, noiseless settings, model sellsspeciiiationn asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect on the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-tex, system might not be used reliably to provide closed captions for online lectures because it fails handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and faimesss.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community, Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory Assumptions and Proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
32

---

Answer: [Yes]
Justification: All theoretical results have been proved. Most of the formal proofs are part of the appendix.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental Result Reproducibility
Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless whether the code and data are provided or not)?
Question:
Answer: [Yes]
Justification: The experiments demonstrate the results for solving few least-squares problems using the Matrix sketching technique described in the paper. All the experimental results are easily reproducible.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, No answer this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how replicate the results, access a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make clear how reproduce that algorithm. (b) If the contribution is primarily
model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (e.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should possible for other researchers to have some path to reproducing verifying the results.
33

---

5. Open access to data and code
Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Question:
Answer: [Yes]
Justification: A supplementary zipfile containing the code is provided to reproduce the experimental results. One of the experiments requires the file 'YearPrediction- MSD.txt file, which can be obtained at (https://archive ics. suuuucieuuuu edu/dataset. 203/yearpredictiooo ctinosdd).....
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips.. cc/
lttpc//gudees/oooeeeeeee wttpauammisssooooooooooooooo wwwsssoliiiiyyo.ooooociiiilliili olicyy for more details.
While encourage the release of code and data, we understand that this might not be possible, so "No" an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https: / //nips. httppppullllldeeeeee Sudeede//ooeebmiiiiccc iessonPoll onnoliccy for more details.
The authors should provide instructions data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
The paper describes the parameters e.g. sketch size, sparsity of sketch (nnz per row), and number estimatee(nummer of machines) along with dataset details to understand the results.
Justification:
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment Statistical Significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes]
Justification: The shaded areas our plots denote the standard error obtained after averaging the result 100 trials.
Guidelines:
The answer NA means that the paper does not include experiments.
34

---

The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
It should be clear whether the error bar is the standard deviation or the standard error of the mean.
is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out range (c.g. negative error rates).
If error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information the puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: All experiments were run on an 19- 13900k with 128GB of RAM and a NVIDIA RTX 3090 GPU.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental rums as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper).
9. Code Of Ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics https: pentps://neeriae /tt.rriur...p//e http/pbbii//thhcaaaaiaiieieeeeee lnes?
Answer: [Yes]
Justification: This work aligns with the NeurIPS code of Ethics.
Guidelines:
The answer NA means that the authors have reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- eration due to laws or regulations in their jurisdiction).
10. Broader Impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [NA]
35

---

Justification: Our theoretical work provides a Matrix sketching algorithm, recovering nearly unbiased least squares estimator. Solving linear systems upproximately is an important component of a large number of algorithms from varied domains. Our work can thus provide a way to recover approximation solutions with smaller bias. This does not lead any negative societal impact.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: This work does not release any data or model, and carries no risk of misuse.
The answer NA means that the paper poses no such risks.
Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging. and many papers do not require this, but we encourage authors to take this into account and make a best faith effort.
12. Licenses for existing assets
Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification: We test our theoretical results on a few publicly available datasets that have been cited in the paper as well.
Guidelines:
The answer NA means that the paper does not existing assets.
36

---

The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset is used and, if possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset. For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswi thcode. com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's s creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification:
No new assets released in this work.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset is used.
At submission time, remember anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification: This work does not involve crowdsourcing or research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According the NeurIPS Code of Bthies, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: This work does not involve crowdsourcing research with human subjects.
37

---

Guidelines:
The answer NA means that the paper does not involve crowdsoureing nor research with
human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent)
may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such the institution conducting the review.
38


Reviews:
[{"review_id": "XU10Cc4el3", "reviewer": "Reviewer_8xrQ", "summary": "This paper presents new techniques for distributed least squares regression using matrix sketching. The key contributions are:\n\n1. A sparse sketching method that produces a nearly-unbiased least squares estimator in two passes over the data, using optimal space and current matrix multiplication time.\n2. Improved communication-efficient distributed averaging algorithms for least squares and related tasks.\n3. A novel bias analysis for sketched least squares, characterizing its dependence on sketch sparsity. This includes new higher-moment restricted Bai-Silverstein inequalities.\n\nThe theoretical results are backed by experiments on real datasets showing the practical benefits of the approach.", "questions": "1. Are there natural extensions of this work to other loss functions beyond least squares?\n2. Do you expect the techniques developed here to be applicable to other sketching problems beyond least squares? If so, which ones?\n\nMinor Comments:\n\n1. The abstract could more clearly state the key theoretical results/bounds achieved\n2. Some additional discussion of practical implications and potential applications would be valuable\n3. A few typos noted (e.g. line 330)", "limitations": "yes", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1. Provides a sparse sketching method that achieves near-unbiased least squares estimation in optimal space and current matrix multiplication time.\n2. Achieves O(d^2 log(nd)) bits of space, which is optimal. Matches current matrix multiplication time O(d^ω), improving over previous approaches.\n3. Introduces a new bias analysis for sketched least squares that sharply characterizes dependence on sketch sparsity. The techniques developed may be applicable to other sketching and randomized linear algebra problems.", "weaknesses": "1. Experiments are conducted on only a few datasets. And it does not explore a wide range of problem sizes or distributed computing scenarios.\n2. Primarily focused on least squares regression, with limited discussion of extensions to other problems.\n3. Some of the theoretical results rely on assumptions (e.g., about leverage score approximation) that may not always hold in practice.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for comments and feedback. We will revise the abstract, and also expand our discussion of applications beyond least squares as outlined below.\n- **Extensions to other loss functions beyond least squares.** In addition to least squares, we provide a more broadly applicable result in Theorem 3, which in particular applies to minimizing general convex loss functions via the framework of Distributed Newton Sketch (see Corollary 1). Thus, our sketching methods and proof techniques are of broader interest to sketching-based optimization algorithms for a variety of loss functions, including the logistic loss, other generalized linear model losses such as the hinge loss, etc. \n- **Applications to other sketching problems.** Our theoretical analysis, which includes several new ideas (including the higher-moment version of Restricted Bai-Silverstein inequality, as well as the careful use of H\"older's inequality in the analysis) is relevant to many instances of RMT-style analysis (i.e., analysis relying on the Stieltjes transform of the resolvent matrix) for sparse sketching operators. This RMT-style analysis has been used in Newton Sketch [16], Randomized SVD [17], and Sketch-and-Project [21]. We chose to focus on distributed least squares, as this gives the clearest worst-case computational improvements.\n- **Distributed computing scenarios.**  In fact, our methods and results are applicable much more generally than the single-server multiple-machine computation model used in the paper. For example, they can be naturally extended to the multiple-server model [43], where the data is randomly partitioned into multiple chunks stored on separate servers, which is common in the literature (see response to Reviewer h7iG for details). The main reason we focused the paper on the single-server multiple-machine model is because this allowed us to obtain worst-case results that are independent of condition-number type quantities (those are unavoidable in the multiple-server model, given our other computational constraints).\n- **Experiments.** Our main contribution is to provide the first sharp theoretical characterization of the least squares bias for extremely sparse sketches, and this is also where we focused our experiments. Nevertheless, we include additional experiments (see the general response and PDF) with other sketching methods (all of which are more computationally expensive than the fast sparse LESSUniform method we used in the paper): Leverage Score Sampling, Gaussian and Subgaussian sketches, as well as Subsampled Randomized Hadamard Transform.\n- **Assumptions about leverage score approximation.** Our main results, Theorems 2 and 3, do not require any assumptions related to leverage score approximation, as they include leverage score approximation as part of the algorithmic procedure.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "iDK0JGnVCP", "reviewer": "Reviewer_RJVb", "summary": "The paper studies the least squares regression task and improves the space and communication amount in a distributed setting to be independent of $\\epsilon$. The key to achieve this is by sketching the data in blocks that have an $\\epsilon$ dependence which can be reduced to only d-dependencies by aggregating them. For the covariance A^TA this becomes only a d x d matrix and for the other required term A^Tb it can be handled by communicating and aggregating just the solutions which have a dimension of d. Although the idea seems very simple described this way, the analysis seems highly non-trivial and has very interesting aspects and novel techniques (or at least new to me). It analyzes in this setting not only the standard least squares error but also a bias term, which allows for smaller sketch dimensions, and only 1/sqrt(eps) dependencies.", "questions": "* I am slightly confused by the definition of the bias. Is $||AE(\\tilde x)-b||$ standard in some literature? I would say this is the variance of the expected sketched estimator. Why would I be interested in this expected estimator, instead of the actual outcome after sketching?\n* I would rather define the bias as $||\\tilde x -x^*||$ which I think does not allow for any improvements over standard sketching results, right?", "limitations": "f", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "* gives least squares sketching results with lower time, space, and communication complexities\n* interesting techniques\n* some further applications are given, though details are fully in the appendix\n* very good writing", "weaknesses": "* motivation of the model, seems to be a niche\n* very low improvements but there seems no big gap that can be leveraged", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the positive feedback, as well as the comments and questions. We will address them in the final version.\n- **Definition of bias.** At a high level, we rely on the statistical definition of the bias of an estimator, which is (informally) the difference between the expectation of that estimator and its target quantity (the estimand). In the case of a least squares estimator $\\tilde x$, the most standard statistical notion of bias would be $E[\\tilde x]-x^*$. Since the estimator is multivariate, one will typically measure the amount of bias by taking the norm, i.e., $||E[\\tilde x]-x^*||$ (as opposed to $||\\tilde x-x^*||$ which is typically referred to as the estimation error). In our case, since we are interested in the least squares prediction vector, i.e. $A\\tilde x$, as an estimator of the vector $b$, we compute the bias as $||E[A\\tilde x]-b||=||A E[\\tilde x]-b||$. This turns out to be the exact right notion for the purpose of distributed averaging and bounding the least squares loss.\n- **Motivation of the model.** Our model is motivated by the general distributed averaging framework (also known as model averaging, or bagging), which is widely used in many settings beyond least squares. In fact, our methods and results are applicable much more generally than the single-server multiple-machine computation model used in the paper. For example, they can be naturally extended to the multiple-server model [43], where the data is randomly partitioned into multiple chunks stored on separate servers, which is common in the literature (see response to Reviewer h7iG for details). The main reason we focused the paper on the single-server multiple-machine model is because this allowed us to obtain worst-case results that are independent of condition-number type quantities (those are unavoidable in the multiple-server model, given our other computational constraints).", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "5nPpzZxb98", "reviewer": "Reviewer_X2Ks", "summary": "Sketched least squares involve estimating the term $(X^TX)^{-1}$ which has a high bias when the sketch matrix $S$ is not sub-Gaussian. This paper gives a sparse sketching method using a LESS embedding which runs in optimal space and current matrix multiplication time, where $S$ is sparse, and constructed based on the leverage scores of the data matrix $A$ (Definition 2 $(s,\\beta_1,\\beta_2)$-LESS embedding). The paper also improves the sharpness of the probability bounds which is applicable to similar problems in RMT using LESS embeddings. For $s=1$, nothing is different, but when $s > 1$, the bias bound is reduced.", "questions": "1. It would be nice to have experiments with other sketching methods to (empirically) justify some statements on the bias, variance (although not stated) and computational time, e.g. estimators mentioned in Table 1. For example, are there sketching estimators with higher bias, but less variance? \n\n2. Despite being clear to read, I had to go back and forth a bit to find out what the novelty is. I appreciate the clarity and thoroughness of explaining the bounds on the bias and variance, runtime, but I would like it if Definition 1 & 2 came much earlier (or at least maybe an informal Definition 2 after line 85?)", "limitations": "N/A", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 4, "strengths": "Inversion bias for estimates of $(X^TX)^{-1}$ is a challenging problem in least-squares sketching. Sub-gaussian sketches have high computational cost but low bias ; other sketches have low computational cost but high bias. The sketch proposed in this paper minimizes the bias of the estimator only requires 1 parallel pass over the data, and has runtime of nnz($A$) + $\\tilde{O}(d^2/\\epsilon)$ which is much faster. Moreover, computing the estimate can be done in parallel, with only the final result needing to be averaged. \n\nThe authors rigorously prove that their sketch fulfills the above criteria by using techniques from random matrix theory, and this is not trivial at all. I find the key contribution in the paper is proving these statements, and the authors did an excellent job here. The paper is extremely well written and easy to follow. Of particular note is the bias analysis for the least squares estimator (Section 4) which is concisely written, where the main ideas (proof sketch) given in the main paper, and technical details in the appendix. The ideas given in the proof sketch can be used for similar problems, and are straightforward to understand. Moreover, the technical details come with sufficient exposition such that it is straightforward for a reader to understand the direction the proof is going (which is certainly much appreciated).", "weaknesses": "1. It would be nice to have experiments with other sketching methods to (empirically) justify some statements on the bias, variance (although not stated) and computational time, e.g. estimators mentioned in Table 1. For example, are there sketching estimators with higher bias, but less variance? \n\nThere are some minor typos, e.g. line 212 reference missing, lines 282 to Equation (2) at bottom of the page is missing a bracket for the numerous expectations (in contrast to lines 542 onwards in the appendix).\n\nIn Appendix A, notation for concentration inequalities should be looked at and made consistent, e.g. Lemma 6 / H{\\\"o}lder's inequality should have a $\\frac{1}{q}$, Lemma 7 ($\\lambda$ max isn't defined), Lemma 8 / Azuma's inequality ($\\lambda$, $m$ should be consistent), Lemma 10 is missing a bracket for $\\mathbb E[x_i^2]$. \n\nThe presentation of the proof for Theorem 5 was slightly jarring (due to Lemma 11, Lemma 12 appearing in the proof), but there also doesn't seem to be a good way to include them (since referring to the two lemmas requires the upper bounds, and flipping back a page is also inconvenient). Maybe a solution is to indent the Lemmas, or box them up?\n\nThe dot before line 626 (after 72) should be removed.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for a careful read and detailed comments. We will address all of the typos and presentation suggestions in the final version.\n- **Additional experiments.** We include additional experiments with other sketching methods (all of which are more computationally expensive than the fast sparse LESSUniform method we used in the paper): Leverage Score Sampling, Gaussian and Subgaussian sketches, as well as Subsampled Randomized Hadamard Transform (see our general response and the PDF). We note that the fact that sketching methods tend to yield smaller least squares bias than uniform subsampling has been empirically observed in prior works [42], which is why we did not focus on this here. Our main contribution is to provide the first sharp theoretical characterization of this phenomenon for extremely sparse sketches. \n- **Are there sketching estimators with higher bias, but less variance?** By definition, the bias has to be no larger than the variance, but there can be cases where the two quantities are comparable in size (which implies that distributed averaging will not be effective). The main example here is i.i.d. sub-sampling, including leverage score sampling. There are theoretical lower bounds [18] which show that (in some cases) the bias of leverage score sampling may not be much smaller than the variance.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "NNOycvxKIr", "reviewer": "Reviewer_h7iG", "summary": "Sketching is a technique from randomized numerical linear algebra which compresses an input matrix $A$ by multiplication with a random matrix $S$. Recent works [18] have shown how to characterize the bias in the sketched inverse covariance $(\\tilde A^\\top\\tilde A)^{-1}$, where $\\tilde A = SA$ for an input matrix $A$ and a sketching matrix $S$. This work extends this machinery to construct approximately unbiased estimators $\\tilde x$ for the least squares regression problem $\\min_x \\|Ax-b\\|$. As the authors show in Theorem 2, this gives an efficient distributed algorithm by averaging the unbiased estimator over $q = O(1/\\epsilon)$ servers. The techniques introduced by the authors also improve prior results on bias-free estimation of the inverse covariance.", "questions": "- The authors may consider discussion the related work https://arxiv.org/abs/2203.09755 on distributed least squares and https://arxiv.org/abs/2007.01327 on using debiasing approaches\n\nOther comments\n\n- Lemma 2 should have a probability statement (e.g., with probability at least $1-\\delta$)", "limitations": "The authors claim that the work does not have limitations in the checklist. Some of the points I raised in the weaknesses section could be discussed in the work, such as restrictions on the computational model and the number of servers.", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The authors introduce a new way to exploit unbiased estimators, which appears to be different from prior works, which incorporate debiasing in a second order optimization framework (e.g. https://arxiv.org/abs/2007.01327). The new technique simply averages local solutions, and conceptually seems simpler than prior techniques.\n\nTechnically, the authors improve prior results on bias-free estimation of the inverse covariance by showing an improved moment bound on a certain random variable by carefully exploiting better bounds for smaller moments via a Holder’s inequality. Other new tools from random matrix theory are additionally needed to handle this change. These ideas, as well as their improvements to bias-free estimation of the inverse covariance, are interesting and may be useful for future results in this area.", "weaknesses": "The result and techniques seem to be a bit niche and incremental. The notion of distributed computation is different from ones that I see in other works in the sketching literature (e.g. https://arxiv.org/abs/1408.5823, https://arxiv.org/abs/1504.06729) where the input is partitioned across multiple servers. Instead, this work considers a setting where there is just one stream which contains the entire input, and the central server has access to $q$ machines that can each access this stream (see the Computational model section), which seems nonstandard. The number of servers needed is also rather large (e.g. $q = O(1/\\epsilon)$), which also seems restrictive. I encourage the authors to include references that prove results in this setting if there are any others. The main technical novelty (line 275) also seems to be an improvement to a bound in the analysis of LESS embeddings from prior work, while the overall sketching framework is largely unchanged.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for the comments and questions. We are glad that the reviewer appreciates the simplicity of our algorithmic approach, and that the reviewer finds our technical contributions interesting and useful. Below, we provide clarifications on the reviewer's remarks, including on the distributed computation framework and the number of machines, as well as providing a comparison with the four references provided by the reviewer. We will be sure to include all of those details in the final version. If you think the responses adequately address your concerns, we encourage you to consider increasing your score.\n- **Novelty in our techniques.** Our main contribution is in the theoretical analysis, which includes several new ideas (including the higher-moments version of Restricted Bai-Silverstein inequality, as well as the careful use of H\\\"older's inequality in the analysis of the dominant term). In fact, these ideas are not niche, since they are relevant to many instances of RMT-style analysis (i.e., analysis relying on the Stieltjes transform of the resolvent matrix) for sparse sketching operators, which has been used for Newton Sketch [16], Randomized SVD [17], and Sketch-and-Project [21]. We chose to focus on least squares, as this gives the clearest computational improvements.\n- **Number of machines required.** The use of $O(1/\\epsilon)$ machines in our main results comes from the small space complexity, namely $O(d^2\\log nd)$ bits per machine, that we impose in our setting. However, if we relax this constraint, then there is a natural trade-off between the space complexity and number of machines required by our methods: If we allow $((1/\\theta)d^2\\log nd)$ bits per machine for some $\\theta\\in [\\epsilon,1]$, then we only need $O(\\theta/\\epsilon)$ machines. Setting $\\theta=\\epsilon$, we recover the standard sketched least squares on a single machine, while $\\theta=1$ recovers our main results, but of course we can freely interpolate between those two extremes.\n- **Comparison with [BP22, arXiv:2203.09755].** This work indeed considers essentially the same problem setting as we do. In fact their Algorithm 1 is the same basic distributed averaging procedure that we use, and their computational model described in Remark 2.3 (Option 2) matches our single-server multiple-machines setup (except that they allow random access instead of streaming access to the data). Our results can be viewed as a direct improvement over their guarantees for Gaussian sketches (e.g., their Theorem 2.2), since we obtain analogous guarantees for extremely sparse sketches which are far more computationally efficient. In fact, we already mention these Gaussian sketch guarantees as a baseline for our work (see Table 1 for a comparison), and we will certainly add the reference to this work in that context.\n-  **Comparison with [DBPM20, arXiv:2007.01327].** This work also considers a very similar problem setting to us, except they focus on regularized least squares, and how the regularization parameter affects the bias of the sketched least squares estimator, with distributed averaging as a main motivation (their main results also assume that all workers have access to the centralized data for the purpose of sketching). Similarly as for the earlier reference, their theoretical results require expensive sketching methods (called \"surrogate sketches\"), which are based on Determinantal Point Processes. We use these sketches as one of our theoretical baselines (see Table 1 for a comparison), and we will certainly add this work as a reference. It is worth noting that our random matrix theory techniques can likely be used to extend the regularization-based debiasing techniques developed in that work to fast sparse sketching. We leave this as a promising direction for future work.\n- **Comparison with [BKLW14, arXiv:1408.5823] and [BWZ16, arXiv:1504.06729].** The distributed computation framework considered in those works partitions the data across multiple servers. Importantly, note that these works solve the task of Principal Component Analysis, which is different from our task of least squares regression, so these works are not directly comparable. Nevertheless, our methods and results can be naturally extended to the multiple-server setting, as we discuss below. The main reason we focused the paper on a more centralized computational model, with one server and multiple machines, is because this allowed us to obtain worst-case results which are independent of condition-number type quantities (such as the data coherence defined below).\n- **Results for data partitioned into multiple servers.** As mentioned above, our results can be extended to the setting where a dataset is partitioned into multiple servers, so that each one of the $q$ machines is accessing a separate chunk of the data. Suppose that a dataset $A,b$ of size $N\\times d$ is partitioned uniformly at random into smaller chunks of size $n=N/q$, and each machine constructs an estimate $\\tilde x_i$ based on a sketch of its own chunk. Then with the same computational guarantees on each machine as in Theorem 2, the averaged estimator $\\tilde x=\\frac1q\\sum_i\\tilde x_i$ will with high probability enjoy a guarantee of $||A\\tilde x-b||\\leq (1+\\epsilon + \\tilde O(\\mu/n))||Ax^*-b||$, where $\\mu=N\\max l_i(A)$ is the coherence of the dataset (here, $l_i$ denotes the $i$th leverage score). The additional error term $\\tilde O(\\mu/n)$, which arises from the partitioning (e.g., see [43]), is often negligible compared to the sketching error $\\epsilon$ when the chunk size $n$ is sufficiently larger than the sketch size. An analogous result can be obtained for Distributed Newton Sketch. We will add these claims to the final version.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:
Thanks to all reviewers for the positive feedback and comments. We responded to those comments in the individual responses to each reviewer. We also provided additional experimental results on four different sketching methods (included in the PDF), and discussed the implications of our theoretical results beyond least squares and in other distributed models. All of this will be included in the final version of the paper, alongside other reviewer suggestions. Here, we summarize the main takeaways.
- **Theoretical implications beyond least squares.** Our main contribution is a set of new theoretical techniques (e.g., higher-moment Restricted Bai-Silverstein, Lemma 3) for analyzing sparse sketching methods, which goes far beyond least squares, as we showed in Theorem 3 and Corollary 1, with an application to general optimization over convex losses (e.g., logistic loss, GLMs, etc) via a variant of the Distributed Newton Sketch. These techniques have wide implications for the analysis of other sketching problems where Restricted Bai-Silverstein-type inequalities have been used, including Randomized SVD [17] and Sketch-and-Project [21].
- **Extensions to other computation models.**  In fact, our methods and results are applicable much more generally than the single-server multiple-machine computation model used in the paper. For example, they can be naturally extended to the multiple-server model [43], where the data is randomly partitioned into multiple chunks stored on separate servers, which is common in the literature (see response to Reviewer h7iG for details). The main reason we focused the paper on the single-server multiple-machine model is because this allowed us to obtain worst-case results that are independent of condition-number type quantities (those are unavoidable in the multiple-server model, given our other computational constraints).
- **Additional experiments.** The fact that sketching methods tend to yield smaller least squares bias than uniform subsampling has been empirically observed in prior works [42]. Our main contribution is to provide the first sharp theoretical characterization of this phenomenon for extremely sparse sketches, and this is also where we focused our experiments. Nevertheless, we include additional experiments (see the PDF) with other sketching methods (all of which are more computationally expensive than the fast sparse LESSUniform method we used in the paper): Leverage Score Sampling, Gaussian and Subgaussian sketches, as well as Subsampled Randomized Hadamard Transform (SRHT). These numerical results further support the claim that sketching enjoys small least squares bias.

Meta-review: This paper proves an interesting headline theoretical result on faster small space distribution regression, requires non-trivial analysis, and also provides nice initial experiments on the method. There are some limitations — for example, that the distributed regression setup the authors consider is fairly niche. Nevertheless, we recommend accepting this paper to NeurIPS.
Decision: Accept (poster)