From Cross-Validation to SURE: Asymptotic Risk of Tuned Regularized Estimators

Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.

42,701 characters · 7 sections · 28 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

\onehalfspacing

abstractWe derive the asymptotic risk function of regularized empirical risk minimization (ERM) estimators tuned by $n$-fold cross-validation (CV). The out-of-sample prediction loss of such estimators converges in distribution to the squared-error loss (risk function) of shrinkage estimators in the normal means model, tuned by Stein's unbiased risk estimate (SURE). This risk function provides a more fine-grained picture of predictive performance than uniform bounds on worst-case regret, which are common in learning theory: it quantifies how risk varies with the true parameter. As key intermediate steps, we show that (i) $n$-fold CV converges uniformly to SURE, and (ii) while SURE typically has multiple local minima, its global minimum is generically well separated. Well-separation ensures that uniform convergence of CV to SURE translates into convergence of the tuning parameter chosen by CV to that chosen by SURE.

Introduction

\paragraph{Background} The goal of supervised learning is to produce good predictions for new observations.\footnote{We would like to dedicate this paper to Gary Chamberlain, whose conversations provided the original inspiration for this project.} An important class of estimators for supervised learning can be described as {regularized empirical risk minimization (ERM) estimators} that are {tuned using cross-validation (CV)}.

{ERM estimators} minimize in-sample average prediction loss (empirical risk), among a given class of predictors. Examples include ordinary least squares and maximum likelihood. ERM estimators are prone to overfitting if the class of predictors is large. Such estimators achieve low in-sample loss, but perform poorly for new observations. To counter overfitting, {regularization} is used. Regularization adds a penalty term to the ERM objective; common penalties include the $L^{2}$ norm of parameters (in Ridge regression) and the $L^{1}$ norm (in Lasso regression). Adding a penalty avoids overfitting, by reducing the estimator variance, at the cost of introducing some bias, which might result in underfitting.

To achieve good performance, avoiding both overfitting and underfitting, the amount of penalization needs to be carefully {tuned}. This can be done by choosing the weight on the penalty term as the minimizer of a {cross-validation} estimate of predictive loss. We focus on n-fold CV, where predictions are evaluated for one hold-out observation at a time, and predictive loss is estimated by averaging evaluations over each of the $n$ observations.

\paragraph{Risk functions} The present paper characterizes the behavior of estimators of this form by deriving an asymptotic approximation to their {risk function}. A large literature in learning theory characterizes such estimators by proving bounds on their worst-case regret---the supremum over data-generating processes (DGPs) of the difference between an estimator's risk and the expected loss of the best predictor in the given class. Such bounds provide strong robustness guarantees, but they might not be informative about the behavior of predictive algorithms for realistic DGPs. By focusing on the risk function, we obtain a more fine-grained characterization: the risk function tells us how expected predictive performance depends on the DGP, while worst-case regret bounds only characterize performance for the least favorable DGP.

Our asymptotic characterization relates the risk function of regularized ERM estimators tuned using CV to the risk function of the James-Stein (JS) shrinkage estimator james1961estimation, and generalizations thereof. JS shrinkage famously dominates maximum likelihood estimation (MLE) in the normal means setting: The risk function (mean squared error) of JS shrinkage, as characterized in stein1981estimation, is lower than the risk of MLE, for any possible DGP. Figure (ref) illustrates. Our result suggests that this same risk improvement carries over, asymptotically, to CV-tuned penalized estimation in general parametric models. Like CV, SURE provides an unbiased estimate of mean squared error; the JS estimator approximately minimizes SURE over the shrinkage intensity.

figure[figure omitted — 3,183 chars of source]

\paragraph{Main result} Our main theorem states that the distribution of the out-of-sample prediction loss of CV-tuned regularized ERM estimators converges to the distribution of the squared-error loss of the corresponding SURE-tuned shrinkage estimator in the Gaussian limit experiment. In particular, the risk function (expected out-of-sample prediction loss, as a function of the true parameter) converges to the mean squared error of SURE-tuned penalized estimation in the normal means model. Two of our intermediate results are of independent interest: the uniform approximation of $n$-fold CV by SURE, and the generic well-separation of the global minimum of SURE.

\paragraph{Key steps} Let us briefly outline the three main parts of our proof. First, we show that in large samples ERM estimators are approximately normally distributed, and that out-of-sample predictive loss is approximately equal to squared error loss. These are standard results, and we follow van2000asymptotic in proving this step. We furthermore need to show that this approximation carries over to penalized estimators, for fixed tuning parameters. Our asymptotic approximations are based on local-to-0 asymptotics: As sample size $n$ increases, the parameter vector drifts to $0$ (which is the minimizer of the penalty term) at a rate of $1/\sqrt{n}$. This rate is such that both bias and variance remain non-negligible in the limit. The drifting-parameter framework is natural for studying regularized estimators, because it is the regime in which the penalty has a first-order effect: if the true parameter were fixed away from $0$, the penalty would become asymptotically irrelevant, while if it were exactly $0$, no bias--variance tradeoff would arise.

Second, we need to show that n-fold CV, as a random function that maps tuning parameter values to estimates of predictive loss, converges uniformly to SURE. In this step of the proof, we build on the prior work of wilsonkasymackey2018. This step involves an influence function approximation for leave-one-out estimators, and a second-order approximation to predictive risk. Uniformity of convergence in the tuning parameter is key to this step, and we need to carefully specify regularity conditions such that uniformity is guaranteed.

Third, we need to show that convergence of the CV criterion function for tuning is sufficient for convergence of its minimizer, the tuned parameter. This is non-trivial, because both CV and SURE typically have multiple local minima, and might have multiple global minima. We need to show that generically the global minimum is well-separated. We do this using separate arguments for $L^{1}$ and $L^{2}$ penalties, characterizing the shape and behavior of SURE in either case.

We should emphasize some limitations of our analysis: We do not consider penalties beyond $L^1$ and $L^2$, or $k$-fold CV with $k < n$; extending our results in these directions is left for future work. Our asymptotic approximations are furthermore not appropriate in the over-parametrized regime, when the number of parameters is of similar or larger magnitude than the number of observations.

\paragraph{Literature} The analysis in this paper connects several lines of work in statistics, econometrics, and machine learning. Leave-one-out cross-validation was formalized by Stone1974; its asymptotic optimality for model selection was established by Li1987. A closely related family of risk estimators includes Mallows' $C_p$ Mallows1973, generalized cross-validation Golub1979, and Stein's unbiased risk estimate stein1981estimation. efron2004estimation provides a unifying perspective, showing that these criteria all take the form of in-sample error plus a covariance penalty. arlot2010survey provide a comprehensive survey of cross-validation procedures and their theoretical properties. A noteworthy contrast with our results is provided by shao1993linear, who shows that leave-one-out CV is asymptotically inconsistent for model selection---it selects overfitted models with positive probability---while $k$-fold CV with $k = o(n)$ is consistent. Our result is complementary in focus: rather than asking which of a finite list of models has the best predictive ability, we characterize the limiting distribution and risk function of the estimator selected by leave-one-out CV, showing it converges to that of the SURE-tuned normal-means estimator.

Shrinkage estimation originates with james1961estimation and was analyzed in depth by stein1981estimation. SURE-based tuning of shrinkage was extended to wavelet thresholding by Donoho1995. The two penalty families we study---Ridge Hoerl1970 and Lasso tibshirani1996regression---are the most widely used forms of regularization in practice, and both are commonly tuned by cross-validation. The close relationship between leave-one-out CV and covariance-penalty criteria such as SURE is further illuminated by zou2007degrees, who show---using SURE as the analytical tool---that the effective degrees of freedom of the Lasso equals the number of nonzero fitted coefficients.

Using local asymptotic frameworks to characterize decision problems was pioneered by le1972limits; work using this approach is reviewed in hirano2020asymptotic. The use of shrinkage asymptotics for parametric models in econometrics is discussed in hansen2016efficient. We build directly on wilsonkasymackey2018, who provide non-asymptotic deterministic guarantees for approximate cross-validation, showing that leave-one-out estimators can be well approximated by a single Newton step from the full-sample estimator.

\paragraph{Roadmap} The remainder of this paper is structured as follows: In Section (ref), we introduce our model and assumptions, and define all relevant notation. In Section (ref), we first provide a heuristic outline of our proof, and then state a series of intermediate lemmas. Proofs of all lemmas are collected in the appendix. In Appendix (ref), we prove the lemmas corresponding to the first part of our argument, involving influence function approximations and asymptotic normality. In Appendix (ref), we prove the second part of the argument, namely the (uniform) approximation of n-fold CV by SURE. In Appendices (ref) and (ref), we show that the global minimizer of SURE is generically well-separated for both $L^{2}$ and $L^{1}$ penalties, thereby proving the third part. In Appendix (ref), we conclude our derivation, proving the convergence of risk for tuned estimators.

Setup

In the following, we first set up our estimators, and their asymptotic counterparts, in a series of definitions. We then state the assumptions that will be invoked to justify our asymptotic approximations.

Throughout this paper, we consider the problem of estimating a parameter vector $\beta_0$ which is defined as the minimizer of expected loss $E[l(\beta, Z)]$. For prediction problems, typically $Z=(W,Y)$, for predictive features $W$ and outcomes $Y$. Examples include linear OLS regression, where $l(\beta, Z) = (Y- W \cdot \beta)^2$, as well as the use of neural nets\footnote{With a small number of parameters relative to the sample size.} or other parametric models for classification, where $l(\beta, Z) = -\log(f(Y|W, \beta))$, and $f(Y|W,\beta)$ is the probability assigned to outcome $Y$ by the model.

We use local-to-0 coordinates, $\theta = \sqrt{n} \cdot \beta$, for sample size $n$, and correspondingly $\theta_0 = \sqrt{n} \cdot \beta_0$. For each $n$, the random vectors $Z_n^i$ are i.i.d. draws from the distribution $\mu_n$, across $i$. Any finite-sample object will be denoted by a subscript $n$; objects without subscripts correspond to the limiting experiment.

Definitions

Definition (ref) introduces notation for loss functions and for limiting loss functions. We evaluate estimates of $\theta$ in terms of their expected loss $\bar L_n(\theta, \theta_0)$. For supervised learning, $\bar L_n(\theta, \theta_0)$ is the out-of-sample expected prediction error.

definition[Loss, empirical loss, and expected loss]$\;$\\ Given the loss function $l(\beta, z)$, define the following. \bals l_n(\theta, z) &= l(\theta / \sqrt{n} ,z) & Loss function in local parameter\\ L_n(\theta) &= \sum_{i=1}^n l_n(\theta,Z_n^i) & Empirical loss\\ \bar L_n(\theta, \theta_0) &= E\left[ L_n(\theta) - L_n(\theta_0)\right] & Expected loss\\ \bar L(\theta, \theta_0) &= \lim_{n\rightarrow \infty} \bar L_n(\theta, \theta_0). & Limiting expected loss \eals We assume that the sequence of distributions $\mu_n$ is such that the limit in the last definition is well-defined.

\paragraph{Scaling} Note that we could have equivalently defined $$\bar L_n(\theta, \theta_0) = n\cdot E\left[ l_n(\theta,Z_n^{n+1}) - l_n(\theta_0,Z_n^{n+1})\right],$$ that is, $\bar L_n(\theta, \theta_0)$ is the expected regret for out-of-sample predictions, multiplied by the sample size $n$. The multiplication by $n$ is required because in Definition (ref), we do not scale empirical loss $L_n(\theta)$ by a factor $\tfrac1n$, and therefore $L_n(\theta)$ diverges. In the definition of the local parameter vector $\theta$ we have however re-scaled the parameter vector $\beta$ by a factor of $\tfrac1{\sqrt{n}}$. This implies that the Hessian (second derivative) of $L_n(\theta)$ with respect to $\theta$, if it exists, is given by $$ \nabla_\theta^2 L_n(\theta) = \frac1n \sum_{i=1}^n \nabla_\beta^2 l(\theta/\sqrt{n},Z_n^i). $$ We can thus expect that this Hessian converges, by a law of large numbers, under some additional regularity conditions.

\paragraph{Estimators of $\theta$ and their asymptotic counterparts} We next specify a series of estimators, in Definition (ref). We start with the standard ERM estimator $\hat \theta_n = \operatorname*{argmin\;}_\theta L_n(\theta)$, and its limiting counterpart $\hat \theta$, where the latter is normally distributed with mean $\theta_0$ and variance $\Sigma$. We then consider the regularized versions of these estimators, that is, the penalized ERM estimator with penalty $\lambda \cdot \pi(\theta)$, and its limiting counterpart. Our choice of local parametrization ensures that a constant value of $\lambda$ along the sequence indexed by $n$ leads to a non-degenerate limit for the penalized ERM estimator, where neither variance nor bias of this estimator vanish.

In Definition (ref), we furthermore introduce leave-one-out loss, and the corresponding leave-one-out ERM and penalized ERM estimator. These will serve as the building blocks of n-fold cross-validation.

definition[Estimators of $\theta_0$] We will consider the following estimators of $\theta_0$, for finite $n$, and in the limit experiment: \bals \hat \theta_n &= \operatorname*{argmin\;}_\theta L_n(\theta) & ERM estimator\\ \hat \theta &\sim N(\theta_0, \Sigma) & Limiting ERM estimator\\ \hat \theta_n^\lambda &= \operatorname*{argmin\;}_\theta \left[L_n(\theta) + \lambda \cdot \pi(\theta)\right] & Penalized ERM estimator\\ \hat \theta^\lambda &= \operatorname*{argmin\;}_\theta \left[\tfrac12\|\theta - \hat \theta\|^2 + \lambda \cdot \pi(\theta)\right], &Limiting penalized ERM estimator \eals where $\pi( \cdot )$ is convex and attains its minimum at $0$.\\ We furthermore consider the following leave-one-out (LOO) loss and estimators of $\theta_0$: \bals L_n^{-i}(\theta) &= \sum_{j\neq i} l_n(\theta,Z_n^j) & LOO empirical loss\\ \hat \theta_n^{-i} &= \operatorname*{argmin\;}_\theta L_n^{-i}(\theta) & LOO ERM estimator\\ \hat \theta^{\lambda,-i}_n &= \operatorname*{argmin\;}_\theta \left[L_n^{-i}(\theta) + \lambda \cdot \pi(\theta)\right]. &\text{LOO penalized ERM estimator} \eals

We can rewrite the limiting penalized ERM estimator as $$ \hat \theta^\lambda = \hat \theta + g^\lambda(\hat \theta), $$ where \bals g^\lambda(\theta) = \operatorname*{argmin\;}_g \tfrac12\|g\|^2 + \lambda \cdot \pi(\theta + g). \eals Denote $\nabla g^\lambda(\theta)$ the derivative of $g^\lambda(\theta)$ where it exists, and define $\nabla g^\lambda(\theta) = 0$ at points where $g^\lambda(\theta)$ is not differentiable.\footnote{This convention is adopted for convenience; we will use it to handle Lasso ($L^1$) penalties in the proof of Lemma (ref) below.}

\paragraph{Estimators of risk} The preceding definition introduced penalized estimators for given, fixed values of the tuning parameter $\lambda$. We are interested in estimators which choose this tuning parameter in a data-dependent way, where $\lambda$ minimizes an estimator of risk. For finite sample size $n$, we consider the n-fold crossvalidation (CV) criterion as an estimator of the risk of penalized ERM estimation. For the limit experiment, we consider Stein's Unbiased Risk Estimator (SURE) as an estimator of the risk of the penalized limiting ERM estimator.

definition[Estimators of risk] We consider the following estimators of risk for penalized ERM estimators with fixed tuning parameter $\lambda$. \bals CV_n(\lambda) &=\sum_i l_n(\hat \theta^{\lambda,-i}_n,Z_n^i), &n-fold CV\\ SURE(\lambda, \hat \theta, \Sigma) &= \operatorname{trace}(\Sigma) + \|g^\lambda(\hat \theta)\|^2 + 2 \operatorname{trace} \left ( \nabla g^\lambda(\hat \theta)\cdot \Sigma\right ). &SURE \eals

\paragraph{Tuned estimators} We can now formally define our tuned estimators. $\hat \theta_n^{*}$ is the penalized ERM estimator using a tuning parameter $\lambda_n^*$ which minimizes the n-fold CV estimator of risk. $\hat \theta^{*}$ is the penalized limiting ERM estimator using a tuning parameter $\lambda^*$ which minimizes the SURE estimator of risk. The tuning parameter is chosen from a set $\Lambda \subset \mathbb R$. Later, we will consider $\Lambda = \mathbb R$ (for Ridge penalties), and $\Lambda$ arbitrary but finite (for Lasso penalties).

definition[Tuned estimators of $\theta$] \bals \hat \theta_n^{*} & = \hat \theta_n^{\lambda_n^*}, & Penalized ERM tuned using CV\\ &\lambda_n^* = \operatorname*{argmin\;}_{\lambda \in \Lambda} CV_n(\lambda) \\ \hat \theta^{*} &= \hat \theta^{\lambda^*}, & Limiting penalized ERM tuned using SURE\\ &\lambda^* = \operatorname*{argmin\;}_{\lambda \in \Lambda} SURE(\lambda, \hat \theta, \Sigma). \eals

We evaluate estimators based on their expected loss for new data-points. For supervised learning, this corresponds to the out-of-sample expected prediction loss. In this paper, we do not consider global criteria such as worst-case risk (maximizing over $\theta_0$) or Bayes risk (averaging over a prior distribution for $\theta_0$). Instead, we are interested in the dependence of expected loss on the parameter $\theta_0$, which is captured by the risk function. Our notation makes this dependence on $\theta_0$ explicit. The risk function gives a more fine-grained picture of estimator performance, relative to global criteria such as worst-case risk, Bayes risk, or worst-case regret.

The parameter $\theta_0$ enters the following expressions both directly, as an argument of $\bar L_n$ and $\bar L$, and implicitly, via the distribution of $Z$ that the expectations are averaging over.

definition[Risk functions] \bals R_n(\theta_0) &= E\left[\bar L_n(\hat \theta_n^{*}, \theta_0) \right] & Finite sample risk, tuned using CV\\ R(\theta_0) &= E\left[\bar L(\hat \theta^{*}, \theta_0)\right] & Limiting risk, tuned using SURE. \eals

Assumptions

Having defined our estimators and evaluation criteria, we next specify the assumptions invoked in our asymptotic analysis. Assumption (ref) sets up a sequence of experiments, indexed by $n$. We assume that, for each $n$, the minimizer of expected loss $E[l(\beta, Z_n^i)]$ is given by $\theta_0/\sqrt{n}$. Put differently, the minimizer $\beta$ of expected loss drifts towards $0$. We furthermore assume that the variance $\Sigma$ of the score $\nabla_\beta l(\theta_0/\sqrt{n}, Z_n^i)$ remains constant along our sequence.

assumption[Sequence of experiments] For each $n$, the random vectors $Z_n^i$ are i.i.d. draws from the distribution $\mu_n$, across $i$. The distributions $\mu_n$ are such that $\theta_0$ and $\Sigma$ do not depend on $n$, where \bals \theta_0 &= \operatorname*{argmin\;}_\theta E[l(\theta /\sqrt{n}, Z_n^i)],\\ \Sigma &= \operatorname{Var}\left(\nabla_\beta l(\theta_0/\sqrt{n}, Z_n^i)\right). \eals

The limiting Hessian $H = \nabla_\theta^2 \bar L(\theta, \theta_0)$ is typically non-degenerate because of our scaling of $\bar L_n(\theta, \theta_0)$ and of $\theta$. The following Assumption (ref) is made for notational convenience. This assumption states that the Hessian $H $ is equal to the identity $I$. This is a coordinate normalization that can be imposed without loss of generality.\footnote{By suitable choice of coordinates we can normalize either $H$, or the asymptotic variance $\Sigma$ of Assumption (ref), or the Hessian of the penalty function $\pi$ (when the latter exists), {but not more than one of these three matrices}, in general. After normalizing the Hessian, we can however diagonalize one more matrix, without loss of generality.}

assumption[Normalized loss function] $$\nabla_\theta^2 \bar L(\theta, \theta_0)|_{\theta = \theta_0} = I.$$

The last part of our proof requires showing that the global optimum of $SURE$ with respect to $\lambda$ is generically unique and well-separated. We will prove this fact for both Ridge and Lasso penalties, using separate arguments for either case.

assumption[Penalty function and grid for tuning]$\;$\\ The penalty $\pi(\theta)$ and the set $\Lambda$ take one of the following two forms: \begin{enumerate} • Ridge: $\pi(\theta) = \tfrac12 \theta \cdot A^{-1} \cdot \theta$, where $A$ is positive definite,\\ and $\Lambda = \mathbb R^+$. • Lasso: $\pi(\theta) = \|A^{-1} \cdot \theta\|_1$, where $A$ is an invertible matrix,\\ and $\Lambda \subset \mathbb R^+$ is finite. \end{enumerate}

The remaining assumptions state regularity conditions. The first item in Assumption (ref) is a condition on the loss function which allows us to invoke results from empirical process theory. Similar assumptions are invoked in van2000asymptotic when deriving the properties of M-estimators. The second item in Assumption (ref) is a weak high-level condition ruling out divergence of ERM estimators, which ensures the applicability of empirical process results.

assumption[Conditions for convergence of the ERM estimator]$\:$ \begin{enumerate} • Lipschitz loss\\ The loss function $l(\beta, z)$ satisfies $$|l(\beta_1, z) - l(\beta_2, z)| \leq m(z) \cdot \|\beta_1 -\beta_2\|, $$ for all $\beta_1, \beta_2$ in a neighborhood of $0$, where $\sup_n \operatorname{Var}(m(Z_n^i)) \leq \infty$. Furthermore, $l(\beta, Z_n^i)$ is differentiable w.r.t. $\beta$, for all $\beta$ in a neighborhood of $\beta=0$, with probability $1$. • Stochastically bounded ERM estimator\\ The sequence $\hat \theta_n = \operatorname*{argmin\;}_\theta L_n(\theta)$ is bounded in probability. \end{enumerate}

To prove the (uniform) convergence of $CV_n$ to $SURE$, we finally impose the following additional regularity conditions.

assumption[Conditions for convergence of CV] $\;$ \begin{enumerate} • Conditions on loss\\ There exist $\mu>0,\nu<\infty$ independent of $n$ such that $L_{n}(\theta)$ is $\mu$-strongly convex and has $\nu$-smooth Hessians with probability approaching 1 under $\mu_{n}$.\footnote{$L_{n}(\theta)$ is $\mu$-strongly convex if $\nabla^2 L_{n}(\theta) - \mu I$ is positive semi-definite for all $\theta$. A function is $L$-smooth if its gradients are Lipschitz continuous with Lipschitz constant $L$.} • Conditions on scores\\ The function $\sqrt{n}\nabla_{\theta}l_{n}(\theta,Z_{n}^{i})$ is Lipschitz continuous almost everywhere, i.e., there exists $B_{n}(Z_{n}^{i})$ such that \[ \left\Vert \sqrt{n}\nabla_{\theta}l_{n}(\theta,Z_{n}^{i})-\sqrt{n}\nabla_{\theta}l_{n}(\theta^{\prime},Z_{n}^{i})\right\Vert \le B_{n}(Z_{n}^{i})\left\Vert \theta-\theta^{\prime}\right\Vert \ \forall\theta,\theta^{\prime}\in\Theta \] and $E_{\mu_{n}}\left[\left\Vert B_{n}(Z_{n}^{i})\right\Vert ^{2}\right]<\infty$. Additionally, there exists $M<\infty$ independent of $n$ such that \[ \mathbb{E}_{\mu_{n}}\left[\left\Vert \sqrt{n}\nabla_{\theta}l_{n}\left(\theta_{0},Z_{n}^{i}\right)\right\Vert ^{4}\right]\le M. \] • Conditions on Hessians\\ The Hessian $\nabla_{\theta}^{2}l_{n}(\theta,Z_{n}^{i})$ is such that \[ \frac{1}{n}\sum_{i=1}^{n}\left\Vert n\nabla_{\theta}^{2}l_{n}\left(\theta_{0},Z_{n}^{i}\right)\right\Vert ^{2}=O_{\mu_{n}}(1). \] Furthermore, there exists $C_{n}(Z_{n}^{i})$ such that \[ \left\Vert n\nabla_{\theta}^{2}l_{n}(\theta,Z_{n}^{i})-n\nabla_{\theta}^{2}l_{n}(\theta^{\prime},Z_{n}^{i})\right\Vert \le C_{n}(Z_{n}^{i})\left\Vert \theta-\theta^{\prime}\right\Vert \ \forall\theta,\theta^{\prime}\in\Theta \] and $\sup_n E_{\mu_{n}}\left[C_{n}(Z_{n}^{i})^{2}\right]<\infty$. • Conditions on Fourth Derivatives\\ The fourth derivative tensor, $D_\theta^4 l_{n}(\theta,Z_{n}^{i})$, of $l_n(\cdot, Z_n^i)$ is such that \[ \sup_n \mathbb{E}_{\mu_{n}}\left[\sup_{\theta\in\Theta}\left\Vert n^{2}D_\theta^{4}l_{n}\left(\theta,Z_{n}^{i}\right)\right\Vert ^{4}\right]\le M, \] for some $M <\infty$. \end{enumerate}

Main result and intermediate lemmas

Our main goal in this paper is to prove the following result:

theorem$$ \bar L_n(\hat \theta_n^*, \theta_0) \rightarrow_d \tfrac12 \|\hat \theta^* - \theta_0\|^2. $$

In words, the distribution of loss of the penalized ERM estimator tuned using n-fold cross-validation converges to the distribution of squared error of the corresponding shrinkage estimator in the normal means model, tuned by Stein's Unbiased Risk Estimate. A special case of these limiting estimators are James-Stein shrinkage estimators, for which closed form characterizations of the risk function are known stein1981estimation. An immediate corollary of Theorem (ref) is the convergence of risk functions, subject to possible truncation of tail events.\footnote{Truncation is necessary because, even for estimators such as linear OLS regression, risk is typically undefined, since some eigenvalues of the design matrix might be close to 0, so that the moments of $\hat\theta_n$ might not exist.}

corollaryLet $M>0$. Then $$ E\left[ \min\left(\bar L_n(\hat \theta_n^*,\theta_0), M \right) \right] \rightarrow E\left[ \min\left(\|\hat \theta^* - \theta_0\|^2,M\right) \right]. $$

Outline of proof

We will build up our argument that proves Theorem (ref) in a series of lemmas. Before doing so, however, we first provide an intuitive sketch of our argument, while neglecting remainder terms. Subsequently, we will prove that these remainder terms are indeed asymptotically negligible.

\paragraph{Influence function approximation} We start by noting that empirical risk is asymptotically equivalent to quadratic error loss relative to the sample mean $\tilde \theta_n$, \bals L_n(\theta) &\approx const. + \tfrac12 \|\theta-\tilde \theta_n\|^2, &\tilde \theta_n &= \theta_0 + \tfrac{1}{\sqrt{n}} \sum_i X_n^i,& \eals where $$X_n^i = -\nabla_\beta l(\theta_0/\sqrt{n}, Z_n^i)$$ is the influence function. Recall that we have normalized the Hessian in Assumption (ref), which simplifies the expression for $\tilde \theta_n$. This approximation of empirical risk immediately implies an asymptotic linear approximation of the empirical risk minimization (ERM) estimator, $\hat \theta_n \approx \tilde \theta_n$. These are standard approximations that deliver asymptotic normality of ERM estimators; see for instance Theorem 5.21 in van2000asymptotic.

We then get the corresponding approximation for the penalized ERM estimator, for fixed $\lambda$, \bals \hat \theta_n^\lambda &\approx \tilde \theta_n^{\lambda} = \tilde \theta_n + g^\lambda(\tilde \theta_n), \eals where we recall the definition $g^\lambda(\theta) = \operatorname*{argmin\;}_g \tfrac12\|g\|^2 + \lambda \cdot \pi(\theta + g)$.

\paragraph{Convergence of CV to SURE} An analogous approximation holds for leave-one-out (LOO) loss. To obtain the LOO sample mean, the influence function $\tfrac{1}{\sqrt{n}} X_n^i$ is subtracted from the sample mean $\tilde \theta_n$, which gives \bals L_n^{-i}(\theta) &\approx const. + \tfrac12 \|\theta-\tilde \theta_n^{-i}\|^2,& where \tilde \theta_n^{-i} &= \tilde \theta_n - \tfrac{1}{\sqrt{n}} X_n^i. \eals The penalized LOO estimator is then approximately given by $$ \hat \theta^{\lambda,-i}_n \approx \tilde \theta_n^{-i} + g^\lambda(\tilde \theta_n^{-i}) \approx \tilde \theta_n^{\lambda} - \tfrac{1}{\sqrt{n}}(I+\nabla g^{\lambda}(\tilde \theta_n)) \cdot X_n^i. $$ In the last step we have replaced $g^\lambda$ by its first-order Taylor approximation around $\tilde \theta_n$, at points $\tilde \theta_n$ where $g^\lambda$ is differentiable. (This approximation won't hold at kink-points of $g^\lambda$, which exist for Lasso penalties, in particular.)

The n-fold cross-validation estimator of the risk of $\hat \theta^{\lambda}_n$ can be approximated by $$ CV_n(\lambda) =\sum_i l_n(\hat \theta^{\lambda,-i}_n,Z_n^i) \approx const. + \tfrac1n \sum_i \|\hat \theta^{\lambda,-i}_n - \theta_0 - \sqrt{n} X_n^i\|^2. $$ When we take this expression, plug in the approximate form of $\hat \theta^{\lambda,-i}_n$, multiply out the inner products, and omit terms which do not depend on $\lambda$, we obtain \bals CV_n(\lambda) &\approx const. + \tfrac1n \sum_i \|\underbrace{\tilde \theta_n + g^\lambda(\tilde \theta_n) - \tfrac{1}{\sqrt{n}}(I+\nabla g^{\lambda}(\tilde \theta_n)) \cdot X_n^i}_{\approx\hat \theta^{\lambda,-i}_n } - \theta_0 - \sqrt{n} X_n^i\|^2\\ &\approx const.+\tfrac1n \sum_i \|g^\lambda(\tilde \theta_n)\|^2 + \tfrac2n \sum_i \langle \nabla g^{\lambda}(\tilde \theta_n) \cdot X_n^i, X_n^i \rangle\\ &\approx const. + \| g^\lambda(\hat \theta_n)\|^2 + 2 \operatorname{trace} ( \nabla g^\lambda(\hat \theta_n)\cdot \hat\Sigma_n )\\ &= const. + SURE(\lambda, \hat \theta_n, \hat \Sigma_n)\\ &\approx const. + SURE(\lambda, \hat \theta_n, \Sigma), \eals where $\hat\Sigma_n$ is the sample second moment of $X_n$. In the second line, $const.$ subsumes any terms that do not depend on $\lambda$, while the approximation omits terms that depend on $\lambda$ but are of order $1/\sqrt{n}$. This approximation to $CV_n$ has the form of Stein's Unbiased Risk Estimate. The first term in this approximation to $CV_n(\lambda)$ corresponds to the average in-sample error, the second term has the form of a covariance penalty efron2004estimation.

\paragraph{Convergence of tuning parameter and tuned estimators} We will need to show that this approximation is uniformly valid in $\lambda$. We will furthermore need to show that uniform proximity of $CV_n$ to $SURE$ is enough to guarantee proximity of the corresponding optimized tuning parameters, $$ \operatorname*{argmin\;}_\lambda CV_n(\lambda) \approx \operatorname*{argmin\;}_\lambda SURE(\lambda, \hat \theta_n, \hat \Sigma_n). $$

figure[figure omitted — 1,593 chars of source]

This latter step is non-trivial, because the criterion function $SURE(\lambda, \hat \theta_n, \hat \Sigma_n)$ typically has multiple local minima. For certain values of $\hat \theta_n$ this function furthermore has multiple global minima. When $SURE$ has multiple global (near-)minima, uniform closeness of $CV$ to $SURE$ is not enough to ensure closeness of the minimizer of $CV$ to the minimizer of $SURE$. The plots in Figure (ref), which are reproduced from wilsonkasymackey2018,\footnote{The numerical values corresponding to these examples are as follows: (a) SURE for Ridge: $\hat \theta = (1.3893, 1.5)$, $L(\theta) =(\theta - \hat \theta)\operatorname{diag}(1,40)(\theta - \hat \theta)$, $\pi(\theta) = \|\theta\|^2$, (b) SURE for Lasso: $\hat \theta = \frac{1}{\sqrt{n}}(\sqrt{1/8}, \sqrt{9/8}, 2 )$, $\pi(\theta) = \sum |\theta_j|$. } illustrate two numerical examples (realizations of $\hat \theta_n$ and of $\hat \Sigma_n$) for which $SURE$ indeed has multiple global minima.

Using separate arguments for Ridge (Appendix (ref)) and Lasso (Appendix (ref)), we will prove, however, that the global minimum of $SURE$ with respect to $\lambda$ is unique and well separated almost everywhere, in a suitable sense. Put differently, cases such as those represented in Figure (ref) are non-generic, such that they do not lead to a breakdown of convergence for the optimized tuning parameter. The arguments proving that multiple global minima only occur on a set of measure $0$ might be the most non-standard part of our proof

Well-separation ensures that the $\operatorname*{argmin\;}$ functional is continuous at almost every realization of $\hat\theta$. From these results we thus conclude that the mapping from $\hat \theta_n$ and $CV_n$ to the tuned estimate $\hat \theta_n^*$ is almost everywhere continuous. This allows us to invoke the continuous mapping and dominated convergence theorems, and to conclude the proof of Theorem (ref).

Intermediate lemmas

Let us now turn to a more formal exposition of our argument. We will prove Theorem (ref) in a series of Lemmas. The lemmas are stated in this section, their proofs in the appendices. All results impose the assumptions stated in Section (ref).

lemma[Lipschitz $g^\lambda$] For any $\lambda \geq 0$, if $\pi( \cdot )$ is convex then $g^\lambda(\theta)= \operatorname*{argmin\;}_g \tfrac12\|g\|^2 + \lambda \cdot \pi(\theta + g)$ is Lipschitz with Lipschitz constant 1.

lemma[Influence function approximation] \be L_n(\theta) - L_n(\theta_0) = \tfrac12 \|\theta - \tilde \theta_n\|^2 + \epsilon_n(\theta), \ee where $\sup_{\theta: \| \theta \| \le C} \epsilon_n(\theta) = o_{\mu_n}(1)$ and $\sup_{\theta: \| \theta \| \le C} \nabla \epsilon_n(\theta) = o_{\mu_n}(1)$ for any $C < \infty$, and \bals \tilde \theta_n &= \theta_0 + \tfrac{1}{\sqrt{n}} \sum_i X_n^i, & X_n^i &= -\nabla_\beta l(\theta_0/\sqrt{n}, Z_n^i). \eals The ERM and penalized ERM estimators \bals \hat \theta_n &= \operatorname*{argmin\;}_\theta L_n(\theta), & \hat \theta_n^\lambda &= \operatorname*{argmin\;}_\theta \left[L_n(\theta) + \lambda \cdot \pi(\theta)\right] \eals satisfy \bals \hat \theta_n &= \tilde \theta_n+ o_{\mu_n}(1),& \sup_\lambda \Vert \hat \theta_n^\lambda &- \tilde \theta_n - g^\lambda(\tilde \theta_n) \Vert = o_{\mu_n}(1). \eals

lemma[Limiting squared error loss] The limiting expected loss is well defined and given by $$\bar L(\theta, \theta_0) =\tfrac12 \|\theta - \theta_0\|^2.$$ Convergence of $\bar L_n(\theta, \theta_0)$ to this limit is uniform in any bounded neighborhood of $\theta_0$: $\sup_{\theta: \|\theta - \theta_0\| \leq C} |\bar L_n(\theta, \theta_0) - \bar L(\theta, \theta_0)| \to 0$ for all $C<\infty$.

corollary[Asymptotic distribution for fixed tuning parameter] The ERM and penalized ERM estimators satisfy \bals \hat \theta_n &\rightarrow^d \hat \theta \sim N(\theta_0, \Sigma),& \hat \theta_n^\lambda &\rightarrow^d \hat \theta + g^\lambda(\hat \theta). \eals

lemma[Convergence of CV to SURE]$\;$\\ The n-fold crossvalidation criterion satisfies $$ \sup_{\lambda \in \Lambda} \left|CV_n(\lambda) - SURE(\lambda, \hat \theta_n, \Sigma)\right| \rightarrow^{\mu_n} 0. $$

lemma[Joint convergence of tuning parameter and tuned estimators] $$ (\lambda_n^*, \hat \theta_n) \rightarrow^d (\lambda^*, \hat \theta). $$ and $$ \hat \theta_n^{*} \rightarrow^d \hat \theta^{*}. $$

From (ref), we then show (ref). The appendices prove each of these Lemmas in turn.