Faster Algorithms for Structured Linear and Kernel Support Vector Machines

We appreciate the reviewer for their insightful comments, and we would like to address your concerns. Before doing so, we would like to highlight one important contribution that may have been overlooked.

We develop a generic algorithm for convex quadratic programming, with a runtime of $O(n(d+m)^{(\omega+1)/2}\log(1/\epsilon))$, assuming the quadratic objective matrix $Q$ has a rank-$d$ factorization and $m$ is the number of linear constraints. To the best of our knowledge, this is the first algorithm for convex quadratic programming that has linear dependence on $n$ and solves the program to high precision.

As a direct application, we adapt this algorithm to linear SVM: given a size-$n$ dataset in $\mathbb{R}^d$, we can solve the linear SVM in $O(nd^{(\omega+1)/2}\log(1/\epsilon))$ time. Prior fastest algorithms with $\log(1/\epsilon)$ dependence all have super-linear dependence on $n$, making them inefficient for large datasets. Moreover, our approach can be generalized to various kernel SVM schemes where the kernel matrix can be low-rank factored in subquadratic time. While we agree it is important to evaluate the practical performance of the algorithm, we also wish to highlight its theoretical contribution.

Q: The paper states that the dimension $d=O(\log n)$ seems to be very restrictive for real-world data analysis.

A: As proven in our hardness result (Theorem 1.5), there exists some $d=\Omega(\log n)$ such that the corresponding kernel SVM cannot be solved in $O(n^{2-o(1)})$ time unless the Strong Exponential Time Hypothesis is false. In fact, our proof of Theorem 1.5 explicitly constructs an instance of Gaussian kernel SVM with $d = C\log n$ for some large constant $C$, such that solving it efficiently is hard. Essentially, this implies that the best one can hope for in the $d = \Omega(\log n)$ regime is to explicitly form the kernel matrix by filling all $n^2$ entries. Therefore, we respectfully disagree that our assumption is “restrictive,” as it represents the theoretical best-case scenario.

Q: The application to the Gaussian kernel seems narrow. Can it be extended to other kernels? How does the width of the Gaussian kernel affect the runtime?

A: Yes, our technique can indeed be extended to a large family of kernels, as parametrized in Alman, Chu, Schild, and Song (2020). Specifically, their work shows that as long as the kernel function takes the form ${\sf K}(u, v) = f(\|u-v\|_2^2)$ and $f$ can be approximated by a low-degree polynomial of degree $o(\log n)$, one can compute a rank-$n^{o(1)}$ factorization of the kernel matrix in $n^{1+o(1)}$ time. This family includes polynomial kernels and $p$-convergent kernels (Song, Woodruff, Yu, and Zhang, 2021).

Regarding the width parameter of the Gaussian kernel, it directly impacts the squared radius of the dataset. Let $R$ denote the original squared radius, and $\sigma$ the width. Then the new squared radius is $\frac{R^2}{2\sigma^2}$. As long as $\frac{R^2}{2\sigma^2} \leq O(\frac{\log n}{\log\log n})$, our algorithm can be applied.

Q: How do your methods compare to the Nyström approximation, which also provides a low-rank approximation to the kernel matrix?

A: The Nyström method does indeed produce good low-rank approximations to kernel matrices. For instance, Musco and Musco (2017) show that one can construct a rank-$s$ factorization of the kernel matrix using $O(ns)$ kernel function evaluations and an additional $O(ns^{\omega-1})$ time. However, there are two significant challenges in applying the Nyström approximation and then solving the SVM:

1. Rank-Approximation Tradeoff: If one desires a small rank (e.g., $s=n^{o(1)}$), a large regularization parameter $\lambda$ is required, leading to poor approximation.
2. Impact of Regularization: The additional $\lambda I$ term complicates the downstream SVM algorithm, introducing an error term dependent on the $\ell_2$ norm squared of the optimal solution and $\lambda$. This increases runtime by a factor of $\log(\lambda)$, and there is no clear theoretical guideline for selecting $\lambda$.

Q: For large-scale datasets, deep neural networks practically offer better computational complexity. How does your method compare with deep neural networks?

A: Theoretically, deep neural networks are non-convex models, making their convergence analysis challenging. The fastest known algorithms for training over-parameterized networks with $m$ neurons per layer require $m \sim n^4$, and the runtime is $O(m^{1.93})$ (Song, Zhang, and Zhang, 2024). Since $m = n^4$, this leads to $O(n^{7.72})$ time—significantly more expensive than our approach, which runs in $O(n^{1+o(1)})$ time. Additionally, kernel SVMs are much simpler and easier to analyze compared to deep neural networks. We believe simpler, convex methods like kernel SVMs are preferable in many settings and designing theoretically sound and efficient algorithms for them is valuable.

We thank the reviewer for their thoughtful comments and insights. Regarding your suggestion to compare with algorithms such as LIBSVM or SVM-Light, it is worth noting that these algorithms are all first-order methods with runtime dependence polynomially on $\epsilon^{-1}$. As we illustrated in the manuscript, this leads to a significant runtime discrepancy even when $\epsilon = 10^{-4}$. Below, we address your questions in detail:

Q: In Section 1, the authors used $B$ as the squared radius, while in Section 2, it was used as $B = Q + tH$. It would be good to make them consistent.

A: Thank you for pointing this out. We have updated the notation for $Q + tH$ to $M$ to avoid confusion and ensure consistency.

Q: It seems the assumption that the squared radius $B = O(\log n / \log\log n)$ is very strong.

A: Indeed, this assumption might appear restrictive at first glance. However, as our lower bound result shows, assuming the Strong Exponential Time Hypothesis, even for $B = \Omega(\log^2 n)$, there is no algorithm that solves the Gaussian kernel SVM in $O(n^{2-o(1)})$ time. This implies that one would need to explicitly form the kernel matrix before running a downstream algorithm to solve the SVM. We conjecture that the lower bound can be strengthened to $B = \Omega(\log n)$, which would reflect a phase transition for fast algorithms in kernel linear algebra, as indicated by the works of Alman, Chu, Schild, and Song (2020) and Aggarwal and Alman (2022).

Q: Section 2 is very heavy in technical details. Could the authors provide some intuitions behind these techniques?

A: Thank you for this suggestion. The current structure of Section 2 highlights the challenges in solving the primal-dual central path equations of QPs efficiently. We further demonstrate that, if the quadratic objective matrix has certain structures, such as low-rank or low-treewidth, it is possible to overcome these obstacles and design efficient algorithms. These techniques rely on crafting nontrivial maintenance data structures and performing robustness analyses to ensure that the approaches converge to a solution.

We believe that Section 1 provides a comprehensive summary of our results, including comparisons with other works and potential applications. Additionally, Section 2.1 offers a high-level overview of the generic framework for solving the program, while Section 2.4 discusses how these techniques extend to the kernel SVM setting. To address your feedback, we have added a paragraph recommending readers to skip Sections 2.2 and 2.3 on their first read, as they are highly technical. We hope this adjustment will make the manuscript more accessible without sacrificing the necessary details.

Q: There are typos such as on line 1540, where “memebers” should be “members.”

A: We appreciate you pointing this out. We have corrected the typos and conducted a thorough proofread of the manuscript to ensure it is free of errors.

References

• Alman, Chu, Schild, and Song (2020): Algorithms and Hardness for Linear Algebra on Geometric Graphs. FOCS 2020.
• Aggarwal and Alman (2022): Optimal-degree polynomial approximations for exponentials and Gaussian kernel density estimation. CCC 2022.

We thank the reviewer for their insightful comments and for pointing out the reference that establishes an $O(\epsilon^{-1})$ iteration complexity for first-order methods. We have incorporated a discussion of this work into our updated manuscript, with the changes highlighted in red. Regarding your questions:

Q: The lower bound for SVM also gives a lower bound for QP. How does this compare to other existing lower bounds for QPs?

A: We are not aware of any direct lower bounds specifically for QPs. However, in addition to the work of Backurs, Indyk, and Schmidt (2017) and our own, there are lower bounds established for linear programs, such as those in Kyng, Wang, and Zhang (2020) and Ding, Kyng, and Zhang (2022). These works show that any linear program with polynomially bounded integer entries can be converted into a 2-commodity flow instance with $O(\mathrm{nnz}(A))$ edges. Consequently, for any $a > 1$, if one can solve an integer-capacitated 2-commodity flow in $O(|E|^a \mathrm{poly} \log(\epsilon^{-1}))$ time, this would imply the existence of an LP solver with a runtime of $O(\mathrm{nnz}(A)^a \mathrm{poly} \log(\epsilon^{-1}))$. Since LPs are a subclass of QPs, this provides a type of lower bound for QPs, albeit without directly addressing the quadratic objective matrix $Q$.

Q: Is there a reference for the result on $T$ in line 291?

A: Yes, this result is standard in the primal-dual IPM literature. See, for instance, Ye (2020), Dong, Lee, and Song (2021), Gu and Song (2022), and it is also implicit in Cohen, Lee, and Song (2019).

Q: What is the graph structure of the bags mentioned in line 310?

A: A "bag" is simply a collection of vertices satisfying the following properties:

1. The union of all the bags equals $V$.
2. For any edge $(u, v) \in E$, there exists a bag that contains both $u$ and $v$.
3. If each bag is treated as a supervertex and edges are added accordingly (i.e., for two bags $B_1$ and $B_2$, if there exists $u \in B_1$, $v \in B_2$, and $(u, v)$ is an edge, then $B_1$ and $B_2$ are connected in the resulting tree), the induced graph is a tree.

This structure explains the term "tree decomposition." Furthermore, for any vertex $u \in V$, all the bags containing $u$ must form a subtree. For more details, please refer to Definition A.1 in the appendix.

Q: Can the authors provide references on the Cholesky factorization mentioned in line 318?

A: Sparse Cholesky factorization for matrices with sparsity patterns corresponding to graphs of small treewidth was first proposed in Schreiber (1982) and later developed in Davis and Hager (1999). More recent applications of this procedure appear in Ye (2020) for linear programs, Dong, Lee, and Ye (2021), and Gu and Song (2022) for semidefinite programs. For further exposition, see Section A.3.

Q: Regarding $d$ and $\log n$: Must $d$ grow with $n$? Can the authors explain the assumption $d = \Theta(\log n)$?

A: Thank you for highlighting this. Upon review, we recognize this as an oversight. The correct assumptions are as follows:

• For our almost-linear time algorithm, we require $d = O(\log n)$.
• For our hardness result, we require $d = \Omega(\log n)$.

Our lower-bound instance has a dimension $\Omega(\log n)$. The origin of this lower bound traces back to Alman and Williams (2015), which establishes hardness results for Bichromatic Hamming Closest-Pair using Orthogonal Vectors. The precise hardness result we used is from Rubinstein (2018). We have revised all occurrences of $d = \Theta(\log n)$ to $d = O(\log n)$ or $d = \Omega(\log n)$, and these changes have been marked in red.

References

• Backurs, Indyk, and Schmidt (2017): On the Fine-Grained Complexity of Empirical Risk Minimization: Kernel Methods and Neural Networks. NeurIPS 2017.
• Kyng, Wang, and Zhang (2020): Packing LPs are Hard to Solve Accurately, Assuming Linear Equations are Hard. SODA 2020.
• Ding, Kyng, and Zhang (2022): Two-Commodity Flow Is Equivalent to Linear Programming Under Nearly-Linear Time Reductions. ICALP 2022.
• Ye (2020): Fast Algorithm for Solving Structured Convex Programs. 2020.
• Dong, Lee, and Ye (2021): A Nearly-Linear Time Algorithm for Linear Programs with Small Treewidth: A Multiscale Representation of Robust Central Path. STOC 2021.
• Gu and Song (2022): A Faster Small Treewidth SDP Solver. 2022.
• Cohen, Lee, and Song (2019): Solving Linear Programs in the Current Matrix Multiplication Time. STOC 2019.
• Schreiber (1982): A New Implementation of Sparse Gaussian Elimination. TOMS 1982.
• Davis and Hager (1999): Modifying a Sparse Cholesky Factorization. SIAM Journal on Matrix Analysis and Applications 1999.
• Alman and Williams (2015): Probabilistic Polynomials and Hamming Nearest Neighbors. FOCS 2015.
• Rubinstein (2018): Hardness of Approximate Nearest Neighbor Search. STOC 2018.

We appreciate the reviewer for their thoughtful comments. We address your concerns below.

Concern: The result works for linear SVM and Gaussian kernels; however, the analysis on other kernel options, such as polynomial kernels, is missing.

Answer: Our technique can indeed be extended to other kernels. As indicated in Alman, Chu, Schild, and Song (2020), as long as the kernel function takes the form ${\sf K}(u, v) = f(\\|u - v\\|_2^2)$ and $f$ can be approximated by a low-degree polynomial of degree $o(\log n)$, there exists a rank-$n^{o(1)}$ factorization of the kernel matrix that can be computed in $O(n^{1+o(1)})$ time. This family includes a large subclass of polynomial kernels called $p$-convergent kernels (Song, Woodruff, Yu, and Zhang, 2021).

Concern: What if the kernel is a combination of kernels, such as $K_{\rm new} = K_1 + K_2$, where both $K_1$ and $K_2$ are Gaussian kernels?

Answer: In this case, we could run our algorithm separately on $K_1$ and $K_2$. Suppose $K_1$ has width $\sigma_1$, $K_2$ has width $\sigma_2$, and let $R$ denote the radius of the dataset. As long as $\frac{R^2}{2\sigma_1^2}, \frac{R^2}{2\sigma_2^2} \leq o(\frac{\log n}{\log\log n})$, our algorithm can solve each individual kernel in $O(n^{1+o(1)})$ time. The results for $K_1$ and $K_2$ can then be combined.

Question: The inequality constraint should be generalized into $Ax \leq b$ instead of $x \leq 0$.

Answer: We are using the slack form of convex quadratic programming, where the constraint region is the intersection of a linear subspace ($Ax = b$) and the non-negative orthant ($x \geq 0$). The formulation with $Ax \leq b$ that you are referring to is called the standard form, where the constraint region is the intersection of the polytope $Ax \leq b$ and the non-negative orthant $x \geq 0$.

One can reduce a standard form constraint to the slack form by increasing the dimension from $n$ to $2n$ and introducing slack variables: $Ax + s = b$, $x \geq 0$, $s \geq 0$. Slack form is commonly used for designing algorithms, such as the simplex method. Here, we also adopt the slack form and note that this formulation is without loss of generality.

Question: When $\epsilon = 10^{-3}$, the best first-order algorithm takes time $\widetilde{O}(\epsilon^{-1} d)$, so it would take $10^3$ iterations rather than $10^6$.

Answer: Thank you for pointing this out. We have corrected this in our updated manuscript.

Question: It would be better if you could perform some experiments and compare against LIBSVM.

Answer: We agree that it is important to evaluate the performance of a theoretical algorithm on real-world datasets. However, due to limited time and computational resources, we were unable to perform experiments for our algorithm.

We note that this line of works (Cohen, Lee and Song, 2019; Lee, Song and Zhang, 2019; Jiang, Song, Weinstein and Zhang, 2021; Dong, Lee and Ye, 2021; Gu and Song, 2022) are all theoretical studies without experimental evaluations. In fact, a common belief in this area is that the robust interior point method studied in these works, as well as our work, is efficient in practice. However, this remains a significant open question due to the nontrivial interplay of linear algebraic data structures. We leave the practical implementation of the algorithm as an open problem.

References

• Alman, Chu, Schild, and Song (2020): Algorithms and Hardness for Linear Algebra on Geometric Graphs. FOCS 2020.
• Song, Woodruff, Yu, and Zhang (2021): Fast Sketching of Polynomial Kernels of Polynomial Degree. ICML 2021.
• Cohen, Lee and Song (2019): Solving Linear Programs in the Current Matrix Multiplication Time. STOC 2019.
• Lee, Song and Zhang (2019): Solving Empirical Risk Minimization in the Current Matrix Multiplication Time. COLT 2019.
• Jiang, Song, Weinstein and Zhang (2021): Faster Dynamic Matrix Inverse for Faster LPs. STOC 2021.
• Dong, Lee and Ye (2021): A Nearly-Linear Time Algorithm for Linear Programs with Small Treewidth. STOC 2021.
• Gu and Song (2022): A Faster Small Treewidth SDP Solver. 2022.