Classification vs regression in overparameterized regimes: Does the loss function matter?

Introduction

Paradigmatic problems in supervised machine learning (ML) involve predicting an output response from an input, based on patterns extracted from a (training) data set. In classification, the output response is (finitely) discrete and we need to classify input data into one of these discrete categories. In regression, the output is continuous, typically a real number or a vector. Owing to this important distinction in output response, the two tasks are typically treated differently. The differences in treatment manifest in two phases of modern ML: optimization (training), which consists of an algorithmic procedure to extract a predictor from the training data, typically by minimizing the training loss (also called empirical risk); and generalization (testing), which consists of an evaluation of the obtained predictor on a separate test, or validation, data set.

Traditionally, the choice of loss functions for both phases is starkly different across classification and regression tasks. The square loss function is typically used both for the training and the testing phases in regression. In contrast, the hinge or logistic (cross-entropy for multi-class problems) loss functions are typically used in the training phase of classification, even though the 0-1 loss function is used for testing. The use of the logistic and hinge losses can be motivated by their computational and statistical properties. For example, the theory of surrogate losses (Zhang, 2004; Bartlett et al., 2006; Ben-David et al., 2012) gives theoretical arguments favoring the logistic and hinge losses over other convex surrogates, including the square loss.Also see, e.g., Section 8.1.2 in Goodfellow et al. (2016) for a representative informal discussion. Yet, there have been indications that the reality is more complex, both in underparameterized and overparameterized regimes of ML. For example, Rifkin (2002) extensively compared the hard-margin support vector machine (SVM), which minimizes the hinge loss, and regularized least-squares classification (RLSC), which minimizes the square loss — ultimately concluding that “the performance of the RLSC is essentially equivalent to that of the SVM across a wide range of problems, and the choice between the two should be based on computational tractability considerations.” Quite similar resultsIn fact, while the results were generally close, in a majority of classification tasks models trained using the square loss outperformed models trained with cross-entropy. for a comparison between the square loss and cross-entropy loss have recently been obtained in Hui and Belkin (2021) for a range of modern neural architectures and data sets across several application domains. The latter is the current dominant standard for training neural networks.

In an important separate development, we have recently seen compelling evidence that overparameterized deep neural networks, as well as other models trained to interpolate the data (i.e. achieve zero, or near zero, training loss), are capable of good test performance, questioning conventionalFor example, Hastie, Tibshirani and Friedman say, on page $221$ of their popular statistical learning textbook (Hastie et al., 2009), that “a model with zero training error is overfit to the training data and will typically generalize poorly.” statistical wisdom (Neyshabur et al., 2014; Zhang et al., 2016; Belkin et al., 2018; Geiger et al., 2019; Belkin et al., 2019). Since then, the theoretical ML community has identified regimes for regression tasks under which overfitting is benign (to borrow the language from Bartlett et al., 2020), and interpolating noisy data with solutions arising from empirical risk minimization is compatible with good generalization (Bartlett et al., 2020; Belkin et al., 2020; Hastie et al., 2019; Mei and Montanari, 2019; Muthukumar et al., 2020). It is worth noting that the ensuing test loss of interpolating solutions is unrelated to their training loss, which is identically zero for all such solutions. This demonstrates, again, that the relationship between training and test losses — both for regression and classification tasks — is more complex than often assumed.

This paper introduces a direct comparison between the different loss functions used in classification and regression, in both the training and testing phases. We analyze the modern overparameterized regime under the linear model with Gaussian features and uncover a remarkable phenomenon of overparameterized training: in sufficiently overparameterized settings, with high probability, every training data point is a support vector. Consequently, the outcome of optimization (with gradient descent) is the same whether we use the hinge loss, logistic loss, or the square loss.

On the other hand, we show that the choice of test loss function results in a significant asymptotic difference between classification and regression tasks. In particular, we identify truly overparameterized regimes for which predictors will generalize poorly for regression tasks (measured by the square loss), but well for classification tasks (measured by the 0-1 loss). The fact that regression and classification can be different has been understood for some time. For example, Devroye, Gyorgi and Lugosi point out in Chapter 6.7 of their classic textbook (Devroye et al., 1991) that “classification is easier than regression function estimation”, in reference to sample complexity. Artificial examples are also given where accurate regression estimates of the density are impossible given the hypothesis class, but classification succeeds because of the separability of the two classes. In contrast, our paper demonstrates that in overparameterized regimes, this phenomenon can be quite generic. Approximability is not the underlying issue; rather, it is a consequence of learning.

In contrast, we show that the choice of loss function used on test data yields significant differences between classification and regression. Depending on the extent of “effective overparameterization”, the same minimum-norm solution can:

succeed at both regression and classification,

succeed at classification and fail at regression, or

Related work

The phenomenon of overparameterization and interpolation yielding significantly improved empirical performance across a variety of models as well as tasks (Neyshabur et al., 2014; Zhang et al., 2016; Geiger et al., 2019; Belkin et al., 2019) has received significant research attention over the last few years. In this section, we contextualize our results in this research landscape.

At a high level, any solution obtained in an overparameterized regime that generalizes well must have some sort of regularization, i.e. special structural constraints on the values it can take. Thus, we need to understand the influence of the choice of training loss function on the resulting solution and its generalization guarantee. In the overparameterized regime, there are infinitely many solutions that interpolate training data, and indeed even more that separate discretely labeled data. Thus, characterizing the implicit regularization (Ji and Telgarsky, 2019; Soudry et al., 2018; Gunasekar et al., 2018; Woodworth et al., 2019; Nacson et al., 2019; Azizan et al., 2020) induced by the choice of optimization algorithm is important to understand properties of the obtained solutions. For the linear model, we have a concrete understanding of the solutions obtained by the most common choices of training loss functions:

If we minimize the logistic loss using gradient descent on separable training dataThe implicit bias has also been characterized for the more difficult non-separable case (Ji and Telgarsky, 2019), but we focus here on separable training data as this will always be the case for an overparameterized setting., we will converge to the hard-margin SVM (Ji and Telgarsky, 2019; Soudry et al., 2018).

As mentioned in the introduction, conventional wisdom recommends the choice of the logistic loss, or the hinge loss, for classification tasks. It is sometimes implied (without theoretical justification) that instead minimizing the square loss would be suboptimal for generalization. However, our first main result (Theorem 11) shows that with sufficient overparameterization, the SVM itself interpolates the binary labels — as pictured in Figure 2, this implies an equivalence in solutions corresponding to several choices of training loss function. Moreover, our subsequent Theorem 13 shows that the interpolating solution generalizes well in classification tasks, for a wide range of overparameterized regimes. And when it does generalize poorly, so does the SVM! These results add theoretical weight to the empirical evidence (Rifkin, 2002; Que and Belkin, 2016) that the hinge loss (and, by extension, the cross-entropy loss) is not necessarily the superior choice for classification tasks.

The SVM maximizes training data margin in feature space. Theoretical analyses of generalization error as a function of the margin have been proposed to explain the success of models such as boosting and neural networks (Schapire et al., 1998; Bartlett, 1998; Bartlett et al., 2017). Explanations based on margin bounds are sometimes credited, in a heuristic manner, for the empirical success of interpolated models in classification tasks. This is, in fact, a misleading explanation (as also noted in Shah et al., 2018); in Section 6, we provide experimental evidence for the tautology of generalization upper bounds as a function of the feature-space margin, when applied to sufficiently overparameterized models. This evidence corroborates the recent perspectives on modern ML which argue against generalization bounds that tie training loss to the expected loss on test data (Belkin et al., 2018; Nagarajan and Kolter, 2019). We instead favor a first-principles approach to analyzing high-dimensional models for classification, inspired by recent progress in regression.

2 Insights from least-squares regression

In this paper, we build on these insights from regression tasks, to show that overparameterized models can similarly generalize well for classification tasks. In fact, it turns out that the conditions for classification are milder, and there is an intermediate regime of overparameterization where the regression problem is “hard” — but the classification problem is “easy”. Focusing on the well-specified case, we show that the balance between preserving signal and absorbing noise does not have to be as delicate for classification tasks as it is for regression tasks. This conclusion cannot be made directly from the regression analyses, as 0-1 classification error is quite different from the mean-square regression error. To bridge the gaps, we use a signal processing perspective on the overparameterized regime that was first developed in Muthukumar et al. (2020), where the conditions for low test error are linked to notions of survival of the true features and contamination by falsely discovered features. We will see (in Section 5.2) that these same quantities show up explicitly in the analysis of classification test error.

3 Recent work on high-dimensional classification/logistic regression

Concurrent to our work, Chatterji and Long (2020) directly upper bound the generalization error of the max-margin SVM under overparameterized linear discriminant models that are not isotropic and, like us, explore how anisotropic the situation needs to be for good generalization. Both their results and techniques are quite different from ours in that they study the iterates of gradient descent leveraging the implicit regularization perspective of optimization algorithms.

Setup

We begin with some basic notation. Thereafter, we describe the setup for training and test data, evaluation of classification and regression tasks, and choices of featurization (in that order).

We describe basic notation for vectors, matrices, and functions.

Let $\mathbf{e}_{i}$ represent the $i^{th}$ standard basis vector (with the dimension implicit). For a given vector $\mathbf{v}$, the functional $\mathsf{sgn}(\mathbf{v})$ denotes the sign operator applied element-wise. Let $\mu_{i}(\mathbf{M})$ denote the $i^{th}$ largest eigenvalue of positive semidefinite matrix $\mathbf{M}$, and $\mu_{\mathsf{max}}(\mathbf{M})$ and $\mu_{\mathsf{min}}(\mathbf{M})$ denote in particular the maximal and minimal eigenvalue respectively. Further, we use $||\mathbf{M}||_{\mathsf{op}},\mathsf{tr}(\mathbf{M})$ and $||\mathbf{M}||_{\mathsf{F}}$ to denote the operator norm, trace norm, and Frobenius norm respectively.

1.2 Function-specific notation

For two functions $f(n)$ and $g(n)$, we write $f\asymp g$ iff there exist universal positive constants $(c,C,n_{0})$ such that

(In most places where we apply the above inequality, the functions $f$ and $g$ are positive valued and so we automatically drop the absolute value signs.)

2 Data

Here, the feature map $\boldsymbol{\phi}$ is known, but the target parameter $\boldsymbol{\alpha^{*}}$ (which we refer to as the signal) is unknown. The label noise in $Y_{i}$ is assumed to be independent of everything else.

We define shorthand notation for the training data: let

denote the data (feature) matrix; $\mathbf{Z}_{\mathsf{train}}:=\begin{bmatrix}Z_{1}&\ldots&Z_{n}\end{bmatrix}^{\top}\in n$ denote the regression output vector; and $\mathbf{Y}_{\mathsf{train}}:=\begin{bmatrix}Y_{1}&\ldots&Y_{n}\end{bmatrix}^{\top}$ denote the classification output vector. Note that if there is no label noise (i.e. $\nu^{*}=0$), then we have $\mathbf{Y}_{\mathsf{train}}=\mathsf{sgn}(\mathbf{Z}_{\mathsf{train}})$.

3 Classification, regression, and interpolation

The overparameterized regime constitutes the case in which the dimension (or number) of features is greater than the number of samples, i.e. $d\geq n$. We define the two types of solutions that we will primarily consider in this regime, starting with interpolating solutions.

We consider solutions $\boldsymbol{\alpha}$ that satisfy one of the following feasibility conditions for interpolation:

Recall from our discussion in Section 2.1 that these interpolations arise from minimizing the square loss on training data. If we instead minimized the logistic or hinge loss, we would obtain the hard-margin support vector machine (SVM), defined below.

Note that data is defined to be linearly separable iff the constraints in Equation (4) can be feasibly satisfied by some parameter vector $\boldsymbol{\alpha}$.

As long as $d\geq n$, any solution that interpolates the binary labels $\{Y_{i}\}_{i=1}^{n}$ satisfies Equation (4) with equality almost surely for any continuous distribution on the features. Thus, in the overparameterized regime, the training data is trivially linearly separable. Note, however, that the feasibility constraints do not require the SVM solution to interpolate the binary labels.

The standard metrics for test error in regression and classification tasks are, respectively, the mean-square-error (MSE) and classification error, defined as follows. In these definitions, we have ignored the irreducible error terms arising from possible additive noise in real outputs and label noise in binary outputs respectively. This reflects the practical goal of all prediction to get the underlying true output right, as opposed to matching noisy measurements of that underlying true output.

Here, all expectations (and ensuing probabilities) are only over the random sample $X$ of test data. As is standard, we will characterize the regression and classification test errors with high probability over the randomness in the training data $\{X_{i},Y_{i}\}_{i=1}^{n}$.

As a final comment, we will typically construct an empirical estimate of both test error metrics from $n_{test}$ test samples of data drawn without any label noise. This is for ease of empirical evaluation.

4 Featurization

We consider zero-mean Gaussian featurization, i.e. for every $i\in\{1,\ldots,n\}$, we have

We denote the spectrum of the (positive definite) covariance matrix $\boldsymbol{\Sigma}$ by the vector $\boldsymbol{\lambda}:=\begin{bmatrix}\lambda_{1}&\ldots&\lambda_{d}\end{bmatrix}$, where the eigenvalues are sorted in descending order, i.e. we have $\lambda_{1}\geq\lambda_{2}\geq\ldots\geq\lambda_{d}>0$.

Throughout, we will consider various overparameterized ensembles obtained by scaling the covariance parameter $\boldsymbol{\Sigma}$ as a function of both the number of training data points, $n$, and the number of features, $d$. We theoretically characterize the performance of solutions for classification and regression tasks using two representative ensembles, defined below.

The isotropic ensemble, parameterized by $(n,d)$, considers isotropic Gaussian features, $\boldsymbol{\Sigma}=\mathbf{I}_{d}$. For this ensemble, we will fix $n$ and study the evolution of various quantities as a function of $d\geq n$.

Note that the isotropic ensemble constitutes the “maximal” level of effective overparameterization (as defined in the second effective rank in Bartlett et al. (2020)) for a given choice of $(n,d)$.

The bi-level ensemble is parameterized by $(n,p,q,r)$, whereWe restrict $(p,q,r)$ to this range to ensure that a) the regime is truly overparameterized (choice of $p$), b) the eigenvalues of the ensuing covariance matrix are always positive and ordered correctly (choice of $q$), c) the number of “high-energy” directions is sub-linear in $n$ (choice of $r$). $p>1,0\leq r<1$ and $0<q<(p-r)$. Here, parameter $p$ controls the extent of artificial overparameterization), $r$ sets the number of preferred features, and $q$ controls the weights on preferred features and thus effective overparameterization. In particular, this ensemble sets parameters

The covariance matrix of the Gaussian features $\boldsymbol{\Sigma}(p,q,r)$ is set to be a diagonal matrix, whose entries are given by:

For this ensemble, we will fix $(p,q,r)$ and study the evolution of various quantities as a function of $n$.

The bi-level covariance matrix is parameterized by the choice for the top $s$ eigenvalues and the bottom $(d-s)$ eigenvalues, with the sum of eigenvalues being invariant(equal to $d$). The parameters of critical importance are $p$, which determines the extent of overparameterization (i.e. number of features), $r$, which determines the number of larger eigenvalues, and $q$, which determines the relative values of larger and smaller eigenvalues (all as a function of the number of training points $n$). We make a few remarks below on this ensemble.

This bi-level ensemble is inspired by the study of estimation of high-dimensional spiked covariance matrices (e.g. Wang and Fan, 2017; Mahdaviyeh and Naulet, 2019) when the number of samples is much smaller than the dimension. In these spiked matrices, the parameter $s$ is typically set to a constant (that does not grow with $n$), and the top $s$ eigenvalues are highly spiked with respect to the other $(d-s)$ eigenvalues. In fact, it is assumed that there exists a universal positive constant $C$, such that the smaller eigenvalues are bounded and the top (larger) eigenvalues grow with $(d,n)$ in the following way:

Under these conditions, the ratio of the top to the bottom eigenvalues grows as $\Omega\left(\frac{d}{n}\right)$, and Wang and Fan (Wang and Fan, 2017) show that the top $s$ estimated eigenvalues of the high-dimensional covariance matrix can be estimated reliably from samples, even when less than the dimension (i.e. $n<d$). This condition, which is also critical for good generalizationIn particular, avoiding signal shrinkage, as also shown in Bartlett et al. (2020). in regression problems, can be verified to be equivalent to the condition $q\leq(1-r)$ in our bi-level ensemble (see Theorem 13 for a full statement). Our definition of the bi-level ensemble allows further flexibility in the choice of these parameters, and we will later show that classification tasks can generalize well even in the absence of this condition.

The bi-level ensemble can be verified to match the isotropic ensemble (Definition 4) as a special case when the parameters are set as $q+r=p$. This case represents the maximal level of effective overparameterization, and in general we take $q\leq(p-r)$ to ensure correct ordering of the eigenvalues. The smaller the value of $q$, the less the effective overparameterization. The models of Chatterji and Long (2020) are spiritually related in how they also use an exponent like $q$ to control the effective overparameterization.

We know that for “benign overfitting” (Bartlett et al., 2020) of additive noise to occur in regression problems, we need to have sufficiently many (growing super-linearly in $n$) “unimportant” directions, corresponding to the lower level of eigenvalues. The choice of parameters $p>1$ and $r<1$ ensures that the number of such “unimportant” directions is equal to $(d-s)=(n^{p}-n^{r})\gg n$, and so the bi-level ensemble as defined does not admit the regime of harmful overfitting of noise for any choice of parameters $(p,q,r)$. This allows us to isolate signal shrinkage as the principal reason for large regression error, and also study the ramifications of such shrinkage for classification error.

In addition to the above, we empirically study (in Section 6) the behavior of various quantities for two other ensembles defined below, both of which have been previously studied in regression tasks.

This ensemble is a simplification of feature families introduced by Belkin et al. (2020); Hastie et al. (2019); Mei and Montanari (2019), all of which demonstrate an explicit benefit of overparameterization in generalization for regression tasks. The features consist of raw, $1$-dimensional random variable $X_{i}\sim\mathcal{N}(0,\sigma^{2})$, and independent $d$-dimensional random variables $\mathbf{W}_{i}\text{ i.i.d.}\sim\mathcal{N}(\mathbf{0},\mathbf{I}_{d})$ for $i\in\{1,\ldots,n\}$. Then, for a given value of $d>n$, each lifted feature is given by:

We define the real output $Z_{i}$, corresponding to input $X_{i}$, as $Z_{i}=X_{i}$, and similarly the binary output $Y_{i}$ is defined as $Y_{i}=\mathsf{sgn}(X_{i})$. For this ensemble, we will fix $n$ and study the evolution of various quantities as a function of $d$. This model allows us to study model mis-specification because the real output is not exactly representable in the feature space.

This ensemble is inspired by commonly chosen reproducing kernel Hilbert spaces and is parameterized by $m\geq 0$. We set the spectrum of the covariance matrix $\boldsymbol{\Sigma}$ to be

For this ensemble, we will fix $(n,d)$ and study the evolution of various quantities as a function of the parameter $m$.

Approximating the SVM by exact interpolation

From the optimization objective and constraints defined in Equation (4), we can see that there is a continuum of margins, defined by $Y_{i}\cdot\phi(X_{i})^{\top}\alpha$, that is possible for each training point. Thus, unlike in least-squares regression, even obtaining an exact expression for the margin-maximizing SVM solution, $\boldsymbol{\widehat{\alpha}}_{\mathsf{SVM}}$, appears difficult in the overparameterized regime.

Such a phenomenon, when true, explicitly links the concept of support-vector-machines with a positive margin constraint to exact interpolation of the training data labels, and suggests a roadmap to analyzing the generalization error of the SVM by analyzing the latter solution (which we do subsequently in Section 5). We show in the theorem below that this phenomenon manifests with high probability provided there is sufficient effective overparameterization.

Let $\boldsymbol{\Phi}_{\mathsf{train}}$ follow the Gaussian featurization from Equation (7) with covariance matrix $\boldsymbol{\Sigma}$, and let $\boldsymbol{\widehat{\alpha}}_{\mathsf{SVM}}$ be the solution to the optimization problem in Equation (4).

the vector $\boldsymbol{\widehat{\alpha}}_{\mathsf{SVM}}$ satisfies the binary label interpolation constraint (Equation (3a)) simultaneously for every $\mathbf{Y}_{\mathsf{train}}\in\{\pm 1\}^{n}$ with probability at least $\left(1-\frac{2}{n}\right)$.

If $\boldsymbol{\Sigma}=\mathbf{I}_{d}$ (i.e., $\boldsymbol{\Phi}_{\mathsf{train}}$ follows the isotropic ensemble), and

then the vector $\boldsymbol{\widehat{\alpha}}_{\mathsf{SVM}}$ satisfies the binary label interpolation constraint (Equation (3a)) for any fixed $\mathbf{Y}_{\mathsf{train}}\in\{\pm 1\}^{n}$ with probability at least $\left(1-\frac{2}{n}\right)$.

Theorem 11 is proved in Appendix C. Both conditions are proved by showing that a complementary slackness condition on the dual of the SVM optimization problem holds with high probability — where the conditions differ is in the application of concentration bounds. The condition in Equation (11) is proved using a broadly applicable “epsilon-net” argument to bound the operator norm of a random matrix, while the sharper condition in Equation (12) leverages Gaussian isotropy and precise properties of the inverse Wishart distribution. Note that the first result holds for all label vectors $\mathbf{Y}_{\mathsf{train}}\in\{\pm 1\}^{n}$ simultaneously, while the second result holds for any fixed $\mathbf{Y}_{\mathsf{train}}$ (but independent of the features).

We now remark on the result for the isotropic and the bi-level ensemble.

Plugging in the condition in Equation (11) into the bi-level ensemble (Definition 5), the following conditions on $(p,q,r)$ are sufficient for all training points to become support vectors with high probability (see Appendix C.2 for a full calculation):

There is an intuitive interpretation for each of these conditions in light of the second “effective rank” condition that is sufficient for benign overfitting (Bartlett et al., 2020) of noise (although our proof technique is quite different). First, the condition $p>2$ mandates an excessively large number of unimportant directions, i.e. corresponding to lower-level (smaller) eigenvalues ($(n^{p}-n^{r})$ of them). Second, the condition $q>\left(\frac{3}{2}-r\right)$ mandates that the ratio between the important directions, i.e. higher-level eigenvalues, and the unimportant directions, is sufficiently small — thus, the unimportant directions are sufficiently weighted. This second condition appears to be strictly stronger than what is required for benign overfitting of noise.

Equation (13) is quite strong as a sufficient condition, but nevertheless admits non-trivial regimes for which classification can generalize well or poorly (see the text accompanying Theorem 13 for a full discussion). However, there is ample evidence to suggest that this condition is not necessary. Notably, Figure 3(b) shows that with a choice of parameterization $(p=3/2,r=1/2)$, the fraction of support vectors becomes equal to $1$ around when $q\geq 0.7$. This choice of parameters for the bi-level ensemble clearly violates both conditions in Equation (13). Thus, all training points become support vectors more often than our theory predicts. Subsequent work to ours (Hsu et al., 2021) tightened the condition in Equation (11) by providing a new deterministic equivalent to the phenomenon of all training points becoming support vectors.

From this section, we have identified an interesting phenomenon by which all training points become support vectors with sufficient effective overparameterization. Moreover, this phenomenon is even more prevalent empirically than our current theory predicts.

Generalization analysis for interpolating solution with Gaussian features

Even with clean data (i.e. zero label noise), the classification setup admits misspecification noise of the form $Y_{i}-\boldsymbol{\phi}(X_{i})^{\top}\boldsymbol{\alpha^{*}}$. The misspecification noise is clearly non-zero mean, and is non-trivially correlated with the features. This resists a clean decomposition of generalization error into the error arising from signal identifiability (or lack thereof) + error arising from overfitting of noise, as in Bartlett et al. (2020).

For a given interpolation $\boldsymbol{\widehat{\alpha}}$, the expression for classification error is distinctly different from mean-square-error (we will see this explicitly in Theorem 17). In particular, we will see that characterizing this expression sharply requires novel analysis of the individual recovered coefficients as a result of interpolation.

We state our main result for this section in the context of the bi-level ensemble (Definition 5). We fix parameters $p>1$ (which represents the extent of artificial overparameterization), and $r\in[0,1)$ (which sets the number of preferred features), and $q\in[0,p-r]$ (which controls the weights on preferred features, thus effective overparameterization); and study the evolution of regression and classification risk as a function of $n$. For the purpose of this section, we denote the regression and classification test losses under the bi-level ensemble as $\mathcal{R}(\boldsymbol{\widehat{\alpha}}_{2,\mathsf{real}};n)$ and $\mathcal{C}(\boldsymbol{\widehat{\alpha}}_{2,\mathsf{binary}};n)$, to emphasize that these losses vary with $n$.

In addition to this and the broad setup as described in Section 3 we make a $1$-sparse assumption on the unknown parameter vector $\mathbf{\boldsymbol{\alpha^{*}}}$, as described below.

Under Assumption 8, we now show the existence of a regime, corresponding to choice of $(p,q,r)$ above, for which the regression test loss stays prohibitively high, but the classification test loss goes to as $n\to\infty$. (We also derive non-asymptotic versions of these results in Appendix E, but only state the asymptotic results here for brevity.)

In this regime, both regression and classification generalize well.

For $(1-r)<q<(1-r)+\frac{(p-1)}{2}$, we have

In this regime, classification generalizes well but regression does not.

For $(1-r)+\frac{(p-1)}{2}<q\leq(p-r)$, we have

In this regime, the generalization is poor for both classification and regression.

Note that the presence of label noise $\nu^{*}$ does not affect these asymptotic scalings (since $\nu^{*}<0.5$).

The new regime of principal interest that we have identified is values of $q\in(1-r,1-r+\frac{p-1}{2})$ for which classification generalizes, but regression does not. The entire proof of Theorem 13 is deferred to Appendices D and E, but we briefly illustrate the intuition for this discrepancy between classification and regression tasks in Section 5.2. In particular, we will see that good generalization for classification requires a far less stringent condition on coefficient recovery than regression.

We now provide some intuition for the scalings described in Theorem 13 for the bi-level ensemble.

The regime that we have identified that is of principal interest is intermediate values of $q$, i.e. $(1-r)<q<(1-r)+\frac{(p-1)}{2}$. This highlights a fascinating role that overparameterization, in the form of the parameter $p$, plays in allowing the good generalization of interpolating solutions in classification tasks. Recall that the larger the value of $p$, the larger the total number of features $d=n^{p}$. Thus, there are several “unimportant directions” in the bi-level ensemble all corresponding to the smaller eigenvalue — which helps in harmless absorption of effective noise. In the proof of Theorem 13, we will identify an explicit mechanism by which having many unimportant directions helps in good generalization for classification, even though the signal is not preserved. At a high level, this mechanism constitutes the spreading out of attenuated signal across several features in a relatively “harmless” way, to exhibit minimal influence on classification performance. In fact, this influence is quantified by a notion of “contamination” by falsely discovered features (defined in Section 5.2) that can be directly linked to the contribution of noise overfitting to regression error.

For $\left(\frac{3}{2}-r\right)<q<(1-r)+\frac{(p-1)}{2}$, we have

For $(1-r)+\frac{(p-1)}{2}<q\leq(p-r)$, we have

Observe that Corollary 16 directly follows from plugging in the condition required in the bi-level ensemble for all training points usually becoming support vectors (Equation (13)), and noting that for $p>2$, we have

Importantly, we have identified that even highly overparameterized regimes, in which all training points become support vectors, can yield good generalization for classification tasks when the hard-margin SVM is used.

2 Path to analysis: Classification vs regression test error

The first step to proving Theorem 13 is obtaining clean expressions for both classification and regression test error. The $1$-sparsity assumption that we have made on the unknown signal enables us to do this as a function of natural quantities corresponding to the preservation of the true feature (survival) and the pollution due to false features (contamination). If we assume that the real labels are generated by the ${t}^{\text{th}}$ feature, $\alpha^{*}_{t}$, then we can define these quantities for any solution $\boldsymbol{\widehat{\alpha}}$. First, as classically observed in statistical signal processing, the estimated coefficient corresponding to the true feature $\alpha^{*}_{t}$ will experience shrinkage and be attenuated by a factor that we denote as survival. From Assumption 8, we defined $\boldsymbol{\boldsymbol{\alpha^{*}}}:=\frac{1}{\sqrt{\lambda_{t}}}\cdot\mathbf{e}_{t}$, and so we have

Second, we have the false discovery of features. We measure the effect of this false discovery for prediction on a test point $X$ by a contamination term:

Recall that $X$ is random, and the features $\boldsymbol{\phi}(X)$ are zero-mean. Therefore, $B$ is a zero-mean random variable. Accordingly, we can define the standard deviation of the contamination term on a test point as below:

where the last step follows from the orthogonality of the $d$ features. The ideas of survival and contamination can be related to the classical signal processing concept of aliasing; Figure 7 in Appendix A provides an illustration.

We state and prove the following proposition, which directly expresses regression and classification test loss in terms of these terms.

Under the $1$-sparse noiseless linear model, the regression test loss (excess MSE) is given by:

and the classification test loss (excess classification error) is given by:

We can think of the quantity $\mathsf{SU}(\boldsymbol{\widehat{\alpha}},t)/\mathsf{CN}(\boldsymbol{\widehat{\alpha}},t)$ as the effective “signal-to-noise ratio” for classification problems.

Proof We first prove Equation (17). Recall that for any estimator $\boldsymbol{\widehat{\alpha}}$, the excess MSE is given by

and then substituting in the $1$-sparse Assumption 8 gives us Equation (17).

Next, we prove Equation (18). Since $\boldsymbol{\phi}(X)=\boldsymbol{\Sigma}^{1/2}\mathbf{W}$ for $\mathbf{W}=(W_{1},\dotsc,W_{d})\sim\mathcal{N}(\mathbf{0},\mathbf{I}_{d})$, we can write $\boldsymbol{\phi}(X)^{\top}\boldsymbol{\alpha^{*}}=W_{t}$ and $\boldsymbol{\phi}(X)^{\top}\boldsymbol{\widehat{\alpha}}=\sum_{j=1}^{d}\sqrt{\lambda_{j}}W_{j}\widehat{\alpha}_{j}$. Thus, the excess classification error of $\boldsymbol{\widehat{\alpha}}$ is given by

Now, the random sum $\sum_{j\neq t}\sqrt{\lambda_{j}}\widehat{\alpha}_{j}W_{j}$ has a Gaussian distribution with mean zero and variance $\mathsf{CN}(\boldsymbol{\widehat{\alpha}},t)^{2}$. Since the $\{W_{j}\}_{j=1}^{d}$ are independent, the classification test error of $\boldsymbol{\widehat{\alpha}}$ is the probability of the following event:

where $U$ and $V$ are independent standard Gaussian random variables. This event is equivalently written as

Since $V/U$ follows the standard Cauchy distribution with cumulative distribution function $F(t)=\tfrac{1}{2}+\tfrac{1}{\pi}\mathsf{tan}^{-1}\left(t\right)$, the claim follows. Equations (17) and (18) give us an initial clue as to why classification test error can be easier to minimize than regression test error. For the right hand side of Equation (17) to be small, we need $\mathsf{SU}\to 1$ to avoid shrinkage, as well as $\mathsf{CN}\to 0$ to avoid contamination. However, for the right hand side of Equation (18) to be small, we only require the ratio of contamination to survival to be small (i.e. $\mathsf{CN}/\mathsf{SU}\to 0$). Clearly, the former condition directly implies the latter, showing that classification is “easier” than regressionOur decomposition of classification error is reminiscent of the decomposition by Friedman (1997) into the ratio of terms depending on the variance (like contamination) and bias (like survival) respectively. Because our data is Gaussian, Proposition 17 allows an exact decomposition.. Theorem 13 is proved fully in Appendices D and E in the following series of steps:

Matching (non-asymptotic) upper and lower bounds are proved on both survival and contamination for interpolation of both real and binary labels. The full statements for these bounds are contained in Theorems 22 and 23 in Appendix D.1.

These bounds are substituted into the bi-level ensemble to get asymptotic scalings for classification and regression test error (Appendix E).

The bulk of the technical work is involved in proving the matching bounds on survival and contamination, i.e. Theorems 22 and 23. Although these results are inspired by the calculations provided in Appendix A.2 for the Fourier case, we build on the techniques provided in Bartlett et al. (2020) for Gaussian features, particularly making use of fundamental concentration bounds that were proved on “leave-one-out” matrices in that work. We build on these techniques to sharply bound both the “survival” and “contamination” terms, and thus obtain matching upper and lower bounds for the classification test error. Crucially, our analysis needs to circumvent issues that stem from effective misspecification in the linear model that arise from the sign operator. While we do not provide a generic analysis of “misspecification noise,” we exploit the special misspecification induced by the sign operator in a number of technical equivalents of the aforementioned random matrix concentration results.

In fact, the non-asymptotic scalings of survival and contamination terms are unaffected even by non-zero label noise on classification training data, provided that the label noise still preserves non-trivial information about the signal. The survival is further attenuated by a non-zero factor of $(1-2\nu^{*})$, which is strictly positive as long as $\nu^{*}<1/2$. Observe that this is equivalent to a hypothetical scenario where the binary labels take on “shrunk” values $\{-(1-2\nu^{*}),(1-2\nu^{*})\}$ instead of the usual $\{-1,1\}$. As long as $\nu^{*}<1/2$, the magnitude of the labels is strictly non-zero and so the labels still provide useful information for classification.

Finally, it is natural to ask how fundamental our assumptions of Gaussianity on data and bi-level covariance structure are to our main generalization result (Theorem 13). We chose the bi-level ensemble to illustrate the separation between classification and regression in the cleanest possible way. However, Theorems 22 and 23 do provide non-asymptotic expressions for survival and contamination for arbitrary covariance matrices. In principle, these expressions can be plugged into Proposition 17 to get upper and lower bounds on classification error for arbitrary covariance matrices. Further, the analysis of benign overfitting in linear regression (Bartlett et al., 2020; Muthukumar et al., 2020) extends to sub-Gaussian features. In the same spirit, we can show that the results — including the existence of the intermediate regime, in which classification works but regression does not — extend to a weaker assumption of independence and sub-Gaussianity on the underlying features. This extension uses an argument similar to the Fourier-case argument given in Appendix A.2 but requires a more direct treatment of the approximation error arising from misspecification. We provide this argument in a forthcoming note. Our results do not extend to kernel settings, where there can be complex dependencies among the (infinite-dimensional) features. This is an important direction for future work.

Examining margin-based explanations for generalization

In this section, we explore the potential for generalization bounds as a function of training data margin to explain the behavior we have observed for classification tasks in the overparameterized regime. Through simple experiments, we demonstrate that margin-based generalization bounds are uninformative in sufficiently overparameterized settings.

For a particular function class $\mathcal{F}$, uniform convergence bounds conservatively approximate the generalization error of $f\in\mathcal{F}$ by that of the least generalizable function in $\mathcal{F}$. The ensuing generalization bounds typically depend on measures of complexity, such as the Vapnik-Chervonenkis dimension, which increase with the number of parameters in the model. Thus, the uniform convergence approximation is not as good when $\mathcal{F}$ is large, e.g. the model has several parameters. This shortcoming of uniform convergence-based bounds was first brought into focus by the remarkable success of boosting with a very large number of primitive classifiers (Schapire et al., 1998). The main observation was that even after the training 0-1 loss became zero, increasing the number of primitive classifiers in the boosted model still reduced the test error.

An analysis in terms of the training data margin was proposed as a possible explanation for this behavior for classifiers $f(\cdot)$ that make their predictions by discretizing the outputs of a real-valued function $g\in\mathcal{G}$, i.e. $f(X)=\mathsf{sgn}(g(X))$. The training margin, $\gamma:=\min_{i}Y_{i}g(X_{i})$ can be intuitively interpreted as a measure of prediction confidence; for linear classifiers, it is precisely the minimum (over training points) distance to the decision boundary. The worst-case margin is not the only quantity that has been considered: generalization bounds based on a weighted combination of margin on all training data points have also been considered and demonstrated to be sharper in certain settings (Gao and Zhou, 2013). In the settings we investigate, all training data points become support vectors – therefore the margins at each training point are equal, and all such notions of margin become equivalent.

Under certain conditions, margin-based generalization bounds can scale far slower with the number of parameters in the model than uniform convergence bounds; for example, in boosting, the dependence is reduced to $\ln(\text{\# of primitive classifiers})$. Since the margin $\gamma$ could be artificially increased (without changing any of the predictions) simply by rescaling the real-valued function $g(\cdot)$ , the quantity of interest is an appropriately normalized margin, e.g. the margin normalized by the Lipschitz constant of the learned function $g(\cdot)$ or its approximation.

The decrease of generalization error despite increasing complexity in “modern” overparameterized regimes is strongly reminiscent of the observations from boosting with a large number of primitive classifiers. It is of particular interest to examine the ensuing generalization bounds for the hard-margin SVM, which maximizes margin on linearly separable data. For the case of linear classifiers, the normalized margin is defined as $\gamma_{N}=\frac{\gamma}{\left\lVert\boldsymbol{\widehat{\alpha}}\right\rVert_{2}}$. We can now state the ensuing classification test error (i.e. 0-1 test loss) as a function of the normalized margin. Notation in the statement is adapted to be consistent with the notation in this paper — for an elementary verification, see Appendix H.

For a random test point $(X,Y)$ drawn from the same distribution as the training data, the following holds with probability $(1-\delta)$ over the training data $\{X_{i},Y_{i}\}_{i=1}^{n}$:

where the ramp loss function $l_{\gamma}$ is defined by

When the training data are separable, we apply the above bound setting the first (average training loss) term to , and only consider the the second term in the bound, i.e. we ignore the high-probability term. Equation (19) reminds us that there is a critical dependence on the intrinsic data dimension, captured by the term $\|\boldsymbol{\Phi}_{\mathsf{train}}\|_{\mathsf{F}}$. We will shortly see that this dependence is critical to track in the overparameterized regime.

2 Can margin track performance of overparameterized models?

We now investigate whether this generalization bound is effective in tracking the true test classification error for the hard-margin SVM in our setting for a number of choices of featurization. Importantly, we consider the solution $\boldsymbol{\widehat{\alpha}}_{\mathsf{SVM}}$ only in sufficiently overparameterized settings under which all training points become support vectors with high probability; therefore, the un-normalized margin $\gamma=1$ and the normalized margin of the SVM solution is exactly equal to $\gamma_{N}=\frac{1}{\left\lVert\boldsymbol{\widehat{\alpha}}\right\rVert_{2}}$.

We study the evolution of margin, the ensuing upper bound in Equation (19), and the true test classification error as we increase the level of overparameterization for two choices of featurizations: isotropic Gaussian features (Definition 4) which generalize poorly according to Theorem 13 and weak features (Definition 9), which are known to exhibit the double-descent behavior. For the case of isotropic features, we retain our 1-sparse assumption from Section 5. For the case of weak features, we consider $Y_{i}=\mathsf{sgn}(U_{i})$ for $i\in\{1,\ldots,n\}$.

Figure 5 plots the isotropic case, and Figure 6 plots the weak features case. For both figures, we hold the number of training points, $n$ constant and vary the number of features, $d$, to tune the extent of overparameterization. In both Figures 5(a) and 6(a), the normalized margin increases with increasing $d$, since the optimizer can use more features to meet the constraint in Equation (3a). The generalization bounds in Figure 5(b) and Figure 6(b) are consequently very similar as well. However, while the test classification loss increases with $d$ for isotropic features, it decreases with $d$ for weak features.

Figures 5 and 6 together show that the relationship between margin and generalization is more complex than typically assumed in highly overparameterized regimes. We highlight a few observations:

Whether margin is qualitatively predictive of generalization is also unclear, as evidenced by the contrasting examples of weak features and isotropy. Under both featurizations, the normalized margin increases with increased overparameterization; but the actual test error behaves very differently (decreasing for weak features, but increasing for isotropy).

Thus, we see that margin-based bounds are not predictive of the behavior of overparameterized models in our setting. It is still possible that an appropriate sense of large margin implies good generalization in certain cases. In particular, for linear models, maximizing the margin is equivalent to minimizing the norm — which, as we have seen, has important generalization properties. However, evidence of this needs to come from first-principles analysis, not from the existing bounds.

The recent work of Negrea et al. (2019) suggests a way forward in the interpolating regime via the introduction of appropriate surrogates, that implicitly capture the good generalization properties vis-a-vis the underlying patterns and the learning algorithm used. It would be interesting to see how these ideas could be unified with the survival/contamination perspective developed here across all three regimes identified in Theorem 13.

We would like to acknowledge the Simons Institute Summer 2019 program on “Foundations of Deep Learning”, which facilitated the initial collaboration that led to this work. More broadly, we thank the participants of this program for many stimulating research discussions that inspired this collaboration — especially Suriya Gunasekar and Matus Telgarsky.

MB acknowledges federal support from NSF and a Google Research Award. DH acknowledges support from NSF grant CCF-1740833 and a Sloan Research Fellowship. AS acknowledges the support of the ML4Wireless center member companies and NSF grants AST-144078 and ECCS-1343398.

A Fourier features on regularly spaced training data: An “ultra-toy” model

As we illustrate in this section, appropriate weightings of these features under this “ultra-toy” model also helped us conjecture all of the main results of this paper.

We consider $n$ (odd) regularly spaced training points from $(-\pi,+\pi)$ — specifically the sequence $(-\pi+\frac{\pi}{n},-\pi+\frac{3\pi}{n},\ldots,-\frac{2\pi}{n},0,+\frac{2\pi}{n},\ldots,+\pi-\frac{\pi}{n})$, a test distribution of $X$ drawn uniformly at random from $(-\pi,+\pi)$, and the $d$ (odd multiple of $n$) features chosen to be the standard real orthonormal Fourier features:

For doing interpolative inference using a weighted norm minimization, we define the weightsIf $\alpha_{j}$ represents the learned coefficient on the cosine at frequency $f$, and $\beta_{j}$ the learned coefficient on the sine at frequency $f$, the minimization is of $\sum_{j}\frac{\alpha_{j}^{2}+\beta_{j}^{2}}{\lambda_{j}}$. A higher $\lambda_{j}$ means that frequency is favored. corresponding to sines and cosines of frequency $j$ by $\{\lambda_{j}\}_{j=0}^{\frac{(d-1)}{2}}$. Following the convention of the rest of the paper, we take the weights $\{\lambda_{j}\}$ to be a decreasing, strictly positive sequence.

Exact aliases are defined as distinct features that agree with each other (possibly up to a constant multiple) on all the sampled points. The Fourier featurization allows exact aliases to exist. There are three different groups of these exact aliases:

The initial constant feature is essentially aliased by the cosines at every multipleFor ease of exposition, the minor issue of the constant feature having a slightly different scaling vis-a-vis its aliases is going to be ignored in this treatment, but this is simply a matter of keeping track of notation. Alternatively, we could eliminate this by using complex Fourier features. We will finesse this issue here by simply not allowing the true signal to have a constant term in it. of $n$.

Each cosine feature in the first $n$ features (namely corresponding to a frequency $j\in\{1,2,\ldots,\frac{n-1}{2}\}$) picks up $(\frac{d}{n}-1)$ cosine aliases with frequencies $(n-j),(n+j),(2n-j),(2n+j),\ldots$. This is because cosine is an even function and the training samples are symmetrically distributed about 0.

Similarly, each sine feature in the first $n$ features (corresponding to a frequency $j\in\{1,2,\ldots,\frac{n-1}{2}\}$) picks up $(\frac{d}{n}-1)$ sine aliases with frequencies $(n-j),(n+j),(2n-j),(2n+j),\ldots$. However, because sine is an odd function, these aliases have their signs alternating with the $(kn-j)$ ones being multiplied by $(-1)$ and the $(kn+j)$ ones being exact aliases.

A.2 Regression vs Classification

Because of the presence of exact aliases, the training data matrix $\boldsymbol{\Phi}_{\mathsf{train}}$ consists of $n$ distinct columns that repeat again and again. In keeping with the bi-level covariance model in Definition 5, we scale the parameters $(s,\lambda_{H},d)$ with $n$ in a coordinated way. Recall that the number of prioritized features is given by $s:=n^{r}$ for $r\in[0,1)$, and the number of features $d=n+n^{p}$ for $p>1$. (We added an extra term of $n$ to make it easier to count the aliases. This has no asymptotic effect when $p>1$ and $n\rightarrow\infty$.) The $\lambda_{H}$ represents how much we favor the special features and in keeping with the scaling in Definition 5, we set $\lambda_{H}=n^{p-r-q}$ for some $q\in[0,p-r]$.

Because of the known orthogonality of the sine and cosine features on $n$ regularly spaced points, the first $n$ columns of $\boldsymbol{\Phi}_{\mathsf{train}}$ are orthogonal. This means that the solution $\boldsymbol{\widehat{\alpha}}$ will only have nonzero entries in the positions that correspond to the $\frac{d}{n}=1+n^{p-1}$ different columns of $\boldsymbol{\Phi}_{\mathsf{train}}$ that are copies of the column corresponding to the feature $\cos(x)$. Since $s<n$, exactly one of these will be favored and so the optimization problem in Equation (20) turns into the much simpler problem:

where $a$ represents the recovered coefficient corresponding to the true underlying feature $\cos(x)$ and $b$ represents the coefficients on all of its exact aliases.

An elementary calculus calculationSee the appendix of Muthukumar et al. (2020) for more discussion of this calculation and its connection to matched filtering in signal processing. shows that Equation (21) is solved by:

The reader can verify that $a$ represents the survival of the true signal. For large enough $n$, this is approximatedIn the style of the Bode Plot of a one-pole low pass filter. by

This approximation is guaranteed to be good to within a factor of 2 everywhere and is usually much better. Notice that Equation (23) is the Fourier-feature counterpart of the upper and lower bounds on survival in Lemmas 61 (binary labels) and 63 (real-valued output). Now, taking $n\rightarrow\infty$, we get

which shows that the signal only fully survives if $q<(1-r)$.

Let us now measure the contaminating effect of falsely discovered features. Following Equation (15), we denote $B(\mathbf{X})$ as the random variable that represents the contribution of all of the aliases to the predictions. Each of the Fourier features of non-zero frequency is zero-mean and has variance $1$. From the orthonormality (in expectation over test data) of the aliases, we get

where in the last step, we substituted Equation (22b). Notice that $\frac{(p-1)}{2}>\frac{p}{2}+\frac{1}{2}-r-q$ whenever $q>(1-r)$, and so asymptotically we get

This approximation is guaranteed to be good to within a factor of 2 everywhere and is usually much better. Notice that this expression is the Fourier-feature counterpart of the lower bound on contamination established for Gaussian features in Lemma 36.

Thus, provided that $q<(1-r)$, the expression in Equation (26) always decays to zero as $n\rightarrow\infty$, regardless of which case we are in. The combination of Equations (26) and (24) tells us that regression in this problem can workIn fact, the argument also works for noisy training data, i.e. $\mathbf{Z}_{i}=\sin(x_{i})+\mathbf{W}_{i}$ where the noise $\mathbf{W}_{i}$ is iid and has zero mean, bounded support, and finite variance $\sigma^{2}$. The formal argument is in Muthukumar et al. (2020), but is summarized here for the ease of the reader. From the central limit theorem and the theory of wide-sense-stationary random variables, in the limit, the representation of the noise part will look marginally Gaussian in the basis of the first $n$ columns of $\boldsymbol{\Phi}_{\mathsf{train}}$, where each of them will be $~{}\mathcal{N}\left(0,\frac{\sigma^{2}}{n}\right)$. The first $s$ of these will survive and thereby contribute a variance of approximately $\sigma^{2}n^{r-1}a^{2}$ to test points, while the other $(n-s)$ of these will be absorbed across all the aliases and thereby contribute a contamination variance of $\sigma^{2}n^{-p}$. The total contribution will be dominated by the $\sigma^{2}n^{-(1-r)}a^{2}$ term. If $q<(1-r)$ and $a\approx 1$, this term vanishes as $n\rightarrow\infty$ and so additive noise does not contribute to regression error unless the noise variance $\sigma^{2}$ itself grows with $n$ (at a rate faster than $n^{1-r}$). to get mean-square-error approaching zero as long as $q<(1-r)$. On the other hand, when $q>(1-r)$, signal does not asymptotically survive and regression MSE approaches the null risk.

For classification, we only care about predicting $\mathsf{sgn}(\cos(\mathbf{X}))$ correctly with high probability when $\mathbf{X}\sim\text{Unif}\left[-\pi,\pi\right]$. Clearly, classification also works under the conditions for which regression works (i.e. $q<(1-r)$), but, as we showed in Theorem 13, can work even in the absence of these conditions. Recall that when $q>(1-r)$, the survival factor $a\rightarrow 0$ as $n\rightarrow\infty$. However, if the contamination is small enough, i.e. $\sigma_{CN}\ll a$, the probability of classification error is extremely low, as illustrated in Figure 9. We observe from Equations (26) and (23) that $\sigma_{CN}\ll a$ if $q<(1-r)+\frac{(p-1)}{2}$. When that happens, classification will asymptotically work.

To see this more formally, we can upper bound the expression of classification error and show that it goes to zero as $n\to\infty$ under these conditionsThe exact Gaussian-feature expression for classification error in Proposition 17 depends solely on the ratio $a/\sigma_{CN}.$ Characterizing the exact expression for Fourier features is more challenging because the contamination does not have a clean distribution, but we can upper bound the probability of classification error using the standard deviation alone.. We use a union bound together with Chebyshev’s inequality in a manner reminiscent of typicality proofs in information theory (Cover and Thomas, 2012). Let $\epsilon=\frac{(p-1)}{2}-(q-(1-r))$ be the difference between the relevant two exponents of $n$. Then, we define the events $\mathcal{A}:=\{x|\ |\cos(x)|<2n^{-\frac{\epsilon}{2}}\}$ and $\mathcal{B}:=\{x|\ |B(x)|>n^{-\frac{\epsilon}{2}}n^{-(q-(1-r))}\}$. The event $\mathcal{A}$ corresponds to having an atypically weak signal in the true feature, and the event $\mathcal{B}$ corresponds to having an atypically large contamination term. Observe that if neither event $\mathcal{A}$ nor event $\mathcal{B}$ holds, we can substitute Equation (23) to get $|a\cos(x)|\geq 2|B(x)|$, and this implies that a classification error will not occur. Therefore, the probability of classification error is upper bounded by $\Pr\left[\mathcal{A}\cup\mathcal{B}\right]$, and by the union bound it suffices to upper bound the probability of each of these events individually. We start with the “weak signal” event $\mathcal{A}$. Because $\cos(x)$ is a function that is always differentiable in the neighborhood where $\cos(x)=0$, this means that $\cos(\mathbf{X})$ as a random variable has a densityThis is known as a shifted arc-sine distribution. in the neighborhood of . Consequently, we have $\Pr\left[\mathcal{A}\right]=\int_{-n^{-\frac{\epsilon}{2}}}^{+n^{-\frac{\epsilon}{2}}}\frac{1}{\pi\sqrt{1-y^{2}}}dy=\frac{2}{\pi}\sin^{-1}(n^{-\frac{\epsilon}{2}})$ which goes to zero as $n\rightarrow\infty$. We now turn to the “unusually large contamination” event $\mathcal{B}$. Because $q<(1-r)+\frac{(p-1)}{2}$, we have $\Pr\left[\mathcal{B}\right]=\Pr\left[|B(\mathbf{X})|>n^{-\frac{\epsilon}{2}}n^{-(q-(1-r))}\right]\leq\Pr\left[|B(\mathbf{X})|>n^{\frac{\epsilon}{2}}\sigma_{CN}\right]$. By Chebyshev’s inequality, we have $\Pr\left[\mathcal{B}\right]\leq n^{-\frac{\epsilon}{2}}$, which goes to zero as $n\rightarrow\infty$.

Since the probabilities of both events $\mathcal{A}$ and $\mathcal{B}$ have been shown to go to as $n\to\infty$, the limiting classification error will also be zero when $q<(1-r)+\frac{(p-1)}{2}$. In fact, this argument can be extended to the case of interpolation of binary labels by using the Fourier series representation of the underlying true label function. Since there is now misspecification induced by the sign operator, this requires understanding the approximation-theoretic properties of the Fourier series by its first $s$ terms as $s\rightarrow\infty$. Such a treatment is beyond the scope of this paper, and we defer it to a separate forthcoming note.

Finally, it is worth noting that the above calculation only upper bounds the classification error; whether this upper bound is matched by a lower bound remains an open question. This would shed light on whether, in fact, classification generalizes poorly when $q>(1-r)+\frac{(p-1)}{2}$ for the case of Fourier featurization. Figure 10 illustrates the three regimes for Fourier featurization. While the first two regimes display behavior that parallels the Gaussian-feature results in Theorem 13, the third regime is inconclusive with respect to classification performance. This is an important question for future work.

B Additional notation for proofs

We consider zero-mean Gaussian featurization, i.e. $\phi(X_{i})=\mathcal{N}(\mathbf{0},\boldsymbol{\Sigma})$. For ease of exposition, we consider $\boldsymbol{\Sigma}$ to be diagonalThis is without loss of generality: if $\boldsymbol{\Sigma}$ were not diagonal, we could first do a coordinate transformation to the basis of the eigenvectors of $\boldsymbol{\Sigma}$.. Corresponding to a given index $t\in\{1,\ldots,d\}$, we define the “leave-one-out” matrix ${\boldsymbol{\Sigma}}_{-t}$ whose eigenvalues are given by: $\mu_{j}({\boldsymbol{\Sigma}}_{-t})={\widetilde{\lambda}}_{j}$ for $j\in\{1,\ldots,d-1\}$. The relation between the spectrum $\{{\widetilde{\lambda}}_{j}\}_{j=1}^{d-1}$ and $\{\lambda_{j}\}_{j=1}^{d}$ is given by

Consider $\{\mathbf{z}_{i}\}_{i=1}^{d}$ i.i.d. with $\mathbf{z}_{i}\sim\mathcal{N}\left(\mathbf{0},\mathbf{I_{n}}\right)$. Observe that we can write effective Gram matrices corresponding to the full as well as the “leave-one-out” spectrum of the covariance matrix:

Using Equation (27), we can also express the “leave-one-out” Gram matrix $\mathbf{A}_{-t}$ as follows:

We will use both of the above expressions for the leave-one-out matrix $\mathbf{A}_{-t}$ in our analysis.

C Support vector proofs and calculations

Recall that we defined the random Gram matrix $\mathbf{A}$ as

where $\mathbf{z}_{j}\text{ i.i.d.}\sim\mathcal{N}(\mathbf{0},\mathbf{I}_{n})$ reflects all the randomness in the matrix $\mathbf{A}$. Note that the spectrum $\{\lambda_{j}\}_{j=1}^{d}$, and all functionals of it, are deterministic.

The dual to the optimization problem (4) can be expressed as below (Boser et al., 1992):

Note that the unconstrained solution of the above is: $\boldsymbol{\beta}^{*}:=\mathbf{A}^{-1}\mathbf{Y}_{\mathsf{train}}$. By complementary slackness, all of the constraints in the optimization problem (4) will be satisfied with equality at optimum i.e. all training points are support vectors, if we have

Thus, it suffices to establish conditions under which Equation (29) holds with high probability.

We start by showing that this is the case, provided that the condition in Equation (11) holds. To do this, we use the following lemma.

Let $\mathbf{E}:=\mathbf{A}-||\boldsymbol{\lambda}||_{1}\mathbf{I}_{n}$. Then, for any choice of positive constant $0<\epsilon<1$ and $\tau>0$, we have (for large enough $n$),

with probability at least $(1-2e^{-\tau})$ over the randomness in the matrix $\mathbf{A}$.

Lemma 20, which essentially controls the operator norm of the error matrix $\mathbf{E}$ using a union bound with discretization (also known as the “epsilon-net” argument), is proved in Appendix F.1. Now, substituting $\tau:=\ln n$ and $\epsilon:=\frac{1}{36\sqrt{n}}$, all of the following inequalities can be verified to hold for large enough $n$:

with probability at least $(1-2e^{-\ln n})=\left(1-\frac{2}{n}\right)$. Now, observe that as a consequence of Equation (11), we have

Substituting these inequalities above finally gives us

Thus, we have shown that for large enough $n$, we have

with probability at least $\left(1-\frac{2}{n}\right)$.

Now, we denote $\mathbf{E}^{\prime}:=\frac{1}{||\boldsymbol{\lambda}||_{1}}\mathbf{I}_{n}-\mathbf{A}^{-1}$. Observe that when Equation (30) holds, we have

where the last inequality again holds for large enough $n$. Thus, for large enough $n$, we have

Furthermore, since we can write $\mathbf{E}^{\prime}=\frac{1}{||\boldsymbol{\lambda}||_{1}}\cdot\mathbf{A}^{-1}\cdot\mathbf{E}$, we get

where inequality $\mathsf{(i)}$ uses the standard inequality on product of operator norms, inequality $\mathsf{(ii)}$ substitutes Equation (31), and inequality $\mathsf{(iii)}$ substitutes Equation (30). Thus, we get, for every $i\in[n]$,

for large enough $n$ and with probability at least $\left(1-\frac{2}{n}\right)$. Here, inequality $\mathsf{(i)}$ follows from the inequality $\mathbf{a}^{\top}\mathbf{M}\mathbf{b}\geq-||\mathbf{a}||_{2}||\mathbf{b}||_{2}||\mathbf{M}||_{\mathsf{op}}$, and inequality $\mathsf{(ii)}$ follows by substituting the upper bound we just derived on $||\mathbf{E}^{\prime}||_{\mathsf{op}}$. This completes our proof for Part 1 of Theorem 11.

Next, we show that Equation (29) holds with high probability under the condition provided in Equation (12), which is the strictly sharper condition for the isotropic Gaussian case. For every $i\in[n]$, we denote $\mathbf{v}_{i}:=\sqrt{n}\cdot Y_{i}\mathbf{e}_{i}$. We add and subtract terms to get

We use the following technical lemma that shows concentration on quadratic forms of the inverse Wishart matrix $\mathbf{A}^{-1}$. From here on, we denote $d^{\prime}(n):=(d-n+1)$ for shorthand.

Let $\mathbf{A}\sim\text{Wishart}(d,\mathbf{I}_{n})$. For any vector $\mathbf{u}\in S^{n-1}$ and any $t>0$, we have

provided that $d^{\prime}(n)>2\max\{t,1\}$.

Lemma 21 is proved in Appendix F.2. Substituting the lower tail bound of Lemma 21 with $t:=2\ln n$ gives us

with probability at least $\left(1-\frac{1}{n^{2}}\right)$. Similarly, substituting the upper tail bound with $t:=2\ln n$ gives us

with probability at least $\left(1-\frac{1}{n^{2}}\right)$. Noting that $||\mathbf{v}_{i}+\mathbf{Y}_{\mathsf{train}}||_{2}^{2}=2(n+\sqrt{n})$ and $||\mathbf{v}_{i}-\mathbf{Y}_{\mathsf{train}}||_{2}^{2}=2(n-\sqrt{n})$, we then get

which is precisely the condition in Equation (12). Under this condition, we have proved that for any training data point corresponding to $i\in\{1,\ldots,n\}$, we have $Y_{i}\beta^{*}_{i}>0$ with probability at least $\left(1-\frac{2}{n^{2}}\right)$. Finally, applying the union bound on all $n$ training data points gives us

with probability at least $\left(1-\frac{2}{n}\right)$. This completes the proof of Theorem 11. $\hfill\blacksquare$

C.2 Implications of Theorem 11 for the bi-level ensemble

In this section, we provide the calculations that help us understand the ramifications of Theorem 11 — in particular, the condition in Equation (11) — for the bi-level ensemble (Definition 5). We reproduce Equation (11) below:

We substitute the parameters of the bi-level ensemble into the left hand and right hand sides of the inequality. Recall that, by definition, we have $||\boldsymbol{\lambda}||_{1}=d=n^{p}$ for the bi-level ensemble and so the left hand side is equal to $n^{p}$. On the other hand, for the right hand side, a simple calculation shows that

where the scaling in $\mathsf{(i)}$ follows because the bi-level ensemble defines $r<1<p$ and $q>0$ (so $(1-a)\asymp 1$ and $(d-s)\asymp d$). Moreover, we have

Putting these together, the right hand side of Equation (11) scales as

and so, for Equation (11) to hold, we get the following sufficient conditions on the parameters $(p,q,r)$ of the bi-level ensemble for sufficiently largeThe reason for requiring sufficiently large $n$ in these statements is the application of the $\asymp$ relation in multiple places. (Also note that Theorem 11 also required sufficiently large $n$.) Accordingly, we can also omit constants from consideration. $n$:

Now, observe that $1-\frac{r}{2}\leq\frac{3}{2}-r$ for all $0\leq r\leq 1$, and so we get sufficient conditions as follows:

These are precisely the conditions in Equation (13).

D Proof of Theorem 13: Bounds on survival and contamination

In this section, we obtain a general, non-asymptotic characterization of classification (and regression) error by bounding survival and contamination terms. As described in Section 5.2, this is then plugged into the expressions in Proposition 17 to prove Theorem 13.

First, we define shorthand notation that is useful for this section, in addition to the notation already defined in Appendix B. For ease of notation, we denote the survival and contamination factors under the $1$-sparse model for the case where we interpolate binary labels as

and for the case where we interpolate real output as

Finally, for a given index $t\in\{1,\ldots,d\}$, we denote as shorthand $\mathbf{z}_{t}:=\mathbf{Z}_{\mathsf{train}}$. It is easy to verify that $\mathbf{z}_{t}\sim\mathcal{N}(\mathbf{0},\mathbf{I}_{n})$ under the $1$-sparse Assumption 8. We also denote $\mathbf{y}_{t}:=\mathbf{Y}_{\mathsf{train}}$. Recall that we consider the possibility of label noise probability equal to $\nu^{*}$: from the generative model defined in Equation (1), we have

The notions of survival and contamination were first introduced in Muthukumar et al. (2020), and characterized there with equality for Fourier featurization on regularly spaced training data. Here, we characterize these quantities for Gaussian features. We state our upper and lower bounds on survival and contamination respectively for two cases — when the output being interpolated is binary, and when the output being interpolated is real. We start with upper and lower bounds on the survival factor.

There exist universal positive constants $(b,b_{2},c,c_{3},c_{4})$ (that do not depend on parameters $(n,d,k,\boldsymbol{\Sigma})$) such that if $r_{k}(\boldsymbol{\Sigma})\geq bn$ and $r_{k}({\boldsymbol{\Sigma}}_{-t})\geq b_{2}n$, we have the following characterizations of the survival factor for any $k\geq t$:

with probability at least $(1-3e^{-\sqrt{n}}-2e^{-\frac{n}{c}})$ over the randomness in the training data $\{X_{i},Y_{i}\}_{i=1}^{n}$.

with probability at least $(1-2e^{-\sqrt{n}}-2e^{-\frac{n}{c}})$ over the randomness in the training data $\{X_{i},Y_{i}\}_{i=1}^{n}$.

We will see subsequently (in Appendix E) that the survival bounds, whether binary labels or real output are interpolated, are matching in their dependence on $n$ up to constants. We now state our characterization of the contamination factor.

almost surely for any realization of the random quantity $\mathsf{SU}_{b}(t)$, and with probability at least $\left(1-\frac{3}{n}\right)$ and $(1-2e^{-\frac{n}{c_{8}}})$ respectively over the randomness in the training data $\{X_{i},Y_{i}\}_{i=1}^{n}$.

almost surely for any realization of the random quantity $\mathsf{SU}_{b}(t)$, and with probability at least $\left(1-\frac{2}{n}\right)$ and $(1-2e^{-\frac{n}{c_{8}}}-e^{-n\delta^{2}})$ respectively over the randomness in the training data $\{X_{i},Y_{i}\}_{i=1}^{n}$.

Observe that the high-probability characterizations of contamination in Theorem 23 themselves hold almost surely for every realization of the respective survival factors for binary and real interpolation, which are random variables. In Appendix E, these expressions will be used together (with a simple union bound) with the matching high-probability characterization of survival factor in Theorem 22. Unlike for the case of survival, the upper and lower bounds for contamination are not necessarily matching — however, as we will see in Appendix E, they turn out to match for all parameterizations of the bi-level ensemble.

As a final remark, in both theorem statements, the only randomness over which all probabilities are taken is solely in the training data $\{X_{i},Y_{i}\}_{i=1}^{n}$. Further, all universal positive constants are taken to be independent of the parameters $(n,d,k,\boldsymbol{\Sigma})$, which entirely describe the problem. In the proofs of Theorems 22 and 23, we will follow these conventions unless specified otherwise.

D.2 Background lemmas

We begin our proofs of Theorems 22 and 23 by stating lemmas that serve as background for our analysis. The first lemma is from Bartlett et al. (2020).

(Concentration of eigenvalues, Lemmas 9 and 10 in Bartlett et al., 2020) There exist universal positive constants $(b,c)$ such that:

For any $k\geq 0$ such that $r_{k}(\boldsymbol{\Sigma})\geq bn$, we have

with probability at least $(1-2e^{-\frac{n}{c}})$ over the random matrix $\mathbf{A}$.

For any $k\geq t$ such that $r_{k}(\boldsymbol{\Sigma})\geq bn$, we have

with probability at least $(1-2e^{-\frac{n}{c}})$ over the random matrix $\mathbf{A}_{-t}$.

Further, as corollaries to the above, we have the following statements:

For any $k\geq 0$ such that $r_{k}(\boldsymbol{\Sigma})\geq bn$, we have

with probability at least $(1-2e^{-\frac{n}{c}})$ over the random matrix $\mathbf{A}$.

For any $k\geq t$ such that $r_{k}(\boldsymbol{\Sigma})\geq bn$, we have

with probability at least $(1-2e^{-\frac{n}{c}})$ over the random matrix $\mathbf{A}_{-t}$.

Note that using Equation (28) to express $\mathbf{A}_{-t}$, we can rewrite the bounds in the above lemma in terms of the quantities ${\boldsymbol{\Sigma}}_{-t}$ and ${\widetilde{\lambda}}_{j}$. In particular, it follows that each of

holds with probability at least $(1-2e^{-\frac{n}{c}})$. We will also apply Equation (38) with $\mathbf{A}_{-t}$ instead of $\mathbf{A}$, and use the corresponding condition $r_{k}(\boldsymbol{\Sigma})\geq b_{2}n$.

The next lemma is the Hanson-Wright inequality, which shows that the quadratic form of a (sub)-Gaussian random vector concentrates around its expectation.

(Hanson-Wright inequality, Rudelson and Vershynin, 2013) Let $\mathbf{z}$ be a random vector composed of i.i.d. random variables that are zero mean and sub-Gaussian with parameter at most $1$. Then, there exists universal constant $c>0$ such that for any positive semi-definite matrix $\mathbf{M}$ and for every $t\geq 0$, we have

We will apply this inequality in two ways. First, we will note that $||\mathbf{M}||_{\mathsf{F}}^{2}\leq n||\mathbf{M}||_{\mathsf{op}}^{2}$ and substitute $t:=c_{1}||\mathbf{M}||_{\mathsf{op}}\cdot n^{3/4}$ (where $c_{1}^{2}=\frac{1}{c}$) to get

with probability at least $(1-2e^{-\sqrt{n}})$. Second, we will note that $||\mathbf{M}||_{\mathsf{op}}\leq\mathsf{tr}(\mathbf{M})$ and moreover, $||\mathbf{M}||_{\mathsf{F}}^{2}=\mathsf{tr}(\mathbf{M}^{2})\leq(\mathsf{tr}(\mathbf{M}))^{2}$. Then, substituting $t:=\frac{1}{c}\cdot\mathsf{tr}(\mathbf{M})\cdot(\ln n)$, we get

with probability at least $(1-\frac{1}{n})$. Finally, note that all probabilities are only over the random vector $\mathbf{z}$. We will frequently apply Lemma 25 as a high-probability statement conditioned on the realization of a random, almost surely positive semi-definite matrix $\mathbf{M}$ which is independent of $\mathbf{z}$.

Finally, the following lemma bounds the squared norm of a Gaussian random vector by a standard tail bound on chi-squared random variables (for e.g. see Wainwright, 2019, Chapter 2), stated for completeness.

Let $\mathbf{z}\sim\mathcal{N}(\mathbf{0},\mathbf{I}_{n})$. Then, for any $\delta\in(0,1)$, we have

with probability at least $(1-2e^{-n\delta^{2}})$.

D.3 Proof of Theorem 22

We first prove Theorem 22, i.e. upper and lower bounds on survival when binary labels or real output are interpolated. We start with the slightly more difficult case of interpolation of binary labels (Equations (33a) and (33b)).

Recall that, by Assumption 8, we have $\alpha^{*}_{t}=\frac{1}{\sqrt{\lambda_{t}}}$. A standard argument based on Moore-Penrose pseudoinverse calculations shows that $\boldsymbol{\widehat{\alpha}}_{2,\mathsf{binary}}=\boldsymbol{\Phi}_{\mathsf{train}}^{\top}(\boldsymbol{\Phi}_{\mathsf{train}}\boldsymbol{\Phi}_{\mathsf{train}}^{\top})^{-1}\mathbf{Y}_{\mathsf{train}}$. We get

where $\mathbf{z}_{t},\mathbf{y}_{t}$ are as defined at the beginning of Appendix D, and $\mathbf{A}$ is the Gram matrix defined in Appendix B. Next, we use the Sherman-Morrison-Woodbury identity to get

Adding and subtracting terms to the numerator, we get

Because of the “leave-one-out” property, note that $\mathbf{A}_{-t}^{-1}\perp\{\mathbf{z}_{t},\mathbf{y}_{t}\}$. Also note that $\mathbf{A}_{-t}^{-1}$ is almost surely positive semidefinite. Thus, we can upper and lower bound the numerator of Equation (47) around its expectation using the Hanson-Wright inequality. First, we calculate the conditional expectation:

Recalling the expression for $\mathbf{y}_{t}$ from Equation (32), a simple calculation yields that

where the last step follows because $z_{t,1}\sim\mathcal{N}(0,1)$.

Now, we apply Equation (43) (the Hanson-Wright inequality) almost surely for every realization of the random matrix $\mathbf{A}_{-t}^{-1}$, and simultaneously to the quadratic forms $(\mathbf{z}_{t}+\mathbf{y}_{t})^{\top}\mathbf{A}_{-t}^{-1}(\mathbf{z}_{t}+\mathbf{y}_{t})$ and $(\mathbf{z}_{t}-\mathbf{y}_{t})^{\top}\mathbf{A}_{-t}^{-1}(\mathbf{z}_{t}-\mathbf{y}_{t})$. Thus, we have each of

with probability at least $(1-2e^{-\sqrt{n}})$ over the randomness in $\{\mathbf{z}_{t},\mathbf{y}_{t}\}$. Similarly, to bound the the denominator, we have each of

with probability at least $(1-e^{-\sqrt{n}})$ over the randomness in $\{\mathbf{z}_{t},\mathbf{y}_{t}\}$. Substituting these bounds into Equation (47), we get each of

with probability at least $(1-3e^{-\sqrt{n}})$ over the randomness in $\{\mathbf{z}_{t},\mathbf{y}_{t}\}$. It remains to obtain high-probability bounds on the random quantities $\mathsf{tr}(\mathbf{A}_{-t}^{-1})$ and $||\mathbf{A}_{-t}^{-1}||_{\mathsf{op}}$. Note that we need both lower bounds and upper bounds on the quantity $\mathsf{tr}(\mathbf{A}_{-t}^{-1})$, but we only need an upper bound on the quantity $||\mathbf{A}_{-t}^{-1}||_{\mathsf{op}}$.

We assume that we can choose $k\geq t$ such that $r_{k}(\boldsymbol{\Sigma})\geq bn$ and $r_{k}({\boldsymbol{\Sigma}}_{-t})\geq b_{2}n$ for universal positive constants $(b,b_{2})$. Consider any such choice of $k$ (which in general could depend on $(n,d)$). First, we use Equation (41) from Lemma 24 to upper bound the quantity $||\mathbf{A}_{-t}^{-1}||_{\mathsf{op}}$ as

with probability at least $(1-e^{-\frac{n}{c}})$ over the random matrix $\mathbf{A}$. Next, we turn to the quantity $\mathsf{tr}(\mathbf{A}_{-t}^{-1})$. To lower bound this quantity, we notice that

Now, from Equation (38) in Lemma 24 applied with $\mathbf{A}_{-t}$, we have

with probability at least $(1-e^{-\frac{n}{c}})$ provided that $r_{k}({\boldsymbol{\Sigma}}_{-t})\geq b_{2}n$. This gives us:

with probability at least $(1-e^{-\frac{n}{c}})$. On the other hand, the upper bound on the trace follows simply by

where the last inequality substitutes Equation (42a), which again holds with probability at least $(1-e^{-\frac{n}{c}})$. Noting that the upper bound on $\mathsf{SU}_{b}(t)$ is monotonically increasing in both $\mathsf{tr}(\mathbf{A}_{-t}^{-1})$ and $||\mathbf{A}_{-t}^{-1}||_{\mathsf{op}}$, and the lower bound on $\mathsf{SU}_{b}(t)$ is monotonically increasing in $\mathsf{tr}(\mathbf{A}_{-t}^{-1})$ but decreasing in $||\mathbf{A}_{-t}^{-1}||_{\mathsf{op}}$, we can substitute the above bounds on these quantities. This completes our characterization of survival when binary labels are interpolated, with the probability of this characterization lower bounded by taking a union bound over the complement of all the above events. After taking this union bound, the probability of each of the lower bound (Equation (33a)) and upper bound (Equation (33b)) holding is at least $(1-3e^{-\sqrt{n}}-2e^{-\frac{n}{c}})$.

D.3.2 Interpolation of real output

Again, using the Sherman-Morrison-Woodbury identity, we have

From Equations (48a) and (48b) above, the following statements each hold with probability at least $(1-e^{-\sqrt{n}})$ over the randomness in $\mathbf{z}_{t}$ and for every realization of the random matrix $\mathbf{A}_{-t}^{-1}$:

Here, $c_{2}$ is a universal positive constant.

Observe that the right hand side of Equation (52) is increasing in the quantity $\mathbf{z}_{t}^{\top}\mathbf{A}_{-t}^{-1}\mathbf{z}_{t}$. Thus, substituting the lower bound for $\mathsf{tr}(\mathbf{A}_{-t}^{-1})$ from Equation (50) and the upper bound for $||\mathbf{A}_{-t}^{-1}||_{\mathsf{op}}$ from Equation (49) lower bounds the quantity $\mathbf{z}_{t}^{\top}\mathbf{A}_{-t}^{-1}\mathbf{z}_{t}$, yielding the lower bound for $\mathsf{SU}_{r}(t)$. Similarly, substituting the upper bound for $\mathsf{tr}(\mathbf{A}_{-t}^{-1})$ from Equation (51) and the upper bound for $||\mathbf{A}_{-t}^{-1}||_{\mathsf{op}}$ from Equation (49) upper bounds the quantity $\mathbf{z}_{t}^{\top}\mathbf{A}_{-t}^{-1}\mathbf{z}_{t}$, yielding the upper bound for $\mathsf{SU}_{r}(t)$. This completes the proof of Theorem 22. Again, a simple application of the union bound shows that each of the lower bound (Equation (34a)) and the upper bound (Equation (34b)) hold with probability at least $(1-2e^{-\sqrt{n}}-2e^{-\frac{n}{c}})$.

D.4 Proof of Theorem 23

We next prove Theorem 23, i.e. upper and lower bounds on contamination, for the cases of interpolating binary labels and real output. Since the contamination factor is intricately related to the contribution of additive noise to regression test error, the proof primarily consists of refinements of the arguments in Bartlett et al. (2020).

We start with a useful set of expressions for the contamination factor in the following lemma. The proof of this lemma is contained in Appendix F.3.

We will use the expression in Equation (53b) to prove an upper bound on contamination, and the expression in Equation (53a) for the lower bound.

We start with the proof for the upper bound on contamination for interpolation of binary labels (Equation (35a)). From Equation (53b) in Lemma 27, we have $\mathsf{CN}_{b}^{2}(t)=\widetilde{\mathbf{y}_{t}}^{\top}\mathbf{\widetilde{C}}\widetilde{\mathbf{y}_{t}}$. Note that by construction, $\mathbf{\widetilde{C}}$ has no dependence on $\{\mathbf{z}_{t},\mathbf{y}_{t}\}$ and thus $\mathbf{\widetilde{C}}\perp\widetilde{\mathbf{y}_{t}}$. The next lemma upper bounds the term $\widetilde{\mathbf{y}_{t}}^{\top}\mathbf{\widetilde{C}}\widetilde{\mathbf{y}_{t}}$ in terms of $\mathsf{tr}(\mathbf{\widetilde{C}})$ and is proved in Appendix F.4.

There exists universal positive constant $c_{6}$ such that when $n\geq c_{6}$, we have

almost surely for every realization of the random matrix $\mathbf{\widetilde{C}}$, and with probability at least $\left(1-\frac{2}{n}\right)$ over the randomness in $\widetilde{\mathbf{y}_{t}}$.

almost surely for every realization of the random matrix $\mathbf{\widetilde{C}}$, and with probability at least $\left(1-\frac{2}{n}\right)$ over the randomness in $\widetilde{\mathbf{y}_{t}}$. The next lemma, which is taken from Bartlett et al. (2020), provides a high-probability upper bound on the quantity $\mathsf{tr}(\mathbf{\widetilde{C}})$.