Stochastic Newton Proximal Extragradient Method

Q1 Running the accelerated version in the experiment.

A1 We added numerical results for the original SNPE with the extragradient step and compared it with the variant without it (Figure 3 in the shared PDF). We observe that the modified variant outperforms the original, suggesting the extragradient step may not be beneficial. However, the original SNPE still outperforms the stochastic Newton method in [1].

Q2 No guarantees for convex.

A2 Due to space limitations, please see our response to Q1 for Reviewer qDEj.

Q3 The complexity of the proposed method with exact Hessians.

A3 With exact Hessians, our method achieves a complexity of $O((\frac{L_2D}{\mu})^{2/3}+\log\log(1/\epsilon))$, though this requires a different analysis (we have an unpublished note). This matches the state-of-the-art complexity for second-order methods in strongly-convex strongly-concave min-max problems [R1]. However, for minimization problems, using Nesterov's acceleration achieves a better complexity of $O((\frac{L_2D}{\mu})^{2/7}+\log\log(1/\epsilon))$ [R2].

[R1] Jiang, R. and Mokhtari, A. Generalized optimistic methods for convex-concave saddle point problems. (2022).

[R2] Arjevani, Y. et al. Oracle complexity of second-order methods for smooth convex optimization. (2019).

Q4 Choosing $x_{k+1} = \hat{x}_k$ is not covered by the theory.

A4 Good point. Our current theory only applies to SNPE in Algorithm 1 and cannot be easily extended to the modified version. However, drawing an analogy from first-order methods, it is reasonable to believe that both versions would share similar convergence guarantees. Specifically, the first-order instantiation of the HPE framework is the extragradient method, where the first step is $\hat{x}_t=x_t-\eta_t\nabla f(x_t)$. Since both extragradient and gradient descent achieve the same linear rate, it is plausible that the same analogy holds for SNPE and its modification. That said, we do not have concrete proof and leave this for future work.

Q5 Dffierence with stochastic Newton when $x_{k+1} = \hat{x}_k$?

A5 Our modified scheme has a different update rule from stochastic Newton, in addition to the difference in step size selection. Our modified scheme is $x\_{t+1}=x\_t-\eta\_t(I+\eta_t\tilde{H}\_t)^{-1}\nabla f(x_t)$ as shown in Eq. (6). In contrast, stochastic Newton follows $x_{t+1}=x_t-\lambda_t\tilde{H}_t^{-1}\nabla f(x_t)$. Hence, the update directions are different, leading to distinct trajectories.

Q6 Using Pearlmutter's implicit Hessian-vector product for inexact second-order updates?

A6 Good point. First, Pearlmutter's technique is not always applicable as it requires access to the computational graph of the objective. Moreover, in the Hessian-free methods, the number of Hessian-vector products per iteration can be substantial to maintain the superior convergence rate of the second-order update, often scaling with the (square root of the) condition number of the Hessian and the desired accuracy. Thus, the stochastic Hessian approximation can sometimes be more efficient for leveraging second-order information.

Q7 Relaxing the requirement of exact gradients.

A7 Extending our analysis to noisy gradients is challenging since we cannot reliably check Condition (3), which is key in our proof. This limitation arises from the HPE framework and is not specific to our algorithm. Moreover, unless we impose strong assumptions on the gradient noise, it is unlikely to achieve superlinear convergence, which is our focus.

Q8 The choices of hyperparameters.

A8 For our Algorithm 1, the hyperparameters are the line-search parameters $\alpha,\beta\in(0,1)$ and the initial step size $\sigma_0$. We chose $\alpha = 1,\beta = 0.5$, and $\sigma_0=1$ without optimizing their choices. For the stochastic Newton method, we followed the default line search parameters in the official GitHub implementation. Both algorithms are relatively robust to the choice of these hyerparameters. We will highlight this point in the revision.

Q9 The stochastic version v.s. an exact Hessian?

A9 We included the numerical results for NPE (our method with exact Hessian) in the shared PDF; see Figures 1 and 2. Figure 1 shows that NPE achieves a faster superlinear convergence rate and converges in fewer iterations due to the use of the exact Hessians. However, our SNPE method also performs comparably to NPE, demonstrating the effectiveness of the Hessian averaging scheme. Moreover, in terms of runtime, SNPE outperforms NPE due to its lower per-iteration cost.

Q10 Logarithmic factors.

A10 Eliminating the logarithmic factors seems difficult. These factors originate from Lemma 3, where we bound the average Hessian noise. The logarithmic dependence on $t$ arises from the union bound, while the dependence on $d$ is due to the matrix concentration inequality (Theorem 3 in [1]). Hence, these logarithmic factors cannot be eliminated without additional assumptions on the stochastic noise.

Q11 Figure 1: is it a labeling mistake?

A11 It is not a labeling mistake; in this experiment, SNPE-UnifAvg is indeed better than SNPE-WeightAvg. This might be due to the small subsampling size, making the stochastic Hessian noise the limiting factor. In our new experiment in the shared PDF, we set the subsampling size to match the dimension, and we observed that the weighted averaging outperforms the uniform averaging in all cases.

Q12 An addition $\Upsilon$ dependency in stochastic NPE?

A12 The transition points and superlinear rates of stochastic Newton also depend on $\Upsilon$. This dependence is not explicit in Table 1 because we focus on the setting where $\Upsilon = O(\kappa)$, which typically holds in practice. In this case, $O(\kappa^2 + \Upsilon^2) = O(\kappa^2)$, simplifying the expressions in Table 1. We will add a remark in the revision to clarity this.

Q1 The assumptions are very restrictive. e.g., Assumptions 3 and 5.

A1 We note that our assumptions are standard in the study of (stochastic) second-order methods, including Subsampled Newton (Erdogdu & Montanari, 2015; Roosta-Khorasani & Mahoney, 2019), Newton Sketch (Pilanci & Wainwright, 2017; Agarwal et al., 2017), and notably, the recent work by Na et al. (2022). In particular, Assumption 3 (Lipschitz Hessians) provides the necessary regularity condition for achieving superlinear convergence and is commonly used in the study of second-order methods, e.g., in Section 9.5.3 of the textbook by Boyd & Vandenberghe (2004). For instance, it is satisfied by the regularized log-sum-exp function and the loss function of regularized logistic regression. Moreover, we also discussed in Section 2 when Assumption 5 holds for stochastic Hessian approximation. Specifically, for Hessian subsampling, this is satisfied when each component function $f_i$ is convex, while it is automatically satisfied for Hessian sketching.

Q2 How do parts (c) and (d) of Theorem 1 prove superlinear convergence? I expected to see $\rho^{2^t}$ where $\rho<1$.

A2 It appears that the reviewer confuses “superlinear convergence” with “quadratic convergence”. Specifically, the rate of $\rho^{2^t}$ mentioned by the reviewer is “quadratic convergence”, which is a special case of, but not equivalent to, “superlinear convergence”. In the optimization literature, the convergence is said to be Q-superlinear if $\lim_{t\rightarrow\infty}\frac{\\|x_{t+1}-x^\*\\| }{\\|x_t-x^\*\\|}=0$ (see, e.g., Appendix A.2 of Nocedal & Wright (2006)). Note that in Theorem 1 (c) and (d), we showed $\\|x_{t+1}-x^\*\\|=O(\frac{1}{t})\\|x_t-x^\*\\|$ and $\\|x_{t+1}-x^\*\\|=O(\frac{1}{\sqrt{t}})\\|x_t-x^\*\\|$, respectively. Since $\lim_{t \rightarrow \infty}\frac{1}{t}=\lim_{t\rightarrow\infty}\frac{1}{\sqrt{t}}=0$, these are superlinear convergence results by definition.

Q3 In Table 2, the iteration complexity of SNPE is $\log(\epsilon^{-1})$ similar to AGD.

A3 Thank you for raising this point. Since our proposed SNPE method achieves superlinear convergence, it has a strictly better dependence on $\epsilon$ compared to AGD. In fact, in the iteration complexity of SNPE presented in Table 2, the dependence on $\epsilon$ can be replaced by $\frac{\log(\epsilon^{-1})}{\log(\log(\epsilon^{-1}))}$, which is provably better than the complexity of AGD by at least a factor of $\log\log(\epsilon^{-1})$. We note that similar superlinear convergence rates have also been established in the prior work on stochastic Newton methods (Na et al., 2022) and in the literature on quasi-Newton methods (Rodomanov & Nesterov, 2022; Jin & Mokhtari, 2023). In our submission, we chose to use the simpler expression $O(\log(\epsilon^{-1}))$ to save space. However, we will update the table and include a discussion in the revision to more accurately reflect our superlinear convergence rate.

Q4 In lines 342-344, the dominating term in the complexity is the one that depends on $\epsilon$. The Damped Newton has better dependence than SNPE.

A4 We agree with the reviewer that damped Newton's method has a better dependence on $\epsilon$ than SNPE, which we also mentioned in lines 342-344. However, it is important to note that the damped Newton's method is deterministic and requires computing the exact Hessian, resulting in a per-iteration cost of $O(nd^2)$. In contrast, our method only requires a stochastic Hessian approximation and typically incurs a total per-iteration cost of $O(nd+d^3)$, as discussed in Section 6. This distinction is crucial because it allows us to achieve a better runtime compared to the damped Newton's method, especially when $n\gg d$. Thus, while the iteration complexity of damped Newton's method exhibits a better dependence on $\epsilon$, its overall arithmetic complexity can be much worse than ours. This is also demonstrated in our experiment presented in the shared PDF file. From Figures 1 and 2, we observe that while the damped Newton's method requires fewer iterations to converge, it takes more time overall to achieve the same accuracy as our method.

Q5 Empirical comparison between SNPE and damped Newton in terms of run-time.

A5 Thank you for your suggestion. We have included the additional experiment in the shared PDF file; please see Figure 2. We would like to remark that damped Newton also performs a backtracking line search in every iteration, resulting in an overhead similar to ours. As expected, when the number of samples $n$ and the dimension $d$ are large, our method has a significant runtime gain compared to damped Newton, and the gap widens as the number of samples increases.

Q6 The authors ignore the complexity of line-search while computing the complexity of SNPE.

A6 We respectfully disagree with the reviewer. As outlined in Remark 3 and proven in Appendix A.4, we explicitly characterize the cost of the line search in our SNPE method. Specifically, after $t$ iterations when $t = \tilde{\Omega}(\Upsilon^2/\kappa^2)$, the total number of line search steps can be bounded by $2t + \log(3M_1\sigma_0 /(\alpha \beta))$. Also, note that each line search step requires computing one gradient and one matrix inversion. Consequently, our method requires, on average, a constant number of gradient evaluations and matrix inversions per iteration. This leads to a constant overhead in the complexity bound, which is effectively hidden by the big O notation.

Additional References:

Boyd, S. and Vandenberghe L. Convex Optimization. Cambridge University Press, 2004.

Nocedal, J. and Wright, S. J. Numerical Optimization. Springer, 2006.

Rodomanov, A. and Nesterov, Y. Rates of superlinear convergence for classical quasi-Newton methods. Math. Program., 2022.

Jin, Q. and Mokhtari, A. Non-asymptotic superlinear convergence of standard quasi-Newton methods. Math. Program., 2023.

We thank the reviewer for their positive feedback.

Q1 Empirical comparison with AGD.

A1 Following your suggestion, we compared the performance of our SNPE method against AGD in our new experiment; please see Figures 1 and 2 in the shared PDF file. From Figure 1, we observe that our SNPE method, with either uniform or weighted averaging, requires far fewer iterations to converge than AGD due to the use of second-order information. Consequently, while SNPE has a higher per-iteration cost than AGD, it converges faster overall in terms of runtime, as demonstrated in Figure 2.

Q1 Authors consider only a strongly convex setup and this seems to me as a major limitation of this work.

A1 Thank you for raising this point. To begin with, we note that the strong convexity assumption is common in the study of stochastic second-order methods, including Subsampled Newton (Erdogdu & Montanari, 2015; Roosta-Khorasani & Mahoney, 2019), Newton Sketch (Pilanci & Wainwright, 2017; Agarwal et al., 2017), and notably, the recent work by Na et al. (2022). Hence, by establishing our results within the strongly convex setting, we can better position our contribution in relation to prior work.

Moreover, the focus on strongly convex functions in these works, as well as in ours, stems from the clear advantages that stochastic second-order methods offer over first-order methods, such as gradient descent. Specifically, stochastic second-order methods achieve a superlinear convergence rate, as demonstrated in this paper, which is superior to the linear rate attained by first-order methods. We also note that the assumption of strong convexity is necessary for achieving superlinear convergence rates, as only a sublinear rate can be achieved in the convex setting, even with exact Hessian information.

Finally, while we believe it is possible to extend our techniques to the convex setting, developing the necessary theory and discussing the results would require more space than is available in the submission. Therefore, this extension is beyond the scope of this paper.

Q2 To check condition (3) authors employ line-search, which introduces an additional logarithmic factor in the convergence rate. However, this can be the burden of the HPE framework, it seems not very significant, but still an issue.

A2 The reviewer is correct in noting that we need to employ line search to ensure Condition (3). However, we would like to mention that most stochastic second-order methods require some form of line search to ensure global convergence, and this limitation is not unique to our methods. For instance, Pilanci & Wainwright (2017), Roosta-Khorasani & Mahoney (2019), and Na et al. (2022) all used a backtracking line search to ensure a sufficient decrease condition.

Moreover, we wish to clarify that the line search scheme only introduces an additive logarithmic factor, instead of a multiplicative one, in our final complexity bound. Specifically, as mentioned in Remark 3 and proved in Appendix A.4, after $t$ iterations when $t = \tilde{\Omega}(\Upsilon^2/\kappa^2)$, the total number of line search steps can be bounded by $2t + \log(3M_1 \sigma_0 /(\alpha \beta))$. Also, note that each line search step requires computing one gradient and one matrix inversion. Consequently, our method requires, on average, a constant number of gradient evaluations and matrix inversions per iteration, leading to a constant overhead in the complexity bound.

Q3 Why do you use $\sigma_{t+1} = \eta_t/ \beta$?

A3 This is a good question. Our motivation behind the choice $\sigma_{t+1} = \eta_t/ \beta$ is to allow the step size to grow, which is necessary for achieving a superlinear convergence rate. Specifically, our entire convergence analysis builds on Proposition 1, which demonstrates that $\\|x_{t+1}-x^\*\\|^2 \leq \\|x_{t}-x^\*\\|^2 (1+2\eta_t \mu)^{-1}$. Consequently, we require the step size $\eta_t$ to go to infinity to ensure that $\lim_{t \rightarrow \infty} \frac{\\|x_{t+1}-x^\*\\|}{\\|x_{t}-x^\*\\|} = 0$. Note that this would not be possible if we simply set $\sigma_{t+1} = \eta_t$, since it would automatically result in $\eta_{t+1} \leq \sigma_{t+1} \leq \eta_t$. Moreover, this condition $\sigma_{t+1} = \eta_{t}/\beta$ is explicitly utilized in Lemmas 8 and 16 in the appendix, where we demonstrate that $\eta_t$ can be lower bounded by the minimum of $\sigma_0/\beta^t$ and another term.

We should note that this more aggressive choice of the initial step size at each round could potentially increase the number of backtracking steps. However, as mentioned in the response to Q2 above, this does not cause a significant issue, since the average number of backtracking steps per iteration can be bounded by a constant close to 2.

Thank you again for the question and we will include the discussions above in our revision.

Q4 Minor remarks.

A4 Thank you for catching the typos. We will fix them all in the revision.

References:

Erdogdu, M. A. and Montanari, A. Convergence rates of sub-sampled newton methods. Advances in Neural Information Processing Systems, 2015.

Roosta-Khorasani, F. and Mahoney, M. W. Sub-sampled Newton methods. Mathematical Programming, 2019.

Pilanci, M. and Wainwright, M. J. Newton sketch: A near linear-time optimization algorithm with linear-quadratic convergence. SIAM Journal on Optimization, 2017.

Agarwal, N., Bullins, B., and Hazan, E. Second-order stochastic optimization for machine learning in linear time. The Journal of Machine Learning Research, 2017.

Na, S., Derezinski, M., and Mahoney, M. W. Hessian averaging in stochastic Newton methods achieves superlinear convergence. Mathematical Programming, 2022.