Leverage Score Sampling Overview

• Leverage Score Sampling is a method that uses projection metrics to select data points based on their influence in the underlying subspace.
• It reduces problem size for tasks like least squares and kernel regression while preserving essential spectral and geometric structures.
• Extensions to ridge regularization, neural tangent kernels, and dynamic data enable scalable learning across diverse applications.

Searching arXiv for recent and foundational papers on leverage score sampling to ground the article. Leverage score sampling is a family of importance-sampling methods that select rows, columns, data points, features, or neurons according to quantities derived from projection matrices or their regularized analogues. In its classical matrix form, the statistical leverage score of a row measures how strongly that row participates in the relevant subspace of the data matrix, and leverage-based sampling yields smaller sketches that preserve the geometry of the original problem for least squares and related tasks (Lee et al., 2020). Subsequent work extends this principle to ridge-regularized kernel methods, random features, neural tangent kernels, active learning, distributed optimization, generative modeling, streaming algorithms, and explainable machine learning, with the common goal of biasing computation toward directions that are structurally important rather than uniformly frequent (Chen et al., 2021).

1. Classical definition and geometric interpretation

For a data matrix $A \in \mathbb{R}^{n \times d}$, with rows $a_i^\top$, the statistical leverage scores are defined through the orthogonal projector onto the column space of $A$. Writing

$$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$

the $i$-th leverage score is

$$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$

If $A = U \Sigma V^\top$ is a thin singular value decomposition, then $P = UU^\top$ and

$\sigma_i(A) = \|U_{i,:}\|_2^2.$

These scores lie in $[0,1]$, sum to $a_i^\top$0, and the coherence is $a_i^\top$1 (Lee et al., 2020).

This formulation identifies leverage as a measure of alignment with the principal subspace. High-leverage rows contain directions that are poorly represented by other rows; low-leverage rows are near-redundant. The same interpretation appears in matrix and kernel settings: ALEVS defines leverage scores for a symmetric positive semidefinite kernel matrix $a_i^\top$2 using the top $a_i^\top$3 eigenvectors $a_i^\top$4, with

$a_i^\top$5

so a point has high leverage when it strongly participates in the top $a_i^\top$6 eigendirections of the kernel matrix (Orhan et al., 2015).

Leverage score sampling converts these scores into a probability distribution. In the matrix case,

$a_i^\top$7

Sampling rows according to $a_i^\top$8 and rescaling each selected row by $a_i^\top$9 produces a sketch $A$0 that preferentially retains high-leverage structure while reducing the problem size (Lee et al., 2020). A closely related deterministic variant selects the columns with the largest leverage scores until their cumulative mass exceeds a threshold, and can achieve relative-error guarantees under a power-law decay assumption on leverage scores (Papailiopoulos et al., 2014).

2. Matrix sketching, regression, and deterministic variants

Leverage score sampling originated in numerical linear algebra and theoretical computer science, where it underlies randomized algorithms for least squares, low-rank approximation and CUR decompositions, linear programming, semidefinite programming, graph sparsification, maximum matching, max-flow, and random spanning trees (Lee et al., 2020). The common pattern is to replace a large matrix or graph by a sampled and rescaled subset whose spectral or geometric properties approximate the original.

For overdetermined least squares,

$A$1

one forms a sampling matrix $A$2, obtains $A$3 and $A$4, and solves

$A$5

With high probability, leverage-score-based sampling yields a solution $A$6 satisfying

$A$7

when $A$8 or $A$9 (Lee et al., 2020). These sketches also preserve quadratic forms and spectral structure through subspace embedding bounds of the form $$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$0.

Deterministic leverage-score selection replaces randomized sampling by thresholding the largest rank-$$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$1 leverage scores $$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$2. If the selected columns capture all but $$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$3 of the total leverage, then

$$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$4

and under power-law decay of the leverage scores, the number of selected columns is polynomial in $$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$5 and $$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$6 (Papailiopoulos et al., 2014). This provides relative-error guarantees in both Frobenius and spectral norm for a purely deterministic column subset selection rule.

Ridge leverage scores regularize the classical notion by replacing $$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$7 with $$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$8. For a matrix $$P \;=\; A (A^\top A)^\dagger A^\top \in \mathbb{R}^{n \times n},$$9, the ridge leverage score of column $i$0 is

$i$1

A deterministic ridge leverage score sampler selects columns until the leftover ridge leverage mass is below $i$2, producing a column subset $i$3 that satisfies additive-multiplicative spectral bounds, $i$4 error column subset selection, $i$5 error projection-cost preservation, and a $i$6 bound on statistical risk for ridge regression under the prescribed regularization (McCurdy, 2018). This suggests a direct connection between leverage-based sampling and sparse, interpretable feature selection in regularized regression.

A separate line of work uses leverage-aware sampling to accelerate iterative solvers rather than merely reduce matrix dimension. For regression $i$7, leverage-score sampling combined with proximal point methods and accelerated coordinate descent yields a running time

$i$8

improving over the previous $i$9, where $$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$0 and $$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$1 is the maximum number of non-zero entries in a row of $$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$2 (Agarwal et al., 2017). In that framework, leverage scores act as a preconditioning device that balances coordinate-wise smoothness constants.

3. Ridge leverage scores in kernel methods and random features

In kernel ridge regression, with kernel matrix $$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$3, labels $$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$4, and regularization $$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$5, one solves

$$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$6

or equivalently predicts at a new point $$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$7 via

$$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$8

The cubic cost of forming $$\sigma_i(A) \;:=\; P_{ii} \;=\; a_i^\top (A^\top A)^\dagger a_i.$$9 motivates random feature and Nyström approximations (Lee et al., 2020).

For kernels admitting a random feature representation

$A = U \Sigma V^\top$0

ridge leverage scores quantify the importance of random features relative to the regularized kernel matrix. Defining

$A = U \Sigma V^\top$1

the ridge leverage function is

$A = U \Sigma V^\top$2

Its total mass

$A = U \Sigma V^\top$3

is the statistical dimension of $A = U \Sigma V^\top$4 at regularization $A = U \Sigma V^\top$5 (Lee et al., 2020).

Sampling random features according to

$A = U \Sigma V^\top$6

and reweighting by $A = U \Sigma V^\top$7 yields a feature matrix $A = U \Sigma V^\top$8 satisfying the spectral approximation

$A = U \Sigma V^\top$9

with $P = UU^\top$0 samples, rather than a feature count tied directly to $P = UU^\top$1 (Lee et al., 2020). This extends earlier Fourier-feature results to any kernel representable as an expectation over finite-dimensional features, including random Fourier features, dot-product kernels, and neural tangent kernels (Lee et al., 2020).

A separate line investigates approximation of the diagonal entries of

$P = UU^\top$2

which are the statistical leverage scores for KRR in Nyström methods. For stationary kernels, an analytic approximation uses the input density $P = UU^\top$3 and the spectral density $P = UU^\top$4 of the kernel. The central approximation is

$P = UU^\top$5

which approximates the rescaled leverage score $P = UU^\top$6 in large samples (Chen et al., 2021). This gives a linear-time, modulo poly-log terms, algorithm for stationary-kernel-based KRR, with asymptotic relative error guarantees and Nyström risk guarantees matching full KRR up to constants (Chen et al., 2021). Conceptually, it shows that leverage-score non-uniformity is governed by local input density and kernel smoothness.

Random-feature leverage scores can also be approximated empirically. Given $P = UU^\top$7 candidate features and $P = UU^\top$8 data points, the empirical matrix

$P = UU^\top$9

induces a discrete distribution

$\sigma_i(A) = \|U_{i,:}\|_2^2.$0

Sampling features from this empirical leverage distribution yields improved out-of-sample performance over vanilla Monte Carlo sampling, with a generalization bound that separates a kernel approximation term from a spectral term involving the smallest eigenvalue of the feature Gram matrix $\sigma_i(A) = \|U_{i,:}\|_2^2.$1 (Shahrampour et al., 2019). This suggests a trade-off: increasing the initial feature pool improves approximation of the target kernel but can worsen conditioning in feature space.

4. Neural tangent kernels and neural-network initialization

Leverage score sampling has been extended to deep learning theory through the neural tangent kernel (NTK). For a two-layer ReLU network

$\sigma_i(A) = \|U_{i,:}\|_2^2.$2

with fixed second-layer weights $\sigma_i(A) = \|U_{i,:}\|_2^2.$3, the NTK with respect to $\sigma_i(A) = \|U_{i,:}\|_2^2.$4 can be written for two-layer ReLU as

$\sigma_i(A) = \|U_{i,:}\|_2^2.$5

Thus the NTK is a kernel of the same random-feature form addressed by generalized leverage score sampling (Lee et al., 2020).

From this viewpoint, finite-width network initialization selects a finite set of random features that approximates the continuous NTK. The initial kernel $\sigma_i(A) = \|U_{i,:}\|_2^2.$6 is a Monte Carlo approximation of the continuous NTK matrix $\sigma_i(A) = \|U_{i,:}\|_2^2.$7. Replacing Gaussian initialization by leverage-score-based initialization samples first-layer weights from

$\sigma_i(A) = \|U_{i,:}\|_2^2.$8

and reweights neurons by $\sigma_i(A) = \|U_{i,:}\|_2^2.$9 (Lee et al., 2020). Under this initialization,

$[0,1]$0

with $[0,1]$1 neurons (Lee et al., 2020).

The same work establishes equivalence between $[0,1]$2-regularized training of a two-layer ReLU network and NTK kernel ridge regression under both classical Gaussian initialization and leverage-score-based initialization (Lee et al., 2020). For Gaussian initialization, the neural network predictions converge to the NTK KRR solution in the infinite-width limit, on both training data and test points, under appropriate width and gradient-flow assumptions. With leverage-score-based initialization, an analogous training-data equivalence theorem holds, with the theoretical advantage that the required width can scale with the statistical dimension $[0,1]$3 rather than a worst-case dimension (Lee et al., 2020).

This suggests a reinterpretation of initialization in the NTK regime: initialization is not merely a random seed but a feature-sampling mechanism. Under that interpretation, leverage score sampling identifies neurons that are important for approximating the regularized NTK, thereby reducing the width needed to emulate NTK ridge-regression dynamics (Lee et al., 2020).

5. Sampling data points, coalitions, and partitions

Leverage score sampling has also been adapted to settings where the sampled objects are data points, active-learning queries, Shapley coalitions, or data partitions rather than matrix rows in a classical least-squares sketch.

In active learning, ALEVS computes leverage scores on class-specific kernel matrices. At each round, a classifier $[0,1]$4 partitions data into predicted-positive and predicted-negative subsets; for each class-specific kernel matrix, the top $[0,1]$5 eigenvectors are computed, and the leverage score of point $[0,1]$6 is

$[0,1]$7

The query rule is deterministic: $[0,1]$8 On ten binary classification tasks, ALEVS outperforms the baselines on most datasets for a substantial portion of the query range, particularly in early iterations, while class-wise leverage scores outperform leverage scores computed on all data at once (Orhan et al., 2015). The method is structure-based rather than margin-based: it queries points that are influential in the top eigenspaces of their predicted class manifolds.

In explainable machine learning, Leverage SHAP formulates Shapley value estimation as a linearly constrained weighted least-squares problem and applies leverage score sampling to the projected design matrix arising in Kernel SHAP. After reducing the constraint, the projected matrix $[0,1]$9 has analytically tractable leverage scores: $a_i^\top$00 for any coalition indicator $a_i^\top$01 with $a_i^\top$02 (Musco et al., 2024). Thus the leverage score depends only on coalition size, and the total leverage mass is identical across each size level. This induces a sampling rule that first chooses a subset size uniformly from $a_i^\top$03 and then chooses a coalition uniformly within that size class. Combined with paired sampling and without-replacement selection, the resulting Leverage SHAP algorithm estimates Shapley values with

$a_i^\top$04

model evaluations in expectation while providing non-asymptotic regression-error guarantees and an $a_i^\top$05-error bound for the Shapley vector (Musco et al., 2024).

In distributed optimization, weighted leverage score sampling is performed on data partitions rather than individual rows. Given row-level leverage probabilities $a_i^\top$06, partition-level probabilities are defined by aggregation: $a_i^\top$07 The algorithm samples $a_i^\top$08 partitions with replacement according to $a_i^\top$09, records multiplicities in a weight vector $a_i^\top$10, rescales selected partitions by $a_i^\top$11, and integrates these weights into a gradient coding scheme via a modified coding matrix $a_i^\top$12 (Charalambides et al., 2020). A central theorem shows that the resulting weighted coded gradient coincides with the gradient of the leverage-score-sketched objective, thereby preserving approximation guarantees while maintaining straggler resilience (Charalambides et al., 2020).

These examples illustrate that leverage score sampling is not restricted to classical row-sketching. The sampled objects may be data instances, feature functions, neurons, coalitions, or partitions, provided there is a projection-based or ridge-regularized notion of influence that can be translated into a sampling distribution.

6. Generative modeling, continual learning, and streaming settings

Several later developments use leverage scores to bias learning toward rare, informative, or dynamically updated structure.

For complete mode coverage in generative adversarial networks, ridge leverage scores are used to reweight the empirical data distribution so that underrepresented modes are sampled more frequently during GAN training. Given a feature map $a_i^\top$13, kernel matrix $a_i^\top$14, covariance $a_i^\top$15, and ridge parameter $a_i^\top$16, the ridge leverage scores are

$a_i^\top$17

Training samples are then drawn with probability

$a_i^\top$18

The paper studies Gaussian-kernel feature maps, fixed explicit feature maps from a pre-trained network, and discriminator-based feature maps recomputed online (Schreurs et al., 2021). On Ring, Grid, unbalanced MNIST, and unbalanced CIFAR10 setups, leverage-score sampling improves mode coverage and often improves FID and Inception Score relative to standard baselines; discriminator-based or classifier-based feature maps produce especially strong rebalancing of minority modes (Schreurs et al., 2021). A plausible implication is that leverage scores identify minority-mode points because those points occupy sparse or structurally distinct regions in feature space.

In continual learning, online leverage score sampling maintains a fixed-size buffer of examples across tasks. At each task, the current sketch $a_i^\top$19 is concatenated with the new task data $a_i^\top$20, the SVD of the enlarged matrix is computed, leverage scores

$a_i^\top$21

are formed, and $a_i^\top$22 rows are sampled without replacement with probabilities proportional to the normalized scores. The retained sketch then serves as the training set for the next stage (Teng et al., 2019). This uses leverage scores as a criterion for which examples to remember and which to forget, with the rationale that high-leverage samples preserve the covariance structure of all data seen so far. On Rotated MNIST, Permuted MNIST, and Incremental CIFAR-100, the method is competitive with GEM and more efficient in wall-clock time than replay methods with more expensive gradient constraints (Teng et al., 2019).

In the turnstile streaming model, leverage-score ideas are generalized to $a_i^\top$23 settings. For $a_i^\top$24, the $a_i^\top$25 leverage score of row $a_i^\top$26 is

$a_i^\top$27

The paper first develops an $a_i^\top$28 sampler for row norms under arbitrary turnstile updates, based on randomized scaling and CountSketch, and then converts it into an $a_i^\top$29 leverage score sampler through preconditioning with an $a_i^\top$30 well-conditioned basis $a_i^\top$31 (Munteanu et al., 2024). The resulting sketch returns sampled indices, slightly perturbed sampled rows $a_i^\top$32, and approximate inclusion probabilities, enabling simulation of leverage-score-based coreset constructions for $a_i^\top$33, logistic, ReLU, and probit-like regression in turnstile streams (Munteanu et al., 2024). Notably, for logistic regression the framework yields the first turnstile-stream algorithm achieving a $a_i^\top$34-approximation with polynomial sketch or subsample size (Munteanu et al., 2024).

These developments shift the role of leverage scores from static matrix compression to dynamic data management: rare modes are oversampled, important past examples are retained, and high-sensitivity rows are extracted from fully dynamic streams.

7. Dependence structures, efficiency, and broader significance

A recurring theme in recent work is that leverage-based marginals can be combined with nontrivial dependence structures or domain-specific approximations without losing theoretical control.

In active learning for linear regression, dependent leverage score sampling via pivotal sampling preserves leverage-based marginals while enforcing negative dependence and spatial coverage. If the sampling distribution is $a_i^\top$35-homogeneous, has marginals proportional to leverage scores, and satisfies a one-sided $a_i^\top$36-independence condition, then linear functions can still be actively learned with $a_i^\top$37 samples, matching independent sampling (Shimizu et al., 2023). For polynomial regression on an interval, the pivotal method improves this to $a_i^\top$38 samples (Shimizu et al., 2023). Empirically, the method reduces the number of samples needed to reach a target accuracy by up to $a_i^\top$39, which the authors attribute to improved spatial coverage relative to independent leverage sampling (Shimizu et al., 2023). This suggests that leverage score sampling need not be independent to retain subspace-embedding behavior.

In high-dimensional geometric optimization, leverage score sampling has been used for the minimum volume covering ellipsoid problem. A deterministic top-leverage subset of rows of the data matrix $a_i^\top$40 is selected until the cumulative leverage mass exceeds $a_i^\top$41; the resulting sampled matrix $a_i^\top$42 satisfies

$a_i^\top$43

which yields explicit bounds on the D-optimal design objective and allows the overall complexity to drop from $a_i^\top$44 to $a_i^\top$45 under power-law decay of leverage scores (Harris et al., 2024). Here leverage scores play a geometric role: points with high leverage are also those most important for the enclosing ellipsoid.

Across these domains, several limitations recur. Exact ridge leverage functions are often expensive to compute because they require access to $a_i^\top$46 or its inverse; in NTK settings, the true leverage distribution is primarily of theoretical interest because approximating it remains costly (Lee et al., 2020). In active learning, ALEVS provides no formal label-complexity guarantee and can be computationally heavy because eigendecompositions are repeated every round (Orhan et al., 2015). In GAN training, effectiveness depends strongly on the chosen feature map and incurs additional cost from repeated leverage-score computation or dimension reduction (Schreurs et al., 2021). In kernel methods, analytic approximations are currently specialized to stationary kernels and moderate-dimensional regimes (Chen et al., 2021). In turnstile $a_i^\top$47 streaming, polynomial dependence on $a_i^\top$48 and $a_i^\top$49 remains substantial, and extending small-sketch guarantees to $a_i^\top$50 faces lower bounds (Munteanu et al., 2024).

Even with these limitations, a unifying pattern is evident. Leverage scores, whether classical, ridge-regularized, empirical, continuous, or $a_i^\top$51-generalized, measure how much a sampled object contributes to the effective subspace or function class relevant to the task. Sampling according to those scores yields approximation guarantees that depend on effective rank, statistical dimension, or leverage decay rather than on ambient problem size. In this sense, leverage score sampling functions as a common principle linking matrix sketching, kernel approximation, neural tangent theory, active data acquisition, robust distributed learning, mode-balanced generative modeling, continual learning, and dynamic-stream coresets (Lee et al., 2020).