Fundamental Bias in Inverting Random Sampling Matrices with Application to Sub-sampled Newton

We thank the reviewer for his/her constructive and thoughtful comments. Please find our point-by-point responses below.

Methods And Evaluation Criteria:

The datasets are standard benchmark datasets (MNIST and CIFAR-10).
I do have a question about lines 400-402: "Figures 2 and 3 do not include the time for input data pre-processing......" What is the relative magnitude of the time for input data pre-processing? That is to say, if the relative magnitudes of pre-processing is greater than the reduction in wall clock time for Newton-LESS / SSN-ARLev, then I don't think wall-clock time would be a meaningful metric.

Reply: We thank the reviewer for the constructive question (which is in fact shared by Reviewer GJpF). It is worth noting that various second-order baselines in Figures 2 and 3, including Newton-LESS, ARLev, ALev, SLev, and SRHT, all require preprocessing of the data.

Specifically, Newton-LESS requires a Hadamard transform preprocessing with a complexity of $O(nd \log n)$, see [Derezinski2021a]. To ensure a fair comparison, we treat SRHT as a uniform sampling method that undergoes the same random Hadamard transform preprocessing as Newton-LESS. The leverage-based sampling methods, ARLev, ALev, and SLev, require computing approximate leverage scores before constructing the sketches $\check{\mathbf{S}}\mathbf{A}$ via random sampling, with complexity $O(\text{nnz}(\mathbf{A}) \log n + d^3 \log d)$.

As such, the preprocessing time are in fact comparable across these methods, and it makes sense to exclude the time for the Hadamard transform and leverage score computation ensures a fair comparison. Considering the computational overhead associated with the sketching process itself (rather than the preprocessing steps) allows for a clearer illustration of the "convergence-complexity" trade-off across different stochastic optimization methods.

We have reproduced Figures 1 and 2 by adding (for all methods) the preprocessing time, and observed a similar trends. In particular, SSN-ARLev achieves a significantly lower overall runtime compared to Newton-LESS and establishes better complexity-convergence trade-off than other baselines. These additional numerical results and discussions will be included in the revision.

Additionally, for certain random Newton-type solvers, preprocessing is done only once at the beginning of the iteration and then reused in subsequent iterations. For instance, in [Li2024], distributed least squares regression based on federated Newton iterations only performs data preprocessing (such as the Hadamard transform and computing leverage scores) once, and then reuse the obtained results alongside the optimization procedure. In such cases, the preprocessing time becomes negligible compared to the iteration computation and communication costs. Thus, we believe wall-clock time remains a meaningful metric.

[Derezinski2021a] Derezinski M, Lacotte J, Pilanci M, et al. Newton-LESS: Sparsification without trade-offs for the sketched Newton update. Advances in Neural Information Processing Systems, 2021, 34: 2835-2847.

[Li2024] Li J, Liu Y, Wang W. Fedns: A fast sketching newton-type algorithm for federated learning. Proceedings of the AAAI Conference on Artificial Intelligence. 2024, 38(12): 13509-13517.

Experimental Designs Or Analyses:

I would be interested to see the boxplots of the relative errors over the 30 independent runs in Fig 2, and 10 independent runs in Fig 3 for both Newton-LESS and SSN-ARLev, as the relative error for Newton-LESS and SSN-ARLev seem close to each other as the sketch size m increases. This is not to say that the experimental design / analyses is flawed, but I'm curious if the relative errors are skewed in one direction or equally symmetric about the solid lines in both Figures.

Reply: We thank the reviewer for the insightful comment. Our theoretical results say that the convergence rate of SSN-ARLev can approximate, but not exceed, the convergence rate of Newton-LESS, see Theorem 4.3. Moreover, as $m$ increases, the convergence rate of the de-biased sub-sampled Newton method approaches that of Newton-LESS. Therefore, it is reasonable to see in the experimental results that the relative errors for both Newton-LESS and SSN-ARLev becoming increasingly similar as $m$ increases.

Our contribution in fact relies in showing that SSN-type methods can indeed achieve (local) problem-independent convergence rates (which, to the best of our knowledge, has never been established for SSN-type methods) that

(1) can be significantly faster than standard SSN; and

(2) matches that of Newton-LESS, but with significantly reduced computation complexity, and thus better complexity-convergence trade-off, as confirmed by Figures 1 and 2.

As suggested by the reviewer, we will include boxplots of the relative errors for both methods also converge as $m$ increases.

We thank the reviewer for his/her constructive and thoughtful comments. Please find our point-by-point responses below.

Theoretical Claims:

I did not have time to check all the proofs. The following are a few questions:

1. Page 4, left column, Line 199: what does it mean to say “subspace embeddings … ensure unbiasedness up to the same level of epsilon but are not guaranteed to remain unbiased for a smaller epsilon”?

Reply: Sorry for this discussion that may be misleading. Here, we meant that sketching matrices satisfying the subspace embedding property, with an error $\epsilon$ (say of order $\epsilon = O\left(\sqrt{\rho_{\max} d_{\mathrm{eff}} \log(d_{\mathrm{eff}}/\delta)/m}\right)$ ), can guarantee the same level $\epsilon$ of the inversion bias, see Lemma B.1 in the supplementary material for a proof. It is, however, unclear whether this (upper bound of) inversion bias is tight and/or can be made smaller with refined analysis. This is in fact the very motivation of the proposed de-biased sampling matrix $\check{S}$, which, under the same sampling size, achieves a smaller inversion bias of order $\epsilon = O\left(\sqrt{\rho_{\max}^3 d_{\mathrm{eff}}^3 / m^3}\right)$, see Proposition 3.2, offering an improvement over the standard subspace embedding-type analysis.

2. Page 21, inequality (a) (in the middle of the page), can you explain how you get the constant factor of 38?

Reply: We thank the reviewer for the constructive comments. Regarding the constant factor of 38 in inequality (a) on Page 21, we will provide a detailed explanation. By combining $H^{1/2} Q_{-s} H^{1/2} \preceq 6 I_d$ and the bound $|| H^{\frac12} \widetilde{H}^{-1} H^{\frac{1}{2}}\ || \leq 2$, along with $||v || = 1$ and $||u || = 1$, we can derive the following two inequalities:

$$\sum^n_{j=1} u^\top \widetilde{H}^{\frac12} Q_{-s} a_{j} a_{j}^\top Q_{-s} \widetilde{H}^{\frac12} u=u^\top \widetilde{H}^{\frac12} Q_{-s} A^\top A Q_{-s} \widetilde{H}^{\frac12} u \leq 36$$

and

$$\sum^n_{j=1} v^\top \widetilde{H}^{-\frac12} a_{j} a_{j}^\top \widetilde{H}^{-\frac{1}{2}} v = v^\top \widetilde{H}^{-\frac12} A^\top A \widetilde{H}^{-\frac12} v \leq 2.$$

Adding these two results gives the constant factor of 38. We will provide these proof details in Section D.2 of the supplementary material in the revised version for clarification.

Experimental Designs Or Analyses:

I do not see issues in general but I do not think one can claim “improved convergence” because the convergence is not improved. The main thing seems to be improving the runtime at the same convergence. The numerical experiments also show this: the convergence rate of SSN-ARLev looks about the same (at least the rate looks similar) as Newton-LESS in Figure 1 on both datasets.

Reply: We greatly appreciate the reviewer's helpful feedback and agree that the claim of "improved convergence" is indeed a bit misleading. As discussed in Section 4 of the submission, the proposed de-biased SSN improves the local convergence rate over standard sub-sampled Newton method, to achieve problem-independent convergence rates (that, to the best of our knowledge, has never been established for SSN). It is worth noting that the theoretical convergence rate of de-biasd SSN can approximately match, but not exceed, that of Newton-LESS (that is a closely-related, but different approach replying on random projection), but will significantly reduced computational complexity.

In particular, the proposed SSN-ARLev (that replies on random sampling) offers significant computational efficiency advantages over Newton-LESS, and thus outperforms Newton-LESS in terms of the complexity-convergence trade-off, see confirmed by Figures 1 and 2.

While the numerical experiments in Figure 1 suggest that SSN-ARLev may offer even improved convergence in certain practical scenarios over the Newton-LESS approach, we are unable yet provide a theoretical proof of such superiority in convergence rate.

Additionally, we believe that combining the proposed de-biased approach with more efficient sub-sampled newton algorithms, such as the adaptive SSN-type algorithm proposed by [Lacotte2021], could further accelerate convergence. However, exploring this is beyond the scope of the present work.

[Lacotte2021] Lacotte J, Wang Y, Pilanci M. Adaptive newton sketch: Linear-time optimization with quadratic convergence and effective hessian dimensionality. International Conference on Machine Learning. PMLR, 2021: 5926-5936.

Other Comments Or Suggestions:

Reply: We thank the reviewer again for the thorough and thoughtful review. All minor issues raised will be addressed in the revised version of the paper.

We thank the reviewer for his/her constructive and thoughtful comments. Please find our point-by-point responses below.

Essential References Not Discussed:

Newton Meets Marchenko-Pastur: Massively Parallel Second-Order Optimization with Hessian Sketching and Debiasing (ICLR 2025, https://arxiv.org/abs/2410.01374)

Optimal Shrinkage for Distributed Second-Order Optimization (ICML 2023, https://arxiv.org/abs/2402.01956)

Reply: We thank the reviewer for pointing to the two helpful references. They will be included and discussed in the revised version of the paper.

Questions For Authors:

1. In Definition 2.3,$ \rho_{\min} \equiv \min_{1\leq i \leq n}\ell_i^{C}/(\pi_i d_{\mathrm{eff}})$, isn't $\pi_i d_{\mathrm{eff}}=\ell_i^{C}$ and thus cancels the numerator? Also, how does this sampling factor related to  $S$, I think it's just how much we weigh each row of $A$ when we sample, but the Definition sounds like $\rho_{\min} $ is also related to  $S$?

Reply: We thank the reviewer for this question. In Definition 2.3, we define $\rho_{\min} \equiv \min_{1 \leq i \leq n} \frac{\ell_i^{\mathbf{C}}}{\pi_i d_{\mathrm{eff}}}$. While for exact (ridge) leverage score (when the random sampling probability $\pi_i$ is chosen according to $\ell_i^{\mathbf{C}}/d_{\mathrm{eff}}$), one has $\pi_i d_{\mathrm{eff}} = \ell_i^{\mathbf{C}}$ and thus cancels the numerator, this is not true for generic random sampling methods.

The parameter $\rho_{\min}$ is introduced to capture the interaction between the choice of sampling probabilities $\pi_i$ and the data structure of $\mathbf{A}$ (e.g., its leverage scores) and how it affects the associated random sampling inversion bias.

2. In Figure 2, the authors present numerical results showing the proposed method is much more efficient compare to Newton-LESS, I feel it would be more helpful if the authors can discuss more about the baseline and what's the main computation overhead there.
3. Also in the experiments, the authors mention they exclude the time for computation of exact or approximate leverage score, does is mean the main workload of computing  $\check S$ is removed? Is there any justification for such exclusion, i.e., is it a fair thing to do?

Reply: We thank the reviewer for the comments. It is worth noting that various second-order baselines in Figures 2 and 3, including Newton-LESS, ARLev, ALev, SLev, and SRHT, all require preprocessing of the data.

Specifically, as performed in [Derezinski2021a], Newton-LESS first performs the Hadamard transform preprocessing on the data, and then generates LESS-uniform sketching matrix $\mathbf{S}$ to get the sketch $\tilde{\mathbf{A}} = \mathbf{S} \mathbf{A}$. The computational costs of these two steps are $O(nd \log n)$ and $O(msd)$, where $s$ represents the sparsity level of (each row in) $\mathbf{S}$. To ensure a fair comparison, we treat SRHT as a uniform sampling method that undergoes the same random Hadamard transform preprocessing as Newton-LESS.

The three leverage-based sampling methods, ARLev, ALev, and SLev, all involve the computation of approximate leverage scores and the sampling process, and has complexities of $O(\text{nnz}(\mathbf{A})\log n + d^3 \log d)$ and $O(m + md)$, respectively.

Thus, the preprocessing time is in fact comparable across all baseline methods, and by excluding the time for the Hadamard transform and leverage score computation, we ensure a fair comparison when evaluating these methods. Considering the computational overhead associated with the sketching process itself (rather than the preprocessing steps) allows for a clearer illustration of the "convergence-complexity" trade-off across different stochastic optimization methods.

We have reproduced Figures 1 and 2 by adding (for all methods) the preprocessing time, and observed a similar trends. In particular, SSN-ARLev achieves a significantly lower overall runtime compared to Newton-LESS and establishes better complexity-convergence trade-off than other baselines.

These additional numerical results and discussions will be included in the revision.

[Derezinski2021a] Derezinski M, Lacotte J, Pilanci M, et al. Newton-LESS: Sparsification without trade-offs for the sketched Newton update. Advances in Neural Information Processing Systems, 2021, 34: 2835-2847.

We thank the reviewer for his/her constructive and thoughtful comments. Please find our point-by-point responses below.

Other Comments Or Suggestions:

1. Are the constants C in the theorems, propositions, and lemmas in Section 2 and 3 the same, or at least on the same order? It would be helpful to clarify this in the paper.

Reply: While the constants $C$ in the theorems, propositions, and lemmas in Sections 2 and 3 are not identical, they are all absolute constants independent of the dimensions $n$, $d_{\mathrm{eff}}$, and $m$.

2. The preliminary in Section 2 seems to be a bit long, especially given the 8-page limit of the main text. Compared to the preliminaries, some proof sketches and intuitions in the appendix seem to be more insightful and may deserve a place in the main text.

Reply: We thank the reviewer for this constructive comment. We agree that the current section of preliminaries are a bit too long and that the proof sketches and intuitions (currently in the appendix) should be moved to the main text to better highlight the connection between RMT and RandNLA. Due to the 8-page limit for the submission, we were unable to include them in the current version.

If the paper gets accepted, we will revise the content accordingly (with the one additional page) in the camera-ready version by incorporating the proof sketches and intuitions from Appendices D.1 and E.1 into the main text, placing them directly after Theorem 3.1 and Proposition 3.2, respectively.

Questions For Authors:

Is the log factors in m in most of the theorems, propositions, and lemmas in Sections 2 and 3 a consequence due to independent sampling of the m rows (uniformly or according to leverage scores)? If so, would it be possible to remove the log factors via dependent sampling, e.g., volume sampling, or some faster alternatives like adaptive leverage score sampling (see, e.g., [1])? It would be insightful to provide some discussions on this, even just as future directions.

[1] Cortinovis, Alice, and Daniel Kressner. "Adaptive randomized pivoting for column subset selection, DEIM, and low-rank approximation." arXiv preprint arXiv:2412.13992 (2024).

Reply: We thank the reviewer for the constructive comment and the reference. We agree with the reviewer's perspective. As shown in Section B of the supplementary materials, the subspace embedding property in Lemma 2.7 relies on the independence of the random row sampling procedure. This implies that the logarithmic factors in $m$ present in most of the theorems, propositions, and lemmas in Sections 2 and 3 arise from the independent sampling of the $m$ rows (including, but not limited to, uniform sampling or leverage score-based sampling).

As such, it is very possible that dependent sampling methods, such as volume sampling or more efficient techniques like adaptive leverage score sampling (see, e.g., [1]), may provide tighter bounds and potentially eliminate the logarithmic factors. For now, we are not able to include these efficient methods in the proposed analysis framework, since our current proof approach strongly relies on the independence of random sampling.

We also agree with the reviewer that these dependent sampling methods could yield better computational efficiency and/or convergence rate for SSN. We leave them for future work. We will include a discussion of these methods in the conclusion section (Section 6) of the revised version, addressing the potential benefits and challenges of using volume sampling, adaptive leverage score sampling, or other dependent sampling approaches.