Improved Active Learning via Dependent Leverage Score Sampling
Authors: Atsushi Shimizu
Xiaoou Cheng
Christopher Musco
Jonathan Weare
Venue: ICLR 2024 oral (2024)
Keywords: leverage score sampling
active learning
polynomial regression
differential equations
pivotal sampling
Primary area: general machine learning (i.e., none of the above)
TL;DR: Better active learning (in theory and practice) in the presence of adversarial noise via non-independent leverage score sampling.

Abstract: We show how to obtain improved active learning methods in the agnostic (adversarial noise) setting by combining marginal leverage score sampling with non-independent sampling strategies that promote spatial coverage. In particular, we propose an easily implemented method based on the \emph{pivotal sampling algorithm}, which we test on problems motivated by learning-based methods for parametric PDEs and uncertainty quantification. In comparison to independent sampling, our method reduces the number of samples needed to reach a given target accuracy by up to $50\%$.

We support our findings with two theoretical results. First, we show that any non-independent leverage score sampling method that obeys a weak \emph{one-sided $\ell_{\infty}$ independence condition} (which includes pivotal sampling) can actively learn $d$ dimensional linear functions with $O(d\log d)$ samples, matching independent sampling. This result extends recent work on matrix Chernoff bounds under $\ell_{\infty}$ independence, and may be of interest for analyzing other sampling strategies beyond pivotal sampling. Second, we show that, for the important case of polynomial regression, our pivotal method obtains an improved bound of $O(d)$ samples.

IMPROVED ACTIVE LEARNING VIA DEPENDENT LEVERAGE SCORE SAMPLING
Atsushi Shimizu, Xiaoou Cheng, Christopher Musco, Jonathan Weare New York University
(as151 06, chengxo, cmusco, wea re @nnyu.ed
ABSTRACT
We show how to obtain improved active learning methods in the agnostic (ad- versarial noise) setting by combining marginal leverage score sampling with non- independent sampling strategies that promote spatial coverage. In particular, we propose an easily implemented method based on the pivotal sampling algorithm. which we test on problems motivated by learning- based methods for parametric PDEs and uncertainty quantification. method reduces the number of samples needed to reach given target accuracy by up 50%. We support our findings with two theoretical results. First, we show that any non-independent leverage score sampling method that obeys a weak one- sided too independence condition (which includes pivotal sampling) can actively learn d dimensional linear functions with Odd lol d) samples, matching indepen- dent sampling. This result extends recent work on matrix Chernoff bounds under f. independence, and may be of interest for analyzing other sampling strategies beyond pivotal sampling. Second, we show that, for the important case of polyno- mial regression, our pivotal method obtains an improved bound on O(d) samples.
comparison to independent sampling, our
INTRODUCTION
In the active linear regression problem, we are given a data matrix A E Raxd with n >> d rows and query access to a target vector R". The goal is learn parameters x E Rd such that Ax zs while observing as few entries in b as possible. We study this problem in the challenging agnostic learning or "adversarial noise" setting, where we do not assume any underlying relationship between A and b. Instead, our goal is to find parameters competitive with the best possible fit, good or bad. Specifically, considering & loss, let x" arg min Ax bli be optimal model parameters. We want to find x using small number of queried target values in b such that
AA - bl2 S (1 e)]lAx* - bll..
(1.1)
for eтror parameter O. Beyond being fundamental learning problem, active regression has emerged as a fundamental tool in learning based methods for the solution and uncertainty anal- ysis of parametric partial differential equations (PDEs) (Chkifa al. [015 Guo al., 2020). For such applications, the agnostic setting is crucial, as a potentially complex quantity of interest is ap- proximated by a simple surrogate model (c.g. polynomials, sparse polynomials, single layer neural networks, etc.) (Lüthen et al. 2021 Hokanson and Constantinc) 2018)). Additionally, reducing the number of labels used for learning is crucial, as each label usually requires the computationally expensive numerical solution of a PDE for new set of parameters (Cohen and De Vore 2015).
1.1 LEVERAGE SCORE SAMPLING
Active linear regression has been studied for decades in the statistical model where b is assumed to equal Ax" plus i.i.d. random noise. In this case, the problem can be addressed using tools from opti- mal experimental design (Pukelsheim) 2006). In the agnostic case, near-optimal sample complexity results were only obtained relatively recently using tools from non-asymptotic matrix concentration Tropp 2012). In particular, it was shown independently in several papers that collecting entries from b randomly with probability proportional to the statistical leverage scores of rows in A can achieve (I.I) with Od + d/e) samples (Sarlós) 2006 Rauhut and Ward 2012) Hampton and
Doostan 2015 Cohen and Migliorati] 2017). The leverage scores are defined as follows:

---

frequete
dorune free Imqquece
spring constant
S0r constant
spring constant
Target Function
(b) Berouuli Sampling
(c) Pivotal Sampling (our method)
Figure I: Polynomial approximations the maximum displacement of damped harmonic oscilla- tor, as a function of driving frequency and spring constant. (a) is the target value, and samples can be obtained through the numerical solution of differential equation governing the oscillator. Both (b) and (c) draw 250 samples using leverage score sampling and perform polynomial regression of degree 20. (b) uses Bernoulli sampling while (c) uses our pivotal sampling method. Our method gives better approximation, avoiding artifacts that result from gaps between the Bernoulli samples.
Definition 1.1 (Leverage Score). Let U E Raxr be any orthogonal basis for the column span of matrix A nnx Let a, and u, be the i-th rows of A and U, respectively. The leverage score т₃ of the i-th row in A can be equivalently written as:
u. 1 = a? (A Notice that T; = -| Σm=
max (a{x))/////))11 /
(1.2)
xER
when A is full-ankk and thus has T = columns.
See (Avron et al. 2017) for short proof of the final equality (.2).. This last definition, based on maximization problem, gives intuitive understanding of the leverage scores. The score of row a, higher it is more "exceptional", meaning that we can find a vector x that has large inner product with a; relative its average inner product (captured by ||x]) with all other rows in the matrix. Based on leverage score, rows that are more exceptional are sampled with higher probability.
sampling, either with without replacement. The typical approach for sampling without replacement, which we call "Bernoulli sampling" is as fol- lows: Each row a; is assigned a probability p; min(1, e. for oversampling parameter e 1. Then each row is sampled independently with probability p;- We construct a subsampled data matrix Á and subsampled target vector b by adding a /Pr to and b:/ vPi to b for any index i that is sampled. To solve the active regression problem, return x arg minx || - bll2-
Prior work considers independent leverage
1.2 OUR CONTRIBUTIONS
In applications to PDEs, the goal is often to approximate a function over a low dimensional distribu- tion X. E.g. X' might be uniform over an interval [-1, C R or over box [-1, x...X [-, C Ry, In this setting, the length d rows of A correspond to feature transformations of samples from X. For example, the ubiquitous task of polynomial regression, we start with x X and add to
2q A a row containing all combinations of entries in with total degree i.c., e1 2/ ..... for all non-negative integers f,....,4 such that Σ 4, S p. For such problems, "grid" based interpola- tion is often used in place of randomized methods like leverage scores sampling, i.c., the target is queried on a deterministic grid tailored X'. For example, when is uniform on a box, the stan- dard approach is to use grid based on the Chebyshev nodes (Xiu) 2016). Pictured in Figure 2 the Cheybshev grid concentrates samples near the boundaries the box, avoiding well known issue of Runge's phenomenon for uniform grids. Leverage score sampling does the same. In fact, the methods are closely related: in the high degree limit, the leverage scores for polynomial regression over the box match the asymptotic density of the Chebyshev nodes (Lüthen al. (021).
So how do the deterministic and randomized methods compare? The advantage of randomized sampling is that they yield strong provable approximation guar-
methods based leverage antees, and easily generalize to any distribution X' Deterministic methods are less flexible on the
no deterministic method can provably solve the agnostic regression problem with few samples. Since we make no assumptions on b, all error in Ax' b could he concentrated only the deterministic indices selected. Randomization is needed avoid high error outliers in b.
As discussed in Chen and Price, 2019), 2019))
2

---

(a) Bernoulli Sampling
(b) Pivotal Sampling (our method)
(c) Chebyshev Grid
The results of three different active learning methods used to collect samples to fit polynomial over [-1,1] x [-1,1]. The image on the left was obtained by collecting points inde- pendently at random with probability according to their statistical leverage scores. The image on the right was obtained by collecting samples at the 2-dimensional Chebyshev nodes. The image in the middle shows our method, which collects samples according to leverage scores, but using non-independent pivotal sampling strategy that ensures samples are evenly spread in spatially.
Figure 2:
other hand, and do not yield provable guarantees. However, the advantage of grid based methods is that they more "evenly" distribute samples over the original data domain, which can lead to better performance in practice. Randomized methods are prone to "missing" larger regions of X's support, as shown in Figure 2 The driving question behind our work is:
Is it possible to obtain the "best of both worlds" for fitting functions over low- -indeiinsiooal domains? Le., can we match or improve on the strong theoretical guarantees of leverage score sampling with method that produces spatially well-distributed samples?
We answer this question in the affirmative. Instead of sampling rows from A independently with probability proportional to the leverage scores, we adopt tool from survey sampling known as pivotal sampling (Deville and Tille| 1998). Our specific version of pivotal sampling is spatially- aware, meaning that it covers the domain in a well-balanced manner, while the marginal probabilities remain proportional to the leverage score. At a high-level, the pivotal method a "competition" based sampling approach, where candidate rows compete in a binary tree tournament. By structuring the tournament so that spatially close points compete at lower levels (we use a novel recursive PCA procedure to build the tree), we ensure better spatial spread than Bernoulli leverage score sampling.
We show that our pivotal sampling method matches or beats the complexity of independent leverage score sampling in theory, while performing significantly better in practice. On the practice side, offer Figure [ as an example from a PDE test problem. In comparison to independent sampling, our spatially-aware method obtains much better approximation to the target for fixed number of samples (more details in Section 4. On the theory side, we prove two results. The first is general: we show that, as long as it samples rows from A with marginal probabilities proportional to the leverage scores, any sampling strategy that obeys weak "one-sided foe independence" condition (which includes pivotal sampling) matches the complexity of independent leverage score sampling:
Let E TRuxd be a data matrix and b E Th be a target vector. Consider any algorithm which samples exactly k rows from A (and observes the corresponding entries in b) from a distribution that 1) satisfies one-sided foo independence (Defit
Theorem 1.1.
3.7) with parameter Dig and 2) the murginal probability of sampling any row is proportional to Tę. Let A and bo be the scaled sampled data and target, as defined in Section 1.1 and let x* arg minxER=
bl1.. As long as k 2 c. (dlog d. D2 + d. Din) for fixed positive constant then with probability 99/100,
A b1B (1 +]]lxXX bb2
(1.3)
One-sided Eoo independence was introduced in recent paper by Kaufman et al. (2022) matrix Chernoff bounds (we provide a formal definition in Sec. (3.. It is the weakest condition under which tight matrix Chernoff bound is known to hold. For example, the condition is implied with stant D = 0(1) by almost all existing notions of negative dependence between random variables, including conditional negative association (CNA) and the strongly Rayleigh property (Pemantle) 2000). As discussed in Sec. B tight matrix Chernoff bound is prerequisite for proving relative
Formally, we assume that the marginal probability of sampling a: equals min (1, eT)) for fixed constant 2 1. Our proof easily generalizes to the case when some probabilities exceed this hound (since sampling hurts) although the total sample complexity will depend on the sum of the marginal prohabilities.
3

---

1.3 however does not imply such a result alone. Our proof of Theorem requirrs adapting an approximate matrix- moatiix-uulipliaaaiinn method from (Drineas) iet al., 2006 to non-independent sampling. It can be viewed extending the work of Kaufman et al. 2 to show that essentially all sampling distributions known to yield tight matrix Chernoff bounds also yield near optimal active regression bounds in the agnostic setting. Importantly, this includes binary-tree-basee pivotal sampling methods like those introduced in this work. Such meth-
error active learning results like Equation
ods are known to satisty the strongly Rayleigh property (Brändén and Jonasson, 2012), and thus one-sided &o independence with Dint = O1)). So, as a corollary of Theorem 1.1 we obtain:
Corollary 1.1. The spatially-aware pivotal sampling methods introduced in Section (which use fixed binary tree) return with probability 99/ 19/100 a vector X satisfying ||AX* 2 S (1
c)ll - b while only observing O
d + :) entries in b.
We hope that Theorem 1.1 will be valuable in obtaining similar results for other sampling methods beyond our own. However, the result falls short of justifying why pivotal sampling performs better than independent leverage score sampling in experiments. Towards that end, we prove a second result specific to pivotal sampling, which shows that the method actually improves on the complexity of independent sampling by a log factor in the important special case of polynomial regression:
Theorem 1.2. Consider any function b [e, u] R. defined on an interval [L,u] R. and consider fitting b with a degree d polynomial based on evaluations of the function at Ii,.... 2K E [t. u). If T1,.... Fkk are collected via pivotal sampling with leverage score marginals (see Appendix [oor details), then as long as k 2 c. (d+ + ) a fixed positive constant e, there is a procedure that uses these samples to construct a degree dl polynomial p which, with probability 19/100, satisfies:
1p - b12 S (1+0) c)
min degree d polynomial
|Ip - -112-
Here |1S denotes the average squared magnitude Je [(x) of function f-
The problem of finding a polynomial approximation to a real-valued function that minimizes the average square error ||p - 6112 can be modeled as an active regression problem involving a matrix A with d + 1 columns and an infinite number of rows. In fact, polynomial approximation is one of the primary applications of prior work on leverage score-based active learning methods (Cohen and Migliorati 2017 Avron et al. 2019). Such methods require O log d + 9) samples, so The- orem 1.2 Is better by a log d factor. Since polynomial approximation is a representative problem where spatially- distributed samples are important, Theorem 1.2 provides theoretical justification for the strong performance of pivotal sampling in experiments. Our proof is inspired by a result Kane jet al. (2017), and relies on showing tight relation between the leverage scores of the polynomial regression problem and the orthogonality measure of the Chebyshev polynomials on [e, u).
1.3 RELATED WORK
The application of leverage score sampling to the agnostic active regression problem has received significant recent attention. Beyond the results discussed above, extensions of leverage score sam- pling have been studied for norms beyond lg (Chen and Derezinski] 2021 Musco et al. Meyer et al.) 2023 Parulekar et al. (0221) ) in the context where the sample space is infinite (i.e. A is operator with infinite rows) (Erdélyi et al., 2020) Avron et al. 2019), and for functions that involve non-linear
transformations (Gajjar al. 2023 Mai al. 202
Munteanu al. 2018).
Theoretical improvements on leverage score sampling have also been studied. Notable is a recent result that improves on the OOd + d/e) bound by log d factor, showing that the active least squares regression problem can be solved with (dde) samples (Chen and Price 2019). This is provably optimal. However, the algorithm in (Chen and Price 2019) is complex. and appears to involve large constant factors: in our initial experiments. it did not empirically improve on indepen- dent leverage score sampling. In contrast, by Theorem 1.2 our pivotal sampling method matches the theoretical sample complexity of (Chen and Price 2019) for the special case of polynomial regression, but performs well in experiments (significantly better than independent leverage score sampling). There have been few other efforts to develop practical improvements on leverage score sampling. Similar to our work, Derezinski et al.] (2018) study a variant of volume sampling that matches the theoretical guarantees of leverage score sampling, but performs better experimentally. However, this method does not explicitly take into account statial- in the underlying regres- sion problem. While the method from (Dereziński et al. 2018) does not quite fit our Theorem 1.1
4

---

(e.g., it samples indices with replacement) we expect similar methods could be analyzed as a special case of our result, as volume sampling induces a strongly Rayleigh distribution.
While pivotal sampling has not been studied in the context of agnostic active regression, it is widely used in other applications, and its negative dependence properties have been studied extensively. (Dubhash et al. 2007) proves that pivotal sampling satisfies the negative association (NA) property. Borcea et al.] 2009) introduced the notion of a strongly Rayleigh distribution and proved that implies stronger notion of conditional negative association (CNA), and (Brändén and Jonasson showed that pivotal sampling run with an arbitrary binary tree is strongly Rayleigh. It follows that the method satisfies CNA. (Greene et al., 2022) discusses an efficient algorithm for pivotal sampling by parallelization and careful manipulation of inclusion probabilities. Another variant of pivotal sampling that is spatiall--wwar
proposed in (Grafstrrmm et al. 2011). Though their approach is out of the scope of our analysis as it involves randomness in the competition order used during sampling, our method is inspired by their work.
1.4 NOTATION AND PRELIMINARIES
Notation. We let [n] denote {1. .... E E denotes the expectation of random variable X. We use bold lower-case letters for vectors and bold upper-case letters for matrices. For vector E IR with entries 21,-.. Zn. ||z]|2
(Ei=1 z 2/21/2 denotes the Euclidean norm of z. Given a matrix A E Raxd, we let a, denote the i-th row, and aay denote the entry in the i-th row and j-th column.
Importance sampling. All of the methods studied in this paper solve the active regression problem by collecting a single random sample of rows in A and corresponding entries in b. We introduce a vector of binary random variables& & = (L,, ). where & is if a; (and thus b;) is selected, and O otherwise. f1.... En will not necessarily be independent depending on our sampling method. Given a sampling method, let p; = [[[,] denote the marginal probability that row is selected. We return an approximate regression solution as follows: let A E Rixd contain a./ /P: for all é such that & = 1, and similarly let b E g* contain bs/ /P5 for the same values of i. This scaling ensures that, for any fixed x, [P| Ax bl = || blß2. To solve the active regression problem, we return X = arg min |Ax b 1 Computing x only requires querying target values b.
Leverage Score Sampling. We consider methods that choose the marginal probabilities propor- tional to A's leverage scores. Specifically, our methods sample row a, with marginal probability p; min(1, Cp -Ta,, where eą is chosen so that E p:- k. Details of how to find e are discussed in Appendix A
We note that we always have C₁ 2 k/d since E Pi S Zim Ti r 2/2 d k.
2
OUR METHODS
In this section, we present our sampling scheme which consists of two steps; deterministically con- structing binary tree, and choosing samples by running the pivotal method this tree. The pivotal method is described in Algorithm [ It takes as input a binary tree with n leaf nodes, each corre- sponding single index i to be sampled. For each index, also have associated probability P;- The algorithm collects a set of exactly l samples S where = Σ/= Pi:- does so by perco- lating up the tree and performing repeated head-to-head comparisons of the indices sibling nodes in the tree, After each comparison, one node promotes to the parent node with updated inclusion probability, and the other node is determined to be sampled or not to be sampled.
It can be checked that, after running Algorithm | index is always sampled with probability p.: regardless of the choice of T. However, the samples collected by the pivotal method are not dependent, but rather negatively correlated: siblings in T' are unlikely to both be sampled, and in general, the events that close neighbors the tree are both sampled are negatively correlated. In particular, if index i could at some point compete with an index j in the pivotal process, the chance of selecting j decreases if we condition on i being selected. We take advantage of this property to generate spatially distributed samples by constructing binary tree that matches the underlying geometry of our data. In particular, assume ally be used to construct a regression matrix A E gaxd via feature transformation (c.g. by adding polynomial features). However, we construct the sampling tree using X alone.
given set of points X Rnxd', X will eventu-
5

---

Algorithm Binary Tree Based Pivotal Sampling (Deville and Tille 1998)
Input: Depth t full binary tree with nt leaves, inclusion probabilities {pi, . . for each leaf. Output: Set of k sampled indices S.
I: Initialize S = 6.
2: while T has least two remaining children nodes do
3: Select any pair of sibling nodes Si, S2 with parent P. Let i, be the indices stored at S1,S2.
4: if px + Pj then
5: With probability p.tpy"' set ps +- P: + Ps. Pj + 0. Store at P.
6: Otherwise, set Pj < Pi + Pý. Pz +- 0. Store at P.
7: else if Pr + Pj > 1 then
8: With probability 1-9. 2-p:-p; set p; t p. 1. P3 = 1. Store at Pand set + SU()).
9: Otherwise, set β3- < Pi + P3 1, P; < Store j at P and set S < U {i}.
Remove S2 from T.
11: return S
Algorithm Binary Tree Construction by Coordinate or PCA Splitting
Input: Matrix X E Raxd', split method E [PCA, coordinate). inclusion probabilities #1. ........................... Pn.
Output: Binary tree T where each leaf corresponds to a row in X.
1: Create a tree T with single root node. Assign set R to the root where R = [m]: < 1)..
2: while There exists node in T that holds set K. such that [X]] > 1 do
3: Select any such node N and let be its level in the tree. Construct X X(K) R K
4: if split method = PCA then
5: Sort X (x) according the direction of the maximum variance.
6: else if split method = coordinate then
7:
Sort X(x) according to values in its ((t mod d') + 1)-th column.
8: Create left child of N. Assign to it all indices associated with the first [ rows of X (K))
9: Create a right child of N. Assign the all remaining indices in X(x)- Delete K. from N.
10:
return T
Our tree construction method is given Algorithm 2 Xx x denotes the subset of rows of X with indices in the set K. First, the algorithm climinates all data points with inclusion probability p. = 1. Next, it recursively partitions the remaining data points into two subgroups of the same size until all the subgroups have only one data point. Our two methods, PCA-based and coordinate-wise only differ in how to partition. The PCA-based method performs principal component analysis to find the direction of the maximum variance and splits the space by hyperplane orthogonal to the direction so that the numbers of data points on both sides are equal. The coordinate-wise version takes a coordinate (corresponding to column in X) in cyclic order and divides the space by hyperplane orthogonal to the chosen coordinate. An illustration of the PCA-based binary tree construction run on fine uniform grid of data points in R² is shown in Figure 3 Note that tree construction method ensures the the number of indices assigned to each subgroup (color) at each level is equal to with 1 point. As such, we end of with an even partition of data points into spatially correlated sets. Two indices will be more negatively correlated if they lie in the same set at a higher depth number.
3 THEORETICAL ANALYSIS
As will be shown experimentally, when using probabilities pI.....P, proportional to the statistical leverage scores of A, our tree-based pivotal methods significantly outperform its Bernoulli coun- terpart for active regression. We provide two results theoretically justifying this operation. We first show that, no matter what our original data matrix X is, and what feature transformation is used to construct A, our methods never perform worse than Bernoulli sampling. In particular, they match independent leverage score sampling. This result is stated as Theorem Its proof is given to Appendix B but we outline main approach here.
the (d log dl d/e/)) sample complexity of
Following existing proofs for independent random sampling (e.g. (Woodruff) 2014) Theorem I.]
requires two main ingredients: subspace embedding result, and an approximate matrix-vectes multiplication result. In particular, let U e Radd be any orthogonal span for the columns of A. Let S E k be subsampling matrix that contains a row for every index selected by our sampling

---

depth-1
depth-2
depth-3
depthm4
depth-5
depth-6
depth-7
depth=8
Figure 3: Visualization binary tree constructed via our Algorithm using the PCA method for matrix X E 1xxx containing points a uniform square grid. For each depth, data points are given the same color they compose a subtree with root at that depth. As we can see, the method produces uniform recursive spatial partitions, which encourage spatially separated samples.
scheme, which has value 1/V at entry i, and is 0 everywhere else. So, in the notation of Theorem 1.1 = SA and Sb.. To prove the theorem, suffices to show that with high probability,
1. Subspace Embedding: For all x E Rd, |||||| S
|Sx S 1.5//x||e.
2. Approximate Matrix- Multiplication: |UUSS' S(b - Ax Ax)) 5 c Ax" l1..
The first property is equivalent to |U"S""U-III I 1/2. Le., after subsampling, U should remain nearly orthogonal. The second property requires that after subsampling with S, the optimal residual, b Ax", should have small product with U. Note that without subsampling, ||U (b- Ax" )13 = 0.
We show that both of the above bounds can be established for any sampling method that 1) samples index i with marginal probability proportional its leverage score 2) homogeneous, meaning that it takes fixed number of samples k, and 3) produces distribution over binary vectors satisfying the following property:
Definition 3.1 (One-sided fec-independence), Let Êi..... .n {0, 13"h be random variables with joint distribution p. Let S C [n]] and let i, E [n] Define the one-sided influence matrix I as:
(1, = Pr(s; 1/5 = AE = TVE E) - Pr(c, 11
= 11 14l E)
Let | moo = maxie(m) Eje(m) (i,))) )- ps is one-sided Lo -esiindepenllnnt with param. Dins if. for all S C [n]. || Diff Note that if ,I,,,,, En are truly independent, we have Dinf
To prove Theorem 1.1 our required subspace embedding property follows immediately from cently established matrix Chernoff-type bounds for sums random matrices involving one-sided lo -iiieeppendent random variables in (Kaufman et al.) 2022). In fact, this work is what inspired us to consider this property, as it is the minimal condition under which such bounds are known to hold. Our main theoretical contribution is thus to prove the matrix-vector multiplication property. We do this by generalizing the approach of (Drineas et al. 2006) (which holds for independent random samples) to any distribution that satisfies one-sided f lox independence.
3.1 IMPROVED BOUNDS FOR POLYNOMIAL REGRESSION
We do not believe that Theorem I.I can be strengthened in general, as spatially well-spread samples may not be valuable all settings. example, in the case when = X. the points in X live dimensional space, and we are only collecting ((d log d) samples, so intuitively any set of points is well-spread. However, we can show that better spatial distribution does offer provably better bounds in some settings where X is low-dimensional and A is high dimensional feature transformation. In particular, consider the case when X is a (quasi) matrix containing just one column, with an infinite number of rows corresponding to every point t in an interval [e, u). Every row of A is polynomial
feature transformation, meaning of the form a, = [1, t, t ..... td1 for degree d. Now, consider target vector b that can be indexed by real numbers in [e, u). Solving the infinite dimensional

---

regression problem minx |Ax blž minx Sz (aax" bt) is equivalent to finding the best polynomial approximation to b in the l2 norm, a well studied application of leverage score sampling Cohen and Mighorati (2017). For this problem, we show that pivotal sampling obtains sample complexity of O(d/e), improving on independent leverage score sampling by a log d factor.
This result is stated Theorem 1.2 and proven in Appendix C Importantly, we note that the de- pendence on dllog in the general active regression analysis comes from the proof of the subspace embedding guarantee the required approximate matrix-muliplication guarantee already follows with O(d/e) samples. Our proof eliminates the log d by avoiding the use of a matrix Chernoff bound entirely. Instead, by taking advantage of connections between leverage scores of the polyno- mial regression problem and the orthogonality measure of the Chebyshev polynomials, we directly use tools from polynomial approximation theory to prove a subspace embedding bound that holds deterministically with just ((d) samples. Our approach is similar to recently result of Kane et al. 2017), which also obtains O(d/c) sample complexity for active degree-d polynomial regression, albeit with sampling method not based on leverage scores.
4 EXPERIMENTS
We experimentally evaluate our pivotal sampling methods on active regression problems with low- dimensional structure. The benefits of leverage score sampling over uniform sampling for such problems has already been established (see e.g. (Hampton and Doostan 2015) or (Gajjar et al.] 2023) and provide additional evidence in Appendix D So, we focus on comparing our pivotal methods to the widely used baseline of independent Bernoulli leverage score sampling.
Test Problems. We consider several test problems inspired by applications parametric PDEs. In these problems, we are given a differential equation with parameters, and seek to compute some observable quantity of interest (QoI) for different choices of parameters. This can be done directly by solving the PDE, but doing so is computationally expensive. Instead, the goal is to use a small number of solutions to the PDE to fit a surrogate model (in our case, a low-degree polynomial) that approximates the Qol well, either over given parameter range, or on average for parameters drawn from a distribution (e.g., Gaussian). The problem is naturally an active regression problem because we can choose exactly what parameters to solve the PDE for. The first equation we consider models the displacement x of a damped harmonic oscillator with a sinusoidal force over time
dx
dr
dx (0) F1. dt
dta (1)
+ dt
(1) + kx()) cos(wt), x00)= To.
The equation has four parameters; damping coefficient c, spring constant k. forcing amplitude f, and frequency W. As a Qol, we consider the maximum oscillator displacement after 20 seconds. We fix
c, f 0.5, and seek to approximate this displacement over the range domain x w = [1. 3] x [0. 2].
We also consider the heat equation for values of æ E [0. 1] with a time-dependent boundary equation and sinusoidal initial condition parameterized by frequency W. The heat equation that we consider describes the temperature fxx, t) by the partial differential equation.
0ʃ 8² Ot Oz2'
f(0,t) = 0, f(x, 0) = sin(war), xe + 0f(1,t) Ot
0.
Leviegee the Netol
Poote Sestering Mehood
Leateri Our Mehoo
(a) Damped Harmonic Os- (b) Heat Equation, degree (c) Surface Reaction, de- (d) Surface Reaction, de- cillator. degree p 12. p 12.
gree p 15.
gree p 20.
Figure 4: Results for active polynomial regression for the damped harmonic oscillator QoL, the heat equation Qol, and the surface reaction model with polynomials varying degree. Our leverage- score based pivotal method outperforms standard Bernoulli leverage score sampling, suggesting the benefits of spatially-aware sampling.
8

---

1.00
(a) Target Function.
(b) Bernoulli Sampling.
(c) Pivotal Sampling (our method).
Figure 5: Polynomial approximation the maximum temperature heat diffusion problem, as function of time and starting condition.
(a) is the target value and both (b) and (c) draw 240 samples using the leverage score and perform polynomial regression of degree 20. However, (b) uses Bernoulli sampling while (c) employs our PCA-based pivotal sampling.
40
(a) Leverage Score.
(b) Target Function.
(c) Bernoulli Sampling.
(d) Pivotal Sampling.
Figure 6: Approximation of the surface reaction model when 425 samples used to fit degree 25 polynomial. To help readers focus on the area near the origin, we draw [-10, 10]² box.
As Qol, we estimate the maximum temperature over all values of æ for t E [0,3] and E [0, 5].
Data Matrix. For both problems, construct A by uniformly selecting n 10 data points in the 2 2-dimensional parameter range of interest. We then add all polynomial features of degree
12
as discussed in Section 1.2 We compute sampled entries from the target vector using standard MATLAB routines. Results comparing our PCA-based pivotal method and Bernoulli leverage score sampling are show in Figure 4 We report median normalized error ||A b [[]]] after 1000 trials. By drawing more samples from b, the errors of all methods eventually converge to the optimal error ||Ax* b / /|bl. but clearly the pivotal method requires less samples to achieve given level of accuracy, confirming the benefits of spatially-awaid
sampling. We also visualize results for the damped harmonic oscillator in Figure [ showing approximations obtained with 250 samples. Visualizations for the heat equation given in Figure For both targets, can directly that pivotal sampling improves the performance over Bernoulli sampling.
We also consider chemical surface coverage problem from (Hampton and Doostan, 2015). Details are relegated to Appendix[D] Again, we vary two parameters and seek to surrogate model using small number of example pairs of parameters. Instead of uniform distribution, for this problem, the rows of X are drawn from a Gaussian distribution with O-mean and 7.5 standard deviation. We construct using polynomial features of varying degrees. Convergence results are shown in Figure 4 and the fit visualized for data near the origin in Figure ( This is a challenging problem since the target function has sharp threshold behavior that difficult approximate with polynomial. How- ever, our pivotal based leverage score sampling performs well, providing better fit than Bernoulli sampling. Additional experiments, including on 3D problems are reported in Appendix D
5 CONCLUSION AND FUTURE WORK
In this paper, we introduce spatially-aware pivotal sampling method and empirically demonstrate its effectiveness for active linear regression for spatially-structure problems. We prove a general theorem that can be used to analyze the number of samples required in the agnostic setting for any similar method that samples with marginal probabilities proportional to the leverage scores using a distribution that satisfies one-sided fox independence. We provide stronger bound for the important special case of polynomial regression, showing that our method can obtain a sample complexity of (d/e), removing log(d) factor from independent leverage score sampling. This result provides initial theoretical evidence for the strong performance pivotal sampling in practice. Extending other spatially structured function classes would be an interesting direction for future work.

---

ACKNOWLEDGEMENTS
Xiaoou Cheng and Jonathan Weare were supported by NSF award 2054306. Christopher Musco was supported by DOE award DE-SC002226 and NSF award 2045590. We would like to thank NYU IT for the use of the Greene computing cluster.
REFERENCES
Haim Avron, Michael Kapralov, Cameron Musco, Christopher Musco, Ameya Velingker, and Amir Zandieh. Random Fourier features for kernel ridge regression: Approximation bounds and sta- tistical guarantees. In Proceedings of the 34th International Conference Machine Learning (ICML), 2017.
Haim Avron, Michael Kapralov, Cameron Musco, Christopher Musco, Ameya Velingker, and Amir Zandich. A universal sampling method for reconstructing signals with simple fourier transforms. In Proceedings of the 51st Anmual ACM Symposium Theory of Computing (STOC), 2019.
Julius Borcea, Petter Brändén, and Thomas M. Liggett. Negative dependence and the geometry of polynomials. Journal of the American Mathematical Society, 22(2):521-567, 2009.
Peter Borwein and Tamás Erdélyi. Polynomials and Polynomial Inequalities. Springer New York, NY, 1995.
Petter Brändén and Johan Jonasson. Negative dependence in sampling. Scandinavian Journal of Statistics, (9(4) 838, 2012.
Xue Chen and Michal Derezinski. Query complexity of least absolute deviation regression via robust uniform convergence. In Proceedings of the 34th Annual Conference on Computational Learning Theory (COLT), 2021.
Xue Chen and Eric Price. Active regression via samele-sample sparsification. In Proceedings of the 32nd Anmual Conference on Computational Learning Theory (COLT), 2019.
Abdellah Chkifa, Albert Cohen, Giovanni Migliorati, Fabio Nobile, and Raul Tempone. Discrete least squares polynomial approximation with random evaluations application to parametric and stochastic elliptic PDEs. ESAIM: M2AN, 49(3):815-837, 2015.
Albert Cohen and Ronald Vore. Approximation of high-dimensional parametric pdes. Acta Numerica, 24:1-159 2015.
Albert Cohen and Giovanni Migliorati. Optimal weighted least-squares methods. SMAI Journal of Computational Mathematics, 3:181-203, 2017.
Michal Dereziński. Manfred K. Warmuth, and Daniel Hsu. Leveraged volume sampling for línear regression. In Advances in Neural Information Processing Systems 31 (NeurIPS), pages 2510- 2519, 2018.
Jean-Claude Deville and Yves Tille. Unequal probability sampling without replacement through a splitting method. Biometrika, 85(1):89-101, 1998.
Petros Drineas, Ravi Kannan, and Michael W. Mahoney. Fast Monte Carlo algorithms for matrices I: Approximating matrix multiplication. SIAM Journal Computing, 36(1):132-157) 2006.
Devdatt Dubhashi, Johan Jonasson, and Desh Ranjan. Positive influence and negative dependence. Combinatorics, Probability and Computing, 16(1):29-41, 2007.
Tamás Erdélyi and Paul Nevai. Generalized Jacobi weights, Christoffel functions, and zeros of orthogonal polynomials. Journal of Approximation Theory, 69(2):111-133 1992.
Tamás Erdélyi, Cameron Musco, and Christopher Musco. Fourier sparse leverage scores and ap- proximate kernel learning. In Advances in Neural Information Processing Systems 33 (NeurIPS), 2020.
10

---

Aarshvi Gajjar, Chinmay Hegde, and Christopher Musco. Active learning for single neuron models with Lipschitz non- oon-linarities. In Proceedings of the 26th International Conference on Artificial Intelligence and Statistics (AISTATS), 2023.
Anton Grafström, Niklas L P Lundström, and Lina Schelin. Spatially balanced sampling through the pivotal method. Biometrics, 68(2):514 54--520, 2011.
Samuel M. Greene, Robert Webber, Timothy C. Berkelbach, and Jonathan Weare. Approximating matrix eigenvalues by subspace iteration with repeated random sparsification. SIAM Journal Scientific Computing, (4((5):3367-A3000 2022.
Ling Guo, Akil Narayan, and Tao Zhou. Constructing least-squares polynomial approximations.
SIAM Review, 62(2):483-508, 2020.
Jerrad Hampton and Alireza Doostan. Coherence motivated sampling and convergence analysis of least squares polynomial chaos regression. Computer Methods in Applied Mechanics and Engineering, 290:73-97, 2015.
Jeffrey M. Hokanson and Paul G. Constantine. Data-driven polynomial ridge approximation using variable projection. SIAM Journal on Scientific Computing, 1110:/116666111151 2018.
Daniel Kane, Sushrut Karmalkar, and Eric Price. Robust polynomial regression up to the information theoretic limit. In Proceedings of the 58th Annual IEEE Symposium on Foundations of Computer Science (FOCS), 2017.
Tali Kaufman,
Rasmus Kyng, and Federico Soldá. Scalar and matrix Chernoff bounds from Egp independence. In Proceedings of the 33rd Annual CCMM Symposium Discrete Algorithms (SODA), 2022.
Nora Lüthen, Stefano Marelli, and Bruno Sudret. Sparse polynomial chaos expansions: Literature survey and benchmark. SIAM/ASA Journal Uncertainty Quantification, 9(2):593-649, 2021.
Tung Mai, Anup B. Rao, and Cameron Musco. Coresets for classification simplified and strength- ened. In Advances in Neural Information Processing Systems 34 (NeurIPS), 2021.
Raphael A. Meyer, Cameron Musco, Christopher Musco, David P. Woodruff, and Samson Zhou. Near-linear sample complexity for Lp polynomial regression. In Proceedings of the 34th Annual ACM-SIAM Symposium on Discrete Algorithms (SODA), 2023.
Alexander Munteanu, Chris Schwiegelshohn, Christian Sohler, and David Woodruff. On coresets for logistic regression. In Advances in Neural Information Processing Systems 31 (NeurIPS), 2018.
Cameron Musco, Christopher Musco, David P. Woodruff, and Taisuke Yasuda. Active linear re- gression for lp norms and beyond. In Proceedings of the 63rd Annual IEEE Symposium on Foundations of Computer Science (FOCS).
Paul Nevai. Bernstein's inequality in 1p for 0 < < 1. Journal of Approximation Theory, 27(3): 239-243, 1979.
Aditya Parulekar, Advait Parulekar, and Eric Price. LI Regression with Lewis Weights Subsampling. In Approximation, Randomization, and Combinatorial Optimization. Algorithms and Techniques (APPROX/RANDOM 2021), pages 49:1-49:21, 2021.
Robin Pemantle. Towards theory of negative dependence. Journal of Mathematical Physics, 41, 06 2000.
Robin Pemantle and Yuval Peres. Concentration of Lipschitz functionals of determinantal and other strong rayleigh measures. Combinatorics, Probability and Computing, 23(1):140-160, 2014.
Friedrich Pukelsheim. Optimal Design of Experiments. Society for Industrial and Applied Mathe- matics, 2006.
Holger Rauhut and Rachel Ward. Sparse Legendre expansions via él-minimization. Journal of Approximation Theory, 164(5):517 1 533, 2012.

---

Tamás Sarlós. Improved approximation algorithms for large matrices via random projections. In Proceedings of the 47th Anmual IEEE Symposium on Foundations of Computer Science (FOCS), 2006.
Keeve M. Siegel. Bounds of the Legendre function. Journal of Mathematics and Physics, 34(1-4) 43-49, 1955.
Lloyd N. Trefethen. Houscholder triangularization of a quasimatrix. IMA Journal of Numerical Analysis, 30(4))887-877, 2009.
Joel A. Tropp. User-friendly tail bounds for sums of random matrices. Foundations of Computa- tional Mathematics, (244):3994334, 2012.
David P. Woodruff. Sketching as a tool for numerical linear algebra. Foundations and Trends in Theoretical Computer Science, 10(1-2):1-177, 2014.
Dongbin Xiu. Stochastic Collocation Methods: A Survey, pages 1-18. Springer International Pub- lishing, 2016.
A PROBABILITY PRE-PROCESSING
As discussed in Section I.4] our sampling methods require computing probabilities pi....,pp where 5i min(1, T;) for some fixed constant eę chosen that > Pr k. We can find such probabilities using simple iterative method, which takes as input initial "probabilities" p; & that are proportional to the leverage scores Ti...., Ty of A. Note that p; could be larger than k d. Pseudocode for how to adjust these probabilities is included in Algorithm
Algorithm Probability Ceiling
Input: Number of samples to choose k; inclusion probabilities (pi,"" .Pa).
Set p; = P; for E [n]].
2: while S {i [n]] > 1} is not empty do
3:
pi +- for all i S.
4: D {i [n]
1). R Eie [nli < 1}.
5:
ps - Eee p. p. for all R.
6: return (. Pn
To see that the method returns Ds with the desired properties, first note that ΣieR P. is always greater than 1. In particular, at the beginning of the while loop, we always have the invariant that Σm P. = k. Accordingly, ΣIER Pi = k ΣiED P. < k - [D). As a result, at the end of the 1. And fact, must be that
algorithm, P. certainly equals min(1,c T²) for some constant by the invariant that ET k.
The main loop in Algorithm 3 terminates as soon as S is empty, which happens after at most steps. In practice however, the method usually converges in 2 or 3 iterations. Nevertheless, our primary objective is to minimize the number of samples required for an accurate regression solution we do not consider runtime costs in detail.
B PROOF OF THBOREM
In this section, we prove our first theoretical result, which holds for any row sampling distribution whose marginal probabilities are proportional to the leverage scores of A, and which satisfies few additional conditions. In particular, we require that the distribution is both k-homogeneous and loo de-independent with constant parameter D. We define the first requirement below, and recall the definition of loc -e--indeenndenc from Section For both definitions, we view subsampling process as a method that randomly selects a binary vector & = {E..... ,.,,, E. from {0, 10". Entries of 1 in the vector correspond to rows that sampled from A (and thus labels that observed
12

---

in b). Under this notation, requirement marginal probabilities is equivalent to requiring that E|G] = min(1, c. Ti) = p; for some fixed constant c, where F is the i-th leverage score of A.
distribution over binary vectors homogeneous if all possible realizations of & contain exactly k ones.
Definition B.1 (Homogeneity).
(-1,--- ,nn is k-
Definition 3.1 (One-sided fe ---independence). Let E1..... .n {0, 1)" be random variables with joint distribution p. Let S - [n] and let i.j [n]]sS. Define the one-sided influence matrix I as:
(i,5) = Pr(Ei 1|:: 1 L IVLE Pr|5, = llee 141 S] 1
Pr(c;
Lef | Im = max Eje(n) IZ (i, is one-sided lox exc-indepenteet with param. Ding if. for all S C [n]. | loo Dift Note that ,,, En are truly independent, we have Dinf 1.
With these definitions in place, we are ready for our proof:
As outlined in Section 3 we need to establish two results: a subspace embed- ding guarantee and approximate matrix-vectur multiplication result. We prove these separately below, then combine the results to complete the proof of the theorem.
Proof of Theorem 1.1
Subspace Embedding. Let S be subsampling matrix corresponding
vector & selected from a distribution with the properties above. That is, S E R contains a row for every index i where 1. That has value 1/ entry i, and is 0 everywhere else. Let U E Rexd be an
orthonormal span for A's columns. We will show that, if k O dloota/0) Dint), with prob.
1 1-6, 6,
JUSS 11e S xx" (U 1) x s SIUSSS |||2 this is
(B.1)
Note that, since ||SUx|| equivalent to showing that, for all Rdd.
||Ux|
(1 a Ux|
||UX||] S []+a)|
We establish BI.I using the following matrix Chernoff bound from (Kaufman et al. 2022):
Lemma B.1 (Matrix Chernoff for lo -independent Distributions). Let E1.-.. En be binary ran- dom variables with joint distribution μ that 2-homogeneous for any positive integer and too- independent with parameter Ding Let Y...... Y, Rauxd be positive semidefinite (PSD) matrices with spectral norm Y I R for some
Let pmax Amax (angseBeg.. (Bę EE Y.Y,)) and
Mmin Ami (Bęup E EL Y.)).. Then, for any (0, 1) and a fixed constant G,
a² Pmax cRD inf
Pr
Amax
Σ Y.Y,
2 (1 am/max
$ dexp
a² Mmin cRD22 inf
Pr
Xmin
Σ Y,Y,
< (1 - a)
S ddcpp
Let D fi : 1} and R = {i E [n] p. < 1}. Letting denote the i-th row of U, UTSTSU can be decomposed as:
y P:
UTS
Eu.u?
VED
Recall that úú = min(1,c.m) c. T;) for some constant c, and in fact, since E p. = k and Σ T; S d, must be that c 2 #. Define p: fT. For each E D, let m, be an integer such that m; pi- For E D, P: = 5T. 2 1, so m, 2 1. Additionally let { " em,
& be variables that are deterministically 1. Then we can trivially split up the second sum above so that:
nnu.
UTS Σ & P:
(B.2)
ER
j=
13

---

It is easy to that {& : R : D.j E [m,]] are fe--independent random variables with parameter Dut and form a distribution that is k homogeneous for k = k + Eie - 1). So, noting that u.u. is PSD,
apply Lemma B.1 to (..2. We are left to bound R and Mmax-
u.ud
To bound R, for i E R, let Y;
For such i, we have that i; 2 P;. So,
P.
|..u
uu
Ti
Ti d
Pi
Pi
k
Le.., for all i E R. || S R for R 1- Additionally, for i E D, let Y;
u.u." m,
Since we chose
my 2 Pr. by the same argument, we have that || S R for R K Next, we bound Mmux by simply noting that
E
=
Σ u.u
+ Σ Lu,u =I
Pi
iER
SED
ER
sED
where I denotes the d x d identity matrix. It follows that Mmaa
Mmin 1.
Plugging into Lemma B. and applying a union bound, we obtain
Pr [|UTSSUU 1| 2 a) S dexp (-
ka DD in
where e' is a fixed constant. Setting O dlog(d/8)) Darr), we have probability at least 1 - 5. This establishes (B.1. a²
| - -I|2 S a with
With our subspace embedding result in place, we move onto proving the necessary approximate matrix-vector multiplication result introduced in Section 3 In particular, we wish to show that, with good probability, || Sb - Ax*) S <llb Ax where x" arg ming | blβ. Since the residual of the optimal solution, Ax", is the same under any column transformation of A, we have b Axx = b Uy* where y" arg miny || bll2- We can then equivalently show that, if we take k O ( Dmn)) samples, with probability - 5,
Approximate Matrix- Multiplication.
S Uy < </bb Uy
(B.3)
As the proof for independent random samples (Drineas et al.] 2006), we will prove (B.3) by bounding the expected squared error and applying Markov's inequality. In particular, we have
[[[UUS S S(b Uy
2 el|b Uy Pr [ 11
(B.4)
b Uy 12
Let z b Uy and note that U z = 0. The numerator on the right side can be transformed as
| Sz| E U Sz U 2 E [| I)2[[]]
Note that above S I is a diagonal matrix with i-th diagonal entry equal to 1 1 if &; = 1 and -1 if = O. Expanding the f--norm and using that & = and i. for i & R, we have
2
27
E S
Σ E
Σ 1=1
E+
E
S 1) miyzi
P
/=1
VER
ΣΣΣ
Co P:pi
JERRRRE
where ea Cov(E,, Et) = BB((s - p.)(E: p.) -E(GEE - p.P..
We will show that
d
22
ΣΣΣ
Ca
Σ
uyzz & Dint
P.
j=1 FERIER
j=1 ER
14

---

where Dutt is the lo -ndeeeddeeeeence parameter. It suffices to show that for every j. we have
" z2
Gi
Dint
(B.5)
VERIER PiP
Pt
EER
Define a symmetric matrix M [
with entries mar
Get VDLvi
and a vector v E IRR with
entries U;
The desired result in B.5 can be expressed as
vpi
Mv I Dirllv/l2-
With S = 0, the one-sided fxc -a-ineependence condition implies that, for all i E [n].
E|
Pt
mal S Dinr.
Pi
IER
ER
Equivalently, if we define diagonal matrix A E IR
[R || such that Aui vpc. have shown:
|AA
"MAl S
Dint.
Bxll= It follows that the largest where for a matrix B, |B|] denotes max; Σ, |aal maxx x
eigenvalue of A MM is at most Dint and thus, the largest eigenvalue of M is also at most Dint Therefore, we have v Mv I Dinr|lv]2.. Considering that p; 2 p. = fT, for E R, we have
27
Wsy
inf
Dinf
Dine
P
Pi
Pt
j=l
VER
SER
SER
iER
d Dint 12ll I.
Recalling that z Uy b, we can plug into (B.4) with 0 ( Dmr)) samples, which proves that (B.3) holds with probability at least - 5.
and B.3 in place, we can prove our main result, Theorem 1.1] We follow similar approach to Woodruff 2014). By reparameterization, proving this inequality equivalent to showing that
Putting it all together. With
1.3) of
|Uý bll < (1 e) Uyy -bll..
(B.6)
where yr = arg miny |SUy Sbll3 and y arg miny ||y - -ll.. Since y" is the minimizer of |[y bll2. have | -b1í 2U (Uy b) y² This indicates that Uy* b is
orthogonal to any vector in the column span of U. Particularly, Uy b orthogonal to Uý* -Uy"
Therefore, by the Pythagorean theorem, have
|Uý* bli = ||yy* - bll2 + ||Uy* - Uy = ||y - b1l2 + 11y
(B.7)
So to prove B.6). it suffices to show that
1y - y" 113 S Uy - -b1
Applying 99.5/100, U S7 SU - I|2 S 2- Then, by triangle inequality.
B.
with k: = Od d- D2 +d Dint/c) ) and 8 = 1/200, we have that, with probability
1y yyll2 < s -y" + I -y") -
(f" yy"lll
SU yy')la | - - yy
< SU'S'U((" - y') + 11 y-
(B.8)
Rearranging, we conclude that y y" <4 of |[Uyy Sb I2 we have Vy SUy Sb|B 2 (SU Sb) Thus,
> Since ÿ is the minimizer
| -y (UUs" Sb Sb SUy*)
= sT - "|2
Applying B33 with 5 1/200 and combining with B.S) using union bound, we thus have that with probability 99/100,
115" - 3 < 4 ST - y
- = U - < 4e| - -bl (B.9)
Plugging into B.7 and adjusting by constant factor completes the proof of Theorem 1.1
15

---

B.1 PROOF OF COROLLARY
We briefly comment on how to derive Corollary 1.1 from our main result. By definition of the method, immediately have that binary-tree-basee pivotal sampling is k-homogeneous, we just need to show that it produces a distribution over samples that is fo--indepeddent with constant parameter Dint- This fact can be derived directly from a line of prior work. In particular, (Brändén and Jonasson, 2012) proves that binary-tree seeaaryttreebasee breeebased pivotal sampling satisfies negative association, (Pemantle and Peres, 2014) proves that negative association implies a stochastic covering property, and (Kaufman et al., 2022) shows that any distribution satisfying the stochastic covering property has -e-iindeendeecc parameter at most D = 2. We also give an arguably more direct alternative proof below based on natural conditional variant of negative correlation.
Proof of Corollary [.7 Brändén and Jonasson 2012) proves that binary-tree-basee pivotal sam-
pling is conditionally negatively associated (CNA). Given a set C S [n] and a vector с E {0. 1)(c] denote the condition & c; for all i C by C. Conditional negative association asserts that, for all C, any disjoint subsets S and T of (Ei, ; ). and any non-decreasing functions f and g.
...).
[[[(S)CC) [[G((T)C) > ElS(S)9(T)IC)
When S and T are singletons and f and g are the identity functions, we have
E[E/C]] Elc, [C 2 E[G,
(B.10)
Since R[G[C]- Pr[g: = 1|C], we also have
Prig; = 1[C1 Pr(c, = 1/C] 2 Pr(c, = 1 A E, = 1|]]
Pr[s. = 1 1AE; 1[ = Pr(si = 15, = LAC
(B.11)
Pr[s: = 1|C] 2
Pr(k, = 1|C]
In words, the entries of our vector & are negatively correlated, even conditioned on fixing any subset of entries in the vector. We will use this fact to show that EjE(n] [ (i, < 2 for all i E [r]]. where I is as defined in Definition 3.1 For a fixed i, let qi = PIę-z & 1 1AL E S). Then, we have |T (i.5) = for j E S, I (1.j) = 9 for j = i, and Σ _] Suf{) [T (t. j) = - di- The last fact follows from k-homogeneity, i.c. that Σi= 9i = k. and (B.11). which implies that I (i.)) S 0 for all j in [n]] {i}, so Σjennouu (i,)) = Σ ) IZ (i.o)) Thus, have Σje(s) (1,5) 2q 2 2.
C PROOF OF THEOREM
In this section prove Theorem which shows that pivotal obtains better sample complexity for polynomial regression on an interval than independent leverage score sampling. Since we can always shift and scale our target function, without loss of generality we take [e, u] to be the interval [-1, 1j. We will sample from the infinite set of points on this interval. Each point / E [-1, will correspond a row in a regression matrix A with d+1 columns and infinite number of rows. Such an object is sometimes referred to as a quasimatrix (Trefethen) 2009 The row with index t equals ag- 1,t, t2 t ). The leverage score for the row with index is defined analogously to
Equation 1.2 as:
(xTa,)2
(t)) = max XER
(C.1)
(xTaa) )ds
Note that S4 (()dtt = (d + 1) since A has d + 1 columns. We will sample l points from [-1,1] interval with the probability of sampling point & to be proportional to r(t). To do so, we renormalize T(!) and consider sampling point with probability proportional kr(() This is the analog of sampling with probability in the discrete case, since [_
d+1
kr(()) di k. Then apply pivotal
We refer the reader to Avron al. 2011 Frdélyi et al. (020)) or Chen and Price 2019 for more in depth treatment of leverage score sampling and active linear regression for quasimatrices.
16

---

sampling in the infinite point limit, where we choose the pivotal competition order that points in [-1,1] compete with each other from left to right across the interval. In this limit, the probability carried along in the competition will, at one point, accumulate to exactly 1. This defines the interval I whose left endpoint is -1 and its right endpoint is defined so that J1, d+1
kT(f) dt 1. This means
that the winner in the last competition in I will be sampled. A new competition then starts for the next point, until the renormalized leverage score accumulates to again and the winner is sampled. We can define [i 1, 1,2,..... k) to be adjacent intervals with J1, d+I that pivotal sampling in the infinite point limit will always sample exactly one point from each of J,,..., Tx- This is actually the only property of pivotal sampling we will need to prove Theorem 1.2 although we note that, within an interval 4, point selected with probability proportional to (t)/ Jr, (ss)ds.
kr(t)
It can be seen
To prove sample complexity bound on (dd/e) for polynomial regression, it suffices to prove that a constant factor subspace embedding guarantee holds when selecting ((d) samples from [-1,1] using the pivotal sampling process above. Observe that the required e-error approximate matrix multiplication guarantee from DectionBBalreae holds with O(d/e) samples collected via the pivotal method with leverage score marginals, and this analysis extends to the infinite quasi-matrix setting (see Appendix C.T). is only the subspace embedding guarantee that adds an extra log d factor to the sample complexity. Since this factor is inherent to our use of matrix Chernoff bound in analyzing the general case, to prove a tighter bound for polynomial regression we use a direct analysis that avoids matrix Chernoff entirely. In particular, our main result of this section is as follows:
Theorem C.I. Let k O (+) points L,,..., l be selected from [-1, 1] via pivotal sampling. Specifically, exactly one point t; is sampled from each interval II..... Ig with probability propor-
tional to its leverage score. We have (deterministically) that for any degree d polynomial p. with r() w()) d+
[-peffde p(t)
p( w(i))
S S p(t)dtt.
(C.2)
To translate from the notation above to the notation used in B.I note that 1_ ptt = ||Ux|] | where is an orthogonal span for the quasi-matrix A defined above, and x is a + 1 dimensional vector containing coefficients of the polynomial p the basis U. Correspondingly, & Σ p(t.)²
w((.)
|S
Proof of Theorem 1.2
With Theorem providing the required subspace embedding guarantee from (B.1). the proof of Theorem 1.2 is essentially identical to the proof for the discrete case (Theo- rem L.I. The only remaining requirement is the required approximate matrix-vector multiplication
guarantee from (B.3. For completeness, we provide direct proof in the polynomial regression case
in Appendix C
L
We build up to the proof of Theorem C.1 by introducing several intermediate results. Our approach is inspired by a result of Kane et al. (2017) which also obtains O(dfe) sample complexity for active degree-d polynomial regression, albeit with a sampling distribution not based on leverage scores. In particular, they prove the same guarantees Thcorem C.I but where points are selected uniformly from k intervals that evenly partition the Chebyshev polynomial weight function 1/vI- t2. As we will see, the leverage scores for polynomials closely approximate this weight function. This connection is well known, as the leverage score function from (C.I) is exactly proportional to the inverse of the polynomial "Christoffel function" " under the uniform measure, a well-studied function in approximation theory.
For simplicity, from now on we denote f(t) = p(t)² and describe the main structure of the proof of Theorem C.I The left hand side of (C.2) can be decomposed into a sum of errors in individual intervals as
[ s())t
S(4)) ke wl(4))
1(1) s(t) w(() w(k))
Σ
<
w(f)dt.
(C.3)
17

---

In each interval I, noting that by the definition of 4, Jr, w(t)dt f.c So we have
h SQ f(4:)
w(t)dt S. (()))(()--((0 (()) ds w(t)dt
w()) w(t,)
w²(s)
- S1. (J,, (F)(s))(ss) + f(s)
10 (s)) do) w(1)dt
(C.4)
wx(s)
= /. (0((0)
1 [(1))d(() dt.
dt + 1: J1. (²)
w(t)
Therefore,
J f(t)dt
1 f(t)) k wtt))
S /
[((t] wt)
{ f(1))u (t) w²(t) dt.
Σ
dt
(C.5)
Then, proving the upper bound Theorem C.1 boils down to establishing lower bound for w()), an upper bound for Iw'(t), and the connection between S4 to the lower bound for wt)).
[f(()) dt and f_, f())dt, which is related
((f)
Both [(() and w()) are polynomials. This is key to providing both pointwise and (weighted) integral upper bounds of their derivatives by their function values and integrals. polynomial approximation theory, Bernsseenn inequalities address exactly this point.
We begin with the integral form of the Markov-Bernster inequalities from Borwein and Erdélyi 1995) and Nevai (1779).
Proposition C.2 (L² Bemsteen's inequality). For any degree d polynomial p(t), with universal constant Co.
( (p((t dt < Cod S lp(t)ldt.
(C.6)
Proposition C.3 (L* Markov's inequality). For any degree d polynomial p(1), with a universal constant C1.
S [woaass sci sc d2 S", IV())dd.
(C.7)
The integral on the left hand side of (C.6) is weighted by VI- This weight is negligible near the boundary of the interval [-1, 1]. which, to some extent, explains the milder O(d) dependence on the right hand side of (C.6) compared to the (d²) dependence in (C7).
C.6) and (C.7) also indicate that a lower bound on w(t) of the form VI- or a constant can cooperate with [t)) to yield an integral bound. This is indeed realizable, and we will prove that the lower bound on w(t) is of the two different forms in different parts of the interval [-1, It roughly the minimum of these two, as vI-zz blows up near the boundary. We define the middle
region be a centered subinterval in [-1, that is z away from the boundary, and the boundary region to the two subintervals of length zz near the boundary. The exact expression of these two regions will be clear when we state different lower bounds of (t)). There we note that these two regions are overlapping so that we get a lower bound on w(f) on the whole interval. Matching this pattern of lower bounds, we will also apply different upper bounds of w' (t) in these two regions.
Lower bound on (t) in the middle region. Erdélyi and Nevai 1992) states the following bound on Ttt), which shows the relation between the leverage score and the Chebyshev polynomial weight function. See Section 4.3 in Meyer et al. (2023) for further discussion of this result.
Proposition C.4. With a constant C2 > 0.
T(t) 2 C2(d+1 1)
C2
w(t)
for |t| S 11- (d 1)2
(C.8)
1
18

---

Lower bound w(t) in the boundary region. This requires bit more effort than the case of the middle region. Inspired by (C.7). we aim for the following bound on r(t) and w(t).
Proposition C.5. With constants C3 0, 2'
(tt 2 C3d(d+ 1), w(t) 2 Czd, for [t] >>> 1
(C.9)
d²
We specify
9 because this implies that, at least for d large enough, the middle region, defined as and the boundary region, defined as |-1, - +
[-1
9 (d-1)²
7 U - 1/
have overlap.
To prove Proposition C.5 we need employ the explicit expression for 7(l) in terms of Legendre polynomials and apply properties of Legendre polynomials. We denote the unnormalized degree Legendre polynomial by Pa(t) and we fix Pa(1) 1. We use La(t) to denote the normalized Legendre polynomials satisfying _-, () = 1. There are classical, explicit bounds for La(t)
and its derivative on the interval [-1, 1].
Lemma C.6 (See e.g. Siegel 1955)). For It| < 1, the normalized degree d Legendre polynomials L.d, satisfy
d(d+1)
[L((t)] S [a(1) = - d + 2' [E(001 S
d
(C.10)
2
2
Proof of Proposition C
Since rtt) and wo()) are related as w((t) =
a+1 we only need to prove the
-(t)
bound for ((). The proof follows from the fact that, (() can be equivalently written the squared norm the t-th row in any orthonormal basis for A,
d
T(t) = Σ I2())
(C.11)
See e.g. Equation (4) in Meyer et al. (2023) for this equivalency. We will then lower bound r(t) by deriving lower bounds for the individual summands, ultimately arguing that for [4] 1 z.. r(t) 2 d². Lower bounds on individual (t) can be derived with an upper bound on the derivative,
since we know at the boundary (1) = / exactly and [t( is an even function. In particular,
by Lemma C.6
((2) (2T. ili (6
Therefore, for any > and [t]
-e. (t (i + } i(i )t + 1) c
For given e, (L) is thus positive and of order 2 as long as i S d/c' for a sufficiently large constant ç' > A careful algebraic manipulation will show that we can choose 3. With the classical formula for summing a linear growth sequence, we have for some constant C3.
d
a/d
7(t) = Σ 2,02) Σ 17(1) Cod(d+ 1).
Now we move on to upper bounds on w() in the two regions matching the form of w(t) so that some cancellation can occur to leave clean integral Sr, [()ddt in (C.5). The classical Markov- Bernstein inequalities in Borwein and Erdelyi (1955) bound the pointwise value of the derivative of a polynomial by its maximum function value.
Lemma C.7 (Bernstein's 's inequality). For any degree d real polynomial p(t).
d sup (p(t),, for < <.
1p((t) < V1 1
(C.12)
tE[-1.1]
19

---

Lemma C.8 (Markov's inequality). For any degree d real potynomial p(t),
In' (1) S d² sup Ip((), for - - 6E|-1,11
(C.13)
Both Lemma C.7 and Lemma C.8 are valid on the whole interval [-1, 1]. but one of them is tighter than the other depending where 1 sits. Although w (t) itself is degree 2d - 1 polynomial, it turns out that the bounds directly given by Lemma C.7 Lemma Lemma[.8] are too crude. However, we can unwrap more structure in (t) by applying the recurrence relation of Legendre polynomials.
Proposition C.9 (Neat expression of r())..
7 (1) (1))(())(())
(C.14)
Proof. We write (1) in terms of the unnormalized Legendre polynomials Pn as
1 PP?
(())
Σ Lİ Σ (i
O
k=0
Pa and 1a are related by the recurrence relation (2n + 1)Pp Pi+ Ph-1 with Po so
7'(1) Σ 22 + 1)P.P
Pi(P P--1) P Potx
Upper bound on w/ (t) in the middle region. Proposition C.10.
d(d + 1) 1- 2
d
7()1 S
[w)(t))
1 2 for - 1<t<l
(C.15)
Proof. This follows from Proposition C.9 and Bernstein's inequality in Lemma C.7 for Legendre
polynomials Pa and Pat1 by noting that the maximum of these polynomials are P. 1, Vi.
Upper bound w/ (l) in the boundary region. The bound in Proposition C.10 will blow up near
the boundary. This is caused by the blowing up of Bernstein's 's inequality in Lemma C.7 We should
use the tighter upper bound offered by Markov's inequality in Lemma C.8 and a similar proof will give the following bound.
Proposition C.11.
(+())) < ddd+1) [w'(t)) < dd(dddd)))) 1), for -
(C.16)
As mentioned, Proposition C.10 and Proposition C.I1 cannot be obtained by treating 7(t) as general 2d polynomial. get instead
In fact, directly applying Lemma C.7 and Lemma C.8 r(),, we will
I0(1 (2d - (dd-1)(d=)) 1)a 2 f
1 / (d 1)3,
and |7(0)) S 1
which is significantly worse in the middle region.
With the previous intermediate results in place, we are now ready to prove our main claim, following the proof of Lemma 2.1 in Kane et al. (2017) for intervals defined by the Chebyshev measure.
Proof of Theorem
C.I
Recall we start with CC.5) as
L f())tt
1 f(ta)
(f) wl
dt + 1 1.5² k. w=(())
(())lw((()
<
dt.
(C.17)
t: v((l))
20

---

We will use different pointwise bounds of w(t) and w') in different parts of the interval. For the middle region, by Proposition C.4 and Proposition we have
(f(t) v(l)
2xCod
[ f(t)dt.
dt
VI 2 (t)ldt -
(C.18)
C2k
C2k
Also, by Proposition C.4 and Proposition
C.10
[w)( w²(
1-
=dd Therefore, 0
c3
22
[())w((1) (t)
di S {
72d C2
f(t)dt.
(C.19)
k
For the boundary region, by Proposition C.5 and Proposition
C.3 we have
1
[f) dt Caki S (f((t)dtt -
4C1d Czk
S" [f(t))dt.
(C.20)
(t)
Also, by Proposition C.5 and Proposition C.11 we have
(w'(t)) w=
ld+1)d
C3 I 2d
Therefore,
cg.
1
[()) (1)
2d
dt
f(tddt.
(C.21)
1: 1-1.-1+ 101- 2.11
(t(t)
Cgk
Finally, noticing that the middle region
and the boundary region
(d-1y= (d id-
1 + d2 U
1 are overlapping, we can further upper bound the right hand side of (C.20), and (.2T). so that
d2
by adding up (C.18).
CC.99.
1 (())dt
f(4) 1e (t:)
$ Co
n²
4C
{ f(t)dt. k -1
d
+
(C.22)
C2
C
C2
C
1=1
Then, by taking O()), we prove the statement.
C.I DIRECT PROOF OF APPROXIMATE MATRIX- VECTOR MULTIPLICATION BOUND
Approximate matrix-vector multiplication bound for Theorem 7.2 Let U be orthogonal basis of the column span of the polynomial regression quasi-matrix A, as defined at the beginning of Ap- pendix C As in the finite dimensional setting, our goal is to prove that, samples via pivotal sampling, using those samples to construct a sampling matrix S, then with prob- ability 1-6. 6,
take k O {
| - Uy") < </bb Uy
(C.23)
We start with the key step
B.4 in the matrix case as
| Uy
Pr ST Uy
1 2 / Uy
(C.24)
[b 1y
Our goal to prove ST - U = O 1.. We have that U's i-th column is the normalized Legendre polynomial L. Let p' = arg min polynomial p ||p bl2
Then Uy* = p", where y* a length + vector holding the coefficients of p when expanded in the basis of Legendre polynomial I,ş So the expectation be written
L,(4)(1 (1:) b(1;)) 27 kw((t)
ST Uy F
(C.25)
21

---

Note that each t, is chosen independently from 4 with probability density function in I as T(t)/ Jr, (t))dt = kt)(t). So we can calculate the expectation (C.25) to get
d
&
(1)(4)) (4) b(()) 2
-
j=0 kw((ts)
(L) b(t)))
kw((t dt
(C.26)
(()((=(t) b(t))²
dt
kw(()
d + k
[] (p'(t) b(t))
where we use (1.11 in the third equality. The integral in the last line is exactly ||b Uy" 1|, so we prove that E|] U S(b Uy
d+1 | Uy 112 Plugging this back into CC244.
we get
+1
Pr [| Uy >1 2 c Uy
1
(C.27)
ck
This shows that if we take k 0 ) samples via pivotal sampling, then with probability 1 6, ||UTS JUTSTSSO U < clbb Uy* 113.
D COMPLEMENTARY EXPERIMENTS
In Section 4 we conduct experiments using three targets; damped harmonic oscillator, a heat equation, and a chemical surface reaction. Here, we show the results of additional simulations. Thus far, the original domain of the experiments is 2D for visualization purposes. Here, consider the 3D original domain by freeing one more parameter in the damped harmonic oscillator model. We also summarize the number of samples required to achieve given target eгror for all four test problems. To corroborate the discussion in Section I.1 we provide a simulation result showing that leverage score sampling outperforms uniform sampling. Lastly, the performance of the randomized BSS algorithm from Chen and Price (2019) on the 2D damped harmonic oscillator target is displayed which supports our discussion in Section 1.3 Before showing the simulation results, we begin this section with deferred detail of the chemical surface coverage target in Section 4
Chemical Surface Coverage. As explained in (Hampton and Doostan 2015). the target function models the surface coverage of certain chemical species and considers the uncertainty which is pa- rameterized by absorption a, desorption y, the reaction rate constant K, and time Given position (Ir, y) in 2D space, our quantity of interest is modeled by the non-linee non-linear evolution equation:
dp
a(1 - 7p - x(1 - p)) p(t 0) 0.9,
dt
(D.I)
a 0.1 exp((.05xz) 7 0.001 + 0.01 exp(0.05y)
In this experiment, we set k 10 and focus after = 4 seconds.
3D Damped Harmonic Oscillator. The damped harmonic oscillator model is given as (restated),
dx. dt (t) + dt
dx (t) (t) + kr(() = f cos(wtt, x=0)= = Tu.
dr
dt (0)
(D.2)
In Section we define the domain as l x [1,3] x [0. 2) while f was fixed as f =05. This time, we extend to 3D space by setting the domain as A x f w = [1.3] x [0.2] [0, 2]. 3D plot of the target is given in Figure (a). Note that if we slice the cube at f - 0.5, we obtain the target in Section 4 This time, create the base data matrix A' ma by constructing
513 fine
grid, and the polynomial degree is set to 12. Figure 7 (b) and (c) give examples of the leverage score of this data matrix.
22

---

10
Azuee
2001 SAyrAD
112
Kiting seepltate
<pring CONMONT:
10
26
30
10
28 25
sortig constant
Soriing constant
(a) Target Function.
(h) Leverage Score at f = 0.D.
(c) Leverage Score at f = 1.0.
Figure 7: (a): The target value of the 3D damped harmonic oscillator model. (b) and (c): Its leverage score. We have f E [0, 2), and the grid is sliced at 0.5 in (b) and at 1.0 in (c).
(a) 3D Damped Harmonic Oscilla- (b) Leverage Score V.s. tor.
Uniform
(c) Randomized BSS.
Sampling.
Figure 8: (a): Results for degree 12 active polynomial regression for the damped harmonic oscil- lator Qol in 3D space. (b). (c): Degree 12 active polynomial regression for the damped harmonic oscillator Qol 2D space. (b) includes the performance of the Bernoulli sampling and our pivotal sampling both using uniform inclusion probability which is given with the thin lines. (c) demon- strates the approximation power of the randomized BSS algorithm. We run the sampling 2000 times with different parameters, round the number of samples to the tens place, and take the median error.
Figure (a) shows the relative error of Bernoulli sampling and our pivotal sampling both using the leverage score. expected, our method shows a better fit than Bernoulli sampling again.
Samples Needed for Certain Error. Table I and summarize the number of samples required to achieve x OPT error and 1.1 OPT error where OPT the error we could obtain when all the data are labeled. In all four test problems, pivotal sampling achieves the target error with fewer samples than the existing method (independent leverlage scores sampling), which denoted by Bernoulli in the tables. Our method is especially efficient when we aim at the target error close to OPT. For instance, to achieve the 1.1 x OPT error in the 2D damped harmonic oscillator model, our method requires less than half samples that Bernoulli sampling requires, showing a significant reduction in terms of the number of samples needed.
Leverage Score Sampling Uniform Sampling. We also conduct complementary experiment to show empirically that leverage score sampling is much more powerful than uniform sampling. We extend the simulation in Fig a uniform inclusion probability setting. This time, draw 350 samples, repeat the simulation 100 times, and report the approximation with a median error in Figure [ The result shows poor performance of the uniform sampling. As they draw fewer samples near the boundaries compared to leverage score sampling, they are not able to pin down the
Oscillator 2D
Heat Eq.
Surface Reaction
Oscillator 3D
10000
10000
10000
51d
poly. deg.
20
20
20
10
Bernoulli (a)
574
554
545
671
Pivotal (b)
398
395
390
533
Efficiency (b/ a)
0.693
0.713
0.716
0.794
Table Number of samples needed to achieve 2 x OPT error.
23

---

Oscillator 2D
Heat Eq.
Surface Reaction
Oscillator 3D
10000
10000
10000
513t
poly. deg.
12
12
Bernoulli (a)
924
814
903
3121
Pivotal (b)
450
442
492
1943
Efficiency (b/ a)
0.487
0.523
0.545
0.623
Table 2: Number of samples needed achieve 1.1 x OPT error.
polynomial function near the edges, resulting in suffering large errors in these areas. The relative error plot is given in Figure 8 (b). We point to three observations. 1) Comparing the thick lines and the thin lines, one can tell that the use of leverage score significantly improves performance. 2) that spatially-aware pivotal sampling outperforms the Bernoulli sampling. 3) By combining the leverage score and spatially aware pivotal sampling (our method), one can attain the best approximation among the four sampling strategies.
Comparing the orange lines and the blue lines,
lpring COEETT
aprng Streart
10
sprig connart
sprig commar
(a) Bernoulli Uniform.
(b) Bernouli Leverage.
(c) Pivotal Uniform.
(d) Pivotal Leverage.
polynomial approximation to the maximum displacement of damped harmonic oscillator, as function of driving frequency and spring constant. X' is uniform dis- tribution over a box. All four draw 350 samples but by different sampling methods; (a) and (b) Bernoulli sampling but (c) and (d) pivotal sampling. Also, (a) and (c) select samples with uniform probability while (b) and (d) employs the leverage score. Clearly, the leverage score suc- cessfully pins down the polynomial function near the boundaries, resulting in better approximation.
Figure 9: Visualizations
Randomized BSS Algorithm. Finally, we run the randomized BSS algorithm (Chen and Price 2019) the 2D damped harmonic oscillator target function, and report the relative error together with Bernoulli leverage score sampling and our pivotal sampling in Figure 8 (c). The setting is the same as the Section 4 that the initial data points are drawn uniformly at random from [1, 3] x [o, 2). and the polynomial degree is set to 12. Even though this algorithm also has a theoretical sample complexity of (dd/c) for the regression guarantee (13) our sampling method achieves a certain relative error with fewer samples.
24


Reviews:
[{"review_id": "jZnYWpyM0G", "reviewer": "Reviewer_Q8k3", "summary": "The authors suggest a modification of the pivotal sampling algorithm which allows for better domain coverage in such an applications as approximating PDE solutions.", "questions": "I would like the authors to clarify the novelty of the suggested theoretical analysis. In particular, I do not understand why the auxiliary results for the proof of Theorem 1.1 do not follows directly from the literature? In particular, why the result of Lemma B.1 does not follow from the existing matrix Bernstein inequalities, see e.g. [Vershynin, Section 5.4]. The results on the approximate matrix multiplication (see the same appendix) also seems to be known [Tropp, 2012, Section 6.4]. It would be great if the authors could better elaborate on the novely, and technical novelty, of Theorems~1.1 and 1.2.\n\nReferences:\n[Vershynin, 2018] Vershynin, Roman. High-dimensional probability: An introduction with applications in data science. Vol. 47. Cambridge university press, 2018. https://www.math.uci.edu/~rvershyn/papers/HDP-book/HDP-book.pdf\n[Tropp, 2015] Tropp, Joel A. \"An introduction to matrix concentration inequalities.\" Foundations and Trends® in Machine Learning 8.1-2 (2015): 1-230. https://arxiv.org/pdf/1501.01571.pdf", "limitations": "", "rating": 8, "confidence": 2, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "The paper is well-written with clear motivation for the main theoretical results and nice experimental illustrations. The modified pivotal sampling algorithm (Algorithm 2) is appealing and seems to be quite easy to implement.", "weaknesses": "My main concerns are stated below and are related with the technical novelty of the results of the analysis of Theorems~1.1 and 1.2. The experimental section is sufficient for the demonstration of the approach, yet it could be elaborated to the higher-dimensional problems.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "WG3KqJ9rWe", "reviewer": "Reviewer_hDW2", "summary": "In this work, the authors improve on active learning methods for linear and polynomial regressions under an agnostic noise setting. Here, the goal is to learn $x \\in \\mathbb{R}^d$ such that $Ax ≈ b$ while observing only a few entries of $b$. Prior works on active linear regression, studied for decades, assume $b$ to be equal to $Ax^*$ plus i.i.d. random noise. In this case, the problem can be addressed using tools from optimal experimental design. In the agnostic case, where such a noise model is not specified, near-optimal sample complexity results and use statistical leverage score-based sampling. In this work, the authors argue that \"any\" sampling strategy that has marginals proportional to leverage scores, satisfying a weak independence condition, will only need samples equal to those required using independent leverage score-based sampling up to constants (theoretical), given by $O(d\\log d)$. They then propose a pivotal sampling-based approach (satisfying the criteria listed) which empirically performs much better, and also show theoretical improvements for the special case of polynomial regression, and uses $O(d)$ samples.", "questions": "Pg 3. Theoretic -> theoretical \n\nPg 6. Bernouilli -> Bernoulli\n\nAs constants play a big role in the motivation of this work and their empirical improvements, a suggestion would be to note them down explicitly in the current work if possible.\n\nThe $\\ell_\\infty$ based independence condition. Can the main theorem result be stated in terms of the parameter $D_{inf}$? It would also be interesting to note the value in the case of independent uniform sampling. It suggests that the $D_{inf}$ is the same for any independent sampling approach.\n\nIn the experiments, it would be interesting to compare with [Chen and Price 2019] and also add a random permutation of rows as the leaves of the tree, in addition to PCA-based ones.", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "— The paper is extremely well-written, providing a thorough discussion on the advantages and limitations of the authors’ work alongside prior research.\n\n— The result on \"any\" sampling strategy is very interesting, even if it doesn’t improve upon previous theoretical guarantees for the given scenario.\n\n— Most sampling-based results obscure significant constants, as noted by the authors regarding the theoretically optimal result due to [Chen and Price 2019], which, in practice, performs poorly. Therefore, the authors' provision of approaches with improved empirical performance, even with similar theoretical guarantees is an important contribution.", "weaknesses": "— The empirical evaluation lacks clarity regarding the improvements. The authors should explicitly state the claimed 50% reduction in samples in the technical sections.\n\n— The theoretical improvements are solely for polynomial regression, which is essentially linear regression with an infinite number of rows. It would be interesting to explore whether the pivotal sampling method has other implications for different norms. For instance, can this \"tournament\" style technique be applied to all norms and their respective optimal sampling strategies?\n\n— The selection of the tree in Algorithm 1 significantly impacts the performance of the sampling method. If the running time is a major concern in large-scale systems or datasets, it isn’t clear if a PCA-based approach is always feasible.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "nXpoyQFJRq", "reviewer": "Reviewer_BstN", "summary": "This paper investigates the active regression problem in the adversarial setting. The authors introduce a novel sampling methodology that promotes spatial coverage, combining it with leverage score sampling. The paper offers two theoretical contributions: the first establishes a sample complexity bound under the one-sided $\\ell_{\\infty}$ independence condition, while the second presents a bound applicable to a specific case of polynomial regression. Empirical results further validate the effectiveness of the proposed approach.", "questions": "a. Is there any theoretical analysis or empirical results regarding the computational complexity of pivotal sampling, particularly concerning the construction of the binary tree?\n\nb. How does the proposed pivotal sampling method compare with maximum leverage sampling strategy?", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "a. The proposed pivotal sampling method using a binary tree tournament is interesting and innovative.\n\nb. The theoretical results are solid. The sample complexity bound for polynomial regression showcases an improvement of a logarithmic factor over the general case. The techniques employed in the analysis have potential implications beyond this specific research.\n\nc. The authors provide empirical results to validate the effectiveness of the proposed pivotal sampling. These results clearly indicate that pivotal sampling significantly outperforms Bernoulli leverage sampling.\n\nd. The presentation is clear and easy to follow.", "weaknesses": "a. My primary concern pertains to the computational complexity of the proposed method. The paper lacks both theoretical analysis and empirical results in this regard, leaving an important aspect unaddressed.\n\nb. The experimental setup in the paper includes only a single baseline, and it would be beneficial to include additional baseline methods, such as maximum leverage score sampling.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "AKZ6O4Dcde", "reviewer": "Reviewer_gr3D", "summary": "This paper proposes a way to perform leverage score sampling that promotes spatial coverage for the problem of active linear regression. Motivated by the empirical successes of deterministic grid-based approaches to PDEs, the presented approach retains the strong theoretical guarantees of leverage score sampling with the spatially-well distributed samples of grid-based approaches. The key insight is to use a pivotal sampling algorithm where a binary tree that matches the geometry of the data is constructed deterministically and the samples are percolated up via head-to-head comparisons (probability weighted coin flip) at sibling nodes. This leads to the desirable behavior of spatial coverage since close neighbors in the tree are less likely to be both included. The authors conduct empirical evaluations that show the improved effectiveness of the proposed method.", "questions": "The authors note that the algorithm of Chen and Price, 2019 is provably optimal, but that “in our initial experiments, it was not competitive with leverage score sampling in practice.” Why were these results not included in the main body of the paper? My understanding is that the Randomized BSS algorithm is only included as a result in Figure 8c of the appendix. It would be compelling to have comparisons to Chen and Price’s method in Fig. 4.", "limitations": "", "rating": 8, "confidence": 4, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "* The problem of active linear regression is interesting and has applications to uncertainty quantification and parametric partial differential equations.\n* The paper is very well-written and organized. There is a great overview of prior work as well as the contextualization of the results (e.g., the fact that the work builds on top of the recent result in Chernoff bounds under $\\ell_\\infty$ independence).\n* The claims of the paper are adequately supported with rigorous theoretical analysis and empirical evaluations.\n* The algorithm and its analysis is novel to the best of my knowledge. The presented analysis is general enough that it could be applied to any non-independent leverage score sampling method that obeys a weak one-sided  $\\ell_\\infty$ independence condition. This might be of independent interest to other researchers in the area. Leveraging \n* The empirical evaluations are compelling and clearly show the improved effectiveness of the dependent leverage score sampling method.", "weaknesses": "* The presented analysis does not provide justification for the improved empirical effectiveness of the method (besides for $\\ell_2$ polynomial regression). So, it is not clear theoretically why the method performs much better than independent leverage score sampling. Nevertheless, this is a limitation that the authors clearly concede and mention a few times throughout.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "HRHSUtM7QW", "reviewer": "Reviewer_gm2B", "summary": "This paper proposes a novel leverage score sampling method based on spatial pivoting and therefore the samples are dependent. The main results involve a nearly optimal sample complexity of $O(d\\log d+d/\\epsilon)$ for active regression which matches the standard independent leverage score sampling, and a optimal bound of $O(d/\\epsilon)$ for polynomial regression. The key proof ingredient for the first result is an adaptation of the $\\ell_\\infty$ independence result due to [KKS22] and how to generalize this property to approximate matrix product, and for the second result is similar to [KKP17]. They also perform extensive experiments for active regressions on PDEs and the experiments are convincing enough that their algorithm has good performance.", "questions": "A few comments regarding citation: the paper [KKS22] was published in SODA 2022. It's better to cite the proceeding version.\n\nA question: it seems the binary tree-based approach to pre-partition the space based on the structure of the problem can also be extended to BSS-type sampling, e.g., see the data structures and algorithms in [SXZ22]. Do you think your argument and techniques can be generalized to BSS, to achieve an optimal sample complexity while leveraging particular structures?\n\n[SXZ22]: Z. Song, Z. Xu and L. Zhang. Speeding up sparsification using inner product search data structures, 2022.", "limitations": "", "rating": 8, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "While the techniques are not completely novel, I do like the theoretical implications of this paper --- they essentially show that leverage score sampling with limited independence also provides theoretically matching bound compared to independent leverage score sampling. This inspires one to design leverage score sampling algorithm that can adaptively choose the samples based on particular structure of the problem. They also show a BSS-type result for polynomial regression without resorting to the BSS barrier function argument. This is very interesting because as far as I know, all sparsification results with $O(d/\\epsilon)$ sample complexity are more or less variants of BSS. I believe results and techniques in this paper have further applications in other sparsification problems.\n\nIn addition, they perform many numerical experiments to show that their algorithm works well in practice. This is nice as it aligns with the idea that one can modify leverage score sampling based on problems to get better (practical) algorithms.", "weaknesses": "As noted in the strengths part, the techniques in this paper are not very novel. To get their first $d\\log d+d/\\epsilon$ sample complexity for active regression, they mainly utilize [KKS22]. The approximate matrix product result requires a bit more work, but it's also standard and not surprising. For the second result regarding polynomial regression, it also mainly follows procedures from [KKP17]. Nevertheless, I still think these results are nice, and are interesting enough to be published in conferences like ICLR.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (oral)