PaperlayerView exact revision on arXiv

Source-linked AI summary

Divide and Conquer Kernel Ridge Regression: A Distributed Algorithm with Minimax Optimal Rates

Yuchen Zhang, John C. Duchi, Martin J. Wainwright

arXiv:1305.5029v2math.STcs.LGstat.ML

TL;DR

Kernel ridge regression is statistically attractive but expensive on large datasets. This paper partitions the data, fits and averages local KRR estimators, and proves that the resulting predictor retains minimax-optimal convergence under bounded splitting. The analysis also identifies computational savings and boundaries on how finely the data can be divided.

  • Problem

    Standard kernel ridge regression has prohibitive O(N^3) time and O(N^2) memory costs for large datasets.

  • Method

    The method randomly partitions N samples into m equal subsets, computes independent KRR estimates with a carefully chosen regularization parameter, and averages them.

  • Results

    The averaged estimator achieves minimax-optimal convergence rates when the number of splits m is not too large, while computational cost scales nearly linearly in N.

  • Takeaways & Limitations

    Variance from local under-regularization is reduced by m-fold averaging, allowing distributed computation without sacrificing the stated optimal convergence rates.

  • Takeaways & Limitations

    When f* does not belong to the RKHS, the current analysis prevents dividing the data into subsets finer than the stated N/γ2(λ) boundary.

Abstract

from arXiv · show

We establish optimal convergence rates for a decomposition-based scalable approach to kernel ridge regression. The method is simple to describe: it randomly partitions a dataset of size N into m subsets of equal size, computes an independent kernel ridge regression estimator for each subset, then averages the local solutions into a global predictor. This partitioning leads to a substantial reduction in computation time versus the standard approach of performing kernel ridge regression on all N samples. Our two main theorems establish that despite the computational speed-up, statistical optimality is retained: as long as m is not too large, the partition-based estimator achieves the statistical minimax rate over all estimators using the set of N samples. As concrete examples, our theory guarantees that the number of processors m may grow nearly linearly for finite-rank kernels and Gaussian kernels and polynomially in N for Sobolev spaces, which in turn allows for substantial reductions in computational cost. We conclude with experiments on both simulated data and a music-prediction task that complement our theoretical results, exhibiting the computational and statistical benefits of our approach.

1 Introduction

Kernel ridge regression is statistically well understood but computationally expensive at large N, motivating a distributed estimator that partitions data, fits local models, and averages them while retaining optimal rates.

  • O(N^3) time and O(N^2) memory make standard kernel ridge regression prohibitive for large datasets.
  • Low-rank kernel approximations and early-stopped optimization reduce computation, but their statistical guarantees require additional analysis.
  • The proposed method randomly partitions N samples into m equal subsets, fits independent KRR estimators, and averages them into a global predictor.
  • O(N^3/m^2) time and O(N^2/m^2) memory result from naive local KRR implementations, with superlinear speedup available using m parallel processors.The processors must still communicate their function estimates.
  • Regularizing each local problem as though it used all N samples makes local fits under-regularized, but averaging reduces their inflated variance enough to preserve optimal convergence.
  • Theoretical results cover RKHS-contained and general regression functions, with minimax-optimal rates illustrated for finite-rank kernels and Sobolev spaces.
  • Simulations and a music-prediction experiment complement the theoretical analysis by examining the method’s computational and statistical benefits.

2 Background and problem formulation

The paper formulates nonparametric regression in an RKHS and introduces KRR as a regularized least-squares estimator whose computation reduces to a kernel-matrix quadratic program.

  • 2.1 Reproducing kernels: A positive-semidefinite kernel defines an RKHS, which is strictly contained within L2(P) under the stated setup.
  • 2.1 Reproducing kernels: Mercer’s theorem represents the kernel through nonnegative eigenvalues and an orthonormal basis for L2(P).
  • 2.1 Reproducing kernels: RKHS functions expand in the eigenbasis, and the space can be viewed as an elliptical subset of the sequence space ℓ2(N).
  • 2.2 Kernel ridge regression: Nonparametric regression estimates the conditional mean f*(x)=E[Y|X=x] from N i.i.d. covariate-response samples by minimizing mean-squared prediction error.
  • 2.2 Kernel ridge regression: KRR minimizes empirical least-squares loss plus a weighted squared-Hilbert-norm penalty controlled by λ>0.
  • 2.2 Kernel ridge regression: The representer theorem restricts KRR solutions to the span of kernel functions evaluated at the training covariates.
  • 2.2 Kernel ridge regression: This restriction converts KRR computation into an N-dimensional quadratic program involving N^2 kernel-matrix entries.

3 Main results and their consequences

Fast-KRR randomly partitions the data, fits local KRR estimators, and averages them, yielding scalable computation with minimax-optimal rates under suitable conditions on the kernel, regularization, and number of partitions.

  • Algorithm: Fast-KRR randomly divides N samples into m subsets, computes one local KRR estimate per subset, and averages the estimates.The algorithm forms a family indexed by λ; its execution is independent of the auxiliary parameter d.
  • General guarantees: Theorem 1 gives an upper bound for the averaged estimator for any trace-class kernel when f∗ belongs to the associated Hilbert space.The theorem provides sharper estimation-error bounds in this setting and yields concrete guarantees for specialized kernel classes.
  • General guarantees: Theorem 2 extends the analysis to general f∗ through an oracle inequality combining approximation error with estimation and residual-error terms.Compared with Theorem 1, its bounds are somewhat looser and impose stronger constraints on the number of partitions.
  • Limitations: The analysis does not establish how to choose λ adaptively in the distributed setting, leaving cross-validation or unbiased risk estimation as an open question.When f∗ does not belong to the Hilbert space, the current analysis also restricts m more strongly than in the f∗∈H case.
  • Kernel-specific consequences: Finite-rank, polynomial-eigendecay, and exponential-eigendecay kernels achieve minimax-optimal rates when m remains below class-specific thresholds.The upper bounds are unimprovable up to constants for Sobolev and exponential kernel classes, and the finite-rank rate is minimax-optimal over all estimators using N samples.

4 Proofs of Theorem 1 and related results

The proof decomposes the averaged estimator’s error into bias and variance terms, then bounds each using auxiliary kernel ridge regression lemmas. These bounds yield conditions on the number of partitions and establish convergence rates for specific kernel classes.

  • Proof strategy: The proof bounds the averaged estimator through bias and variance components of one subsample estimator.The subsample size is n = N/m, and the proof separately controls estimation and approximation contributions.
  • Proof strategy: The proofs of the bias and variance lemmas constitute one main technical contribution.The theorem follows by combining these lemmas with the preceding error decomposition.
  • Partition threshold: m must remain below a kernel- and regularization-dependent threshold for the decomposition-based estimator to retain the target convergence rate.The analysis identifies this threshold by controlling the term T3(d).
  • Kernel-specific rates: For exponentially decaying eigenvalues, choosing λ = N^-2ν yields convergence rate λ = N^-2ν.The truncation choice makes the tail terms negligible under the stated eigenvalue-decay analysis.
  • Kernel-specific rates: For Gaussian-type exponential decay, the convergence rate is approximately √(log N/N) under the stated partition condition.A boundedness assumption gives a corresponding second corollary statement.

5 Proof of Theorem 2 and related results

Theorem 2 is proved using an error decomposition around a regularized population solution and auxiliary bounds for subsample estimation and approximation. The resulting oracle analysis applies when the regression function need not belong to the RKHS, but imposes stricter partition conditions.

  • Technical ingredients: The proof uses conditional mean-zero noise, RKHS minimization, Hölder’s inequality, Jensen’s inequality, and fourth-moment control.These steps establish the moment bound needed for the theorem’s residual analysis.
  • Proof strategy: The proof analyzes subsample error relative to the regularized population solution f*_{λ} using parallel bias and variance lemmas.The subsample estimator is standard KRR with n = N/m samples and regularization λ.
  • Proof strategy: Theorem 2’s proof combines auxiliary bounds under Assumptions A and B′ to control the oracle-style residual error.The argument applies the lemmas after establishing the required conditions and summing the resulting bounds.

6 Experimental results

The experiments test Fast-KRR’s accuracy, partition threshold, runtime, and performance on music year prediction. Across simulations and real data, partitioning preserves accuracy in supported regimes while reducing computation.

  • Simulation studies: m ≤ 16 gives Fast-KRR accuracy comparable to baseline KRR, while m = 64 closely matches the full estimator for N ≥ 2^11.Using λ = N^-2/3 avoids the performance gap caused by over-regularizing each subproblem with λ = n^-2/3.
  • Simulation studies: m ≲ N^1/3 is the theoretical Sobolev threshold, while simulations show error begins increasing around log(m) ≈ 0.45 log(N).The empirical threshold confirms polynomial growth in m can retain minimax optimality and is somewhat better than the predicted threshold.
  • Simulation studies: Fast-KRR yields substantial computational improvements without sacrificing accuracy in the timing experiment.The experiment records mean error and runtime across partition counts and dataset sizes; Fail denotes out-of-memory matrix inversion, while N/A denotes significantly non-optimal error.
  • Real data experiments: Fast-KRR has the best fixed-time performance on held-out song-year prediction, although its margin over Nyström sampling is not substantial.The comparison includes Nyström subsampling and random-feature kernel ridge regression, with error bars from ten experiments.
  • Real data experiments: Fast-KRR achieves much lower error and negligible performance variability than standard KRR trained on only a 1/m fraction of the data.The experiment uses regularization λ ≈ (N/m)^-1 for the fraction-data KRR baseline.

7 Discussion

The discussion emphasizes that Fast-KRR retains minimax convergence when the number of splits is sufficiently controlled, while reducing computational cost. It also identifies unresolved comparisons between Fast-KRR and Nyström kernel evaluations.

  • Theoretical guarantees: Fast-KRR achieves minimax optimal convergence rates when the number of data splits m is not too large.The error guarantees depend on the kernel’s effective dimensionality γ(λ).
  • Theoretical guarantees: m ≲ N/γ(λ)^2 is sufficient for the stated estimation-error guarantee.The threshold is expressed using the kernel effective dimensionality γ(λ).
  • Computational consequences: The convergence rate is minimax optimal while computational cost scales nearly linearly in N.This combines the statistical guarantee with the computational benefit of subsampling schemes.
  • Computational consequences: Fast-KRR may require at most Nγ(λ)^2 kernel evaluations when using m ≍ N/γ(λ)^2 machines.The comparison with Nyström uses N eγ(λ) evaluations, but the two scalings are not currently comparable.
  • Open questions: The optimal computational scaling of Fast-KRR versus Nyström remains unresolved across data-support settings.The discussion notes that eγ(λ) may scale linearly with N or resemble the true effective dimensionality, depending on the setting.

A Proof of Lemma 1

The proof of Lemma 1 bounds the conditional estimation error by decomposing the regression-function coefficients into finite-dimensional and tail components. It controls these terms through matrix concentration, moment assumptions, and norm inequalities.

  • Error decomposition: The proof defines the estimation error Δ = bf − f* and reduces the analysis to bounding its conditional expectation.Jensen’s inequality and the gradient optimality condition provide the initial reduction.
  • Finite-dimensional control: The proof establishes ∥E[Δ | X]∥_H ≤ ∥f*∥_H using positivity of the empirical covariance operator.This bound controls the Hilbert-space norm of the conditional error.
  • Error decomposition: The coefficient sequence of E[Δ | X] is split into a first d-dimensional component δ↓ and an infinite tail δ↑.Orthonormality of the basis separates the squared norm into finite and tail contributions.
  • Finite-dimensional control: The finite-dimensional analysis introduces the design matrix Φ, eigenvalue matrix M, tail vector v, and diagonal scaling Q.These objects rewrite the error terms into matrix expressions suitable for bounding.
  • Concentration bounds: Moment bounds and concentration of the normalized empirical Gram matrix control the finite-dimensional terms.The proof invokes Lemma 5, an event ensuring a lower spectral bound, and elementary norm inequalities.

A.1 Proof of Lemma 5

The proof of Lemma 5 bounds matrix and vector terms arising in the finite-dimensional error analysis. It combines moment estimates, symmetrization, and a matrix concentration inequality.

  • Term-by-term bounds: The proof bounds Q^-1M^-1θ↓ as the first component of the finite-dimensional error expression.Q is diagonal, allowing direct operator-norm control of this factor.
  • Term-by-term bounds: The second component, Q^-1Φ^T v, is rewritten as a product and bounded using Cauchy–Schwarz and moment estimates.The analysis separately controls the design-matrix column and tail-error factors.
  • Matrix concentration: The empirical covariance fluctuation Q^-1((1/n)Φ^TΦ − I)Q^-1 is analyzed through centered rank-one matrix summands.The proof defines the summands and applies symmetrization with Rademacher variables.
  • Matrix concentration: A matrix moment inequality for independent symmetric Hermitian matrices supplies the concentration step.It is applied with Xi = εiAi/n, followed by separate bounds on the resulting terms.

B Proof of Lemma 2

The proof controls the error by decomposing it into finite-dimensional and tail components, then bounding the resulting terms using concentration, moment, and norm inequalities.

  • Error decomposition: The error is expanded in the eigenfunction basis and split into truncated and tail components.Orthonormality yields E[∥∆∥^2] = E[∥∆↓∥^2] + E[∥∆↑∥^2].
  • Termwise bounds: The proof bounds the decomposition terms with Lemma 8, previously established inequalities, conditional variance control, and elementary quadratic inequalities.The argument separately controls the event complement, tail contribution, and transformed noise term before combining them.
  • Finite-dimensional representation: The finite-dimensional analysis uses the basis matrix Φ, eigenvalue matrix M, tail vector v, and Q = (I + λM^-1)^1/2.These objects rewrite the optimality condition into matrix expressions that can be bounded term by term.
  • Noise control: The transformed noise contribution is controlled through the diagonal entries of Q^-1 and the effective-dimension quantity γ(λ).The resulting sum is bounded by (σ^2/n)γ(λ).
  • Conclusion: Combining the component bounds with probability control and expectations completes the lemma after normalization by λn^2.The proof applies inequality (a + b + c)^2 ≤ 3a^2 + 3b^2 + 3c^2 and then divides by λn^2.

C Proof of Lemma 3

The proof of Lemma 3 uses the same truncation framework as earlier lemmas, derives bounds for the principal terms, and controls the remaining term under the stated assumptions.

  • Event-based control: The first term is rewritten using the design points and an event E, whose complement is controlled by bound (45).The resulting expression is then split into terms handled by Lemmas 10–12 and the event-probability bound.
  • Preliminary bounds: Lemma 9 supplies preliminary bounds used to control the tail and higher-order error quantities.The proof invokes Lemma 9 before bounding the remaining terms in the decomposition.
  • Error decomposition: The error is decomposed into finite-dimensional and tail components in the eigenfunction basis.The corresponding squared norms add by orthonormality.
  • Lemma-based bounds: Lemmas 10 and 11 provide bounds for the decomposed terms, with Lemma 11 requiring Assumptions A and B′.The proof combines these bounds through a quadratic inequality before applying the event-complement control.
  • Moment arguments: The second inequality follows from the first by Jensen’s inequality, while the remaining estimates use Hölder’s inequality and moment assumptions.Assumption A is applied with moment 2k ≥ 8, equivalently k ≥ 4.

D Proof of Lemma 4

The proof of Lemma 4 handles non-centered errors through a first-order optimality condition, truncation, and bounds on empirical covariance and moment terms.

  • Optimality condition: Because the errors are not conditionally mean-zero, the proof begins with the empirical objective’s first-order optimality condition.Fréchet differentiability yields the condition used to analyze the estimator error.
  • Population comparison: The optimal population solution cancels the corresponding expected gradient term, producing an analogue of the earlier covariance-based identity.This step uses the optimality of f*_{λ̄} and the sample covariance operator.
  • Truncation: The error is expanded in the eigenfunction basis and split into truncated and tail components before taking expectations.The tail vector, basis matrix Φ, eigenvalue matrix M, and Q encode the resulting finite-dimensional representation.
  • Termwise bounds: Lemma 10 bounds the first two terms in the expanded expression, while Lemma 12 controls the third term.The proof combines these estimates with (a + b + c)^2 ≤ 4a^2 + 4b^2 + 2c^2.
  • Moment control: The final moment bound uses Jensen’s and Cauchy–Schwarz inequalities together with the matrix moment bound from inequality (52).The argument applies (52) with k = 2 to establish the remaining claim.
Loading 1305.5029v2…