Predicting kernel regression learning curves from only raw data statistics

Dhruva Karkada; Joseph Turnbull; Yuxi Liu; James B. Simon

cs.LGOct 16, 2025

Predicting kernel regression learning curves from only raw data statistics

Authors: Dhruva Karkada, Joseph Turnbull, Yuxi Liu, James B. Simon

Organizations: UC Berkeley · Imbue

Abstract

We study kernel regression with common rotation-invariant kernels on real datasets including CIFAR-5m, SVHN, and ImageNet. We give a theoretical framework that predicts learning curves (test risk vs. sample size) from only two measurements: the empirical data covariance matrix and an empirical polynomial decomposition of the target function $f_*$ . The key new idea is an analytical approximation of a kernel's eigenvalues and eigenfunctions with respect to an anisotropic data distribution. The eigenfunctions resemble Hermite polynomials of the data, so we call this approximation the Hermite eigenstructure ansatz (HEA). We prove the HEA for Gaussian data, but we find that real image data is often "Gaussian enough" for the HEA to hold well in practice, enabling us to predict learning curves by applying prior results relating kernel eigenstructure to test risk. Extending beyond kernel regression, we empirically find that MLPs in the feature-learning regime learn Hermite polynomials in the order predicted by the HEA. Our HEA framework is a proof of concept that an end-to-end theory of learning which maps dataset structure all the way to model performance is possible for nontrivial learning algorithms on real datasets.

Figures & tables

Panel of Figure 1 — Figure 1: We provide an end-to-end theory of learning for kernel ridge regression (KRR) that maps minimal statistics of the data distribution to test-time performance. (Top left.) KRR implicitly consists of two steps: (1) the kernel maps the data to high-dimensional nonlinear features ${\bm{x}}\mapsto{\bm{\psi}}({\bm{x}})$ , then (2) it fits a linear estimator to these features. Let the covariance in data space be $\mathbb{E}\!\left[{\bm{x}}{\bm{x}}^{\top}\right]={\bm{U}}{\bm{\Gamma}}{\bm{U}}^{\top}$ and let the covariance in feature space be $\mathbb{E}\!\left[{\bm{\psi}}({\bm{x}}){\bm{\psi}}({\bm{x}})^{\top}\right]={\bm{\Xi}}{\bm{\Lambda}}{\bm{\Xi}}^{\top}$ . (Lower left.) We introduce an ansatz (blue box) that predicts the feature covariance statistics ( ${\bm{\Xi}},{\bm{\Lambda}}$ ) from only the data covariance ( ${\bm{U}},{\bm{\Gamma}}$ ) and the functional form of ${\bm{\psi}}(\cdot)$ . This is sufficient to predict average-case test error using known results (green box). (Top right.) We are able to predict learning curves for KRR on image tasks without requiring omniscient knowledge of the feature statistics (i.e., without ever constructing or diagonalizing a kernel matrix). (Bottom right.) We are able to accurately predict the KRR sample complexity, including constant prefactors, for learning polynomials of the ImageNet dataset. See Figure 3 for additional plots and Appendix D for experimental details.

Figure 2: The Hermite eigenstructure ansatz ( HEA ) accurately predicts the eigenvalues (top) and eigenfunctions (bottom) of various kernel/dataset combinations. For four kernel/dataset settings (columns), we compute the empirical kernel eigensystem $\{(\lambda^{(\text{emp})}_{i},\phi^{(\text{emp})}_{i})\}$ and compare to the theoretical eigenpairs $\{(\lambda^{(\text{th})}_{i},\phi^{(\text{th})}_{i})\}$ obtained from Definition 4 and indexed in order of decreasing $\lambda^{(\text{th})}_{i}$ . In the top plot in each column, the $i$ -th point from the top right has coordinates $(\lambda^{\text{(emp})}_{i},\lambda^{\text{(th})}_{i})$ and its color indicates the polynomial degree of $\phi_{i}^{(\text{th})}$ . In the bottom plot in each column, we bin both the predicted and empirical eigenfunctions into logarithmic spectral bins and visualize the pairwise subspace overlap ( Equation 44 ), with axes matching the top plot. Grey pixels indicate bins with no eigenvalues. In all plots, concentration along the diagonal indicates theory-experiment match. See Section D.2 for further explanation and experimental details.

Figure 3: Using only the empirical data covariance and a polynomial decomposition of the target function, we predict learning curves across a variety of kernels, datasets, and targets. (Top row.) We predict test error on synthetic target functions of real datasets. Each target is the multi-Hermite polynomial ( Equation 3 ) whose leading term is indicated in the plots. (Recall that

z_{i}={\bm{u}}_{i}^{\top}{\bm{x}}/\sqrt{\gamma_{i}}

are the rescaled PCA coordinates.) For each of these targets, our ansatz predicts both learning curves (top left) and the sample complexity required to achieve

\texttt{MSE}\leq 0.5

(top right). (Bottom left and center.) We train on binarized true target functions. See Section D.1 for details. (Bottom right.) We construct synthetic targets by drawing from a Gaussian process. The source exponent controls the difficulty of the target. See Section D.1 for details. For all learning curve predictions, we estimate the coefficients of the target function in the predicted eigenbasis using the Gram-Schmidt process described in Section D.3 . See Section D.4 for full experimental details.

Figure 4: The HEA accurately predicts polynomial learning order in MLPs in the feature-learning regime. We measure the amount of time it takes to train three-layer ReLU MLPs, each dot being a trained MLP on one multi-Hermite polynomial target once a test error of

\texttt{MSE}\leq 0.1

is reached. We find that the optimization time is well-predicted by the inverse square root of the HEA eigenvalue

\lambda_{\bm{\alpha}}^{-1/2}

Appendix figures & tables20 assets

Supplementary material from the paper’s appendix.

Appendix

Figure 5: We compute the Laplace kernel on-sphere coefficients

c_{\ell}

with

\sigma=r=1

and plot

c_{\ell}\gamma_{1}^{\ell}

for various

\ell

. Dashed lines show

r^{2}\gamma_{1}^{-1}

, the value of the predicted minimum.

Figure 6: We show an alternative visualization of the top row of Figure 2 , comparing the predicted and empirical spectra of various kernel/dataset combinations. We see that the HEA accurately predicts the minute details of the kernel spectrum. Furthermore, tail deviations are indeed caused by finite-kernel effects in the empirical spectrum rather than a failure of the HEA. We plot the empirical spectrum with a curve and the predicted spectrum with dots, instead of vice versa, so the colors of the dots can denote the polynomial order of the predicted eigenmode.

Figure 7: We compare the original plot from Figure 3 with kernel width 8 (left) to the same experimental setup except for a kernel width 2 (right). The true eigenfunctions of the narrower kernel deviate from the HEA prediction, especially for the high-degree modes. For a similar comparison, see Figure 13 .

Figure 8: We compare the original experimental setup from Figure 3 (sample normalization, kernel width 10) on the left, to two similar setups: no sample normalization (center) and a narrower kernel width of 4 (right). Interestingly, the narrow kernel width does not change the agreement much; removing sample normalization worsens agreement slightly.

Panel of Figure 9 — Figure 9: Here we show four plots in which the experimental setup does not change; instead, we compare different techniques for estimating the required theoretical quantities. (Upper left.) HEA prediction, identical to Figure 1 . (Upper right.) We numerically diagonalize the kernel matrix (size $25000\times 25000$ ) and use the obtained eigenstructure to make predictions. The accuracy of the prediction degrades as the number of training samples approaches the size of the kernel matrix used to estimate the true eigenstructure. (Lower left.) We use the greedy algorithm described in Algorithm 1 to estimate the target coefficients. We see good initial agreement, but it degrades due to accumulating non-orthogonality. (Lower right.) For comparison, we include the predictions one would obtain if one modeled the data distribution as an isotropic Gaussian. Clearly, a major contribution of our work is to handle the anisotropy in natural data, since it strongly affects the resulting learning curves.

Figure 10: We compare the plot from Figure 3 on the left with an identical experimental setup, except without sample normalization. We see that normalization is necessary for obtaining agreement, since the HEA is derived using the kernel’s on-sphere level coefficients.

Figure 11: For “dense” synthetic targets (i.e., targets whose power is split among many eigenfunctions), the HEA predicts the performance of the Laplace kernel very well. Contrasting this with the predictions for Laplace kernel on real targets ( Figure 10 ), we conclude that dense targets are more forgiving of errors in target coefficient recovery. For this experiment, we use ZCA strength

\omega^{2}=10^{-2}

and sample normalization.

Figure 12: Additional sample complexity plot, this time for Laplace kernel on CIFAR-5m. We use sample normalization and a ZCA strength of

3\times 10^{-3}.

Figure 13: The HEA breaks when kernel width is too narrow. We repeat the experiment of Figure 2 with the Gaussian kernel function $K({\bm{x}},{\bm{x}}^{\prime})=e^{-\frac{1}{2\sigma^{2}}\left|\!\left|{\bm{x}}-{\bm{x}}^{\prime}\right|\!\right|^{2}}$ and synthetic Gaussian data. The data is generated as ${\bm{x}}\sim\mathcal{N}(0,{\bm{\Gamma}})$ , where ${\bm{\Gamma}}$ has eigenvalues $\gamma_{i}\propto i^{-3}$ for $i=1\ldots 30$ , with $\gamma_{1}\approx 0.83$ after normalization. As the kernel width $\sigma$ shrinks to $\sigma\lesssim 1$ , the HEA breaks: eigenvalue prediction fails with higher-order modes dominating predictions, and the clear eigenspace overlap present at large width evaporates.

Figure 14: With the Laplace kernel, the HEA requires high effective data dimension… We repeat the experiment of Figure 13 with two modifications: first, the kernel function is Laplace $K({\bm{x}},{\bm{x}}^{\prime})=e^{-\frac{1}{\sigma}\left|\!\left|{\bm{x}}-{\bm{x}}^{\prime}\right|\!\right|}$ with width $\sigma=8\sqrt{2}$ , and second, the eigenvalues are $\gamma_{i}\propto(i+3)^{-\alpha}$ with variable exponent $\alpha$ . As $\alpha$ increases (moving left to right), the effective dimension $d_{\text{eff}}=\frac{\left(\sum_{i}\gamma_{i}\right)^{2}}{\sum_{i}\gamma_{i}^{2}}$ decreases, and agreement with the HEA degrades. Empirically, we find that $d_{\text{eff}}\approx 20$ is usually high enough to see good agreement.

Figure 15: …but with the Gaussian kernel, it does not need high effective dimension. We repeat the experiment of Figure 14 exactly, but use a Gaussian kernel with width $\sigma=3$ . Thanks to the smoothness of the Gaussian kernel, we see little degradation of the theory-experiment match of the HEA.

Figure 16: Different datasets have different degrees of Gaussianity in their early principal components. For four natural datasets — CIFAR-10, SVHN, MNIST, and the UCI Mushrooms tabular dataset — we compute the first four normalized principal coordinates

(z_{i})_{i=1}^{4}

and plot their distributions as histograms. We find that the first few principal coordinates for CIFAR-10 are fairly close to Gaussian, this is less true for SVHN, yet less true for MNIST, and not at all true for the Mushrooms dataset.

Figure 17: The HEA works better with complex high-dimensional datasets and worse with simpler datasets. We repeat the experiment of Figure 2 with a Gaussian kernel of width $\sigma=3$ on $8000$ samples from four datasets of increasing non-Gaussianity: CIFAR-10, SVHN, MNIST, and the UCI Mushrooms tabular dataset. As the dataset becomes simpler, the HEA works worse and worse, though the predicted trend is still roughly present, if quantitatively wrong, even for the tabular dataset.

Figure 18: Visualization of the first four PCA directions for the CIFAR-10, SVHN, and MNIST datasets.

Figure 19: Validation of MLPs being feature-learning. We investigate the Gram matrix of the first weight matrix, $\text{W}_{1}^{\top}\text{W}_{1}$ . We train MLPs with $\zeta=1$ on $d=10$ Gaussian synthetic data for $1000$ GD steps. Yellow denotes higher entries, while blue denotes near-zero entries. The Gram matrices are heavily weighted to modes present in the target function, confirming the feature learning regime.

Figure 20: Our primary MLP results are largely invariant to hyperparameter choices. We train MLPs on synthetic Gaussian data, with one hyperparameter being motified. The base case is shown in the (top left) for reference. Apart from the (top center) , all modifications produce no notable changes to our observation of

\lambda_{\bm{\alpha}}^{-1/2}\sim\eta\cdot n_{\text{iter}}

scaling. (Top center). The 2 hidden layer case is only different than the rest by a change in eigenvalue-independent prefactor, the scaling law still holding.

Figure 21: Sweeping the richness parameter changes the slope. All. The effective optimization time uses the richness-independent learning rate

\eta_{0}

instead of

\eta

. (Left) We validate taking

\zeta\ll 1

to bring us into the lazy/NTK regime, where we know the eigenvalue to be inversely proportional to the effective optimization time. To speed up training for the lazy runs, we set

\eta_{0}=2\cdot 10^{1}

. This likely is the cause of the constant mode being notably high. (Right) We investigate the so-called ”ultra-rich” regime

\zeta\gg 1

and find there to be a

\lambda_{\bm{\alpha}}^{-0.4}\sim\eta\cdot n_{\text{iter}}

relationship. We suspect the

0.4

slope in the ultra-rich regime would change with variation of hyperparameters and is thus not fundamental, unlike the slopes of

1

and

1/2

in the lazy and rich regimes.

Panel of Figure 22 — Figure 22: Diagonalizing the kernel in 2 dimensions at the $a_{1}/a_{0}\to 0$ limit.

Panel of Figure 23 — Figure 23: Diagonalizing the kernel in 3 dimensions at two limits.

Figure 24: Structure of the kernel spectrum, zooming into one particular

N

-th segment. They converge to a segment of the target spectrum according to multiplicity of the target spectrum.

Explore similar work

Apr 25, 2026stat.ML

Learning Curves and Benign Overfitting of Spectral Algorithms in Large Dimensions

Existing large-dimensional theory for spectral algorithms resolves either the optimally tuned point or the interpolation limit, but leaves the under-regularized regime unexplored. We study the learning curve and benign overfitting of spectral algorithms in the large-dimensional setting where the sample size and dimension are of comparable order, i.e.,

n \asymp d^γ

for some

γ>0

. We first consider inner-product kernels on the sphere

\mathbb{S}^{d-1}

and establish a sharp asymptotic characterization of the excess risk across the full regularization path under various source conditions

s \geq 0

, where

s

measures the relative smoothness of the regression function. Our results reveal that the learning curve is not simply U-shaped but instead consists of three distinct regimes: over-regularized, under-regularized, and interpolation regimes. This characterization allows us to fully capture the benign overfitting phenomenon, demonstrating that benign overfitting arises consistently across both the under-regularized and interpolation regimes whenever

s

is positive but no larger than a critical threshold. We further show that, in the sufficiently regularized regime, the kernel learning curve is recovered by an associated sequence model. Finally, we extend the learning-curve analysis to large-dimensional KRR for a class of kernels on general domains in

\mathbb{R}^d

whose low-degree eigenspaces satisfy spectral-scaling and hyper-contractivity conditions.

Weihao Lu, Qian Lin, Yingcun Xia +1

Department of Statistics and Data Science, National University of Singapore. · Department of Statistics and Data Science, Tsinghua University.

Kernel Ridge Regression Spectral Representation Method Overfitting Regularization Optimal Sample Complexity Finite Time Analysis

Ask PDF

May 14, 2026stat.ML

Average Gradient Outer Product in kernel regression provably recovers the central subspace for multi-index models

We study a prototypical situation when a learned predictor can discover useful low-dimensional structure in data, while using fewer samples than are needed for accurate prediction. Specifically, we consider the problem of recovering a multi-index polynomial

f^*(x)=h(Ux)

, with

U\in\mathbb{R}^{r\times d}

and

r\ll d

, from finitely many data/label pairs. Importantly, the target function depends on input

x

only through the projection onto an unknown

r

-dimensional central subspace. The algorithm we analyze is appealingly simple: fit kernel ridge regression (KRR) to the data and compute the Average Gradient Outer Product (AGOP) from the fitted predictor. Our main results show that under reasonable assumptions the top

r

-dimensional eigenspace of AGOP provably recovers the central subspace, even in regimes when the prediction error remains large. Specifically, if the target function

f^*

has degree

p^*

, it is known that

n\asymp d^{p^*}

samples are necessary for KRR to achieve accurate prediction. In contrast, we show that if a low degree

p

component of

f^*

already carries all relevant directions for prediction, subspace recovery occurs in the much lower sample regime

n\asymp d^{p+δ}

for any

δ\in(0,1)

. Our results thus demonstrate a separation between prediction and representation, and provide an explanation for why iterative kernel methods such as Recursive Feature Machines (RFM) can be sample-efficient in practice.

Libin Zhu, Damek Davis, Dmitriy Drusvyatskiy +1

Department of Mathematics, University of Washington, Seattle, WA 98195 · Wharton Department of Statistics and Data Science, University of Pennsylvania, Philadelphia, PA 19104, USA · Department of Mathematics, U. Washington, Seattle, WA 98195 +1

Kernel Ridge Regression Kernel Method Subspace Low-Dimensional Structure Optimal Sample Complexity Mean

Ask PDF

May 14, 2026stat.ML

On Kernel Eigen-alignments of KRR: Reconstruction and Generalization

This paper investigates the critical role of eigenalignments between the kernel matrix and learning targets in achieving robust generalization in learning problems. We establish a direct connection between generalization performance in kernel methods and the estimation of eigenvectors and eigenvalues of matrices, offering a more intuitive understanding compared to prior work with minimal assumptions. We also show that, since the prediction task in KRR is essentially the weighted sum of eigenvectors/singular vectors, by analyzing how much error can be caused by perturbations to the kernel matrix, we can then derive a bound on this generalization error using the estimation stability of matrix eigenvalues and eigenvectors. Compared with previous work, our analysis concentrates on finite-sample settings and on the generalization error arising from having a suboptimal finite training set. Our findings reveal that in kernel methods, as long as the kernel is of high rank, the near-zero reconstruction error can be trivially obtained, implying that the reconstruction error will have limited predictive power for generalization. Finally, we establish a generalization bound from an eigenvalues/eigenvectors estimation perspective, showing that strong generalization requires increasing eigenvector alignment, eigenvalue magnitude, or gaps between consecutive eigenvalues.

Yang Liu, Ernest Fokoue, Richard Lange +1

Galisano College of Coumputing and Information Science · Rochester Institute of Technology · Rochester, NY +1

Kernel Ridge Regression Spectral Decomposition Generalization Bounds Reconstruction Error Brain Alignment

Ask PDF